R语言实现模型的评估

一粒沙

发布于 2021-07-05 15:30:08

1.6K0

发布于 2021-07-05 15:30:08

在R语言中构建模型，有很多包进行了模型的封装。那么模型的评估在R中也有对应的包ipred。此包利用了bagging和boosting算法进行对模型的的评估。在这里我们介绍下这两个算法的区别：1）样本选择上：Bagging：训练集是在原始集中有放回选取的，从原始集中选出的各轮训练集之间是独立的；Boosting：每一轮的训练集不变，只是训练集中每个样例在分类器中的权重发生变化。而权值是根据上一轮的分类结果进行调整。2）样例权重：Bagging：使用均匀取样，每个样例的权重相等；Boosting：根据错误率不断调整样例的权值，错误率越大则权重越大。3）预测函数：Bagging：所有预测函数的权重相等；Boosting：每个弱分类器都有相应的权重，对于分类误差小的分类器会有更大的权重。4）并行计算：Bagging：各个预测函数可以并行生成；Boosting：各个预测函数只能顺序生成，因为后一个模型参数需要前一轮模型的结果。当然，bagging在学习算法模型不稳定和受训练模型影响很大的模型有更好的效果。接下来我们看下在这个ipred包中如何运行的。首先看下所需要的包：

install.packages(“ipred”)
install.packages(“MASS”)
install.packages(“rpart”)

接下来通过实例看下运行：

library("ipred")
library("rpart")
library("MASS")
library("survival")
data("GlaucomaM",package="TH.data")
gbag <- bagging(Class ~ ., data = GlaucomaM, coob=TRUE)

在这里的erro打分利用的是Brier分数。它是衡量概率校准的一个参数。简单来说，Brier分数可以被认为是对一组概率预测的“校准”的量度，或者称为“ 成本函数 ”，这一组概率对应的情况必须互斥，并且概率之和必须为1。Brier分数对于一组预测值越低，预测校准越好。例如：如果预测为70％（P = 0.70），并且下雨，则Brier评分为（0.70-1）2= 0.09。接下来看下结果的详细信息：

 #comb进行对模型进行自定义。
comb.lda <-list(list(model=lda, predict=function(obj, newdata) predict(obj, newdata)$x))
mod <-bagging(Class ~ ., data=GlaucomaM, comb=comb.lda)
predict(mod,newdata=GlaucomaM[c(1:3, 99:102), ])

接下来是间接分类模型的构建。所谓间接分类模型，就是将数据集分为三种类型的变量:用于预测类的变量(解释变量)、用于定义类的变量(中间变量)和类成员变量本身(响应变量)。中间变量是在解释变量的基础上建模的，响应变量是在中间变量上定义的。我们直接通过一个实例来看下间接分类模型的评估过程：

data("Smoking",package = "ipred")
classify <-function(data){
  data <- data[,c("TVPS","BPNL", "COHB")]
  res <- t(t(data) > c(4438, 232.5, 58))
  res <- as.factor(ifelse(apply(res, 1, sum)> 2, 1, 0))
  res
}
 
response <-classify(Smoking[ ,c("TVPS", "BPNL", "COHB")])
smoking <-data.frame(Smoking, response)
 
formula <-response~TVPS+BPNL+COHB~TarY+NicY+COY+Sex+Age
 
inclass(formula,data = smoking, pFUN = list(list(model = lm, predict =
mypredict.lm)),cFUN = classify)

最后是错误率评估功能，可以通过预测误差估计值的重新抽样评估生存数据的误分类误差、均方根误差或Brier评分。我们首先看下此功能最重要的参数estimator=c("cv","boot", "632plus")指的是进行评估的几种方法：cv交叉验证；boot基于bootstrap；632plus指的是632+ bootstrap基于偏差校正的bootstrap。接下来直接看实例：

data("iris")
library("MASS")
 
# force predictto return class labels only
mypredict.lda<- function(object, newdata)
  predict(object, newdata = newdata)$class
 
# 10-fold cv ofLDA for Iris data
errorest(Species~ ., data=iris, model=lda,
         estimator = "cv", predict=mypredict.lda)