[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.


> Thank you,
> Paul

David Winsemius, MD
West Hartford, CT

More information about the R-help mailing list