从glmnet中提取系数变量名称并放入数据框中

40

我希望从glmnet生成的模型系数中提取数据,并从中创建一个SQL查询。函数 coef(cv.glmnet.fit) 产生了一个 'dgCMatrix' 对象。当我使用as.matrix 将其转换为矩阵时,变量名丢失,只剩下系数值。

我知道可以将系数打印在屏幕上,但是否可能将名称写入数据框?

有人能帮助提取这些名称吗?


你需要发布一个可重现的例子。发布一些数据,制作一个例子告诉我们问题所在以及您希望它成为什么样子。使用glmnet时,我总是默认获取变量名称,因此我不理解这个问题。 - LyzandeR
9个回答

35

更新:我的回答的前两个评论都是正确的。我将以下回答保留下来,仅供后人参考。

以下回答简洁有效,不需要任何其他包:

tmp_coeffs <- coef(cv.glmnet.fit, s = "lambda.min")
data.frame(name = tmp_coeffs@Dimnames[[1]][tmp_coeffs@i + 1], coefficient = tmp_coeffs@x)

加1的原因是@i方法从0开始为截距做索引,但@Dimnames[[1]]从1开始。


OLD ANSWER: (only kept for posterity) 试试这些句子:

非零系数:

coef(cv.glmnet.fit, s = "lambda.min")[which(coef(cv.glmnet.fit, s = "lambda.min") != 0)]

所选的功能:

colnames(regression_data)[which(coef(cv.glmnet.fit, s = "lambda.min") != 0)]

将它们作为数据框放在一起很简单,但如果您也需要代码的这部分,请告诉我。


4
请注意,colnames(regression_data)[which(coef(cv.glmnet.fit, s = "lambda.min") != 0)] 不考虑截距(第一列),因此显示的名称是错误的。 - RUser4512
2
@x 对象方法将为您提供非零系数。 - Davor Josipovic
感谢您的输入。我现在提供了更好的解决方案。 - Mehrad Mahmoudian
@jruf003 我明白你的意思。那我把旧答案和更新后的答案位置交换一下怎么样? - Mehrad Mahmoudian
1
@joel.wilson 也许你还没有使用 glmnet::cv.glmnet() 来拟合模型。是这种情况吗? - Mehrad Mahmoudian
显示剩余7条评论

9

请查看broom包。它有一个tidy函数,将不同的R对象(包括glmnet)的输出转换为数据框。


这是最简洁的答案。 - Paul Bendevis

6

这些名称可以通过dimnames(coef(cv.glmnet.fit))[[1]]访问,因此以下代码可以将系数名称和值放入data.frame中:

data.frame(coef.name = dimnames(coef(GLMNET))[[1]], coef.value = matrix(coef(GLMNET)))

4

在Mehrad的解决方案基础上,这里提供了一个简单的函数来打印仅包含非零系数的表格:

print_glmnet_coefs <- function(cvfit, s="lambda.min") {
    ind <- which(coef(cvfit, s=s) != 0)
    df <- data.frame(
        feature=rownames(coef(cvfit, s=s))[ind],
        coeficient=coef(cvfit, s=s)[ind]
    )
    kable(df)
}

上述函数使用knitr中的kable()函数生成适用于Markdown的表格。

1
s = 'lambda.min 给我一个错误。Error in lambda[1] - s : non-numeric argument to binary operator 有什么想法吗? - joel.wilson

2

这里,我写了一个可重现的例子,并使用 cv.glmnet 拟合了一个二元(逻辑)模型。也可以使用 glmnet 模型拟合。在这个例子的结尾处,我将非零系数和相关特征组装成了一个名为 myResults 的数据框:

library(glmnet)
X <- matrix(rnorm(100*10), 100, 10);
X[51:100, ] <- X[51:100, ] + 0.5; #artificially introduce difference in control cases
rownames(X) <- paste0("observation", 1:nrow(X));
colnames(X) <- paste0("feature",     1:ncol(X));

y <- factor( c(rep(1,50), rep(0,50)) ); #binary outcome class label
y
## [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [51] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## Levels: 0 1

## Perform logistic model fit:
fit1 <- cv.glmnet(X, y, family="binomial", nfolds=5, type.measure="auc"); #with K-fold cross validation
# fit1 <- glmnet(X, y, family="binomial") #without cross validation also works

