I want to do a weighted Cox regression analysis (coxphw command in R) since there is violation of proportional hazards assumption in standard Cox model shown by Schoenfeld residuals. What would be the test/metric for goodness of fit in a weighted Cox regression? Thank you!
As AdamO commented, the
coxphw() function is essentially the
coxph() function with a particular set of weights. The weights, however, are a function of event times rather than fixed for each case. For the default
AHR (average hazard ratio) setting, the weight for all cases at risk at an event time
t is set to the product of the estimated survival probability and the inverse of the estimated probability of censoring at that time. That gives a particular type of "average hazard ratio" that has a reasonable interpretation even if proportional hazards (PH) doesn't hold.
The first implication is that there is no need to test for PH if you use
coxphw() with the
AHR setting. The PH assumption is irrelevant to interpretation of the coefficients.
Second, the "average hazard ratio" for
coxphw() is closely related to the concordance metric for model performance; from page 8 of Dunkler et al:
in a two-group comparison, average hazard ratios approximate the odds of concordance $OC = c/(1 − c) = P(T_A < T_B)/P(T_B < T<A)$ very well, $T_A$ and $T_B$ being the survival times of two randomly chosen subjects of groups $A$ and $B$. Intuitively, this can be explained as each pairwise comparison of a subject that fails with all subjects still at risk being given equal importance.
It is not, however, suitable for predictions as functions of time (Table 4 of Dunkler et al). That would suggest a measure of concordance as most appropriate for testing goodness of fit. The
coxphw package does provide a
concord() function to provide a generalized measure of concordance for each predictor in the model, based on the cases that went into the model.
What would be more satisfying would be a way to validate the model in a way similar to what's done in the
rms package: for example, fitting a model to multiple bootstrapped samples of the data, and testing performance of each of those models on the whole data set. That gives a reasonable idea of how well the results would apply to a new sample from the underlying population.
coxphw package doesn't make that easy, as it seems to have only limited ability to predict outcomes on cases not included in the model. A work-around would be to use the regression coefficients returned by
coxphw() to calculate the linear predictor from the covariate values of cases not in the model, and see how frequently the order of linear predictor values agrees with the order of event times between comparable pairs of cases (all pairs of cases with event times; pairs of censored/event cases with censored time greater than event time).
If times to events are important, you could consider spending more time working around the PH problem. For example, if your model has the wrong functional form for a continuous covariate, that could show up as non-PH. Depending on your goals and those of your audience, understanding the nature of any non-PH behavior might end up being more important than just reporting an average hazard ratio.