1. 交叉验证基本介绍
通常在建立模型后需要使用外部进行验证，以评估模型的外部可用性。然而，获取外部数据并不容易，这时交叉验证（Cross Validation）则是一种较好的可替代方案。交叉验证的方法很多，这里我们介绍最常用的k折交叉验证。
简单解释一下：
假如原始数据为100例患者，建模后我们使用K折交叉验证（K一般取3-10折）。
若取K=10，则：将原始数据分为10份（K值）：
9份做训练，1份做验证，这样循环10次（因为这样的组合有10种）。
取10次评价指标的平均值（如AUC值）作为模型的验证结果。即【数据分为K个子集，1个子集做验证，K-1个子集做训练，这样的取法有K种，所以会产生K个新数据集（训练集+训练集）。我们进行K次训练和验证得到的平均评价指标可能较为准确
为了保证数据分割的影响，目前的k折交叉验证一般会进行多次重复（200-1000次）。即进行200次的10折交叉验证。这样做出的结果可能会更加准确。
如文献所示
2.交叉验证的原理图
3.R Scripts
3.1数据示例
library("caret")
library(pROC)
bc<-read.csv("E:/r/test/buyunzheng.csv",sep=',',header=TRUE)
数据有8个指标，最后两个是PSM匹配结果，我们不用理他，其余六个为：
Education：教育程度，age:年龄，parity产次，induced：人流次数，case：是否不孕，这是结局指标，spontaneous：自然流产次数。有一些变量是分类变量，我们需要把它转换一下
bc$education<-ifelse(bc$education=="0-5yrs",0,ifelse(bc$education=="6-11yrs",1,2))
bc$spontaneous<-as.factor(bc$spontaneous)
bc$case<-as.factor(bc$case)
bc$induced<-as.factor(bc$induced)
bc$education<-as.factor(bc$education)
使用caret包的createFolds函数根据结果变量把数据分成10份，因为是随机分配，为了可重复性，我们先设立一个种子
set.seed(666)
folds <- createFolds(y=bc$case,k=10)###分成10份
接下来我们可以根据每一份的数据建立方程，并求出AUC，95%CI和截点，并且画出ROC曲线。
我们先来做第一个数据的，要提取列表的数据，需要做成[[1]]这种形式
fold_test <- bc[folds[[1]],]#取fold 1数据，建立测试集和验证集
fold_train <- bc[-folds[[1]],]#
建立预测模型
fold_pre <- glm(case ~ age + parity +spontaneous,
family = binomial(link = logit), data =fold_train )###建立模型
fold_predict <- predict(fold_pre,type='response',newdata=fold_test)##生成预测值
roc1<-roc((fold_test[,5]),fold_predict)
round(auc(roc1),3)##AUC
round(ci(roc1),3)##95%CI
绘制ROC曲线
plot(roc1, print.auc=T, auc.polygon=T, grid=c(0.1, 0.2),
grid.col=c("green", "red"), max.auc.polygon=T,
auc.polygon.col="skyblue",
print.thres=T)
嫌一个一个做比较麻烦的话我们也可以做成循环，一次跑完结果
先建立一个auc的空值，不然跑不了
auc_value<-as.numeric()###建立空值
for(i in 1:10){
fold_test <- bc[folds[[i]],] #取folds[[i]]作为测试集
fold_train <- bc[-folds[[i]],] # 剩下的数据作为训练集
fold_pre <- glm(case ~ age + parity +spontaneous,
family = binomial(link = logit), data =fold_train )
fold_predict <- predict(fold_pre,type='response',newdata=fold_test)
auc_value<- append(auc_value,as.numeric(auc(as.numeric(fold_test[,5]),fold_predict)))
}
#找到最大AUC的索引，用该索引对应的数据作为测试集，其它数据作为训练集即可训练出AUC最高的模型
max.index = which(auc_value==max(auc_value),arr.ind=TRUE)
fold_test <- bc[folds[[max.index]],]
fold_train <- bc[-folds[[max.index]],]
网友评论