R语言机器学习与临床预测模型26--综合判别改善指数
本内容为【科研私家菜】R语言机器学习与临床预测模型系列课程
R小盐准备介绍R语言机器学习与预测模型的学习笔记
你想要的R语言学习资料都在这里， 快来收藏关注【科研私家菜】
01 什么是综合判别改善指数？
综合判别改善指数(Integrated Discrimination Improvement, IDI)。这个指标也用于判断预测模型改善情况，与上一讲的净重新分类指数(NRI)有类似也有不同。
IDI是由Pencina等人于2008年提出的一个非常新的判别指标。由于它考虑了不同切点的情况，可以用来反映模型的整体改善状况，在一定程度上补齐了NRI的短板。同时，虽然AUC也考虑到了不同切点，但是AUC的改善情况在临床中不易解释，IDI也因此弥补了AUC的缺陷，可以形象地展示研究对象被准确重新判别的比例。
因此，当我们在进行2个疾病模型比较，或者2个指标诊断效能比较时，除了传统的ROC曲线及其AUC，也可以同时给出NRI和IDI，更加全面多层次的展示模型的改善情况。
02 IDI计算公式
IDI的计算其实也比较简单，它反映的是两个模型预测概率差距上的变化，因此是基于疾病模型对每个个体的预测概率计算所得。它的计算方法为：
其中Pnew,events、Pold,events表示在患者组中，新模型和旧模型对于每个个体预测疾病发生概率的平均值，两者相减表示预测概率提高的变化量，对于患者来说，预测患病的概率越高，模型越准确，因此差值越大则提示新模型越好。
在模型比较时，将两部分相减即可得到IDI，总体来说IDI越大，则提示新模型预测能力越好。与NRI类似，若IDI＞0，则为正改善，说明新模型比旧模型的预测能力有所改善，若IDI＜0，则为负改善，新模型预测能力下降，若IDI=0，则认为新模型没有改善。
我们可以通过计算Z统计量，来判断IDI与0相比是否具有统计学显著性，统计量Z近似服从正态分布，公式如下：
02 二分类变量R语言实现
示例数据来自于survival包里自带的一份梅奥诊所的数据，记录了418位患者的临床指标与原发性胆汁性肝硬化（PBC）的关系。其中前312个患者来自于RCT研究，其他患者来自于队列研究。我们用前312例患者的数据来预测2000天时间点上是否发生死亡。此处需要说明的是原始数据是一个生存数据，我们重新定义一个二分类的结局，死亡or 存活，不考虑时间因素。先载入这份数据，如图4.所示。这个表中的结局变量是status，0 = 截尾（删失），1 = 接受肝移植，2 = 死亡。胆我们的研究目的“死亡与否”是个二分类变量，所以要做变量变换。再看time一栏，有的不够2000天，这些样本要么是没到2000天就死亡了，要么是删失了。我们要删掉2000天内删失的数据。
##hereconsider pbc dataset in survival package as an example
library(survival)
dat=pbc[1:312,]
dat$sex=ifelse(dat$sex==''f'',1,0)
##subjectscensored before 2000 days are excluded
dat=dat[dat$time>2000|(dat$time<2000&dat$status==2),]
##predcitingthe event of ''death'' before 2000 days
event=ifelse(dat$time<2000&dat$status==2,1,0)
##standardprediction model: age, bilirubin, and albumin
z.std=as.matrix(subset(dat,select=c(age,bili,albumin)))
##newprediction model: age, bilirubin, albumin, and protime
z.new=as.matrix(subset(dat,select=c(age,bili,albumin,protime)))
##glmfit (logistic model)
mstd=glm(event~.,binomial(logit),data.frame(event,z.std),x=TRUE)
mnew=glm(event~.,binomial(logit),data.frame(event,z.new),x=TRUE)
##UsingPredictABEL package
library(PredictABEL)
pstd<-mstd$fitted.values
pnew<-mnew$fitted.values
##用cbind函数把前面定义的event变量加入数据集，并定义为dat_new
dat_new=cbind(dat,event)
##计算NRI,同时报告了IDI,IDI计算与cutoff点设置无关。
##cOutcome指定结局变量的列序号
##predrisk1,predrisk2为新旧logistic回归模型
reclassification(data=dat_new,cOutcome=21,
predrisk1=pstd,predrisk2=pnew,
cutoff=c(0,0.2,0.4,1))
# 计算的IDI为0.44%，说明新模型较旧模型预测能力仅改善0.44%。
03 生存资料模型
生存资料的NRI与分类结局的NRI差别在于前者需要构建Cox回归模型，所以我们首先构建新旧Cox回归模型，计算这两个模型的NRI.
##here consider pbc dataset in survival package as an example
library(survival)
dat=pbc[1:312,]
dat$time=as.numeric(dat$time)
##定义生存结局
dat$status=ifelse(dat$status==2,1, 0)
##定义时间点
t0=365*5
##基础回归模型
indata0=as.matrix(subset(dat,select=c(time,status,age,bili,albumin)))
##增加1个预测变量新模型
indata1=as.matrix(subset(dat,select=c(time,status,age,bili,albumin,protime)))
##旧模型中预测变量矩阵
covs0<-as.matrix(indata0[,c(-1,-2)])
关注R小盐，关注科研私家菜(VX_GZH: SciPrivate)，有问题请联系R小盐。让我们一起来学习 R语言机器学习与临床预测模型
