Next, we fit Cox proportional hazards models using the coxph() function. To begin, we consider a model that uses sex as the only predictor.

> fit.cox <- coxph(Surv(time, status) ~ sex)
> summary(fit.cox)
Call:
coxph(formula = Surv(time, status) ~ sex)

  n= 88, number of events= 35 

          coef exp(coef) se(coef)     z Pr(>|z|)
sexMale 0.4077    1.5033   0.3420 1.192    0.233

        exp(coef) exp(-coef) lower .95 upper .95
sexMale     1.503     0.6652     0.769     2.939

Concordance= 0.565  (se = 0.045 )
Likelihood ratio test= 1.44  on 1 df,   p=0.23
Wald test            = 1.42  on 1 df,   p=0.233
Score (logrank) test = 1.44  on 1 df,   p=0.23

Note that the values of the likelihood ratio, Wald, and score tests have been rounded. It is possible to display additional digits.

> summary(fit.cox)$logtest[1]
    test 
1.438822 
> summary(fit.cox)$waldtest[1]
test 
1.420000 
> summary(fit.cox)$sctest[1]
    test 
1.440495 

Regardless of which test we use, we see that there is no clear evidence for a difference in survival between males and females.

> logrank.test$chisq
[1] 1.440495

As we learned in this chapter, the score test from the Cox model is exactly equal to the log rank test statistic!