Comparing AUC between joint model and Cox model?

[RandallJEllis]

Thank you for an incredible package. My goal is to compare time-varying AUROC for a joint model and a Cox model. tvROC takes a joint model, but I don't think it can take as input a Cox model, and I've been trying different functions like timeROC to calculate the AUROC for the Cox model, but it's hard to determine whether the time-varying AUROC is being calculated in the exact same way. Below is my current code, any help is appreciated!

#### JOINT MODEL

lmeFit <- lme(
    ORRES_boxcox ~ AGEYR_z + AGEYR_z_squared + AGEYR_z_cubed + SEX + 
      EDCCNTU_z + APOEGN + AGEYR_z*APOEGN + AGEYR_z_squared*APOEGN + 
      AGEYR_z_cubed*APOEGN,
    random = ~ COLLECTION_DATE_DAYS_CONSENT_yr | BID,
    data = train_df,
    na.action=na.exclude,
    control = lmeControl(opt = 'optim')
  )

coxFit <- coxph(
  Surv(time_to_event_yr, label) ~ AGEYR_z + AGEYR_z_squared + 
    AGEYR_z_cubed + SEX + EDCCNTU_z + APOEGN + AGEYR_z*APOEGN +
    AGEYR_z_squared*APOEGN + AGEYR_z_cubed*APOEGN,
  data = train_df,
  id = train_df$BID,
  x = TRUE
)

jointFit <- jm(
  Surv_object = coxFit,
  Mixed_objects = list(lmeFit),
  time_var = "COLLECTION_DATE_DAYS_CONSENT_yr",
  functional_forms = list("ORRES_boxcox" = ~ value(ORRES_boxcox) + slope(ORRES_boxcox)),
  data_Surv = train_df,
  id_var = "BID"
)

dt = 1
first_year <- as.integer(min(train_df$COLLECTION_DATE_DAYS_CONSENT_yr)) + 2
last_year <- as.integer(max(train_df$COLLECTION_DATE_DAYS_CONSENT_yr)) - 1
print(first_year)
print(last_year)

res_roc_iter <- list()
for (tstart in seq(first_year,last_year)){
  
  thoriz <- tstart + dt
  
  # This is a strange quirk of tvROC (function to calculate ROC curves). The 
  # function subsets newdata before the start of the interval (Tstart) by
  # the longitudinal time variable (COLLECTION_DATE_DAYS_CONSENT_yr), but
  # then subsets the subset for rows where the time to event 
  # (time_to_event_yr) is less than the time horizon (Tstart + Dt) but where
  # the label is 1. With the way the data is set up for the Cox model, where
  # time points before the time to event are labeled 0, and time points at or
  # after the time to event are labeled 1, tvROC will never find events. 
  # The solution is, for all patients who ever experience an event, fix all
  # values in the label column at 1. 
  val_df_case_label_allT <- val_df %>%
    group_by(BID) %>%
    mutate(label = ifelse(any(label == 1 & time_to_event_yr <= thoriz),
                                   1,
                                   label))
  
  print(
    paste0(
      'Interval: Year ',tstart, ' to ',thoriz
    )
  )
  
  roc <- tvROC(jointFit, newdata = val_df_case_label_allT, 
                   Tstart = tstart, 
                   Dt = dt,
                   type_weights = c("IPCW")
      )
  res_roc_iter[[paste0("year_",tstart)]] <- roc
}

#### COX MODEL
risk_scores <- predict(coxFit,
                             newdata = val_df,
                             process = 'event')


# Compute time-dependent ROC and AUROC
print('Starting Cox model ROC analysis')
baseline_roc <- timeROC(
  T = val_df$time_to_event_yr,
  delta = val_df$label,
  marker = risk_scores,
  cause = 1,  # Specify the event of interest
  times = seq(first_year+1,last_year+1),
  # iid = TRUE  # Optional: Compute influence functions for confidence intervals
)

[drizopoulos]

To my knowledge, JMbayes2 and timeROC do not calculate the same AUCs. I plan to implement a method for tvROC() that will also work with Cox models. When this is ready, I'll post a message here.