百度360必应搜狗淘宝本站头条
当前位置:网站首页 > 文章教程 > 正文

文献解析:生存数据和分类结局列线图的做法,史上最全

xsobi 2024-11-24 00:29 1 浏览

今天给大家分析一篇文献,作者是浙大邵逸夫医院的章仲恒副主任,论文的主要内容就是教大家怎么画列线图。写的比较详细,拜读之后顺手写个粗浅的解析分享给大家,依然还是建议大家去看原文哦,论文题目见文章最后。

全文包括2分类逻辑斯蒂回归,多分类逻辑斯蒂回归,生存分析的列线图的画法和解析,基本上涵盖了临床上常用的预测方法,所以很值得大家收藏

接下来我一个一个给大家解析

作者先模拟了一个数据集,大概长这样:

这个数据集有1000个个案,6个变量,其中age和lac都是正态分布的连续自变量,sex和shock为因子,y是二分类结局,Y是多分类结局(3分类),所以我们如果用y做结局,就是拟合一个二分类逻辑斯蒂回归,用Y做结局就是多分类逻辑斯蒂回归。接下来我们一个一个看:

二分类结局的列线图画法

画列线图用到的包是rms。

第一步我们需要用datadist()函数来得到预测变量的分布,以此来规定我们的图形的尺度,紧接着拟合我们的逻辑回归模型(用到的函数为lrm),代码如下:

ddist <- datadist(data)
options(datadist='ddist')
mod.bi<-lrm(y~shock+lac*sex+age,data)

模型拟合好了之后,我们把这个模型直接喂给nomogram函数就行,代码如下:

nom.bi <- nomogram(mod.bi,          
											lp.at=seq(-3,4,by=0.5),
                     fun=function(x)1/(1+exp(-x)),
                     fun.at=c(.001,.01,.05,seq(.1,.9,by=.1),.95,.99,.999),
                     funlabel="Risk of Death",
                     conf.int=c(0.1,0.7),
                     abbrev=TRUE,
                     minlength=1,lp=F)

plot(nom.bi,lplabel="Linear Predictor",
     fun.side=c(3,3,1,1,3,1,3,1,1,1,1,1,3),
     col.conf=c('red','green'),
     conf.space=c(0.1,0.5),
     label.every=3,
     col.grid = gray(c(0.8, 0.95)),
     which="shock")

运行上面的代码即可得到下面的列线图:

上面的代码中我们自行设定了部分参数,需要给大家写写都是什么意思,首先是lp.at=seq(-3,4,by=0.5)表示我们的自变量的线性组合只在-3到4的范围内以步长为0.5进行显示,不过奇怪的是作者设置了这个参数,又把lp设置为F(就是不显示预测变量的线性组合值),所以作者在上面的代码中设置这个完全是没有必要的,如果我们将lp参数改为T则可看到预测变量的线性组合结果,效果如下:

然后是fun这个参数,因为我们的结局是二分类的,这个分类是通过连接函数得到的,我们的预测变量的线性组合直接得到的其实是log odds(看不明白的同学请看之前关于逻辑回归机器学习的文章),我们这儿其实想要的是某一类的概率,所以就用了一个反logit转换,目的是为了得到每一类的概率p,具体,代码中我们将参数fun设定为function(x)1/(1+exp(-x)),其实写成fun=plogis也是可以的,这个plogis就是逻辑斯蒂的累计分布函数(Logistic Cumulative Distribution Function (plogis Function))

为了更好地给大家说明,这儿给大家做个演示,运行如下代码大家就可以看到逻辑斯蒂的累计分布函数长啥样:

x_plogis <- seq(- 10, 10, by = 0.1) 
y_plogis <- plogis(x_plogis) 
plot(y_plogis) 

看到没,纵轴就是我们想要的p,这个就是通过反logit转换得到p的原理,好好体会哈。

