如何在 glm 中指定特定的优化方法答案

【问题标题】：How to specify a particular optimization method in glm如何在 glm 中指定特定的优化方法
【发布时间】：2021-08-15 16:09:21
【问题描述】：

在 glm 的文档中，有一个 method 选项，在文档中它说

用户提供的拟合函数可以作为函数提供或命名函数的字符串，其中函数采用与 glm.fit 相同的论点。如果指定为字符串从 stats 命名空间中查找。

我想使用二项式非负最小二乘优化器，使系数为非负且总和为 1。使用此优化器的一个示例是在带有选项 method = "method.NNLS" 的 SuperLearner 包中.下面是一个可重现的例子：

library(SuperLearner)
# binary outcome
set.seed(1)
N <- 200
X <- matrix(rnorm(N*10), N, 10)
X <- as.data.frame(X)
Y <- rbinom(N, 1, plogis(.2*X[, 1] + .1*X[, 2] - .2*X[, 3] +
                           .1*X[, 3]*X[, 4] - .2*abs(X[, 4])))

SL.library <- c("SL.glmnet", "SL.glm", "SL.knn", "SL.mean")

# least squares loss function
test.NNLS <- SuperLearner(Y = Y, X = X, SL.library = SL.library,
                          verbose = TRUE, method = "method.NNLS", family = binomial())
> test.NNLS

Call:  
SuperLearner(Y = Y, X = X, family = binomial(), SL.library = SL.library, method = "method.NNLS",  
    verbose = TRUE) 


                   Risk      Coef
SL.glmnet_All 0.2460486 0.0000000
SL.glm_All    0.2507033 0.2423697
SL.knn_All    0.2508500 0.3493301
SL.mean_All   0.2475494 0.4083002

请注意，Coefs 是非负数且总和为 1。有没有办法使用与 glm 相同的优化器？还是使用相同约束的不同优化器？我试过了：

glm(Y ~ as.matrix(X), family = "binomial", method = "method.NNLS")

但它不起作用。

【问题讨论】：

glmnet 包中的glmnet 接受边界。
@G.Grothendieck，它可以添加一个约束，使系数总和为 1？
我正在努力解决这个问题，但这并不重要。 (1) nnls 优化器本身仅将参数约束为非负数。 method.NNLS 代码通过简单地将系数除以它们的总和以确保它们总和为 1 来处理 sum-to-1 约束。这与拟合 GLM 时不同 (AFAICT)约束。你能确认你对这种类型的解决方案没问题，或者你真的想用限制为总和为 1 的参数来拟合模型吗？
(2) 您对 SuperLearner 优化器和插入 glm() 的方法的承诺/必要性如何？换句话说，解决问题的替代解决方案（拟合参数限制为非负且总和为 1 的二项式 GLM）可以解决您的问题吗？
(1) 如果可能，我想拟合模型，其参数总和为 1 (2) 我不致力于 SuperLearner 优化器。替代解决方案将解决我的问题。

标签： r optimization glm

【解决方案1】：

一种方法（不是唯一的方法，可能不是最有效或最紧凑的方法）是对变换空间中的参数进行一般最大似然估计；可能最常见的这种转换是加法对数比率转换（例如compositions::alr()）。该解决方案使用bbmle 包，它是optim() 的包装器；在这种情况下，它实际上并没有比直接使用 optim() 提供那么多优势

library(bbmle)
library(compositions)
linkfun <- function(mu) as.numeric(compositions::alr(mu))
invlink <- function(eta) as.numeric(compositions::alrInv(eta))

## negative log-likelihood function: first convert the parameter vector
## from the constrained (alr) space to the unconstrained space, then
## convert the linear predictor from an unconstrained space to the
## constrained (probability, (0,1)) space
lfun <- function(beta) {
    -sum(dbinom(Y, size=1, 
                prob = plogis(as.matrix(X) %*% invlink(beta)),
                log=TRUE))
}
npar <- ncol(X)

使用mle2：

## extra setup for using mle2 with a score function that takes a vector
## parameter
pnm <- paste0("beta", seq(npar))
parnames(lfun) <- pnm
    
m1 <- mle2(lfun, start = setNames(rep(0, npar-1), pnm[-1]), vecpar = TRUE)
invlink(coef(m1))
##  [1] 4.897617e-01 2.221793e-01 1.274830e-05 2.198596e-05 7.541545e-02
##  [6] 1.175966e-01 1.590407e-05 4.662563e-02 4.825376e-02 1.170077e-04

或

opt <- optim(par = rep(0, npar-1), fn = lfun, method = "BFGS")
invlink(opt$par)

注 1

我的模型拟合与您用于生成数据的方法不匹配（这将在公式中编码为~ 0 + x1 + x2 + x3 + x3:x4 + abs(x4)），但您的SuperLearner 示例也不匹配，所以我并不担心。

注2

我的猜测（没有仔细考虑）是上面的方法不会给出与使用 SuperLearner 相同的结果；查看代码，SuperLearner 所做的是找到（加权）最小二乘问题的非负解，然后通过将系数向量除以其总和将系数转换为合成空间。这与解决最小偏差/最大似然问题（参数限制为非负且总和为 1）不同。如果您想要 SuperLearner 解决方案，它可能确实是可以通过将适当的method= 参数传递给glm 来做到这一点，但我看不到如何使用glm 进行受限MLE 解决方案。

【讨论】：