glmnet:对于glmnet()和cv.glmnet()，相同的lambda但是不同的系数

Question

尽管我使用的是相同的lambda，但cv.glmnet()生成的系数似乎与glmnet()生成的系数不同。为什么会这样呢？他们不应该是一样的吗？

library(glmnet)

# Data dimensions
num.samples <- 30
num.genes <- 17000

# Data objects - note that both X and Y are scaled
set.seed(123)
Y <- matrix(rnorm(num.samples), ncol=1)
set.seed(1234)
X <- matrix(rnorm(num.samples*num.genes), ncol=num.genes)

# Run cv.glmnet: get lambda.min and coef
fit.cv <- cv.glmnet(X, Y, nfolds=num.samples, intercept=FALSE)
fit.cv.lambda <- fit.cv$lambda.min
fit.cv.coef <- coef(fit.cv, s = fit.cv.lambda)[,1][2:(num.genes+1)]

# Run glmnet with lambda.min from cv.glmnet: get coef
second.lambda=fit.cv.lambda-0.0001 ## second.lambda included because glmnet manual recommends using >1 lambda for glmnet()
fit <- glmnet(X, Y, lambda=c(fit.cv.lambda,second.lambda), intercept=FALSE) 
fit.lambda <- fit$lambda[1]
fit.coef <- coef(fit, s = fit.cv.lambda)[,1][2:(num.genes+1)]

# Lambda is the same, but coefficients are not
fit.cv.lambda==fit.lambda ## TRUE
not.equal = which(fit.cv.coef != fit.coef)
length(not.equal) ## 18
mean(abs(fit.cv.coef[not.equal] - fit.coef[not.equal])) ## 0.0004038209

(我还注意到，glmnet()和cv.glmnet()的系数在alpha的某个值下并没有不同，但似乎没有明显的模式)

Answer 1

简短答案：这是一个数值精度问题。您所遇到的差异并不是由于cv.glmnet和glmnet之间的差异。更确切地说，这是由以下因素综合造成的：

• 惩罚路径，lambda，在这两个对象之间是不同的，我指的是整个惩罚路径，而不仅仅是利益的惩罚是否同时存在。
• 默认收敛阈值，thresh =

如果希望从具有不同惩罚路径的两个glmnet或cv.glmnet对象获得的估计相等(或至少非常接近)，则在两个函数中使用thresh =选项来降低收敛阈值。此外，我建议在exact = TRUE中设置coef()。

扩展答案：下面我们在几个例子中说明了这一点。在这样做之前，了解coef()函数(它用type = "coefficients"调用predict.glmnet()函数)的逻辑也很重要。

• 如果为对象的原始惩罚路径中已经计算的惩罚请求系数估计，coef()将简单地从原始对象返回估计数。
• 如果为不位于原始路径和exact = FALSE (这是默认的)的惩罚请求系数估计，则根据原始路径中最近的惩罚的估计值进行插值估计系数。
• 如果您请求对未在原始路径和exact = TRUE中的惩罚进行系数估计，则将新的惩罚添加到原始路径中，并对整个模型进行修改以获得估计值。

示例1:相同的惩罚路径，相同的估计数

如果使用默认参数，glmnet()和cv.glmnet()将为给定的数据集计算相同的惩罚路径(但由于停止准则，路径中计算的惩罚数可能有所不同)。我们如下所示：

library(glmnet)

# Data dimensions
num.samples <- 30
num.genes <- 17000

# Data objects - note that both X and Y are scaled
set.seed(123)
Y <- matrix(rnorm(num.samples), ncol=1)
X <- matrix(rnorm(num.samples*num.genes), ncol=num.genes)

# Run cv.glmnet and glmnet, obtain same penalty path up to min(num penalty)
fit <- glmnet(X, Y, intercept=FALSE)
cvfit <- cv.glmnet(X, Y, intercept= FALSE)
min_num_lambdas = min(length(fit$lambda), length(cvfit$lambda))
all.equal(fit$lambda[1:min_num_lambdas], cvfit$lambda[1:min_num_lambdas]) #TRUE

然后，我们可以使用coef()从两个对象获得相同的估计值，而不管新的惩罚是否在原始路径中。

# Requested penalty in original path
coef1 <- coef(fit, s = fit$lambda[10])
coef2 <- coef(cvfit, s = fit$lambda[10])
all.equal(coef1, coef2) #TRUE

# Requested penalty not in original path -- uses interpolation
coef1 <- coef(fit, s = 0.40)
coef2 <- coef(cvfit, s = 0.40)
all.equal(coef1, coef2) #TRUE

# Force glmnet to refit the model with s added to the penalty path
coef1 <- coef(fit, s = 0.40, exact = TRUE, x = X, y = Y)
coef2 <- coef(cvfit, s = 0.40, exact = TRUE, x = X, y = Y)
all.equal(coef1, coef2) #TRUE

示例2:不同的惩罚路径，不同的估计

我们对glmnet() fit做了一个小小的修改，要求在路径上而不是在100次中执行99次惩罚。这使得惩罚路径在fit和cvfit之间存在差异。现在，即使我们使用exact = TRUE选项，估计值也是不同的。

fit <- glmnet(X, Y, intercept=FALSE, nlambda = 99)
coef1 <- coef(fit, s = fit$lambda[10], exact = TRUE, x = X, y = Y)
coef2 <- coef(cvfit, s = fit$lambda[10], exact = TRUE, x = X, y = Y)
all.equal(coef1, coef2) # Mean relative difference: 0.002006215

The Fix

为了确保估计值匹配，我们可以用较低的阈值对两个对象进行重构。

fit <- glmnet(X, Y, intercept=FALSE, nlambda = 99, thresh = 1e-20)
cvfit <- cv.glmnet(X, Y, intercept= FALSE, thresh = 1e-20)
coef1 <- coef(fit, s = fit$lambda[10], exact = TRUE, x = X, y = Y)
coef2 <- coef(cvfit, s = fit$lambda[10], exact = TRUE, x = X, y = Y)
all.equal(coef1, coef2) #TRUE

请注意，您也可以在thresh =中使用coef()选项，但这只有在新的代价尚未出现在两个对象的原始路径时才能起作用。

fit <- glmnet(X, Y, intercept=FALSE, nlambda = 99)
cvfit <- cv.glmnet(X, Y, intercept= FALSE)
coef1 <- coef(fit, s = 0.40, exact = TRUE, x = X, y = Y, thresh = 1e-20)
coef2 <- coef(cvfit, s = 0.40, exact = TRUE, x = X, y = Y, thresh = 1e-20)
all.equal(coef1, coef2) #TRUE

thresh =需要减少的值将完全依赖于数据。