Software 8.6.1 (Factors and their overall test). Example 8.1.1: the numbers of each cell type, the hazard ratios with their intervals, and the likelihood ratio test for cell type by comparing the models with and without it:

> library(survival); options(show.signif.stars = FALSE)
> vet <- veteran; vet$trt <- factor(vet$trt, labels = c("standard", "test"))
> table(vet$celltype)                    # reference level: squamous

 squamous smallcell     adeno     large
       35        48        27        27
> m1 <- coxph(Surv(time, status) ~ trt + celltype + karno, data = vet)
> round(cbind(HR = exp(coef(m1)), exp(confint(m1))), 3)
                     HR 2.5 % 97.5 %
trttest           1.299 0.876  1.926
celltypesmallcell 2.282 1.347  3.865
celltypeadeno     3.171 1.778  5.653
celltypelarge     1.484 0.853  2.580
karno             0.969 0.959  0.979
> m0 <- coxph(Surv(time, status) ~ trt + karno, data = vet)
> anova(m0, m1)  # LR test of the three cell-type terms
Analysis of Deviance Table
 Cox model: response is  Surv(time, status)
 Model 1: ~ trt + karno
 Model 2: ~ trt + celltype + karno
   loglik  Chisq Df Pr(>|Chi|)
1 -483.97
2 -474.91 18.102  3   0.000419

Reading the output. Adenocarcinoma and small-cell tumours carry about three and two times the hazard of squamous tumours; the overall test of cell type is \(\chi ^2_3 = 18.1\), \(p = 0.0004\). The test regimen’s hazard ratio, \(1.30\), is not significant.

Python: CoxPHFitter with formula=’trt + C(celltype) + karno’.

Software 8.6.2 (Ties and concordance). The coefficients of treatment and Karnofsky score under the Breslow, Efron and exact methods for ties; the number of tied death times; and the concordance of the model:

> f <- Surv(time, status) ~ trt + celltype + karno
> ties <- sapply(c("breslow", "efron", "exact"),
+                function(m) coef(coxph(f, data = vet, ties = m))[c(1, 5)])
> round(ties, 5)
         breslow    efron    exact
trttest  0.25731  0.26174  0.26258
karno   -0.03111 -0.03127 -0.03156
> sum(duplicated(vet$time[vet$status == 1]))   # tied death times
[1] 31
> round(c(concordance = m1$concordance[["concordance"]]), 4)
concordance
     0.7371

Reading the output. The three methods differ only in the third decimal place. With \(31\) tied death times out of \(128\), Efron’s method is close to the exact one, Breslow’s slightly further. The concordance is \(0.737\).

Python: lifelines uses Efron’s method; cph.concordance_index_ gives \(C\).

Software 8.6.3 (Interactions, strata and predictions). The likelihood ratio test of a treatment-by-cell-type interaction; the stratified model; and predicted survival at \(30\), \(90\) and \(180\) days (rows) and medians for Karnofsky scores \(40\), \(60\), \(80\) (columns), squamous tumour, standard regimen:

> m2 <- coxph(Surv(time, status) ~ trt * celltype + karno, data = vet)
> anova(m1, m2)                          # does the treatment effect vary?
Analysis of Deviance Table
 Cox model: response is  Surv(time, status)
 Model 1: ~ trt + celltype + karno
 Model 2: ~ trt * celltype + karno
   loglik  Chisq Df Pr(>|Chi|)
1 -474.91
2 -471.59 6.6458  3    0.08409
> ms <- coxph(Surv(time, status) ~ trt + karno + strata(celltype), data = vet)
> round(exp(coef(ms)), 4)                # stratified by cell type
trttest   karno
 1.2622  0.9648
> nd <- data.frame(trt = "standard", celltype = "squamous",
+                  karno = c(40, 60, 80))
> sf <- survfit(m1, newdata = nd)
> round(summary(sf, times = c(30, 90, 180))$surv, 3)  # rows: times
         1     2     3
[1,] 0.774 0.872 0.929
[2,] 0.506 0.695 0.823
[3,] 0.191 0.412 0.623
> summary(sf)$table[, "median"]
  1   2   3
 92 144 283

Reading the output. The interaction test gives \(\chi ^2_3 = 6.6\), \(p = 0.08\). Stratified by cell type, the treatment hazard ratio is \(1.26\). The predictions show the strength of the performance score: six-month survival rises from \(19\%\) at score \(40\) to \(62\%\) at \(80\), and the median from \(92\) to \(283\) days.

Python: cph.predict_survival_function(new_patients).