函数内的update()只会搜索全局环境吗？

Question

函数内的update()只会搜索全局环境吗？

9

我试图编写一个包装函数来批量进行似然比检验。我尝试使用update()函数来更新初始模型。但是，它似乎不是在函数内查找对象，而是在全局环境中搜索对象。

fake <- data.frame(subj= rep(1:5, 4), 
                   factor1 = rep(LETTERS[c(1,2,1,2)], each=5), 
                   factor2 = rep(letters[1:2], each=10), 
                   data=sort(rlnorm(20)))

foo <- function(){
                  temp <- fake
                  model1 <- lmer(data~factor1*factor2 + (1 |subj), temp)
                  model1a <- update(model1, ~.-factor1:factor2)
                  model1a}

并且它会显示以下错误信息：

Error in eval(expr, envir, enclos) : object 'factor1' not found

有没有办法让update（）在函数内进行搜索？谢谢！

编辑：

我犯了一个错误。我想将“temp”传递给lmer，而不是“fake”。

编辑2：一种方便的解决方案是简单地指定数据对象。虽然update（）现在没有问题，但anova（）似乎认为我尝试比较的模型基于不同的数据对象。

 foo <- function(){
                  temp <- fake
                  model1 <- lmer(data~factor1*factor2 + (1 |subj), data=temp)
                  model1a <- update(model1, ~.-factor1:factor2, data=temp)
                  anova(model1, model1a)
            }
 foo()

我收到一个错误信息：

 Error in anova(model1, model1b) : 
   all models must be fit to the same data object

我想这个错误不仅限于update()函数。但我想知道有没有人知道如何解决这个问题。请注意，如果我不使用update()函数而是将模型拼写出来（见下文），上述错误就会消失：

 foo <- function(){
                  temp <- fake
                  model1 <- lmer(data~factor1*factor2 + (1 |subj), data=temp)
                  model1a <- lmer(data~factor1 + factor2 + (1 |subj), data=temp)
                  anova(model1, model1a)
            }
 foo()

 Data: temp
 Models:
 model1a: data ~ factor1 + factor2 + (1 | subj)
 model1: data ~ factor1 * factor2 + (1 | subj)
         Df     AIC    BIC  logLik  Chisq Chi Df Pr(>Chisq)  
 model1a  5 -4.6909 3.7535  7.3454                           
 model1   6 -8.8005 1.3327 10.4003 6.1097      1    0.01344 *
 ---
 Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

编辑3：看起来问题出在anova()函数上。我也尝试了@hadley的建议。

foo2 <- function(){
  my_update <- function(mod, formula = NULL, data = NULL) {
  call <- getCall(mod)
  if (is.null(call)) {
    stop("Model object does not support updating (no call)", call. = FALSE)
  }
  term <- terms(mod)
  if (is.null(term)) {
    stop("Model object does not support updating (no terms)", call. = FALSE)
  }
  if (!is.null(data)) call$data <- data
  if (!is.null(formula)) call$formula <- update.formula(call$formula, formula)
  env <- attr(term, ".Environment")
  eval(call, env, parent.frame())}

      model1 <- lmer(data~factor1*factor2 + (1 |subj), temp)
      model1a <- my_update(model1, ~.-factor1:factor2)
      anova(model1, model1a)
 }
 foo2()

我收到了以下所示的错误提示信息：

 Error in as.data.frame.default(data) : 
   cannot coerce class 'structure("mer", package = "lme4")' into a data.frame

- Alex

除了需要加载lme4包外，在R 2.15.1中对我来说没有错误。 - MattBagg

您是想将 temp 传递给 lmer 而不是 fake 吗？ - BenBarnes

抱歉，是的，我想将temp传递给lmer，而不是fake。 - Alex

就我个人而言，使用CRAN（稳定版本）的lme4出现错误（显示找不到“temp”对象），但使用开发版本（在GitHub上：https://github.com/lme4/lme4/）则可以成功。在开发版本中存在一些GLMM稳定性问题，但就功能和健壮性而言，据我所知，开发版本优于稳定版本... - Ben Bolker

检查编辑 -- 就此而言，anova() 示例似乎也适用于开发的lme4。 - Ben Bolker

2个回答

4

虽然我非常喜欢@Hadley的回答（我自己也可能会使用该函数），但您还可以在update函数中指定data参数。（在这里，我假设您想将temp传递给您的模型。）

model1a <- update(model1, ~.-factor1:factor2, data = temp)

编辑

如果你想用anova比较模型，使用update将会混淆data参数的名称，并且“欺骗”anova，使其认为这两个模型是使用不同数据集拟合的。只更新公式并创建一个新模型将避免这种情况:

foo <- function(){
                  temp <- fake
                  model1 <- lmer(data~factor1*factor2 + (1 |subj), data=temp)
                  newForm <- update.formula(formula(model1), ~.-factor1:factor2)
                  model1a <- lmer(newForm, data=temp)
                  anova(model1, model1a)
            }

- BenBarnes

但是...我认为这实际上在稳定的lme4中不起作用（我同意它应该）...? - Ben Bolker

@BenBolker，有趣的问题。在发布之前，我已经在OSX上使用R 2.15.2成功测试了答案，使用的是（我非常确定）lme4的常规CRAN Mac二进制版本。它还可以在WinXP上使用lme4_0.999999-0和Matrix_1.0-9与R 2.15.2一起工作。 - BenBarnes

抱歉，我误解了你在做什么。 - Ben Bolker

看起来anova()将使用update()更新的对象视为基于不同数据集（请参见我在原始帖子中添加的示例）。您是否知道解决此问题的方法？谢谢！ - Alex

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- hadley · Accepted Answer

我以前也被这种行为所困扰，所以我写了自己的update版本。它会评估公式环境中的所有内容，因此应该相当健壮。

my_update <- function(mod, formula = NULL, data = NULL) {
  call <- getCall(mod)
  if (is.null(call)) {
    stop("Model object does not support updating (no call)", call. = FALSE)
  }
  term <- terms(mod)
  if (is.null(term)) {
    stop("Model object does not support updating (no terms)", call. = FALSE)
  }

  if (!is.null(data)) call$data <- data
  if (!is.null(formula)) call$formula <- update.formula(call$formula, formula)
  env <- attr(term, ".Environment")

  eval(call, env, parent.frame())
}

library(nlme4)

fake <- data.frame(
  subj = rep(1:5, 4), 
  factor1 = rep(LETTERS[c(1,2,1,2)], each = 5), 
  factor2 = rep(letters[1:2], each = 10), 
  data = sort(rlnorm(20)))

foo <- function() {
  temp <- fake
  model1 <- lmer(data ~ factor1 * factor2 + (1 | subj), fake)
  model1a <- my_update(model1, ~ . - factor1:factor2)
  model1a
}
foo()