[R] Cox proportional hazards confidence intervals
David Winsemius
dwinsemius at comcast.net
Mon Nov 21 05:50:29 CET 2011
On Nov 20, 2011, at 6:34 PM, Paul Johnston wrote:
> 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).
I don't understand two things: Why would your report the inferior
result, and I suppose I also wonder why does it make that much
difference? The estimate is what it is and a p-value of .04 is not
that different from one of .06. Or are we dealing with religious
beliefs here?
>
> 2. Should I calculate HR and CIs that "agree" with the LR or logrank
> P? How do I do that?
Therneau and Grambsch show how to calculate profile likelihood curves
that can be used to generate confidence intervals on pages 57-59 of
"Modeling Survival Data". This "survival analysis jockey" considers
that book an essential reference. They basically use the offset
capacity to construct 50 likelihoods around the estimate for one
particular variable in a more complete model and then show where the
97.5th and 0.025th percentile points are for an beta estimate based on
a chi-square distribution for these log-likelihoods.
Further code not possible in the absence of the complete formula.
--
David.
>
> Thank you,
> Paul
David Winsemius, MD
West Hartford, CT
More information about the R-help
mailing list