企业官网建设流程全解析

校准曲线这东西，凡是做过临床预测模型的人应该都不陌生。不管是用Cox回归做生存预测，还是用Logistic回归做二分类诊断，模型建完之后，审稿人几乎必问一句：你们的模型校准度怎么样？尤其如果你准备按TRIPOD声明来报告模型，校准曲线和区分度指标就是绕不开的两块内容。我在之前的系列文章里讲过Cox模型的构建、验证、列线图绘制，这篇单独把校准曲线拎出来，专门聊清楚它在Cox模型里到底怎么画、怎么看、怎么解读，顺便把能直接跑的R语言代码全部贴出来。这篇文章适合正在做临床预测模型、被校准曲线折腾过的医生和科研人员，也适合R语言刚入门但想搞懂校准概念的朋友。

1. 校准曲线到底在做什么

1.1 模型评估的两把尺子：区分度与校准度

很多人一提预测模型性能，第一反应就是C-index或者AUC。C-index确实重要，它回答的问题是：在随机抽取的一对患者中，模型能不能把更早发生事件的那位排在风险更高的位置。换句话说，它衡量的是模型的排序能力，叫区分度（discrimination）。

但区分度高不代表模型的预测概率就对。举个例子，某个模型对所有患者预测的生存率都系统性偏高10个百分点，它依然可能拥有很高的C-index，因为患者之间的相对排序没有变，高的还是高，低的还是低。可如果你用这个模型去告诉一位患者“你两年生存率有70%”，而真实情况可能只有60%，这就不是排序的问题了，这是预测绝对值准不准的问题。

校准度（calibration）衡量的正是这个：模型预测出的概率，和实际观察到的概率，到底有多大差距。校准曲线就是把这种差距画出来给人看。理想状态下，模型预测生存率是70%，那么在实际随访数据中，同样一批预测值为70%的患者，两年后应该大约有70%的人还活着。

所以，建立预测模型，区分度和校准度是两条腿，缺一条都站不稳。C-index看的是“谁能排前面”，校准曲线看的是“预测的概率靠不靠谱”。

1.2 Cox模型校准曲线的独特之处

同样是校准曲线，在Logistic回归和Cox回归里画的逻辑完全不一样，这是很多人第一次上手就懵掉的地方。

Logistic回归预测的是一个固定时间点的结局概率，比如“患者是否患病”。你只要把预测概率和实际患病比例放在一起比较就行，操作简单直观。可Cox模型输出的是生存函数，也就是随时间变化的生存概率。同一名患者，1年生存率是一个值，3年生存率又是另一个值。因此Cox模型的校准曲线必须指定一个时间点，比如“2年生存率校准曲线”，然后在这个时间点上比较预测生存率和观测生存率。

这里还有个隐藏的坑：生存数据里有删失（censoring）。随访结束时还没发生事件的患者，不代表他永远不会发生事件，只是还没观察到而已。所以“观测生存率”不能简单地用“某组内还活着的人数除以总人数”来计算，必须用Kaplan-Meier估计或者伪观测值、逆删失加权这类方法来处理删失。这也是为什么Cox校准曲线比Logistic校准曲线难画、也更容易画错的原因。

搞懂了这两点，再看下面的代码和分析，思路就会清晰很多。

2. 画校准曲线的两套主流方案选型

2.1 实战中最常用的rms路线：建模校准一条龙

在R语言里面画Cox校准曲线，医学统计领域用得最多的方案是rms包。Frank Harrell开发的这个包几乎是临床预测模型的标准工具，从数据描述、模型拟合、列线图绘制到模型验证，全部在同一个生态里完成。

rms方案的核心思路是：先用cph()函数拟合Cox模型，然后直接调用calibrate()函数对模型做校准。calibrate()内部采用Bootstrap重抽样，得到三条曲线：理想曲线（Ideal）、表观曲线（Apparent）和偏差校正曲线（Bias-corrected）。这三条曲线的含义后面会详细讲，简单说，偏差校正曲线才是你论文里最该展示的结果，因为它校正了模型自身过拟合带来的乐观偏差。

