预测模型构建利器——基于logistic的列线图(R语言)

简介: 预测模型构建利器——基于logistic的列线图(R语言)

模型知识准备

什么是列线图?


列线图(Alignment Diagram),又称诺莫图(Nomogram图),它是建立在多因素回归分析的基础上,将多个预测指标进行整合,然后采用带有刻度的线段,按照一定的比例绘制在同一平面上,从而用以表达预测模型中各个变量之间的相互关系。


列线图的基本原理,简单的说,就是通过构建多因素回归模型(常用的回归模型,例如Cox回归、Logistic回归等),根据模型中各个影响因素对结局变量的贡献程度(回归系数的大小),给每个影响因素的每个取值水平进行赋分,然后再将各个评分相加得到总评分,最后通过总评分与结局事件发生概率之间的函数转换关系,从而计算出该个体结局事件的预测值。


列线图将复杂的回归方程,转变为了可视化的图形,使预测模型的结果更具有可读性,方便对患者进行评估。正是由于列线图这种直观便于理解的特点,使它在医学研究和临床实践中也逐渐得到了越来越多的关注和应用。


如何看懂列线图?

image.png


我们可以看到列线图主要由左边的名称以及右边对应的带有刻度的线段所组成。


列线图的名称主要包括三类:


1. 预测模型中的变量名称:例如图中的年龄(Age)、ftv(孕前检查次数),ht(是否高血压),smoke(是否抽烟)等信息,每一个变量对应的线段上都标注了刻度,代表了该变量的可取值范围,而线段的长度则反映了该因素对结局事件的贡献大小。


2. 得分,包括单项得分,即图中的Point(第一行),表示每个变量在不同取值下所对应的单项分数,以及总得分,即Total Point(最后一行),表示所有变量取值后对应的单项分数加起来合计的总得分。


3. 预测概率:例如图中的Low weight rate,表示出生低体重儿的概率。


实例分析:如果产妇35岁,根据列线图,得分10分;无高血压,得分0分;抽烟,得分37分;再加上其他分数,如果总分100分,那么对应的生产低体重儿的概率就约为0.12,那么你就可以拿着这个图表给患者看,告诉她产出低体重儿的概率比较低。这样就到预测的效果了,懂了吗,那么最后就是如何搭建模型,最后达到这些效果。


什么是校准曲线?

校准曲线的解释:事实上,校准曲线是实际发生概率与预测概率的散点图。校准曲线可视化了霍斯默-勒梅休(Hosmer-Lemeshow)拟合优度检验的结果,所以除了校准曲线之外,我们还应该检查霍斯默-勒梅休拟合优度检验的结果。预测率和实际发生率越接近Y = X,Hosmer-Lemeshow拟合优度检验的p值大于0.05,模型校准得越好。在这种情况下,校准曲线几乎Y = X线,表示模型校准良好。


模型搭建

此次模型是基于R语言的多功能程序包,进行logistic回归预测,绘制列线图,绘制校准曲线,包括C指数等计算


为了更好的搭建我们的模型,我们首先对数据进行预处理,对数据进行了清洗和规整化,使得我们的每一个指标都是数值类型


#查看当前工作路径
getwd()
#设置当前工作路径
setwd("E:/个人文件/学术科研/身体约束预测模型")
getwd()
#导入第三方程序包
library(rms)
library(readr)
#导入数据
mydata=read.csv("数据.csv",encoding = 'UTF-8')

哑变量的设置

1、如果是二分类变量不用设为哑变量


2、如果是三分类以上的有序变量尽量不设哑变量,这样可以节省自由度,准确度会更高一些


3、非有序变量的三分类以上的变量必须设为哑变量,否则会影响结果


将哑变量引入回归模型,虽然使模型变得较为复杂,但可以更直观地反映出该自变量的不同属性对于因变量的影响,提高了模型的精度和准确度。但是我们也需要根据我们的实际的场景来对我们的数据进行处理。


如果要处理,我们就可以利用下面的代码进行处理,为什么此处不需要,就是因为元数据已经是一个有序的数值变量,而且没有规律的,不会很大的影响我们的模型,反而设置之后会影响。


处理分类变量

#此处用于哑变量的设置,此数据已经设置好了
mydata$leibie=as.factor(mydata$leibie)
mydata$n0=as.factor(mydata$n0)
mydata$shuxing=as.factor(mydata$shuxing)
mydata$age=as.factor(mydata$age)
mydata$time=as.factor(mydata$time)
mydata$n4=as.factor(mydata$n4)
mydata$n5=as.factor(mydata$n5)
mydata$n8=as.factor(mydata$n8)
mydata$n9=as.factor(mydata$n9)
mydata$n10new=as.factor(mydata$n10new)
mydata$MOBILITY=as.factor(mydata$MOBILITY)
mydata$n13=as.factor(mydata$n13)
mydata$n14=as.factor(mydata$n14)
mydata$n15=as.factor(mydata$n15)
mydata$n16=as.factor(mydata$n16)
mydata$n18=as.factor(mydata$n18)
mydata$n32=as.factor(mydata$n32)
mydata$n40new=as.factor(mydata$n40new)
#查看


