上篇文章依照jama surgery的一篇文章给大家写了生存数据预测模型评价的C指数、校准曲线和模型验证结果的做法，其实生存数据预测模型的评价方法还有很多，本期接着往下看。

Time-dependent ROC

当结局是一个二分类变量的时候，考虑模型性能的两个指标一个叫灵敏度和特异度，我们希望两个都大，模型在不同的cutoff点时这两个指标如何变化的图示就是ROC曲线。详细参考之前的逻辑回归分类器的文章。

上面的说法是没有考虑时间因素的，认为结局是固定的，这个在大多数情况下都是适用的，就是说一般情况下我们的结局都是测一次或者说结局是不随着时间变化的，比如结局只有一个点，可以直接用ROC曲线评估模型。

The classical (standard) approach of ROC curve analysis considers event (disease) status and marker value for an individual as fixed over time, however in practice, both the disease status and marker value change over time

但是遇到生存数据，比如随访的病人这个时间点没死，随访时间延长可能就死了，这个时候模型的灵敏度和特异度在各个时间都不一样，光说ROC就没有意义，需要加上时间，就叫做time-dependent ROC。

捋一捋生存分析模型判断生存情况的原理：

模型会给每一个个案预测一个风险x(文献中叫做个体的marker)，t时刻归类到死亡或者存活需要一个标准c，模型根据这个风险x和标准c的大小比较结果判定个案是否死亡（阳性阴性）。

模型判断的是否死亡与个案的实际状况进行对比就有了混淆矩阵和灵敏度特异度。

D i (t)是t时间的实际状态，t时刻的敏感度和特异度的计算公式就是下图：

就是说t时刻个案确实是阳性D i (t)=1的同时模型也说我是“阳性”、“死亡”（marker比标准c大）叫做模型的真阳性率也就是敏感度，模型能将阳性个案识别为阳性的能力。

同理sp(c,t)就是特异度，是真阴率，模型能将无病识别为无病的能力。

本质上ROC曲线可以根据灵敏度和特异度两个指标来绘制的，我们知道了时间依赖的灵敏度和特异度，曲线也就可以做出来了。

上面的过程简化下就是：时间依赖的ROC就是不同时刻的结局计算出来的灵敏度和1-特异度构成的ROC。理论上生存结局可以有无数个ROC。

研究目的的不同，Time-dependent ROC又可以进行细分：

(1) cumulative/dynamic (C/D),

(2) incident/dynamic (I/D)

(3) incident/static (I/S)

看下图，ABCDEF6个个案进行随访，模型的判定是否死亡的标准是c:

这三种ROC的计算就在于如何界定某时刻t的阳性个案，比如(C/D)就是t时刻所有累计的阳性都算数，所以此时算灵敏度的分母按照下图就是ADE，随着时间的增大，阳性数量会增大，在计算特异度的时候只考虑DF的概率，所以叫做累计动态cumulative/dynamic (C/D)，这种方法也是临床使用最多的，它可以评价到t时刻为止模型的表现。

第二种叫做发病动态incident/dynamic (I/D)，在这种情况下只有A作为t时刻的实际阳性，BE都不是（在计算的时候BE直接扔掉了），然后CDF当作阴性，在算特异度的时候只考虑DF的概率；第三种叫做发病静态incident/static (I/S)，这种情况下依然是看t时刻发病，依然只有A作为阳性，超过某个区间（0，t*）都没有发病的作为阴性，即图中的DF作为阴性。这两种情况文章应用很少，不做重点记忆。

实操

画时间依赖ROC用到的包叫做timeROC，主要参数如下所示：

timeROC函数接受的第一个参数T是时间，delta参数是生存状态，删失一定要设定为0。还有marker参数，这个是模型的预测的自变量，我现在有数据如下，想做一个以

bili，chol，albumin为自变量形成的COX模型的时间依赖ROC，并且将时间设定在时间的每个5分位数上，我就可以写出如下代码：

timeROC(T=data$time,                      delta=data$status,marker=data$bili,                      other_markers=as.matrix(data[,c("chol","albumin")]),                      cause=1,weighting="marginal",                      times=quantile(data$time,probs=seq(0.2,0.8,0.1)),                      iid=TRUE)

输出结果如下，并且可以用confint得到AUC的置信区间：

然后就可以出图了，我们用plot函数设定好time参数（time一定要和timeROC的参数对应起来）就可以直接出图了如下：

