尽管常规代码有效,但在自定义函数中包装 nnet::multinom() + ggeffects::ggemmeans() 失败:'symbol' 类型的对象不可子集化

Wrapping nnet::multinom() + ggeffects::ggemmeans() in a custom function fails though regular code works: object of type 'symbol' is not subsettable

我想用 nnet::multinom() 拟合多项式模型并用 ggeffects::ggemmeans() 获得预测。虽然这样的过程在常规代码中有效,但我无法将其包装在函数中。

例子

数据

library(dplyr)

my_mtcars <- 
  mtcars %>%
  mutate(across(c(vs, carb), as.factor)) %>%
  as_tibble()

拟合和预测的工作方式如下

library(nnet) # 7.3-15
library(emmeans) # 1.5.4
library(ggeffects) # 1.0.2

m <- multinom(carb ~ vs, data = my_mtcars)
ggemmeans(model = m, terms = "vs")

## # Predicted probabilities of carb
## # x = vs

## # Response Level = 1

## x | Predicted |       95% CI
## ----------------------------
## 0 |      0.00 | [0.00, 0.00]
## 1 |      0.50 | [0.43, 0.57]

## # Response Level = 2

## x | Predicted |       95% CI
## ----------------------------
## 0 |      0.28 | [0.24, 0.32]
## 1 |      0.36 | [0.30, 0.42]

## # Response Level = 3

## x | Predicted |       95% CI
## ----------------------------
## 0 |      0.17 | [0.14, 0.19]
## 1 |      0.00 | [0.00, 0.00]

## # Response Level = 4

## x | Predicted |       95% CI
## ----------------------------
## 0 |      0.44 | [0.39, 0.50]
## 1 |      0.14 | [0.12, 0.17]

## # Response Level = 6

## x | Predicted |       95% CI
## ----------------------------
## 0 |      0.06 | [0.05, 0.06]
## 1 |      0.00 | [0.00, 0.00]

## # Response Level = 8

## x | Predicted |       95% CI
## ----------------------------
## 0 |      0.06 | [0.05, 0.06]
## 1 |      0.00 | [0.00, 0.00]

但是当我尝试将此过程包装在自定义函数中时,它失败了

my_multinom <- function(dat, dv, expl) {
  
  frmla <- as.formula(paste0(dv, "~", expl))
  
  model_fit <- nnet::multinom(frmla, data = dat)
  ggemmeans(model = model_fit, terms = expl)
}

my_multinom(dat = my_mtcars, dv = "carb", expl = "vs")

Error in object$call$formula[[2]] :
object of type 'symbol' is not subsettable

值得注意的是,问题似乎出在 multinom()ggemmeans() 之间的交互上。如果我们从 my_multinom() 中省略 ggemmeans() 那么它似乎工作正常:

my_multinom_no_ggemmeans <- function(dat, dv, expl) {
  
  frmla <- as.formula(paste0(dv, "~", expl))
  model_fit <- nnet::multinom(frmla, data = dat)
  model_fit
}

my_multinom_no_ggemmeans(dat = my_mtcars, dv = "carb", expl = "vs")

## # weights:  18 (10 variable)
## initial  value 57.336303 
## iter  10 value 38.192450
## iter  20 value 37.940409
## final  value 37.940164 
## converged
## Call:
## nnet::multinom(formula = frmla, data = dat)

## Coefficients:
##   (Intercept)       vs1
## 2    13.44961 -13.78607
## 3    12.93879 -33.99280
## 4    13.91961 -15.17237
## 6    11.84015 -23.96194
## 8    11.84015 -23.96194

## Residual Deviance: 75.88033 
## AIC: 95.88033 

知道为什么 my_multinom() 包装器失败了吗?


更新


我可能找到了解决方案,但我不明白为什么它有效。基于this github issue(不同的包),我采用了以下解决方案:

my_multinom_with_do.call <- function(dat, dv, expl) {

  frmla <- as.formula(paste0(dv, "~", expl))

  model_fit <- do.call(multinom, args = list(formula = frmla, data = dat))
  ggemmeans(model = model_fit, terms = expl)
}

有效:

my_multinom_with_do.call(dat = my_mtcars, dv = "carb", expl = "vs")

但为什么这行得通,而我原来的 my_multinom() 却不行?

因为惰性求值,它不起作用。 model_fitcall 成员有 formula = frmla,未计算。 emmeans 对该模型的支持需要一个公式。如果您在原始函数中添加一行,它将起作用:

my_multinom <- function(dat, dv, expl) {
    
    frmla <- as.formula(paste0(dv, "~", expl))
    
    model_fit <- nnet::multinom(frmla, data = dat)
    model_fit$call$formula <- frmla
    ggemmeans(model = model_fit, terms = expl)
}

do.call 方法起作用的原因是当您创建传递给 do.call 的列表时,frmla 被评估。