cv.glmnet 与 glmnet 结果；衡量解释力答案

【问题标题】：cv.glmnet vs glmnet results; gauging explanatory powercv.glmnet 与 glmnet 结果；衡量解释力
【发布时间】：2018-05-14 22:11:22
【问题描述】：

当通过 glmnet 包估计套索模型时，我想知道是否更好：(a) 直接从 cv.glmnet 采购的 cv.fit 对象中提取系数/预测/偏差，或者 (b) 使用最小 lambda 从 cv.glmnet 重新运行 glmnet 并从 glmnet 进程中提取这些对象。（请耐心等待——我觉得这已被记录在案，但我在网上看到了这两种方法的示例/教程，并且没有可靠的逻辑可以选择哪种方式。）

也就是说，对于系数，我可以运行（a）：

cvfit = cv.glmnet(x=xtrain, y=ytrain, alpha=1, type.measure = "mse", nfolds = 20)
coef.cv <- coef(cvfit, s = "lambda.min")

或者我之后可以运行 (b)：

fit = glmnet(x=xtrain, y=ytrain, alpha=1, lambda=cvfit$lambda.min)
coef <- coef(fit, s = "lambda.min")

虽然这两个过程选择相同的模型变量，但它们不会产生相同的系数。同样，我可以通过以下两个过程之一进行预测：

prdct <- predict(fit,newx=xtest)
prdct.cv <- predict(cvfit, newx=xtest, s = "lambda.min")

他们预测相似但不相同的向量。

最后，我想我可以通过以下两种方法中的任何一种来解释 % 偏差：

percdev <- fit$dev.ratio
percdev.cv <- cvfit$glmnet.fit$dev.ratio[cvfit$cvm==mse.min.cereal]

但实际上这样拉percdev.cv是不可能的，因为如果cv.glmnet使用的lambda序列少于100个元素，cvfit$glmnet.fit$dev.ratio和cvfit$cvm==mse.min.cereal的长度不匹配.所以我不太确定如何从cvfit$glmnet.fit 中提取最小 lambda dev.ratio。

所以我想我想知道哪个过程最好，为什么，以及人们通常如何提取适当的 dev.ratio 统计数据。谢谢！

【问题讨论】：

来自?glmnet，lambda 的文档说：警告：小心使用。避免为 lambda 提供单个值（对于 CV 之后的预测，请改用 predict()）。而是提供递减的 lambda 值序列。 glmnet 依靠它的暖启动来提高速度，而且它通常更快地适应整个路径而不是计算单个适应。。所以我想这部分回答了你关于使用哪个的问题。
另见：stackoverflow.com/questions/29311323/…
谢谢。是的，我知道通常glmnet 应该与（默认或提供的）lambda 序列一起使用，但是一旦将这样的序列提供给cv.glmnet 和“最佳”lambda（lambda.1se 或lambda.min ) 获得，人们会认为使用该 lambda 会产生与在cv.glmnet 下相同的glmnet 结果。即使计算速度较慢。另外，我有一种预感，从cv.glmnet 获取系数和预测可能会更好，但我不确定如何从cv.glmnet 获取dev.ratio。

标签： r glmnet

【解决方案1】：

正如 cmets 中所指出的，它与提供的 lambda 序列有关，因为如果您查看 cv.glmnet 的源代码，它会在前几行调用 glmnet:::cv.glmnet.raw，运行 glmnet() on定义的 lambda 值。

所以我们可以使用下面的例子：

library(mlbench)
data(BostonHousing)
data = BostonHousing
data$chas=as.numeric(data$chas)

cvfit = cv.glmnet(x=as.matrix(data[,-14]),y=data[,14])
coef.cv <- coef(cvfit, s = "lambda.min")

fit = glmnet(x=as.matrix(data[,-14]), y=data[,14], alpha=1, lambda=cvfit$lambda.min)
coef <- coef(fit, s = "lambda.min")

head(cbind(coef.cv,coef))
6 x 2 sparse Matrix of class "dgCMatrix"
                       1            1
(Intercept)  31.74123706  31.86654225
crim         -0.09834634  -0.09869320
zn            0.04144161   0.04158829
indus         .            .         
chas          2.68518774   2.68163334
nox         -16.30664523 -16.35459059

它们略有不同，如果您提供 cv.glmnet 中使用的 lambda 序列：

fit = glmnet(x=as.matrix(data[,-14]), y=data[,14], alpha=1, lambda=cvfit$lambda)
coef <- coef(fit, s = cvfit$lambda.min)
head(cbind(coef.cv,coef))

6 x 2 sparse Matrix of class "dgCMatrix"
                       1            1
(Intercept)  31.74123706  31.74123706
crim         -0.09834634  -0.09834634
zn            0.04144161   0.04144161
indus         .            .         
chas          2.68518774   2.68518774
nox         -16.30664523 -16.30664523

它们现在是一样的。 dev.ratio 也会匹配：

fit$dev.ratio[fit$lambda==cvfit$lambda.min]
[1] 0.7401482
cvfit$glmnet.fit$dev.ratio[which.min(cvfit$cvm)]
[1] 0.7401482

【讨论】：