I am calculating cox propotional hazards models with the coxph function from the survival package. My data relates to failure of various types of endovascular interventions. I can successfully obtain the LR, Wald, and Score test p-values from the coxph.object, as well as the hazard ratio as follows:
formula.obj = Surv(days, status) ~ type coxph.model = coxph(formula.obj, df) fit = summary(coxph.model) hazard.ratio = fit$conf.int[1] lower95 = fit$conf.int[3] upper95 = fit$conf.int[4] logrank.p.value = fit$sctest[3] wald.p.value = fit$waldtest[3] lr.p.value = fit$logtest[3] I had intended to report logrank P values with the hazard ratio and CI obtained from this function. In one case the P was 0.04 yet the CI crossed one, which confused me, and certainly will raise questions by reviewers. In retrospect I can see that the CI calculated by coxph is intimately related to the Wald p-value (which in this specific case was 0.06), so this would appear to me not a good strategy for reporting my results (mixing the logrank test with the HR and CIs from coxph). I can report the Wald p-values instead, but I have read that the Wald test is inferior to the score test or LR test. My questions for survival analysis jockeys out there / TT: 1. Should I just stop here and use the wald.p.value? This appears to be what Stata does with the stcox function (albeit Breslow method). 2. Should I calculate HR and CIs that "agree" with the LR or logrank P? How do I do that? Thank you, Paul ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
