[R] Multiple comparisons for coxph survival analysis model

Pavel Kúr pavel.kur at prf.jcu.cz
Sat Sep 26 13:43:33 CEST 2009


Hello, all R-users!

I am working on fitting a survival analysis model using the coxph
function for Cox proportional hazards regression model. Data look like
usual:

==========================
group     block    death    censor
Group1     1          4          1
Group1     1         12         1
...
Group2     30        4          1
Group2     30        4          1
...
Group3     57       16        1
Group3     57       16        1
==========================

And I need to compare surviving among the particular groups. Fitting
works normally:

================================================================
> cph.1 <- coxph(Surv(death, censor) ~ group + cluster(block), data = seedlings)
> summary(cph.1)
Call:
coxph(formula = Surv(death, censor) ~ group+ cluster(block),
    data = seedlings)

  n= 27000
                    coef exp(coef) se(coef) robust se     z    p
groupGroup2 0.436      1.55   0.0539     0.296  1.47 0.14
groupGroup3 3.048     21.06   0.0439     0.283 10.77 0.00

                    exp(coef) exp(-coef) lower .95 upper .95
groupGroup2      1.55     0.6467     0.865      2.76
groupGroup3     21.06     0.0475    12.100     36.67

Rsquare= 0.38   (max possible= 0.997 )
Likelihood ratio test= 12892  on 2 df,   p=0
Wald test            = 271  on 2 df,   p=0
Score (logrank) test = 16164  on 2 df,   p=0,   Robust = 84.2  p=0
==================================================================

I have obtained tests of significance for differences between the
second/third group and the first (reference) group, but I want to
compare each group with each other, not only all groups with the first
one!
So I need to use some multiple comparison methods.
I have tried the multcomp library, which I normally use for glm
models. But it hasn't worked:

==============================================================
> summary(glht(cph.1, linfct=mcp(group="Tukey")))
Error in glht.matrix(model = list(coefficients = c(0.435824045783883,  :
  ‘ncol(linfct)’ is not equal to ‘length(coef(model))’
============================================================

So I tried a different approach using the contrast matrix:

======================================================
> summary(glht(cph.1, linfct = contrMat(coef(cph.1),type="Tukey")))

         Simultaneous Tests for General Linear Hypotheses

Multiple Comparisons of Means: Tukey Contrasts


Fit: coxph(formula = Surv(death, censor) ~ group + cluster(block),
    data = seedlings)

Linear Hypotheses:
                                 Estimate Std. Error z value Pr(>|z|)
groupGroup2 - groupGroup3 == 0   2.6117     0.1781   14.66   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
=======================================================

Well, I have achieved comparison between the 2nd and 3rd group, but
this time the 1st group is missing.
In glm the reference group is expressed as the intercept, so comparing
with it is comparing with the intercept. But there is no intercept in
coxph!

Please, is there any way how to accomplish full multiple comparisons in coxph?

Thank you in advance!
Pavel Kur




More information about the R-help mailing list