mydata

#查看数据集(展示数据集里面的各个自变量的信息)

mydata=as.data.frame(mydata)

summary(mydata)

#打包数据

dd=datadist(mydata)

options(datadist='dd')

模型搭建


fit=lrm(n23~age+time+n4+n5+n7+n8+n9+n10new+n13+n14+n15+n16+n18
        +n19f1+n20+n21+n22+n32+n39new+n40new,data = mydata,x=T,y=T)
Logistic Regression Model
 lrm(formula = n23 ~ age + time + n4 + n5 + n7 + n8 + n9 + n10new + 
     n13 + n14 + n15 + n16 + n18 + n19f1 + n20 + n21 + n22 + n32 + 
     n39new + n40new, data = mydata, x = T, y = T)
                        Model Likelihood    Discrimination    Rank Discrim.    
                              Ratio Test           Indexes          Indexes    
 Obs          1026    LR chi2     800.88    R2       0.796    C       0.974    
  1            761    d.f.            20    g        4.249    Dxy     0.948    
  2            265    Pr(> chi2) <0.0001    gr      70.027    gamma   0.948    
 max |deriv| 2e-08                          gp       0.363    tau-a   0.363    
                                            Brier    0.054                     
           Coef    S.E.   Wald Z Pr(>|Z|)
 Intercept -7.7527 1.9490 -3.98  <0.0001 
 age        0.3654 0.1917  1.91  0.0566  
 time      -0.4446 0.1658 -2.68  0.0073  
 n4         0.2440 0.2370  1.03  0.3034  
 n5        -0.9956 0.3846 -2.59  0.0096  
 n7         0.2332 0.5105  0.46  0.6479  
 n8         0.7252 0.1620  4.48  <0.0001 
 n9        -0.6308 0.2890 -2.18  0.0291  
 n10new     0.8700 0.2214  3.93  <0.0001 
 n13        0.7915 0.3066  2.58  0.0098  
 n14       -0.2616 0.3375 -0.77  0.4383  
 n15        0.2198 0.2727  0.81  0.4203  
 n16        0.0793 0.1116  0.71  0.4774  
 n18        0.9823 0.2536  3.87  0.0001  
 n19f1     -0.1460 0.4530 -0.32  0.7473  
 n20        0.4265 0.3507  1.22  0.2240  
 n21       -0.6420 0.4020 -1.60  0.1103  
 n22        2.8574 0.3696  7.73  <0.0001 
 n32       -1.2263 0.1500 -8.17  <0.0001 
 n39new    -1.0837 0.3095 -3.50  0.0005  
 n40new     0.7555 0.3476  2.17  0.0298

数据格式:自变量~因变量,data=源数据,x=T,y=T


我们发现这个模型里面的自变量过多,我们可以通过查看相应的指标,看是否符合统计学意义,也就是我们用于判断模型变量的舍取


image.png


fit=lrm(n23~age+time+n5+n8+n9+n10new+n13+n18
        +n19f1+n20+n21+n22+n32+n39new+n40new,data = mydata,x=T,y=T)

全部满足


image.png


summary(fit)  查看


> summary(fit)
             Effects              Response : n23 
 Factor      Low High Diff. Effect    S.E.    Lower 0.95 Upper 0.95
 age         2   3    1      0.410530 0.18719  0.043632   0.777420 
  Odds Ratio 2   3    1      1.507600      NA  1.044600   2.175900 
 time        1   3    2     -0.876950 0.32704 -1.517900  -0.235960 
  Odds Ratio 1   3    2      0.416050      NA  0.219160   0.789810 
 n5          1   3    2     -2.078500 0.74955 -3.547600  -0.609450 
  Odds Ratio 1   3    2      0.125110      NA  0.028793   0.543650 
 n8          1   3    2      1.420900 0.30820  0.816850   2.025000 
  Odds Ratio 1   3    2      4.140900      NA  2.263400   7.575800 
 n9          1   4    3     -1.707400 0.84598 -3.365500  -0.049303 
  Odds Ratio 1   4    3      0.181340      NA  0.034545   0.951890 
 n10new      1   3    2      1.855200 0.42948  1.013400   2.696900 
  Odds Ratio 1   3    2      6.392800      NA  2.755000  14.834000 
 n13         1   3    2      1.296000 0.48649  0.342520   2.249500 
  Odds Ratio 1   3    2      3.654700      NA  1.408500   9.483300 
 n18         1   3    2      1.987700 0.49964  1.008400   2.967000 
  Odds Ratio 1   3    2      7.298600      NA  2.741200  19.433000 
 n19f1       1   2    1     -0.030359 0.44114 -0.894990   0.834270 
  Odds Ratio 1   2    1      0.970100      NA  0.408610   2.303100 
 n20         1   2    1      0.475350 0.34309 -0.197090   1.147800 
  Odds Ratio 1   2    1      1.608600      NA  0.821120   3.151200 
 n21         1   2    1     -0.537750 0.37185 -1.266600   0.191060 
  Odds Ratio 1   2    1      0.584060      NA  0.281800   1.210500 
 n22         1   2    1      2.751200 0.35705  2.051400   3.451100 
  Odds Ratio 1   2    1     15.662000      NA  7.779000  31.534000 
 n32         2   4    2     -2.455500 0.29487 -3.033400  -1.877500 
  Odds Ratio 2   4    2      0.085824      NA  0.048153   0.152970 
 n39new      1   2    1     -1.096700 0.30151 -1.687600  -0.505700 
  Odds Ratio 1   2    1      0.333990      NA  0.184960   0.603090 
 n40new      1   3    2      1.337900 0.65617  0.051788   2.623900 
  Odds Ratio 1   3    2      3.810900      NA  1.053200  13.790000 