用rms方案有几个明显的好处。第一，代码量少，几行就能出图，适合大多数人快速出结果。第二，和rms生态下的nomogram()列线图、scores()评分表衔接流畅，你前一步刚建好列线图，下一步直接校准，数据格式不用来回转换。第三，calibrate()的输出不仅有图，还有量化的校准统计指标，方便写进论文结果部分。

不过这个方案也有局限。它的参数设置比较讲究，比如必须提前设定好time.inc，模型拟合时必须加上x = TRUE, y = TRUE，很多新手在这些细节上栽跟头。另外它默认的图形样式偏基础，要达到发表级别的美观程度，通常还需要用ggplot2重构。

2.2 更灵活的riskRegression路线：多模型、多时间点一条龙

如果你的需求不是简单画一条曲线，而是想在一个图中比较多个模型的校准情况，或者同时查看多个时间点的校准曲线，那我更推荐riskRegression包。

riskRegression的核心函数是Score()，功能很强大。你传入一个已经拟合好的Cox模型（可以是survival包的coxph()，不要求必须是rms的cph()），它就能计算出Brier分数、时间依赖性AUC、校准曲线等多种评估指标。配合plotCalibration()函数，可以快速画出校准图。

riskRegression还有一个优势，它对时间点非常友好。比如你想同时看1年、3年、5年的校准曲线，一条代码就能在同一个坐标系里叠加显示，这在rms方案里操作起来要麻烦得多。我在做多时间点生存模型验证时，基本都会用riskRegression来出一版对比图。

缺点是riskRegression包的API风格和rms差别较大，它的Hist()对象和Surv()对象经常让初学者头晕，而且它的输出内容相对更“统计风味”一些，需要花点时间理解各项指标的含义。

2.3 方案对比小结

到底选哪套？我的建议很直接：如果你只是为当前这个Cox模型补一张校准曲线，且已经用了rms建模，那就用rms方案；如果你要做多模型对比、多时间点验证，或者用的是survival包建模，我推荐riskRegression。两套代码我下面都会贴出来，你可以都跑一遍，看哪个更顺手。

3. 全套实操：从数据到发表级校准曲线

3.1 准备工作与核心陷阱

在写正式代码前，先把演示环境准备好。下面的代码全部用R语言实现，测试数据是survival包自带的lung数据集。这是一份晚期肺癌随机临床试验数据，包含228名患者的生存时间（time，单位是天）、生存状态（status，1为删失，2为死亡），以及年龄、性别、ECOG评分等变量。用它演示Cox校准曲线非常合适。

如果你电脑上还没装R环境，建议先去R语言官网下载最新版R，再装一个RStudio作为IDE。安装好之后，执行下面的代码安装本次需要用到的包：

install.packages("rms") install.packages("survival") install.packages("riskRegression") install.packages("ggplot2")

数据准备和预处理的时候，我先提醒两个高频错误，几乎所有刚接触这份数据的人都会踩。

第一个坑是lung数据里status变量的编码。在这个数据集里，status = 1表示删失（还没观察到死亡），status = 2表示事件发生（死亡）。而在R的Surv()函数里，内部约定通常是1表示事件发生，0表示删失。因此你需要把原来的status变量做一次映射转换，否则后面的分析会得出完全错误的结果。

第二个坑是缺失值。lung数据里的ph.ecog有几个缺失值，如果不处理，cph()函数在拟合时会默认使用na.action机制删除缺失行，但后续datadist()和calibrate()可能因为变量不同步而报错。稳妥起见，演示代码里直接用na.omit()把缺失行删掉。

library(survival) library(rms) library(riskRegression) library(ggplot2) # 加载数据，并剔除缺失值 data("lung") lung2 <- na.omit(lung) # 转换生存状态编码：原编码2=死亡，这里转换为1=事件 lung2$status2 <- ifelse(lung2$status == 2, 1, 0) head(lung2)

转换完之后，status2就可以直接放进Surv()里了。接下来正式建模。

3.2 rms方案代码：cph + calibrate + plot

rms这套流程有几个必须记住的“仪式”，少一个都可能出问题。第一步是用datadist()设置数据分布，并把它指定为全局选项；第二步是在cph()里设置x = TRUE, y = TRUE，这样calibrate()才能拿到原始数据做重抽样；第三步是设置time.inc，这是校准的时间点，单位必须和time变量一致。

