R语言绘制多时间点生存ROC曲线：从原理到实战

简介：针对医学生物统计与临床试验中的生存预测需求，这份R语言源代码包面向具备一定R基础的研究者，用于绘制SCI科研常见的时间依赖ROC曲线（即多时间点生存ROC），以评估模型在不同随访时刻的区分能力。包内共3个文件，包括一个R主脚本、一份txt格式的输入数据说明或示例、以及一份pdf格式的ROC图输出样例，压缩包整体仅13KB，轻量且聚焦。核心脚本集成了数据读取、时间点设置、ROC及AUC计算和图形修饰等环节，可直接替换自身数据运行，便于对比多个预测时间点的性能变化。资源已有351人学习浏览，适合需要快速产出多时间点生存ROC图的医学、生物统计方向研究者，既能节省调试时间，也有助于理解pROC与survival包在生存预测模型评估中的实际应用。

1. 多时间点生存ROC：一张图回答“这个标记物在第3年还有没有用”

多时间点生存ROC曲线，是肿瘤预后、生信数据挖掘和临床预测模型里最常被问到“能不能画”的一张图。它回答的是单时间点ROC答不了的问题：你手里的标志物在第1年、第3年、第5年分别还有多少区分能力。很多人算了一个时间点AUC=0.78，就直接写进文章，但审稿人追问“随访5年，这个指标到第3年还稳不稳”时，拿不出第二张图。用R语言做这件事，核心不是调一个包，而是把随访数据组织对、把时点选对、把生存删失处理对。这套方案适合做生存分析、生信数据和临床预测模型的人，拿到源码后只需要替换自己的三列数据就能跑通。

2. 为什么单时点ROC不够用：时依ROC原理与R包选型

2.1 生存数据里的“金标准”会过期

普通ROC要求每个样本有一个固定标签：患病或未患病。但生存数据里，标签是随时间变的。同一个病人，第365天没事件，第730天事件发生了；另一个人随访结束也没事件，只留下一个删失状态。如果你图省事，把“第3年是否死亡”当y，删失的病人要么被整行丢掉，要么全算成阴性，两种做法都会扭曲AUC。

时依ROC（time-dependent ROC）把时间轴引进来：在时刻t，病例（case）定义为t之前发生了目标事件的人，对照（control）定义为到t时刻仍未发生事件的人。删失者只要在t之前删失，就不再参与这个时点的计算，但他在更早的时点仍然有贡献。这就是Heagerty在2000年提出的I/D（incident/dynamic）定义。与之相对的I/S定义，则用“t时刻是否存活”当静态标签，对晚删失样本更敏感。

这两种定义没有绝对好坏，但如果你想回答“这个标志物对随访期内累积发生的事件区分能力如何”，I/D更贴近临床直觉。timeROC包实现的就是I/D定义，而且它的方差估计用的是influence function解析法，不需要反复Bootstrap，这是它最省心的地方。

2.2 timeROC、survivalROC、riskRegression怎么选

R语言里做多时间点生存ROC，主流是三个包，它们的分工不太一样：

我一般这样取舍：探索阶段和出图用timeROC，因为它一次调用就能给出一串时间点的AUC和置信区间；如果数据量大、AUC曲线波动厉害，用survivalROC的NNE方法做平滑估计对比一下；最后写文章要报告“内部验证AUC”，再用riskRegression的Score做IPW或交叉验证。标题里的“R语言绘制”，落到实现上，90%的常见源码包都是把timeROC的调用脚本包了一层，加了数据读取和ggplot出图。

2.3 时间点怎么定：从临床问题反推，而不是让数据替你选

这是最容易翻车的一步。多时间点ROC的“时间点”，应该来自临床问题，而不是把随访时间均匀切成20份。肿瘤研究常见的做法是取1年、3年、5年，对应临床上的复查节点；急重症研究可能取28天、90天、180天。定时间点之前，先看数据的随访质量：

