r - 估计命令如何查找 R 公式中的变量名称？

Question

我想在用户定义的函数上使用 R 的 nls() 函数来估计大量模型。由于许多变量在我的规范中是固定的，我想要一种在我的函数中预先设置它们的方法，但我没有正确理解 R 如何在公式中包含的函数中查找变量。

我看过 Hadley Wickham 的高级 R 书中有关元编程的部分，但它并没有启发我。以下是我使用 mtcars 数据集尝试实现的目标的简化示例:

我尝试为跨规范固定的变量设置默认值:

expo <- function(x, theta, weight = wt) {
  x*weight^theta
}

我还尝试过仅使用固定变量的列名作为函数内的变量

expo <- function(x, theta) {
  x*wt^theta
}

如果我只想计算函数，这两种方法都有效，比如

attach(mtcars)
expo(qsec, 1)
detach()

但是，例如，如果我尝试在估计例程中使用我的 expo() 函数

nls(mpg ~ phi + expo(qsec, theta),
    data = mtcars,
    start = c('phi' = -2, 'theta' = 1))

它失败并显示消息expo(qsec, theta) 中的错误:未找到对象“wt”。评论中提出的一种可能性是简单地将数据集(本例中为 mtcars)传递给 expo() 作为参数。但由于我只会在调用 nls() 时调用 expo() ，其中数据集已经是一个参数，所以如果我能找到一种方法来避免这种重复。

在适本地定义或调用 expo() 后，我的最终目标是能够执行以下操作:

depvars <- c('qsec', 'drat', 'dist')
lapply <- (depvars, function(x) {
    formula <- as.formula(paste0('mpg ~ phi + expo(', x, ', theta)'))
    nls(formula,
        data = mtcars,
        start = c('phi' = -2, 'theta' = 1))
}

Answer 1

棘手的是 R 的词法作用域搜索是在封闭环境中进行的， 这在调用过程中可能会令人困惑，因为调用者环境可能都有封闭的环境，并且事情很快就会变得困惑。

我将使用 rlang 包来调试此场景。

首先，如果您在全局环境中定义了 expo， 那么这将是它的封闭环境:

expo <- function(x, theta) {
  x*wt^theta
}

rlang::get_env(expo)
# <environment: R_GlobalEnv>

所以当你调用它时，R会首先在函数调用中搜索变量 (不是来电者!) 环境， 然后在封闭环境(此处为全局环境)中。

我不知道nls到底是做什么的， 但我假设它会根据您提供的数据创建一个环境并评估那里的公式。 然而，它创建的环境似乎只包含它可以在公式中明确看到的变量， 我发现的东西:

expo <- function(x, theta) {
  cat("caller: ")
  print(ls(rlang::caller_env()))
  cat("enclosing: ")
  print(ls(rlang::env_parent(rlang::current_env())))
}

nls(mpg ~ phi + expo(qsec, theta),
    data = mtcars,
    start = c('phi' = -2, 'theta' = 1))
# caller: [1] "mpg"   "phi"   "qsec"  "theta"
# enclosing: [1] "expo"    
# Error ...

正如我们所见，expo 的调用者环境包含我们可以在公式中识别的变量， 其封闭环境仅包含expo的定义 (全局环境)。 不幸的是，这意味着您甚至无法在 expo 中使用类似 eval.parent 的内容， 因为该环境不会包含 data 中的所有变量。

如果您仍然想解决这个问题， 您可以在调用 nls 之前使用您的数据修改 expo 的封闭环境， 像这样:

expo <- function(x, theta) {
  x*wt^theta
}

environment(expo) <- list2env(as.list(mtcars))

nls(mpg ~ phi + expo(qsec, theta),
    data = mtcars,
    start = c('phi' = -2, 'theta' = 1))
# Error ... number of iterations exceeded maximum of 50