1

I've been using Cox for survival analysis with the goal of identifying the effect of 5 different levels of a treatment on survival. My dataset is not small and has 130 covariates after dummy encoding.

About 15 of these - including one level of my treatment - violate the PH assumption.

I think in my situation the best extended Cox option would be adding time interactions to those that violate the PH assumption meaningfully. Another option, of course, could be an AFT, if the assumptions there are also satisfied.

My question is: If I were to try each of these approaches, how could I compare which these are preferred given my objective. That is, how would I know which is better? It's not like I'm making a prediction about survival in which case I could use CI or something more precise. Thus, comparing which is best and knowing how to know this is not obvious to me.

JED HK
  • 409

1 Answers1

2

The most rational and linear approach would be to put Bayesian priors on the amount of non-PH for all variables, before making use of information you get from examining fit to the observed outcome variable. This will give you the correct uncertainties for "final" effects in the model.

Less satisfactory but simpler approaches are (1) use the unitness correlation $\rho$ that goes with the test of PH to gauge the extent of the violation, instead of using P-values, which can detect trivial violations, or (2) fit an overall PH model then fit an overall non-PH model (e.g., log-normal survival time model) and see which one "wins" by some global goodness of fit assessment. I'm less inclined to pick and choose covariate-by-covariate.

Frank Harrell
  • 91,879
  • 6
  • 178
  • 397
  • Thanks Professor Harrell. Can you explain what you mean by ρ here? When I print 'cox.zph' here (on R) I see chi-squared, df, and P-value. – JED HK Aug 10 '22 at 13:15
  • 1
    Hmmm not sure why this was taken out of cox.zph. You can run Spearman's correlation of survival time vs r <- resid(fit, 'scaledsch') to get it. I think you can get the survival time from as.numeric(dimnames(r)[[1]]). $\rho$ gives you a unitness index of how correlated $\beta(t)$ is with $t$, which is a sample-size-less index of non-PH. – Frank Harrell Aug 10 '22 at 14:36
  • It was removed a while ago (https://cran.r-project.org/web/packages/survival/news.html). I followed your advice and now must interpret the output. Most ρ are around 0. The largest in magnitude away from 0 was 0.064 (though this variable was not significant according to cox.zph. This suggests correlation between β(t) and (t) is low, I presume. However according to Schoenfeld plots some hazards do seem to change quite a bit. How can I reconcile this? – JED HK Aug 10 '22 at 15:49
  • I could not find mention of removal of $\rho$ in the link you provided. In terms of reconciling, check the scale on which the residual plot y-axis is drawn. Sometimes the scale exaggerates small effects. 0.064 is not very large. – Frank Harrell Aug 10 '22 at 16:00
  • Ah, sorry about that. I took it from EdM's answer here: https://stats.stackexchange.com/questions/144923/extended-cox-model-and-cox-zph/238964#238964

    Small enough to ignore, then? Feel free to link anything. You have one example in your book in Chapter 20 in which ρ is -0.23 and highly significant but you do not comment on what -0.23 means in terms of magnitude

    – JED HK Aug 11 '22 at 08:01
  • Hard to know the exact impact of -.23 but I think 0.06 is safe to ignore. In general I'm losing favor of goodness of fit tests. Better to look at impact, AIC, and $R^{2}_{adj}$ as in
    https://www.fharrell.com/post/impactpo and even better to embed models inside more general models putting Bayesian priors on parameters you hope you don't need (in your case, parameters allowing for non-PH). We need a linear rather than a dichotomous modeling process.
    – Frank Harrell Aug 11 '22 at 11:47
  • Your help is greatly appreciated Professor. As a non-statistician, I hope I can make sense of it. Final things before I leave you alone: 1) I will need a credible reference for ignoring 0.06. Do you have something suitable? 2) Schoenfeld plots suggest some PHs change quite a bit (e.g., from 2 to 1 for a categorical var) for some variables. How would you reconcile this with a ρ being close to 0? This seems (at least to me) to be conflicting – JED HK Aug 11 '22 at 13:08
  • The only thought I have is to see if the categories for which you see departures are low enough in prevalence (relative frequency) that they don't matter as much as it appears. – Frank Harrell Aug 11 '22 at 13:20