Introduction
生存曲线处理的是**"事件发生了没有"和"什么时候发生"合在一起**的数据,而且必然带删失——有人到研究结束还没发生事件,你只知道他"至少活到了这一天"。
普通的均值比较处理不了这种数据:把删失的人当成"没事件"会低估风险,直接扔掉又引入偏倚。Kaplan-Meier 曲线是专门为此设计的:只有观察到事件时才下降一级,删失只是把这个人从风险集里移走,不让曲线掉。
读的时候按这个顺序:
- 起点是 1。 时间 0 时所有人都还没发生事件。
- 每级台阶是一次事件。 曲线是阶梯状的,不是平滑的——它只在观察到事件的那一刻下降。
- 看两条线分不分得开。 分得越开,两组的生存经验差别越大。
- 看置信区间。 后期变宽是正常的:风险集里人少了。
- 看风险表。 这一行最容易被跳过,却最要紧——某个时间点之后只剩几个人,那之后的曲线就别当真了。
- 最后看 p 值。 log-rank 概括的是整个随访期的差异,不说明差多少。
Example Data
生存曲线要的输入是:随访时间 + 事件状态(0/1)+ 分组变量,一行一个人。时间和状态必须成对,缺一个这个人就进不了模型;分组变量决定画几条曲线。
lung_survival 是 228 例肺癌的随访数据。它由 survival::lung 改写而来,事件列已经固化成明确的 0/1——原始的 status 用的是 1 = 删失、2 = 死亡,而绝大多数函数期待的是 0/1,编码错了曲线照样画得出来,只是上下颠倒。
library(survival)
library(ggsurvfit)
library(dplyr)
library(biopalette)
# 公开地址,和数据集页上「下载 CSV」给的是同一个 —— 不写仓库相对路径:
# 那个目录不进仓库,读者 clone 下来也没有这个文件,这段代码就跑不了
data_km <- read.csv("https://assets.evanzhou.org/tessera/csv/lung_survival.csv") |>
mutate(
event = as.integer(event),
# 显式写 levels:图例和配色都按因子水平的顺序走
sex = factor(sex, levels = c("male", "female"))
) |>
# 三列缺一不可:时间和状态是一对,分组缺了这个人也没法归入任何一条曲线
filter(!is.na(time_days), !is.na(event), !is.na(sex))
# 换自己的数据时保留这一行。事件编码错了不会报错、不会有任何症状,
# 只会得到一张上下颠倒的曲线
stopifnot(all(data_km$event %in% c(0L, 1L)))
# survfit2() 来自 ggsurvfit,和基础的 survfit() 算的是同一件事,
# 但它把调用记在对象里,后面的图层才能自己找到数据
fit_km <- survfit2(Surv(time_days, event) ~ sex, data = data_km)
Palettes
两组用 heat_light 的蓝红配色。两个色相直接区分两组,曲线、置信区间和图例始终使用同一映射。
置信区间沿用曲线本色再降透明度,不额外引入一套浅色——两条带子叠在一起时,多一种颜色就多一层要辨认的东西。
曲线图里曲线比填充重。 置信区间是背景信息,深浅必须和曲线拉开,否则重叠处谁也看不清。
# heat_light 正好两色,按因子水平的顺序绑上去 —— 用名字而不是位置,
# 后面 scale_*_manual() 的 breaks 才不会和图例次序错开
pal_sex <- setNames(
rev(get_palette("heat_light", type = "qualitative")),
levels(data_km$sex)
)
Recipe
| No. | Method | Input Data | Palettes |
|---|---|---|---|
| 1 | ggsurvfit |
fit_km |
pal_sex |
| 2 | survminer |
data_km |
pal_sex |
1 · ggsurvfit
一张完整的 KM 图有五件东西:曲线、置信区间、风险表、中位线、检验 p 值。
ggsurvfit() 返回的是普通 ggplot 对象,所以这五件里除了曲线本身,其余都是标准图层,能和同一篇里的其它图共用主题。
#| fig: km
#| fig-width: 8
#| fig-height: 7
# log-rank:rho = 0 是标准 log-rank,各时间点等权;rho = 1 是 Peto 权重,
# 更看重早期事件。选哪个是分析决定,不是画图参数
sd_logrank <- survdiff(Surv(time_days, event) ~ sex, data = data_km, rho = 0)
p_logrank <- pchisq(sd_logrank$chisq, df = length(sd_logrank$n) - 1,
lower.tail = FALSE)
ggsurvfit(fit_km) +
labs(
title = "Overall Survival by Sex",
subtitle = "Tessera Toy: lung_survival - Kaplan-Meier with 95% CI",
x = "Time (days)", y = "Survival probability"
) +
scale_color_manual(values = pal_sex, breaks = names(pal_sex),
labels = c("Male", "Female")) +
# 置信区间沿用曲线颜色,只降低存在感;不再引入 recipe 外的颜色
scale_fill_manual(values = pal_sex, guide = "none") +
add_confidence_interval(type = "ribbon", alpha = 0.14) +
add_risktable() + # 风险表 —— 判断尾部可不可信全靠它
add_quantile(y_value = 0.5, color = "gray50", linewidth = 0.75) + # 中位生存参考线
scale_ggsurvfit() + # 把 y 轴调成 0–1 的百分比刻度
ggplot2::theme_classic() +
theme(
legend.title = element_blank(),
legend.position = "bottom",
plot.title = element_text(size = 16, face = "bold")
) +
annotate(
"text", x = 500, y = 0.8,
label = sprintf("log-rank p = %.3f", p_logrank),
hjust = 0, size = 4.2
)