假设我们要画两年生存率校准曲线，time是以天为单位，那么两年就是365.25 * 2 = 730.5天。这个换算容易出错，尤其当原始时间是“月”时，千万别习惯性地写24。下面是可以直接跑的完整代码：

# 1. 设置datadist，rms生态必需 dd <- datadist(lung2) options(datadist = "dd") # 2. 设定校准时间点：2年 time_inc <- 365.25 * 2 # 3. 拟合Cox模型 fit_cph <- cph(Surv(time, status2) ~ age + sex + ph.ecog, data = lung2, x = TRUE, y = TRUE, time.inc = time_inc, surv = TRUE) # 4. 查看模型概要 print(fit_cph)

代码里我加了surv = TRUE，这个参数让模型额外保存生存函数估计，保证后续绘图时能顺利计算每个个体的预测生存概率。缺了这个参数，在某些版本的rms里calibrate()会表现出莫名其妙的问题。

接下来就是核心一步：

# 5. Bootstrap法校准，200次重抽样，分20组比较 set.seed(2024) # 设随机种子，保证结果可复现 cal <- calibrate(fit_cph, u = time_inc, # 校准时间点 B = 200, # 重抽样次数 pred.group = 20, # 将预测概率分为20组 cmethod = "hare") # 采用HARE方法估计观测生存率 # 6. 绘图 plot(cal, xlim = c(0, 1), ylim = c(0, 1), xlab = "Predicted 2-Year Survival Probability", ylab = "Observed 2-Year Survival Probability", subtitles = TRUE)

跑完这段代码，你会看到一张经典的校准曲线图。图中通常有三条线：一条是对角虚线，代表理想状态下的完美校准线；一条是表观校准线，显示模型在训练数据上的预测和观测一致性；还有一条是偏差校正后的校准线，它更真实地反映了模型在新患者群体中可能的表现。

B = 200这个参数是Bootstrap重抽样次数。200次是比较常规的默认选择，既能保证结果稳定，速度也不会太慢。如果样本量大、事件数多，可以提高到500次，但代价是计算时间明显增加。我自己的经验是，在228例这种小样本上，200次足够，再往上提升对最终图形的影响肉眼几乎看不出来。

pred.group = 20表示把患者按预测概率分成20个小组，每组内分别比较平均预测生存率和KM估计的观测生存率。分组越多，曲线上的点越密，看起来越细腻，但如果样本量不足，每组内的人数太少，KM估计就会变得很不稳定。小样本数据建议用15到20组，大样本可以适当增加到30组甚至更多。

cmethod = "hare"是观测生存率的估计方法。HARE是一种基于风险回归的平滑估计方法，比直接使用Kaplan-Meier更稳。一般保持默认就行，如果希望结果更传统，也可以改成cmethod = "KM"，两者结果差距不大。

3.3 读懂calibrate的输出与统计表格

图画出来之后，很多人只看一眼曲线形态，其实calibrate()还输出了非常有价值的量化信息。直接打印cal对象可以看到图例中的统计量，比如样本量n、事件数d、分组数p，以及模型在第pred.group百分位点上的C统计量。用summary()还能得到更详细的表格，显示在50%和90%分位点上的预测生存率与观测生存率：

# 查看校准汇总结果 summary(cal)

输出的结果大致长这样：

这个表格很实用。比如在50%分位点，模型预测的中位生存概率是0.321，而实际观察到的同等患者两年生存率是0.332，两者相差0.011，说明在校准中位水平上模型的预测偏差很小。90%分位点类似。我写论文时通常会把这个表整理成文字放进结果部分，比如“模型在50%和90%分位点上的预测生存率与实际观测生存率相差均在0.02以内，校准良好”。

判断校准好坏不要只看整条曲线贴不贴对角线，这两个分位点的数值差其实更直观。差异在0.05以内一般算可接受，超过0.1就要警惕模型在该风险区间的预测能力可能存在问题。

3.4 用ggplot把rms校准曲线重绘成论文风格

rms自带的plot()画出的图准确是准确，但样式比较老旧。很多期刊对图片清晰度、配色、字体都有要求，直接在原图基础上用还是会显得敷衍。我通常的做法是：提取calibrate()对象内部的数据，转成数据框后用ggplot2重新绘制，这样想怎么调就怎么调。