plot(ROC.bili.cox,time=999.2)

输出如下，time改动的话图也会相应改动：

这样出来的画其实是很占空间的，因为你每个时间点都得画个图。常见的在文献中都会将几个时间点的ROC画在一个图上，然后加上图例进行区分，比如下面这篇来自Diabetes Care的时间依赖ROC图就是这样表示的：

那么对于我们的图，我们也可以进行修改成像文献中一样，直接用plot函数将add参数改为TRUE，然后依次画出各个时间点的ROC就可以，参考代码及输出如下：

plot(ROC.bili.cox,time=999.2)plot(ROC.bili.cox,time=1307.4,add = T,col="blue")plot(ROC.bili.cox,time=2555.7,add = T,col="green")legend(.7, .2, legend=c("关注", "公众号","Codewar"),       col=c("red", "blue","green"), lty=1, cex=0.7,bty = "n")

以上就是Time-dependent ROC的所有内容，我们接着往下看决策曲线的做法。

决策曲线

关于预测模型的决策曲线和校准曲线之前我有专门一篇文章做详细的解释，那么具体到生存分析中，我们依然来回忆一下：

Diagnostic and prognostic models are typically evaluated with measures of accuracy that do not address clinical consequences. Decision-analytic techniques allow assessment of clinical outcomes but often require collection of additional information may be cumbersome to apply to models that yield a continuous result.

就是我们知道模型的评估指标，比如对于生存分析来讲我们的模型Brier score和C-index表现的再好，其实也只能告诉你，模型本身不错。

但是指导临床决策我们还要考虑受益的问题，比如对于癌症来讲不必要的活检的危害其实是比没发现癌症的危害小很多的，我们宁愿病人多做一个检查也绝不愿意漏诊一个癌症。那么多少个额外活检的危害抵得上一次漏诊的危害，这个值是可变的，也是我们不知道的，但是我希望的的这个值无论如何变，通过我们的模型做出临床决策都是受益比瞎猜要好的，这个是我们的终极目的。而决策曲线就是看实际情况是不是符合刚刚的终极目的。

内在的逻辑请参看我之前文章。还有JAMA Guide to Statistics and Methods中的相应介绍，贴在下面就不给大家翻译了：

文献中常见的决策曲线图长这样

上图是来自jama子刊的一篇文章，它是将两个模型的dca画在了一个图上，我们依然是在R中进行对照复现。

做决策曲线用到的函数是dca函数，具体到cox模型我们需要指定时间点time：

同时经常出错的地方还在formula参数的写法，形成模型整体的决策曲线需要写模型的marker就是模型的预测风险，范围必须在0-1。这个是很多同学做错的地方。

比如我对做好的cox模型将time设定为15然后做决策曲线的示例代码如下：

dca(Surv(status, time) ~ cancerpredmarker, data = mydata, time = 15)

其中cancerpredmarker是我拟合的cox模型对个案的风险预测值，上图就显示了我的预测模型的决策曲线，多数文章也会对多个预测模型的DCA放一起进行比较，比如我做了模型1得到marker命名为向量“关注”，模型2命名为“公众号Codewar”，下面代码即刻让两个模型的DCA同时显示：

dca(Surv(status, time) ~ 关注+公众号Codewar,     data = mydata,    time = 15,    thresholds = 1:50 / 100) 

形成的图示如下：

具体图的各种风格大家还可根据ggplot的语法进行图里位置和图片主题风格的修饰。此文略过。还要给大家提醒的是如果模型是竞争风险做法依然是一样的，但是一定要将结局编码为因子同时删失水平编码为“censor”，其余做法都一样。

Competing risks endpoints are handled similarly to survival endpoints. The outcome must be defined as a factor, with the lowest level called "censor", and the other levels defining the events of interest. The dca() function will treat the first outcome listed as the outcome of interest

好啦，生存数据的预测模型的做法和评价系列到这儿就告一段落了，整体还是比较全面的，全部吃透发个jama应该在做法上不成问题了，有idea的话尽快联系我实施吧。

小结

今天接着前两期给大家补上了生存数据预测模型的时间依赖ROC和DCA的原理和做法。感谢大家耐心看完，自己的文章都写的很细，重要代码都在原文中，希望大家都可以自己做一做，请转发本文到朋友圈后私信回复“数据链接”获取所有数据和本人收集的学习资料。如果对您有用请先记得收藏，再点赞分享。