library(survival) fit <- survfit(Surv(time, status) ~ 1, data = demo) quantile(fit)

quantile(fit)会给出随访时间的分位数。核心原则是：你选的每个时间点，都要落在随访数据能支撑的范围内，而且到该时点还在风险集里的人数不能太少。我一般把最大时点设在中位随访时间附近，最多不超过最大随访时间的80%。后面第4章会展开讲“时点取得太靠后会看到什么奇怪曲线”。

3. 用R把多时间点生存ROC画出来：数据准备到出图全流程

3.1 数据长什么样：三列是最小骨架

任何多时间点生存ROC，数据都逃不出三列：随访时间、结局状态、标记物值（可以是基因表达、风险评分、影像参数）。先把数据组织成这个结构，后面的坑能少一半。

set.seed(2024) n <- 200 time <- round(rexp(n, rate = 1/36), 1) # 模拟随访月数，均值36个月 status <- rbinom(n, 1, prob = 0.45) # 1=事件，0=删失 marker <- rnorm(n) + 0.6 * status # 让marker与结局相关 demo <- data.frame(id = 1:n, time = time, status = status, marker = marker) head(demo)

这段是构造演示数据，不是分析数据。time是随访时长，单位随意，但后面times参数必须用同一个单位；status必须编码成0/1，1代表目标事件，0代表删失或竞争事件（你要分析的那个结局没发生）；marker是你要评价的预测变量，可以是连续型也可以是风险评分。真实数据替换时，只改这三列的来源就行，列名建议统一成time/status/marker，脚本可复用性会高很多。

如果你手里的zip源码里只有分析脚本没有示例数据，我建议先用上面这段模拟数据把脚本跑通，再换真实数据。别一上来就拿全量数据跑，报错时你分不清是数据问题还是代码问题。

3.2 计算多个时间点的AUC：timeROC核心调用