calibrate()返回的对象里，apparent存储的是表观校准数据，calibrated存储的是偏差校正数据，两者都是两列矩阵，第一列是预测生存率，第二列是观测生存率。提取出来就能画：

# 提取校准数据 cal_app <- as.data.frame(cal$apparent) cal_corrected <- as.data.frame(cal$calibrated) # 添加数据来源标记，方便合并绘图 cal_app$type <- "Apparent" cal_corrected$type <- "Bias-corrected" cal_all <- rbind(cal_app, cal_corrected) colnames(cal_all)[1:2] <- c("Predicted", "Observed") # ggplot2重绘 ggplot(cal_all, aes(x = Predicted, y = Observed, color = type)) + geom_abline(intercept = 0, slope = 1, linetype = "dashed", color = "black") + geom_line(aes(group = type), linewidth = 0.8) + geom_point(size = 2) + coord_equal() + scale_color_manual(values = c("Apparent" = "gray50", "Bias-corrected" = "firebrick")) + labs(x = "Predicted 2-Year Survival Probability", y = "Observed 2-Year Survival Probability", color = "") + theme_bw(base_size = 12) + theme(legend.position = "right")

这种图放在论文里比默认的base R图形美观很多。需要特别说明的是，cal对象里是否总是同时包含apparent和calibrated跨版本可能略有差异，但在我测试过的rms 6.x版本里都是这样的结构。如果你运行时报列名找不到，可以用str(cal)查看一下对象内部结构再提取。

另外提醒一句，保存图片时，务必设置合适的尺寸和分辨率。如果走期刊投稿，一般要求300 dpi以上，可以用ggsave(..., dpi = 300, width = 5, height = 4)保存，字体和线条粗细也要根据期刊要求调整。

3.5 riskRegression方案代码：coxph + Score + plotCalibration

如果你的模型是用survival::coxph()拟合的，不一定非要用rms重写一遍，直接用riskRegression就能完成校准。而且这套方案在绘制多个时间点的校准曲线时优势明显。

先拟合一个普通的Cox模型：

# 用survival包拟合Cox模型 fit_cox <- coxph(Surv(time, status2) ~ age + sex + ph.ecog, data = lung2, x = TRUE, y = TRUE)

然后调用Score()：

# 计算校准指标 score_obj <- Score( list("Cox model" = fit_cox), formula = Hist(time, status2) ~ 1, data = lung2, times = 730.5, # 2年时间点 metrics = "brier", # 同时计算Brier分数 plots = "calibration" ) # 绘制校准曲线 plotCalibration(score_obj, times = 730.5, type = "survival")

这里有两个容易忽略的细节。第一，Score()的formula参数必须用Hist(time, status2) ~ 1，而不是Surv(time, status2) ~ 1，这是riskRegression包的固定风格。第二，times参数指定想要校准的时间点，单位依然要和原始数据一致，用天就填730.5。

plotCalibration()画出来的图是分组散点加理想对角线，每个点代表预测生存率相近的一组患者。如果同时传入多个模型，还可以在图中叠加显示多条校准曲线，用于模型间的横向比较：

# 多模型校准曲线示例 fit_cox2 <- coxph(Surv(time, status2) ~ age + sex, data = lung2, x = TRUE, y = TRUE) score_obj2 <- Score( list("full model" = fit_cox, "reduced model" = fit_cox2), formula = Hist(time, status2) ~ 1, data = lung2, times = 730.5, metrics = "brier", plots = "calibration" ) plotCalibration(score_obj2, times = 730.5, type = "survival")

这样的图特别适合在模型比较章节展示：完整模型和简化模型的校准线都靠近理想线，但完整模型的曲线更贴近对角线，说明增加变量确实改善了校准度。审稿人看到这种图往往会觉得你的验证工作做得很扎实。

如果你想同时看1年、2年、3年三个时间点的校准曲线，times = c(365.25, 730.5, 1095.75)即可，plotCalibration()会按时间点分别出图，非常方便。

4. 校准曲线常见问题与排查实录

4.1 时间点选择的经验与单位转换

校准曲线的“校准时间点”选择，是我见过出错频率最高的地方。绝大多数人的困惑是：我应该选择哪个时间点？