也欢迎大家的意见和建议，大家想了解什么统计方法都可以在文章下留言，说不定我看见了就会给你写教程哦，有疑问欢迎私信，有合作意向请直接滴滴我。

如果你是一个大学本科生或研究生，如果你正在因为你的统计作业、数据分析、模型构建,科研统计设计等发愁，如果你在使用SPSS, R，Mplus中遇到任何问题，都可以联系我。因为我可以给您提供最好的，最详细和耐心的数据分析服务。

如果你对Z检验，t检验，方差分析，多元方差分析，回归，卡方检验，相关，多水平模型，结构方程模型，中介调节，量表信效度等等统计技巧有任何问题，请私信我，获取详细和耐心的指导。

如果你或你的团队需要专业的科研数据清洗，建模服务，教学培训需求等等。请联系我。

If you are a student and you are worried about you statistical #Assignments, #Data #Analysis, #Thesis, #Reports, #Composing, #Quizzes, Exams.. And if you are facing problem in #SPSS, #R-Programming, #Excel, Mplus, then contact me. Because I could provide you the best services for your Data Analysis.

Are you confused with statistical Techniques like z-test, t-test, ANOVA, MANOVA, Regression, Logistic Regression, Chi-Square, Correlation, Association, SEM, multilevel model, mediation and moderation etc. for your Data Analysis...??

Then Contact Me. I will solve your Problem...

If You or Your Research Team Need Professional Scientific Data Cleaning, Model Building Services or Statistical Consulting... Please Contact Me.

往期精彩

R数据分析：生存数据的预测模型建立方法与评价

R数据分析：生存分析的列线图的理解与绘制详细教程

R数据分析：结合APA格式作图大法讲讲ggplot2和ggsci，请收藏

R数据分析：变量间的非线性关系，多项式，样条回归和可加模型

Mplus数据分析：性别差异gendergap的相关研究如何做？

R机器学习：分类算法之logistics回归分类器的原理和实现

R数据分析：PLS结构方程模型介绍，论文报告方法和实际操作

R数据分析：跟随top期刊手把手教你做一个临床预测模型

R数据分析：Lasso回归筛选变量构建Cox模型并绘制列线图

R数据分析：如何用层次聚类分析做“症状群”，实例操练

R数据分析：工具变量回归与孟德尔随机化，实例解析

R数据分析：潜类别轨迹模型LCTM的做法，实例解析

R文本挖掘：中文词云生成，以2021新年贺词为例

R机器学习：分类算法之判别分析LDA,QDA的原理与实现

R可视化：plot函数基础操作，小白教程

R机器学习：重复抽样在机器学习模型建立过程中的地位理解

R数据分析：用lme4包拟合线性和非线性混合效应模型

R数据分析：如何用mice做多重插补，实例解析

R数据分析：孟德尔随机化中介的原理和实操

R数据分析：生存分析的列线图的理解与绘制详细教程

R数据分析：cox模型如何做预测，高分文章复现

R数据分析：广义估计方程式GEE的做法和解释

R数据分析：潜类别轨迹模型LCTM的做法，实例解析

R数据分析:潜变量与降维方法(主成分分析与因子分析)

R数据分析：如何给结构方程画路径图，tidySEM包详解

R数据分析：自我报告的身高数据的离群值探索

R数据分析：生存分析与有竞争事件的生存分析的做法和解释

R机器学习：朴素贝叶斯与支持向量机的原理与实现

R数据分析：混合效应模型的可视化解释，再不懂就真没办法

R数据分析:如何理解模型中的“控制”，图例展示

R数据分析：tableone包的详细使用介绍

R数据分析：如何用lavaan包做结构方程模型，实例解析

R机器学习：分类算法之K最邻进算法（KNN）的原理与实现

R数据分析:潜增长模型LGM的做法和解释，及其与混合模型对比

R数据分析：论文中的轨迹的做法，潜增长模型和增长混合模型

R数据分析：纵向分类结局的分析-马尔可夫多态模型的理解与实操

R数据分析：临床预测模型实操，校准曲线和DCA曲线做法示例

R数据分析：国产新冠口服药比辉瑞好的文章的统计做法分享

R数据分析：再写潜在类别分析LCA的做法与解释

R数据分析：潜在转化分析LTA的做法和解释（一）

R机器学习:分类算法之K最邻进算法(KNN)的原理与实现

R数据分析：交互作用的简单斜率图做法及解释

R数据分析:双连续变量交互作用的简单斜率图作图及解释