library(timeROC) res <- timeROC( T = demo$time, # 随访时间 delta = demo$status, # 0/1结局 marker = demo$marker, # 预测因子 cause = 1, # 指定哪种事件算阳性 times = c(12, 36, 60), # 要评估的时间点（月） iid = TRUE # 计算解析方差，用于置信区间 ) res$AUC # 三个时间点的AUC res$se # 标准误 res$CI_AUC # AUC的95%置信区间

参数说明：cause只在存在竞争风险时要注意，如果有多种结局而你只关心其中一种，cause要指向那一种；times是核心，必须和T的单位一致，且不能有超过最大随访时间的值；iid=TRUE是timeROC最值的参数，它用影响函数直接推出标准误，省掉了Bootstrap的几千次重抽样。

跑完看输出，AUC向量会和times一一对应。常见情况是AUC随时点往后推移而下降，这说明标记物对远期事件的区分能力衰减；也有AUC先升后降的模式，说明标记物在某个窗口期最强，这种信息比单点ROC丰富得多，正是审稿人想看到的。

3.3 画图：AUC随时间变化曲线与单时点ROC

pdf("multi_time_auc.pdf", width = 5, height = 4, family = "Arial") plot(res, conf.int = TRUE, # 画置信区间带 col = "darkred", lwd = 2, xlab = "随访时间（月）", ylab = "AUC") dev.off()

这个图就是多时间点生存ROC最常见的成品形态：横轴是随访时间，纵轴是AUC，带一条随时间变化的曲线和置信区间带。如果想单独看某个时间点的ROC曲线，比如3年时的灵敏度和特异度表现：

plot(res, time = 36, col = "darkred", lwd = 2)

timeROC的plot函数对time参数有两个用法：不指定time画AUC-时间曲线，指定具体time画该时点的ROC曲线。前者适合放正文，后者适合放补充材料。有一点要注意：这是base R的绘图体系，想叠加多个标记物、改主题、拼图，都得在plot后面继续用lines或points，不能像ggplot那样直接加图层。后面3.4会给一个转ggplot的做法，适合要精细排版的人。

3.4 多标记物对比：一条图画两条曲线

实际分析里很少只评价一个指标，常见的是“新标记物 vs 临床传统指标”对比。做法是分别调用timeROC，然后在同一张图上叠加：

res_marker1 <- timeROC(demo$time, demo$status, demo$marker1, cause = 1, times = c(12, 36, 60), iid = TRUE) res_marker2 <- timeROC(demo$time, demo$status, demo$marker2, cause = 1, times = c(12, 36, 60), iid = TRUE) plot(res_marker1, col = "darkred", lwd = 2, xlab = "随访时间（月）", ylab = "AUC") lines(res_marker2$times, res_marker2$AUC, col = "steelblue", lwd = 2) legend("bottomleft", legend = c("新标记物", "临床指标"), col = c("darkred", "steelblue"), lwd = 2, bty = "n")

如果你想把图做得更精细，或者要拼进ggplot体系，直接把timeROC的结果转成数据框：

auc_df <- data.frame( time = res$times, AUC = res$AUC, lower = res$CI_AUC[, 1], upper = res$CI_AUC[, 2] ) library(ggplot2) ggplot(auc_df, aes(time, AUC)) + geom_ribbon(aes(ymin = lower, ymax = upper), alpha = 0.2) + geom_line(linewidth = 1) + coord_cartesian(ylim = c(0.5, 1)) + labs(x = "随访时间（月）", y = "AUC")

把timeROC的结果搬进ggplot，是我自己踩过坑之后固定的写法。因为后续加主题、调字体、拼多图，ggplot比base R顺手得多；但注意CI_AUC是矩阵，取列时要用[, 1]和[, 2]，很多人在这里取错，画出来的置信区间带是歪的。

4. 多时间点生存ROC的常见坑：现象、原因、解决

4.1 时间点超过最大随访时间，AUC曲线尾部“飞起”

现象：AUC曲线在随访后期突然急剧上升或剧烈抖动，置信区间宽到没眼看，甚至出现AUC接近1的“完美表现”。

原因：到随访后期，真正还在风险集里的人只剩十几个甚至几个。timeROC的I/D估计在风险集很小时方差暴涨，少数几个事件就能把AUC拉得很高。这不是你的标记物变强了，是样本量撑不住了。

解决：把times限制在随访时间的合理分位数内，我一般看80%分位数：

quantile(demo$time, c(0.5, 0.8, 0.95))

然后把times的最大值设在不超过P80的位置。比如中位随访36个月，P80是48个月，那就不要选60个月做时点，选12、24、36更稳。如果临床确实需要5年结果，唯一解法是增加随访时间或扩大样本量，而不是硬画。

4.2 用survivalROC想要置信区间，自己写Bootstrap把电脑跑崩了

现象：survivalROC本身不返回置信区间，于是很多人写for循环Bootstrap 500次，每次重新算AUC，最后内存溢出或跑了半小时没结果。

原因：survivalROC的NNE方法每次调用都要做局部平滑和最近邻计算，500次Bootstrap叠加在几千个样本上，计算量是O(B×n²)量级，R的单线程循环扛不住。

解决：如果只是要置信区间，换timeROC并把iid设为TRUE，它用解析方差直接给出SE和CI，不需要Bootstrap。如果你的场景必须用survivalROC（比如要NNE平滑估计），那就减少Bootstrap次数到100~200，并且每次只保存AUC标量，不要保存整个ROC对象：

boot_auc <- numeric(200) for (b in 1:200) { idx <- sample(n, replace = TRUE) fit_b <- survivalROC(demo$time[idx], demo$status[idx], demo$marker[idx], predict.time = 36, method = "NNE", span = 0.25) boot_auc[b] <- fit_b$AUC } quantile(boot_auc, c(0.025, 0.975))

注意这里每次Bootstrap只存一个AUC值，200次循环也就几百KB内存，不会崩。

4.3 删失比例过高时AUC虚高

现象：你的队列删失比例超过70%，画出来的多时间点生存ROC漂亮得不像话，AUC在多个时点都在0.8以上，但换个队列就崩到0.6。

原因：删失比例高意味着大量样本在事件发生前就离开了观察。在I/D定义下，这些人在t之前删失就不进t时刻的分母，导致剩下来参与计算的都是“更容易出事”的人，AUC被系统性拔高。这在生信数据挖掘里尤其常见，很多公共数据库的随访本来就不完整。

解决：先报告删失比例，然后在验证阶段换用IPW加权的AUC估计。riskRegression包的Score函数可以做：

library(riskRegression) m_cox <- coxph(Surv(time, status) ~ marker, data = demo) sc <- Score(list("marker" = m_cox), formula = Surv(time, status) ~ 1, data = demo, cause = 1, times = c(12, 36, 60), plots = "auc")

Score用逆删失概率加权重新估计AUC，对高删失队列更稳健。如果结果和timeROC差很多，说明你的结论依赖删失假设，文章里要如实报告。

4.4 把生存数据硬转成二分类扔给pROC

现象：为了省事，有人先把数据转成“3年是否死亡”，删失样本要么丢掉，要么算成未死亡，然后用pROC包画普通ROC，结果和timeROC算出来的AUC对不上。

原因：这不是哪个包算错了，而是定义完全不同。pROC假设标签固定，删失样本的“3年状态”本来就是未知的，强行归入任何一类都会引入偏差。生存ROC里每个时点用的是累积病例和动态对照，删失者在该时点不参与，但在更早时点仍然贡献过信息。

解决：涉及随访时间的数据，一律用timeROC或survivalROC这类时依ROC工具，不要转成二分类再画。如果审稿人要求你解释两者的差异，你就把I/D定义和删失处理讲清楚，这本身就是方法学亮点。

4.5 时点选太多，图乱得像心电图

现象：有人把times设成seq(6, 60, by = 6)，画出10个点，AUC曲线在0.6和0.9之间来回跳，审稿意见直接是“图不清晰，结论无法评估”。

原因：时点越多，每个点附近的风险集越小，估计噪声越大。而且临床医生没法解释“第30个月的AUC是0.72”这种含义，时点本身没有医学节点支撑。

解决：固定3个临床可解释的时间点，正文放AUC-时间曲线，三个时点的ROC单图放补充材料。时点之间相隔不要太密，比如12、36、60这样。如果随访短，就28、90、180天。核心是每个时点都要能回答“临床上什么时候用这个指标”。

5. 把图调到能投SCI：出图尺寸、字体与验证习惯

5.1 出图尺寸与字体

投稿图不是越大越好，而是符合期刊排版规范。单栏图宽89mm左右，双栏图宽180mm左右，R里直接按毫米指定输出：

png("roc_multitime.png", width = 180, height = 120, units = "mm", res = 300)

res=300对应印刷精度，放到Word里不会虚。字体建议统一Arial或Helvetica，base R的plot里用family参数指定；如果遇到中文标签在PDF里乱码，直接改用英文标签，投稿图用英文本来就是常态。

5.2 给审稿人的验证：时点、风险集、AUC、置信区间一张表

光有图还不够，审稿人通常想看具体数字。我把每个时点的风险集人数、AUC、95%CI整理成补充表，模板如下：

风险集人数可以用survival::survfit的summary拿到。这个表的意义在于让审稿人一眼看到“第60个月风险集只剩41人”，曲线尾巴再好看他也会理解那是噪声。

我早先在一个随访40个月的小队列里，把时间点设到了36个月，结果最后一段置信区间宽到几乎覆盖整张图。现在拿到任何数据，第一件事是看删失比例和中位随访，再决定时间点，而不是直接抄教程里的1/3/5年。这套流程最值钱的地方，不是那个plot函数，而是知道每个数字背后有多少人在支撑。希望帮到你。

本文还有配套的精品资源，点击获取