列线图的绘制
nom=nomogram(fit,fun = plogis,fun.at = c(.001,.01,.05,seq(.1,.9,by=.1),.95,.99,.999),lp=T,
             funlabel = "身体约束预测")
plot(nom)


image.png


到此我们的列线图就绘制好了,至于如何去看,上面的文章也做了详细的解释和介绍!


校准曲线的绘制

cal=calibrate(fit,method = 'boot',B=1000)

plot(cal,xlim = c(0,1.0),ylim = c(0,1.0),xlab = "Predicted Probability",ylab = "observed Probability")


image.png


参考线:说明预测值和理论值百分百重合,这是不可能的,用于我们的参考


校准曲线反映我们的模型的预测的好坏


接下来我看看C指数:0.937


临床预测模型是临床研究的常用方法,通常我们构建出模型,需要对模型进行评估,其中一个重要指标是区分能力,也就是区分有病/无病、有效/无效、死亡/存活等结局的预测能力。区分能力一般用C统计量的大小来表示。


C:表示关联度

一般而言,C统计量的值>0.7,说明关联度较高,模型效果较好


image.png


说明共线性比较的强,综合判断我们的基于logistic机器学习算法,对身体约束的变量筛选和预测的模型搭建还是比较的合理,模型效果也比较的好。


对于结局是二分类的回归模型,ROC曲线下面积的值与C指数是一样的,95%可信区间也是一样的


至此我们的logistic模型就完成了!


相关文章
|
6天前
|
数据采集 机器学习/深度学习 数据可视化
R语言贝叶斯模型预测电影评分数据可视化分析
R语言贝叶斯模型预测电影评分数据可视化分析
|
12天前
|
数据可视化 Python
R语言分析糖尿病数据:多元线性模型、MANOVA、决策树、典型判别分析、HE图、Box's M检验可视化
R语言分析糖尿病数据:多元线性模型、MANOVA、决策树、典型判别分析、HE图、Box's M检验可视化
|
13天前
|
机器学习/深度学习 数据可视化 算法
数据分享|R语言交互可视化分析Zillow房屋市场:arima、VAR时间序列、XGBoost、主成分分析、LASSO报告
数据分享|R语言交互可视化分析Zillow房屋市场:arima、VAR时间序列、XGBoost、主成分分析、LASSO报告
|
13天前
|
机器学习/深度学习 数据可视化 Python
R语言神经网络模型预测多元时间序列数据可视化
R语言神经网络模型预测多元时间序列数据可视化
|
15天前
|
机器学习/深度学习 数据可视化 数据挖掘
PYTHON链家租房数据分析:岭回归、LASSO、随机森林、XGBOOST、KERAS神经网络、KMEANS聚类、地理可视化
PYTHON链家租房数据分析:岭回归、LASSO、随机森林、XGBOOST、KERAS神经网络、KMEANS聚类、地理可视化
|
15天前
|
机器学习/深度学习 存储 算法
【数据分享】R语言SVM和LDA文本挖掘分类开源软件存储库标签数据和词云可视化
【数据分享】R语言SVM和LDA文本挖掘分类开源软件存储库标签数据和词云可视化
|
17天前
|
机器学习/深度学习 自然语言处理 数据可视化
【数据分享】R语言对airbnb数据nlp文本挖掘、地理、词云可视化、回归GAM模型、交叉验证分析
【数据分享】R语言对airbnb数据nlp文本挖掘、地理、词云可视化、回归GAM模型、交叉验证分析
|
18天前
|
机器学习/深度学习 数据可视化 算法
【视频】主成分分析PCA降维方法和R语言分析葡萄酒可视化实例|数据分享
【视频】主成分分析PCA降维方法和R语言分析葡萄酒可视化实例|数据分享
|
19天前
|
自然语言处理 数据可视化 算法
Python主题建模LDA模型、t-SNE 降维聚类、词云可视化文本挖掘新闻组数据集
Python主题建模LDA模型、t-SNE 降维聚类、词云可视化文本挖掘新闻组数据集
|
25天前
|
算法 数据可视化
R语言ARIMA集成模型预测时间序列分析
R语言ARIMA集成模型预测时间序列分析