Jpaete
Jpaete

Reputation: 59

Different interval confidences on the same forest plot (R, survival)

I'm currently working on survival functions and especially on forest plots. I would like to plot the different levels of confidence intervals on the same line, with different colors or shades.

To do so, I use the forest function from the forestploter package. I use this package because it is greatly adaptable. However, I'm not able to have these different CIs. For the moment, I was able to get the different CIs but not on the same line.

Do anyone of you know how I could fix it ?

Here is an example of what I've been doing with the lung data from the survival package.

library(survival)
library(forestploter)

# Fit a Cox proportional hazards regression model
reg <- coxph(Surv(time, status) ~ age + sex + ph.ecog + ph.karno, data = lung)

  # Extract variable names, coefficients, and standard errors
variables <- names(coef(reg))
beta <- coef(reg)
se <- sqrt(diag(reg$var))

  # Calculate hazard ratios (HR)
HR <- exp(beta)

  # Calculate p-values
p <- 1 - pchisq((beta / se)^2, 1)

  # Calculate confidence intervals (CI)
CI_95 <- exp(confint(reg))
CI_90 <- exp(confint(reg, level = 0.90))

  # Format CI string for HR table
CI_hr_95 <- sprintf("%.2f (%.2f to %.2f)", HR, CI_95[, 1], CI_95[, 2])

# Create essai data.table with results
essai <- data.table(
  variables,  # Variable names
  beta,        # Coefficients
  p = round(p, 3),  # Rounded p-values
  HR,          # Hazard ratios
  CI_95_lower = CI_95[, 1],  # Lower CI for 95%
  CI_95_upper = CI_95[, 2],  # Upper CI for 95%
  CI_90_lower = CI_90[, 1],  # Lower CI for 90%
  CI_90_upper = CI_90[, 2],  # Upper CI for 90%
  CI_hr_95     # Formatted HR CI string
)

essai$` `<- paste(rep(" ", 50), collapse = " ")

# Define forest plot theme
tm <- forest_theme(
  base_size = 12,      # Font size
  base_family = "",    # Font family (optional)
  ci_pch = 10,         # Point character for CI lines
  ci_col = "black",    # Color for CI lines
  ci_alpha = 1,        # Transparency of CI lines
  ci_fill = NULL,      # Fill color for CI regions (optional)
  ci_lty = 1,          # Line type for CI lines
  ci_lwd = 2,          # Line width for CI lines
  ci_Theight = NULL,    # Height for CI tick marks (optional)
  legend_name = "CI threshold",  # Legend label for CI lines
  legend_value = c("90%", "95%")  # Legend values for CI lines
)

# Prepare data for forest plot
forest_data <- essai[, c(1, 9, 10, 3)]  # Select relevant columns

# Define forest plot elements
forest <- forest(
  forest_data,
  est = list(essai$HR, essai$HR),           # Estimate values (HR for both lines)
  lower = list(essai$CI_90_lower, essai$CI_95_lower),  # Lower CI for each line
  upper = list(essai$CI_90_upper, essai$CI_95_upper),  # Upper CI for each line
  sizes = 0.5,                               # Point size
  ci_column = 3,                             # Column with CI information
  ref_line = 1,                               # Reference line at HR=1
  ticks_at = c(0, 1, 5.3),                     # X-axis ticks
  xlim = c(0, 5.3),                           # X-axis limits
  nudge_y = 0.3,                              # Vertical nudge for labels
# Even with nudge_y = 0, it is not on the same line
  theme = tm,                                # Forest plot theme
  arrow_lab = c("Longer WT", "Shorter WT")    # Labels for HR lines (optional)
)

# Create the forest plot
plot(forest)

PS : For who might know this handbook, I would like to recreate the Figure 19.10 (p 468) from Regression Modeling Strategies: With Applications to Linear Models, Logistic Regression, and Survival Analysis (2015)

Different CI on the same plot

PS 2 : With my real data, using the forest function is great, because it makes easier to deal with label and to plot HR and CIs for interaction variables.

Thanks.

Upvotes: 0

Views: 91

Answers (0)

Related Questions