4.4 Kaplan–Meier estimation in software
Software 4.4.1 (The estimator by hand). Example 4.1.1: a loop over the distinct event times computing the number at risk, the deaths and the running product; then survfit for comparison:
> library(survival)
> tm <- c(3, 5, 6, 8, 10, 12, 15, 18, 22, 24) # the ten patients of
> ev <- c(1, 0, 1, 1, 0, 1, 0, 1, 0, 0) # Section 3, 0 = censored
> S <- 1; out <- NULL
> for (tj in sort(unique(tm[ev == 1]))) { # event times only
+ n <- sum(tm >= tj); dj <- sum(tm == tj & ev == 1)
+ S <- S * (1 - dj / n)
+ out <- rbind(out, c(t = tj, n_risk = n, deaths = dj, S = round(S, 4)))
+ }
> out
t n_risk deaths S
[1,] 3 10 1 0.9000
[2,] 6 8 1 0.7875
[3,] 8 7 1 0.6750
[4,] 12 5 1 0.5400
[5,] 18 3 1 0.3600
> km <- survfit(Surv(tm, ev) ~ 1)
> summary(km)$surv # the same steps
[1] 0.9000 0.7875 0.6750 0.5400 0.3600
Reading the output. The table reproduces the hand calculation, and survfit gives the same five values: \(0.9\), \(0.7875\), \(0.675\), \(0.54\), \(0.36\).
Python: lifelines.KaplanMeierFitter().fit(tm, ev).survival_function_.
Software 4.4.2 (Survival, intervals and the median). The lung data: the estimated survival at \(180\), \(365\) and \(730\) days with the number at risk and the \(95\%\) interval; the median with its interval; and Greenwood’s standard error at one year by hand, against survfit’s:
> fit <- survfit(Surv(time, status) ~ 1, data = lung) > s <- summary(fit, times = c(180, 365, 730)) > round(cbind(time = s$time, n_risk = s$n.risk, surv = s$surv, + lower = s$lower, upper = s$upper), 3) time n_risk surv lower upper [1,] 180 160 0.722 0.666 0.783 [2,] 365 65 0.409 0.345 0.486 [3,] 730 13 0.116 0.072 0.187 > fit # median and its 95% CI Call: survfit(formula = Surv(time, status) ~ 1, data = lung) n events median 0.95LCL 0.95UCL [1,] 228 165 310 285 363 > use <- fit$time <= 365 # Greenwood by hand, at 1 year > terms <- fit$n.event / (fit$n.risk * (fit$n.risk - fit$n.event)) > v <- sum(terms[use]) > round(c(greenwood_se = s$surv[2] * sqrt(v), survfit_se = s$std.err[2]), 5) greenwood_se survfit_se 0.03582 0.03582
Reading the output. One-year survival is \(40.9\%\) (\(34.5\%\) to \(48.6\%\)). The median is \(310\) days (\(285\) to \(363\)). Greenwood’s sum gives a standard error of \(0.03582\) at one year, exactly the value survfit uses.
Python: kmf.median_survival_time_ and kmf.confidence_interval_.
Software 4.4.3 (Curves by group and median follow-up). Kaplan–Meier estimates for men and women, with the numbers of subjects and events, the medians and their intervals; then the reverse Kaplan–Meier estimate of median follow-up:
> bysex <- survfit(Surv(time, status) ~ sex, data = lung)
> tab <- summary(bysex)$table
> tab[, c("records", "events", "median", "0.95LCL", "0.95UCL")]
records events median 0.95LCL 0.95UCL
sex=1 138 112 270 212 310
sex=2 90 53 426 348 550
> rev_km <- survfit(Surv(time, status == 1) ~ 1, data = lung) # reverse KM
> summary(rev_km)$table["median"] # median follow-up, days
median
588
Reading the output. Men: \(138\) patients, \(112\) deaths, median \(270\) days (\(212\) to \(310\)). Women: \(90\) patients, \(53\) deaths, median \(426\) days (\(348\) to \(550\)). The intervals barely overlap, which suggests a real difference. The median follow-up is \(588\) days.
Python: fit a KaplanMeierFitter to each group and plot both on one axis.
Questions on this section
Stuck on something here? Ask below and it stays attached to this topic.