R数据分析:生存分析的列线图的理解与绘制详细教程
阅读原文时间:2023年07月09日阅读:2

列线图作为一个非常简单明了的临床辅助决策工具,在临床中用的(发文章的)还是比较多的,尤其是肿瘤预后:

Nomograms are widely used for cancer prognosis, primarily because of their ability to reduce statistical predictive models into a single numerical estimate of the probability of an event, such as death or recurrence, that is tailored to the profile of an individual patient.

找个公开数据库做生存分析出个列线图,然后出个文章是很多临床同学可以依赖的较容易的实现路径,之前有给大家介绍过列线图,今天开始再给大家比较详细地写写生存分析列线图系列,希望可以对大家有帮助。

理解列线图

要弄明白生存分析的列线图的出图逻辑。我们首先来回顾下cox模型究竟是拟合是什么东西,看下图:

在基础风险确定后,乘上以e为底数的指数函数(我们关心的协变量的线性部分是在指数上)就可以得到风险函数(为什么能这么做就涉及到比例风险假设)。通过线性部分的指数函数和基础风险我们cox模型最终得到的是hazard function。

通过hazard function我们可以得到hazard rate,但是对于临床应用来讲,临床医师关心的东西更直观,他们关心的是具体协变量条件下个体的生存概率,画出的列线图常常如下面所示:

列线图中的结局常常是某个时间点的生存概率,这就要求我们在统计处理上做出转换。就是将风险函数转换成生存函数进而得到预测的生存概率。

接下来我们就来详细地过一遍实操重点。

本文中涉及到的文献图片和方法描述均来自JAMA Surg杂志的文章,文章引用如下:

Hyder O, Marques H, Pulitano C, et al. A Nomogram to Predict Long-term Survival After Resection for Intrahepatic Cholangiocarcinoma: An Eastern and Western Experience. JAMA Surg. 2014;149(5):432–438. doi:10.1001/jamasurg.2013.5168

变量选择

首先看变量筛选,经常我们用来做模型的数据库中有很多变量,列线图作为一个临床应用工具,变量肯定是越少越好的(让医生算分算半个多小时总分总是不合适推广的嘛,虽然大家做论文都不关心临床转换,但是还要有这个意识比较好),所以必须精选,这篇文章用到的方法叫做Backward stepwise selection:

Backward stepwise selection using the AIC in Cox proportional hazards regression modeling identified 6 variables that were the most associated with survival: age, tumor size, multiple lesions, nodal status, vascular invasion, and presence of cirrhosis of the underlying liver

我看生存分析列线图的文章这个方法用得还是比较普遍的哈,基本都是单变量筛过之后再来个stepwise selection:

整体这个方法在R中操作也是非常方便的,像rms包中专门就有fastbw函数进行倒退法的逐步筛选。

可以用aic为标准,也可以用p值为标准进行筛,很方便的。

用生存分析模型出列线图

首先明确,同学们不要再称呼“列线图模型”了,列线图只是具体模型的可视化、工具化表示,他本身不是模型。模型本身具体要看你到底做的是什么统计模型,比如逻辑模型,比如线性回归模型,再比如今天写的COX模型,这个才叫模型。

A nomogram is a graphical representation of a mathematical model involving several pre- dictors to predict a particular endpoint based on traditional statistical methods such as Cox proportional hazards model for survival data or logistic regression for binary outcome

好多同学问我能不能帮忙做一个列线图模型,其实这种表达我是摸不着头脑的。

出列线图,首先要确定内在的统计模型,比如今天写生存数据的列线图,我就要先做一个COX模型,然后再借助nomogram函数出图,这个函数的参数很多,下图只是部分参数:

可以看到这个函数需要的第一个参数就是一个做好的模型fit。具体到生存分析的列线图,我们就需要先跑一个cox模型出来,然后对跑模型的数据集d进行下面的操作:

ddist <- datadist(d); options(datadist='ddist')

生存分析的列线图需要调节的地方可以有很多。

比如我们做一个生存分析,在nomogram函数中不设定任何参数直接去出列线图的话,出出来的图是这样的:

图中只会有cox模型线性部分的预测值,这个时候我们需要将线性预测值转换为生存概率才符合临床应用实际,就像下图发表中的文献一样:

此时要做的就是进行风险函数和生存概率函数的转换。我们需要设定转换的代码如下:

surv <- Survival(f)

通过上面的代码我们就可以将cph函数拟合出来的风险函数转换成生存函数,从而在列线图的绘制中我们可以规定显示具有临床意义的时间点的某个个体的生存概率。比如我想得到3年和6年的生存概率为结局的列线图,我就可以写出如下代码:

plot(nomogram(f, fun=list(function(x) surv(3, x),
                            function(x) surv(6, x)),
            funlabel=c("age 3 Survival Probability",
                       "age 6 Survival Probability"))

上面的代码中f是cph对象,fun中给定的就是将线性预测值转换成生存概率的函数,运行代码即可出图如下:

并且针对nomogram可以做很多的个性化的修饰,比如lp参数可以控制是否显示线性预测值的打分轴,lp.at和fun.at可以控制线性预测值打分轴和转换函数显示的点。比如对于上图,我希望线性预测值的轴只显示0到4的点,我就可以写参数:

lp.at = c(0,1,2,3,4)

就可以实现。

还有,有时候我们分类变量的水平比较多,名字比较长,我们可能会将abbrev参数设定为TRUE来缩略显示长度,比如你注意下面的图和上面的在sex水平上的显示区别就是因为我们将abbrev参数设定为了T:

我们还可以很方便地改变列线图的轴标签,只需要将变量打上我们想要的标签,比如将两个变量标签分别设定为“关注”和“Codewar”后,将下面的参数在nomogram函数中设定一下:

vnames='labels'

运行后查看效果:

对于上面的列线图,我可能还觉得我们的图轴和标签离得有点远,这个时候我就可以将xfrac设定小一点比如我设定为0.2,这个时候图就会紧凑很多;我们还可以通过tcl参数设定轴的刻度标线的长度,比如我设定为1,这时候图的刻度线就会变长,读图就会更轻松。

xfrac=.2,tcl=1

参数像上面设定后,我们的图就如下所示:

调了一下还是蛮有效果的哈,但是我还是不满意,看人家jama的列线图,背景色都有,淡淡的蓝色显得就很高级,这个操作大家只需要在出图前设定:

par(bg = "aliceblue")

然后再plot效果就有了:

这下一看就是高分杂志的的图,背景色中的aliceblue你也可以任意改成你喜欢的颜色。

nomogram函数还有很多的参数可调,一篇文章肯定是写不完的,其余的调节功能留给大家自己探索了吧。

接下来写读图的部分。

学会读图

为了更加的加深大家对模型本身和列线图这个可视化工具的区别的理解,我们今天带大家结合cox模型来读cox模型的列线图。

首先我们学会读线性预测值,首先再一次回忆模型的表达:

线性部分就是表达式中指数函数的指数部分,比如我现在跑了一个cox模型结果如下:

那么我知道age的线性系数0.0419,sex中male的线性系数是-0.5975 ,所以我们的模型对一个10岁的男性线性部分预测值就应该为10*0.0419-0.5975=-0.178,回到我们这个模型的列线图中

回到列线图:我们可以看到10岁的得分是0分,男性得分为0,总分0分,对应的线性预测值大概也为-0.18(大家可以用尺子比个大概哈),达成一致。

我们再看生存概率的读法,比如对于一个100岁的男人来讲,依照下面的列线图,她的年龄得分应该是100,性别得分是0,总分是100,对应的3年的生存概率应该大约是0.62(大家可以用把尺子比对下哈):

然后我们出列线图的内在模型再一次验证,我们用predictSurvProb函数,将新数据设定为1个100岁的男性,times设定为3,用原来的cox模型预测出来的生存概率确实也是0.627。依然达成一致。

上面就是读图方法与模型结果的相互验证,希望能够进一步加深列线图只是模型的可视化的表示这一概念的理解。

好了,今天的文章重点就放在列线图出图上,文章中还有报告决策曲线和校准曲线,C指数等下一次再给大家详细写。