【问题标题】:merTools predictInterval() for model with nested random effectmerTools predictInterval() 用于具有嵌套随机效应的模型
【发布时间】:2016-02-02 22:58:36
【问题描述】:

merTools 包中的predictInterval() 不喜欢嵌套随机效果吗?例如,使用ggplot2 包中的msleep 数据集:

library("lme4")
library("merTools")
library("ggplot2")
mod <- lmer(sleep_total ~ bodywt + (1|vore/order), data=msleep)
predInt <- predictInterval(merMod=mod, newdata=msleep) 

返回错误:

Error in '[.data.frame'(newdata, , j) : undefined columns selected

这运行正常没问题:

mod <- lmer(sleep_total ~ bodywt + (1|vore) + (1|order), data=msleep)
predInt <- predictInterval(merMod=mod, newdata=msleep)

(实际上它给出了关于随机效应变量中NA 级别的警告,但我并不担心)

更新

正如下面 Ben Bolker 回答的 cmets 中所讨论的,merTools 的新版本考虑了嵌套随机效应。但是,当我尝试预测包含新级别的嵌套随机效应的数据时,会出现错误。

这行得通:

mod <- lmer(sleep_total ~ bodywt + (1|vore/order), data=msleep)
predInt <- predictInterval(merMod=mod, newdata=msleep) 

这很有效,尽管有几个警告(有关警告的其他问题,请参见下文*):

mod <- lmer(sleep_total ~ bodywt + (1|vore) + (1|order), data=msleep)
msleep2 <- msleep %>% mutate(vore = "omni")
predInt <- predictInterval(merMod=mod, newdata=msleep2) 

但这不起作用:

mod <- lmer(sleep_total ~ bodywt + (1|vore/order), data=msleep)
msleep2 <- msleep %>% mutate(vore = "omni")
predInt <- predictInterval(merMod=mod, newdata=msleep2) 

出现以下错误:

Error in `[.data.frame`(tmp, alllvl) : undefined columns selected
In addition: Warning message:
In predictInterval(merMod = mod, newdata = msleep3) :
  newdata is tbl_df or tbl object from dplyr package and has been
              coerced to a data.frame

在这里,"omni" 实际上并不是 vore 的新级别,但是当与 order 结合时,它会创建新的变量嵌套组合。

如果我使用"new" 或其他任何不是vore 观察级别的东西,我会得到类似的结果:它适用于模型的非嵌套版本,但不适用于嵌套版本。


*另外,我是否应该关注上面第二个模型块给出的警告:

> mod <- lmer(sleep_total ~ bodywt + (1|vore) + (1|order), data=msleep)
> msleep2 <- msleep %>% mutate(vore = "omni")
> predInt <- predictInterval(merMod=mod, newdata=msleep2)
Warning messages:
  1: In predictInterval(merMod = mod, newdata = msleep2) :
     newdata is tbl_df or tbl object from dplyr package and has been
       coerced to a data.frame
  2: In chol.default(sigma, pivot = TRUE) :
     the matrix is either rank-deficient or indefinite

我猜第二个是vore 对每个观察值取相同值的结果,但这不应该是预测的问题,不是吗?如果在我拟合模型时变量取相同的值,我可以看到这是一个问题,但在预测新观察时不认为这应该是一个问题?

【问题讨论】:

    标签: r lme4


    【解决方案1】:

    可以(显然)通过明确写出交互术语来解决此问题。 警告:我实际上并没有检查以确保结果预测是正确的,只是看到没有产生错误并且结果对象是近似合理的......

    msleep <- transform(msleep,voreOrder=interaction(vore,order,drop=TRUE))
    mod2 <- lmer(sleep_total ~ bodywt + (1|vore)+(1|voreOrder), data=msleep)
    predInt <- predictInterval(merMod=mod2, newdata=msleep) 
    

    这确实会生成警告消息,但显然它们是由于vore 变量中的&lt;NA&gt; 值(我不知道这个数据集...)

    【讨论】:

    • 这是目前处理此问题的正确方法,我已经检查了预测是否按预期工作。现在我需要处理 merTools 补丁以允许 merTools 解析嵌套规范。
    • @jknowles,你可能会发现一些机器here很有用...
    • 作为一个更新,我认为 merTools 的开发分支非常接近以一种不优雅的方式解决这个问题。如果你有一个可以工作的迫切需要。我正在努力尽快发布新的 CRAN 版本,并希望能解决这个问题。
    • @jknowles 我已经使用 install.packages("merTools", type="source") 安装了 0.2.0 版。然后我在原始问题中运行相同的步骤,但出现新错误: > predInt
    • 是的,我相信这是由于lme4 代码的更改。我刚刚向 CRAN 提交了一个补丁,该补丁需要更新版本的 lme4 以避免此错误。如果您将lme4 安装更新到最新版本,它应该可以工作,或者如果您更新到merTools 的开发版本,它将强制您现在更新lme4
    猜你喜欢
    • 1970-01-01
    • 2021-04-04
    • 2018-08-26
    • 2014-07-24
    • 1970-01-01
    • 2019-02-09
    • 2013-05-03
    • 1970-01-01
    • 1970-01-01
    相关资源
    最近更新 更多