I've created a few Cox regression models and I would like to see how well these models perform and I thought that perhaps a ROC-curve or a c-statistic might be useful similar to this articles use:
J. N. Armitage och J. H. van der Meulen, ”Identifying co‐morbidity in surgical patients using administrative data with the Royal College of Surgeons Charlson Score”, British Journal of Surgery, vol. 97, num. 5, ss. 772-781, Maj 2010.
Armitage used logistic regression but I wonder if it's possible to use a model from the survival package, the survivalROC gives a hint of this being possible but I can't figure out how to get that to work with a regular Cox regression.
I would be grateful if someone would show me how to do a ROC-analysis on this example:
library(survival) data(veteran) attach(veteran) surv <- Surv(time, status) fit <- coxph(surv ~ trt + age + prior, data=veteran) summary(fit)
If possible I would appreciate both the raw c-statics output and a nice graph
Thank you very much for answers. @Dwin: I would just like to be sure that I've understood it right before selecting your answer.
The calculation as I understand it according to DWin's suggestion:
library(survival) library(rms) data(veteran) fit.cph <- cph(surv ~ trt + age + prior, data=veteran, x=TRUE, y=TRUE, surv=TRUE) # Summary fails!? #summary(fit.cph) # Get the Dxy v <- validate(fit.cph, dxy=TRUE, B=100) # Is this the correct value? Dxy = v[rownames(v)=="Dxy", colnames(v)=="index.corrected"] # The c-statistic according to the Dxy=2(c-0.5) Dxy/2+0.5
I'm unfamiliar with the validate function and bootstrapping but after looking at prof. Frank Harrel's answer here on R-help I figured that it's probably the way to get the Dxy. The help for validate states:
... Somers' Dxy rank correlation to be computed at each resample (this takes a bit longer than the likelihood based statistics). The values corresponting to the row Dxy are equal to 2 * (C - 0.5) where C is the C-index or concordance probability.
I guess I'm mostly confused by the columns. I figured that the corrected value is the one I should use but I haven't really understood the validate output:
index.orig training test optimism index.corrected n Dxy -0.0137 -0.0715 -0.0071 -0.0644 0.0507 100 R2 0.0079 0.0278 0.0037 0.0242 -0.0162 100 Slope 1.0000 1.0000 0.2939 0.7061 0.2939 100 ...
In the R-help question I've understood that I should have "surv=TRUE" in the cph if I have strata but I'm uncertain on what the purpose of the "u=60" parameter in the validate function is. I would be grateful if you could help me understand these and check that I haven't made any mistakes.