glm()模型的交叉验证

Question

我想对我之前在R中建立的一些glm模型进行10倍交叉验证。我对boot包中的cv.glm()函数有点困惑，尽管我已经阅读了很多帮助文件。当我提供以下公式时:

library(boot)
cv.glm(data, glmfit, K=10)

这里的 "data" 参数是指整个数据集还是仅指测试集？

到目前为止我看到的例子都将 "data" 参数作为测试集，但那并没有太多意义，比如为什么在同一测试集上进行10次交叉验证？它们都会给出完全相同的结果（我猜测！）。

不幸的是，?cv.glm 以模糊的方式解释了它：

data: 包含数据的矩阵或数据框。行应该是案例，列对应于变量之一为响应

我的另一个问题是关于 $delta[1] 的结果。这是10次试验中的平均预测误差吗？如果我想获得每个折叠的误差怎么办？

这是我的脚本：

##data partitioning
sub <- sample(nrow(data), floor(nrow(x) * 0.9))
training <- data[sub, ]
testing <- data[-sub, ]

##model building
model <- glm(formula = groupcol ~ var1 + var2 + var3,
        family = "binomial", data = training)

##cross-validation
cv.glm(testing, model, K=10)

- Error404

---

请查看 boot:::cv.glm 的示例部分。您应该输入完整的数据、模型和 CV 的折叠。 - Roman Luštrik

感谢您的回复@RomanLuštrik。听起来很不错。但我仍然有一些疑问。这个函数是否使用交叉验证中提供的所有数据？假设我为cv.glm(data, glm, K=10)提供了一个1000行的数据框，它会将数据分成10个大小为100的部分并进行交叉验证吗？抱歉，我已经查看了?cv.glm，但没有找到答案。 - Error404

如果您要进行二倍交叉验证，该函数将使用50％的数据来拟合模型。它将使用另外50％的数据来查看模型描述数据的效果。或者，在留一法交叉验证中，它会对除一个数据“点”外的所有数据点进行模型拟合，并观察这个被单独拿出来的“点”的表现情况。重复N次即可得到结果。 - Roman Luštrik

您好@RomanLuštrik。您说如果我进行2倍交叉验证，该函数将根据50%的数据拟合模型，并使用另外50%作为测试集。如果函数这样做，那么为什么它需要一个名为“glmfit”的参数，该参数是先前拟合的模型？ - Error404

如果您对交叉验证有疑问，我建议您在crossvalidated.com上开一个帖子。 - Roman Luštrik

Answer 1

我对使用各种包中的10倍交叉验证方法总是有点谨慎。 我有自己的简单脚本手动创建任何机器学习包的测试和训练分区：

#Randomly shuffle the data
yourData<-yourData[sample(nrow(yourData)),]

#Create 10 equally size folds
folds <- cut(seq(1,nrow(yourData)),breaks=10,labels=FALSE)

#Perform 10 fold cross validation
for(i in 1:10){
    #Segement your data by fold using the which() function 
    testIndexes <- which(folds==i,arr.ind=TRUE)
    testData <- yourData[testIndexes, ]
    trainData <- yourData[-testIndexes, ]
    #Use test and train data partitions however you desire...
}

- Jake Drew

---

感谢Jake Drew。为了测试目的，我已经将你上面的代码与cv.glm的结果进行了比较，结果是相同的。感谢你的帖子，现在我可以信任cv.glm ；-) - citraL

Answer 2

@Roman在他的评论中提供了一些答案，然而，通过使用cv.glm检查代码可以回答你的问题：

我相信这段代码将数据集随机分成K折，必要时进行四舍五入以使K不能整除n：

if ((K > n) || (K <= 1)) 
    stop("'K' outside allowable range")
K.o <- K
K <- round(K)
kvals <- unique(round(n/(1L:floor(n/2))))
temp <- abs(kvals - K)
if (!any(temp == 0)) 
    K <- kvals[temp == min(temp)][1L]
if (K != K.o) 
    warning(gettextf("'K' has been set to %f", K), domain = NA)
f <- ceiling(n/K)
s <- sample0(rep(1L:K, f), n)

这里的内容表明，delta值不是均方根误差。就像帮助文件所说的那样：默认为平均平方误差函数。 这意味着什么？我们可以通过检查函数声明来了解：

function (data, glmfit, cost = function(y, yhat) mean((y - yhat)^2), 
    K = n) 

这表明，在每个折叠中，我们计算误差平方的平均值，其中误差是指预测响应与实际响应之间的常规差异。

delta[1]只是每个折叠中所有这些项的总和的加权平均值，请参见cv.glm代码中的内联注释：

for (i in seq_len(ms)) {
    j.out <- seq_len(n)[(s == i)]
    j.in <- seq_len(n)[(s != i)]
    Call$data <- data[j.in, , drop = FALSE]
    d.glm <- eval.parent(Call)
    p.alpha <- n.s[i]/n #create weighted average for later
    cost.i <- cost(glm.y[j.out], predict(d.glm, data[j.out, 
        , drop = FALSE], type = "response"))
    CV <- CV + p.alpha * cost.i # add weighted average error to running total
    cost.0 <- cost.0 - p.alpha * cost(glm.y, predict(d.glm, 
        data, type = "response"))
}

- Alex