## Adapted from @Mehrad Mahmoudian:
myCoefs <- coef(fit1, s="lambda.min");
myCoefs[which(myCoefs != 0 ) ]               #coefficients: intercept included
## [1]  1.4945869 -0.6907010 -0.7578129 -1.1451275 -0.7494350 -0.3418030 -0.8012926 -0.6597648 -0.5555719
## [10] -1.1269725 -0.4375461
myCoefs@Dimnames[[1]][which(myCoefs != 0 ) ] #feature names: intercept included
## [1] "(Intercept)" "feature1"    "feature2"    "feature3"    "feature4"    "feature5"    "feature6"   
## [8] "feature7"    "feature8"    "feature9"    "feature10"  

## Asseble into a data.frame
myResults <- data.frame(
  features = myCoefs@Dimnames[[1]][ which(myCoefs != 0 ) ], #intercept included
  coefs    = myCoefs              [ which(myCoefs != 0 ) ]  #intercept included
)
myResults
##       features      coefs
## 1  (Intercept)  1.4945869
## 2     feature1 -0.6907010
## 3     feature2 -0.7578129
## 4     feature3 -1.1451275
## 5     feature4 -0.7494350
## 6     feature5 -0.3418030
## 7     feature6 -0.8012926
## 8     feature7 -0.6597648
## 9     feature8 -0.5555719
## 10    feature9 -1.1269725
## 11   feature10 -0.4375461

s = 'lambda.min 给我一个错误。Error in lambda[1] - s : non-numeric argument to binary operator 有什么想法吗? - joel.wilson
你使用了函数glmnet还是cv.glmnet?它们的结果数据结构不同。 - David C.
"h(simpleError(msg, call)) 函数出错: 在选择 'which' 函数的方法时,评估参数 'x' 时发生错误:无法将 'list' 对象强制转换为 'double' 类型。" - Emmanuel Goldstein

2

使用coef()方法可以对glmnet()对象(即你的模型)进行处理。在下面的例子中,索引[[1]]表示多项逻辑回归中结果类别的数量,对于其他模型,可能需要将其删除。

coef_names_GLMnet <- coef(GLMnet, s = 0)[[1]]
row.names(coef_names_GLMnet)[coef_names_GLMnet@i+1]

row.names()在这种情况下需要增加(+1),因为coef()对象中变量(数据特征)的编号从0开始,但在转换为字符向量后,编号从1开始。


1

如果您知道如何获得lambda,我找到了两种不同的方法来展示特定lambda所需的选定模型中的预测因子。其中之一包括截距项。可以通过“glmnet”库中cv.glmnet的平均交叉验证获得lambda。您可能只想查看每种方法的最后几行:

 myFittedLasso = glmnet(x=myXmatrix, y=myYresponse, family="binomial")
 myCrossValidated = cv.glmnet(x=myXmatrix, y=myYresponse, family="binomial")
 myLambda = myCrossValidated$lambda.1se  # can be simply lambda

 # Method 1 without the intercept
 myBetas = myFittedLasso$beta[, which(myFittedLasso$lambda == myLambda)]
 myBetas[myBetas != 0]
 ## myPredictor1    myPredictor2    myPredictor3
 ##   0.24289802      0.07561533      0.18299284


 # Method 2 with the intercept
 myCoefficients = coef(myFittedLasso, s=myLambda)
 dimnames(myCoefficients)[[1]][which(myCoefficients != 0)]
 ## [1] "(Intercept)"    "myPredictor1"    "M_myPredictor2"    "myPredictor3"

 myCoefficients[which(myCoefficients != 0)]
 ## [1] -4.07805560  0.24289802  0.07561533  0.18299284

请注意,上面的例子暗示了二项分布,但这些步骤可以应用于任何其他类型。

1
# requires tibble.
tidy_coef <- function(x){
    coef(x) %>%
    matrix %>%   # Coerce from sparse matrix to regular matrix.
    data.frame %>%  # Then dataframes.
    rownames_to_column %>%  # Add rownames as explicit variables.
    setNames(c("term","estimate"))
}

没有tibble:

tidy_coef2 <- function(x){
    x <- coef(x)
    data.frame(term=rownames(x),
               estimate=matrix(x)[,1],
               stringsAsFactors = FALSE)
}

0
当我使用tidymodels框架中的glmnet时,遇到了类似的问题,即模型在工作流程中进行训练,但是coef()和上述解决方案都无法解决。然而,对我有用的是glmnet:::coef.glmnet代码的一部分:
# taken from glmnet:::coef.glmnet
coefs <- predict(x, "lambda.min", type = "coefficients", exact = FALSE)

dd <- cbind(
  data.frame(var = rownames(coefs)),
  as.data.table(as.matrix(coefs))
)

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接