其实逻辑回归在做的就是把自变量的线性组合通过连接函数映射到上面这个图上(就是得到p)从而达到对2分类结局分类的目的,在做列线图的时候我们也得把这个东西以fun参数的形式给指出来,不然它就不会对自变量进行转换,出来的概率就是错误的,转换就是通过fun这个参数实现的。

再看fun.at参数,这个参数是设定函数轴的标签的,我们将其设定为fun.at=c(.001,.01,.05,seq(.1,.9,by=.1),.95,.99,.999)就是希望如果自变量的组合会导致这些值的话,就会在相应的地方打上标签,因为我们自变量的线性组合最小为-3其对应的p为0.0474,所以在途中标签是从0.05开始的;funlabel参数是用来规定轴的名字的。conf.int参数用来设定变量的confidence levels,默认不显示。

接下来看看plot函数中的几个参数:lplabel用来设定预测变量线性组合的标签,fun.side用来设定fun刻度值标签的上下位置,防止重叠;col.conf给置信水平上颜色;conf.space规定置信轴的垂直位置;label.every=3意思是轴的每3个值打一个标签;方便我们对齐的垂直轴的颜色可以通过col.grid调节。

有序结局变量的列线图画法

有序逻辑斯蒂回归也是用lrm()函数进行拟合的,代码如下:

mod.ord <- lrm(Y ~ age+rcs(lac,4)*sex)

上面的代码中拟合了加了限制性样条的sex和lac的交互项(之后会给大家写什么是样条,本篇文章略过)。

为了更好地给大家解释有序逻辑回归,我们先看看有序逻辑回归的模型表达:

我们在做有序逻辑回归的时候,比如我们的结局变量有4个水平,这个时候我们会同时拟合3个二分逻辑模型(下图中的1,2,3),我们始终是有一个参考水平的,解释的时候依然和二分逻辑回归一样的解释:

比如我们的例子中,我们就会同时拟合如下3个模型,第1个模型回答了个案是不是大于第1类,第二个模型回答了是不是大于第二类,第3个模型回答了是不是大于第3类,这样就穷尽了所有的可能性。

相应地,我们的概率p也得有三个:

fun2 <- function(x) plogis(x-mod.ord$coef[1]+mod.ord$coef[2])
fun3 <- function(x) plogis(x-mod.ord$coef[1]+mod.ord$coef[3])

以上就规定好了出图时的自变量线性组合的转换方法,接着就可以出图了,代码为:

nom.ord <- nomogram(f, fun=list('Prob Y>=1'=plogis,
                                'Prob Y>=2'=fun2,
                                'Prob Y=3'=fun3),
                    lp=F,
                    fun.at=c(.01,.05,seq(.1,.9,by=.1),.95,.99))
plot(nom.ord, lmgp=.2, cex.axis=.6)

因为我们的模型中是有交互的,是lac和sex交互,其中sex是有两个水平,那么我们出来的列线图就会在sex的两个水平分别画出lac,于是我们做出来的列线图就是下面这个样子的:

我们有看到,在列线图中结局Y不同顺序的概率都有出现,通过这个列线图就可以很方便地得到某个个案到底处于结局变量的哪个水平的概率最大了。

生存分析的列线图

对于生存数据我们可以用半参数模型和参数模型进行建模(见我之前的文章),在生存分析中我们一般会感兴趣特定时点的生存概率和中位生存期,当然啦,这些模型都是可以通过列线图进行可视化的。

首先我们需要拟合生存模型,用到的是psm函数,本例中我们用到的数据长这样:

需要用到的变量是(time,status,ph.ecog,sex和age,其中ph.ecog是病人日常活动功能的一个指标。我们现在就拟合一个模型来研究ph.ecog对病人生存的影响,将sex和age当作协变量,代码如下:

 mod.sur <- psm(Surv(time,status) ~ ph.ecog+sex+age,
 data, dist='weibull')

然后我们可以用下面的代码得到病人的中位生存期和生存概率:

med <- Quantile(mod.sur)
surv <- Survival(mod.sur)
ddist <- datadist(lung)

到这儿,我们画图需要的东西就全部准备好了,接下来进行出图,首先是以生存期为结局的图:

nom.sur1<-nomogram(mod.sur,
                   fun=list(function(x) med(lp=x, q=0.5),
                            function(x) med(lp=x,q=0.25)),
                   funlabel=c("Median Survival Time",
                              "1Q Survival Time"),
                   lp=F)
plot(nom.sur1,
     fun.side=list(c(rep(1,7),3,1,3,1,3),rep(1,7)),
     col.grid = c("red","green"))

运行代码就可以得到下面的图啦:

可以看到我们的结局变量有两,分别是中位生存期和第一四分位生存期,因为我们定义fun的时候是定义了两嘛。

接下来是以生存概率为结局进行出图,代码如下:

nom.sur2 <- nomogram(mod.sur, fun=list(function(x)surv(200, x),
  function(x) surv(400, x)),
  funlabel=c("200-Day Survival Probability",
             "400-Day Survival Probability"),
  lp=F)

上面的代码分别定义了随访200天和400天时病人的生存概率,也就是我们的结局,接下来就用plot函数出图了,代码如下:

plot(nom.sur2,
     fun.side=list(c(rep(c(1,3),5),1,1,1,1),
                   c(1,1,1,rep(c(3,1),6))),
     xfrac=.7,
     col.grid = c("red","green"))

代码中的xfrac参数是用来规定轴标签和图距离的哈,运行代码就可以得到下面的列线图:

半参生存分析

半参生存分析其实就是比例风险模型,为什么叫半参呢,这儿再给大家补习一下:

A parametric survival model is one in which survival time (the outcome) is assumed to follow a known distribution. Examples of distributions that are commonly used for survival time are: the Weibull, the exponential (a special case of the Weibull), the log-logistic, the log-normal, etc.

The Cox proportional hazards model, by contrast, is not a fully parametric model. Rather it is a semi-parametric model because even if the regression parameters (the betas) are known, the distribution of the outcome remains unknown. The baseline survival (or hazard) function is not specified in a Cox model (we do not assume any shape or form).

我们依然是先拟合模型,以ph.ecog为i变量并控制sex和age,比例风险模型的拟合用的是cph()函数,本例中代码如下:

mod.cox <- cph(Surv(time,status) ~ ph.ecog+sex+age,
               lung, surv=TRUE)

接下来是出图前的准备(注意看fun参数):

ddist <- datadist(lung)
options(datadist='ddist')
surv.cox <- Survival(mod.cox)
nom.cox <- nomogram(mod.cox, fun=list(function(x)
  surv.cox(200, x),
  function(x) surv.cox(400, x)),
  funlabel=c("200-Day Sur. Prob.",
             "400-Day Sur. Prob."),
  lp=F)

然后出图:

plot(nom.cox,
fun.side=list(c(rep(c(1,3),5),1,1,1,1),
c(1,1,1,rep(c(3,1),6))))

运行代码即可做出如下的列线图:


通过上图就可以很方便地得到某个病人200天和400天时的生存概率了,好了,以上就是论文的所有内容,就给大家解析到这儿,希望对做列线图的同学有用。

论文原文:Zhang, Z., & Kattan, M. W. (2017). Drawing Nomograms with R: applications to categorical outcome and survival data. Annals of translational medicine, 5(10).

小结

今天给大家写了临床常见预测模型的列线图的做法,感谢大家耐心看完,自己的文章都写的很细,代码都在原文中,希望大家都可以自己做一做,请关注后私信回复“数据链接”获取所有数据和本人收集的学习资料。如果对您有用请先收藏,再点赞转发。

也欢迎大家的意见和建议,大家想了解什么统计方法都可以在文章下留言,说不定我看见了就会给你写教程哦,另,咨询代做请私信。

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

如果你对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...

猜你喜欢

R数据分析:列线图的做法及解释

R数据分析:手把手教你画列线图(Nomogram)及解读结果

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

相关推荐

js向对象中添加元素(对象,数组) js对象里面添加元素

一、添加一个元素对象名["属性名"]=值(值:可以是一个值,可以是一个对象,也可以是一个数组)这样添加进去的元素,就是一个值或对象或数组...

JS小技巧,如何去重对象数组?(一)

大家好,关于数组对象去重的业务场景,想必大家都遇到过类似的需求吧,这对这样的需求你是怎么做的呢。下面我就先和大家分享下如果是基于对象的1个属性是怎么去重实现的。方法一:使用.filter()和....

「C/C++」之数组、vector对象和array对象的比较

数组学习过C语言的,对数组应该都不会陌生,于是这里就不再对数组进行展开介绍。模板类vector模板类vector类似于string,也是一种动态数组。能够在运行阶段设置vector对象的长度,可以在末...

如何用sessionStorage保存对象和数组

背景:在工作中,我将[{},{}]对象数组形式,存储到sessionStorage,然后ta变成了我看不懂的形式,然后我想取之用之,发现不可能了~记录这次深刻的教训。$clickCouponIndex...

JavaScript Array 对象 javascript的array对象

Array对象Array对象用于在变量中存储多个值:varcars=["Saab","Volvo","BMW"];第一个数组元素的索引值为0,第二个索引值为1,以此类推。更多有...

JavaScript中的数组Array(对象) js array数组

1:数组Array:-数组也是一个对象-数组也是用来存储数据的-和object不同,数组中可以存储一组有序的数据,-数组中存储的数据我们称其为元素(element)-数组中的每一个元素都有一...

数组和对象方法&amp;数组去重 数组去重的5种方法前端

列举一下JavaScript数组和对象有哪些原生方法?数组:arr.concat(arr1,arr2,arrn);--合并两个或多个数组。此方法不会修改原有数组,而是返回一个新数组...

C++ 类如何定义对象数组?初始化数组?linux C++第43讲

对象数组学过C语言的读者对数组的概念应该很熟悉了。数组的元素可以是int类型的变量,例如int...

ElasticSearch第六篇:复合数据类型-数组,对象

在ElasticSearch中,使用JSON结构来存储数据,一个Key/Value对是JSON的一个字段,而Value可以是基础数据类型,也可以是数组,文档(也叫对象),或文档数组,因此,每个JSON...

第58条:区分数组对象和类数组对象

示例设想有两个不同类的API。第一个是位向量:有序的位集合varbits=newBitVector;bits.enable(4);bits.enable([1,3,8,17]);b...

八皇后问题解法(Common Lisp实现)

如何才能在一张国际象棋的棋盘上摆上八个皇后而不致使她们互相威胁呢?这个著名的问题可以方便地通过一种树搜索方法来解决。首先,我们需要写一个函数来判断棋盘上的两个皇后是否互相威协。在国际象棋中,皇后可以沿...

visual lisp修改颜色的模板函数 怎么更改visual studio的配色

(defunBF-yansemokuai(tuyuanyanse/ss)...

用中望CAD加载LISP程序技巧 中望cad2015怎么加载燕秀

1、首先请加载lisp程序,加载方法如下:在菜单栏选择工具——加载应用程序——添加,选择lisp程序然后加载,然后选择添加到启动组。2、然后是添加自定义栏以及图标,方法如下(以...

图的深度优先搜索和广度优先搜索(Common Lisp实现)

为了便于描述,本文中的图指的是下图所示的无向图。搜索指:搜索从S到F的一条路径。若存在,则以表的形式返回路径;若不存在,则返回nil。...

两个有助于理解Common Lisp宏的例子

在Lisp中,函数和数据具有相同的形式。这是Lisp语言的一个重大特色。一个Lisp函数可以分析另一个Lisp函数;甚至可以和另一个Lisp函数组成一个整体,并加以利用。Lisp的宏,是实现上述特色的...