我的经验是：优先选择临床上有决策意义的时间点，而不是机械地使用中位随访时间。比如做肺癌术后预测模型，医生和患者最关心的是3年或5年生存率，那校准时间点就选3年或5年，别因为随访时间只够2年就强行画5年节点，会因为没有足够的风险集样本而导致曲线末端剧烈波动。

时间单位这个坑也一样常见。如果你的time变量是以天为单位，时间点就只能写天数；如果是月，就写月份数；如果是年，就写年数。经常有人time.inc默认没设或者设成1，结果模型拟合和校准时间对不上，画出来的曲线完全不知所云。交叉检查的方法很简单：summary(fit_cph)会输出预先设定的time.inc值，先确认这一项是不是你想校准的时间节点再往下走。

如果你的数据是以年为单位，比如time记录为2.3、5.1这种，那time.inc = 2就代表2年。我建议无论原始单位是什么，最好在数据准备阶段把time统一成“年”或者“天”，并在代码注释中标明，免得隔一个月回来看自己都忘了。

4.2 内部校准与外部校准的区别

很多人的校准曲线是在建模数据本身上画的，这在严格意义上是“内部校准”。calibrate()里用Bootstrap重抽样的目的，就是通过反复在重抽样样本上拟合模型、再回到原始数据验证，来估算模型因为过拟合产生的乐观度，最后把这种乐观度从表观校准里减掉，得到偏差校正校准曲线。所以rms画出来的Bias-corrected线，比直接在原始数据上比较预测和观测的Apparent线更有参考价值。

如果你的数据量足够大，更好的做法是预留独立验证集，或者使用外部队列，在验证集上直接绘制校准曲线。这种情况下不需要Bootstrap校正，直接把模型预测值应用到新数据，按预测概率分组，比较各组预测生存率与KM估计的观测生存率即可。

代码实现也不难：把训练集上拟合好的模型拿到验证集上，用predict()算出验证集每位患者的预测生存概率，然后按分位数分组，每组内用survfit()估计观测生存率，最后画散点。这种方式得到的是真正的外部校准证据，临床说服力远高于内部校准。

4.3 常见报错速查表

下面把我这些年跑校准曲线时遇到的典型报错整理成表，每一条都是实测踩过的坑。

4.4 几个容易忽视的小细节

细节一：随机种子的设置。calibrate()依赖Bootstrap重抽样，不设种子的话，每次运行结果会在小数点后有微小波动。虽然不影响整体判断，但投稿时如果审稿人要求复现，不同次运行结果不一致会给人很不严谨的印象。所以在calibrate()之前，一定加set.seed()。

细节二：样本量与分组数的匹配。pred.group不是越大越好。我见过有人把预测概率分成50组，画出来的曲线密密麻麻全是锯齿，反而看不出趋势。样本量在200例左右时，20组是安全的选择；样本量上千的时候，可以考虑30到40组。如果你的模型目的只是看整体趋势，分组少一点反而更稳健。

细节三：一定要结合Brier分数。校准曲线是“看图说话”，Brier分数则把预测误差变成了一个数值，结合使用更利于论文量化报告。riskRegression的Score()函数里metrics = "brier"可以直接算出这个值，Brier分数越低说明预测越准确，生存模型中它同时考虑到了校准度和区分度。

细节四：校准曲线和列线图的关系。如果你之前已经画过列线图，那么这两者是配套的。列线图展示的是“怎么利用模型去预测”，校准曲线回答的是“这个预测到底有多靠谱”。做临床预测模型报告时，先列线图，再校准曲线，再加C-index，这一套组合拳是标准配置。

关于校准曲线，最后再分享一点个人体会。我在做模型验证时发现，校准曲线的形态其实比C-index更容易暴露模型的真实问题。有时候模型C-index看着还行，但校准曲线明显偏离对角线，说明模型在某些风险区间存在系统性偏差。这类问题光靠调阈值、改截断值掩盖不了，得从模型结构、变量筛选或者数据质量上找原因。画校准曲线不是流程上的走过场，它是真正帮你审视模型质量的一盏灯。希望这篇的代码和避坑经验，能让你少走一些我当时走过的弯路。