Emmeans 函数 - 参考网格中没有变量

Posted 2023-02-18

【中文标题】Emmeans 函数 - 参考网格中没有变量

【英文标题】：Emmeans function - no variable in reference grid

【发布时间】：2021-12-02 06:25:14

【问题描述】：

我正在尝试在 lmer 数据集上运行 emmeans 函数，但它不起作用。这是我的数据：

structure(list(Date = structure(c(16578, 16578, 16578, 16578, 
16578, 16578), class = "Date"), Time = c(7, 7, 7, 9, 11, 11), 
    Turtle = c("R3L12", "R3L12", "R3L12", "R3L12", "R3L12", "R3L12"
    ), Tex = c(11.891, 12.008, 12.055, 13.219, 18.727, 18.992
    ), m.Tb = c(12.477, 12.54, 12.54, 12.978, 16.362, 16.612), 
    m.HR = c(7.56457, 6.66759, 17.51107, 9.72277, 19.44553, 13.07674
    ), season = c("beginning", "beginning", "beginning", "beginning", 
    "beginning", "beginning"), year = c(2015L, 2015L, 2015L, 
    2015L, 2015L, 2015L), Mass = c(360L, 360L, 360L, 360L, 360L, 
    360L)), row.names = c(NA, 6L), class = "data.frame") 

型号代码：model1 <- lmer(m.HR ~ season + (1|Time) + (1|Date) + (1|Turtle), turtledata)

em意思是代码：

model1.emmeans <- emmeans(model1, "Turtle")

这些是我得到的错误：

To enable adjustments, add the argument 'pbkrtest.limit = 20608' (or larger)
[or, globally, 'set emm_options(pbkrtest.limit = 20608)' or larger];
but be warned that this may result in large computation time and memory use.
Note: D.f. calculations have been disabled because the number of observations exceeds 3000.
To enable adjustments, add the argument 'lmerTest.limit = 20608' (or larger)
[or, globally, 'set emm_options(lmerTest.limit = 20608)' or larger];
but be warned that this may result in large computation time and memory use.
Error in emmeans(model1, "Turtle") : 
  No variable named Turtle in the reference grid

我不确定为什么它说没有 Turtle，因为它是我数据集中的字符变量。

基本上，我只是想让 emmeans 运行，但我也担心它不会运行，因为完整的数据集有 20,000 行长。

【问题讨论】：

乌龟是随机效应，你希望emmeans做什么？

@GeorgeSavva 我想说明 Turtle 是重复测量的事实。有没有更好的方法来做到这一点？

很明显，您没有收到显示的消息与模型拟合显示的数据集，因为该消息暗示您有超过 20,000 个观察值。请给出可重现的例子。

@RussLenth 显示的数据集只是前 6 行。我做了 dput(head(df))

那你应该在 OP 里这么说。

【参考方案1】：

函数emmeans 测试固定效果（被操纵的东西），而不是随机效果（由于设计而发生的东西）。以下示例显示了这一点，也是创建最小可重现示例的一种方法：

library(emmeans)
library(lme4)

# some artificial data
set.seed(143)
foo <- data.frame(
  m.HR   <- rnorm(100, mean=c(rep(c(5, 6), 25), rep(c(8, 9), 25))),
  season <- rep(c("a", "b"), each=50),
  Turtle <- rep(c("T1", "T2"), 50)
)

# simplified model with one fixed and one raqndom effect
model1 <- lmer(m.HR ~ season +  (1|Turtle), foo)

(model1.emmeans <- emmeans(model1, "Turtle"))
# --> error as Turtle is a random effect

(model1.emmeans <- emmeans(model1, "season"))
# --> works as season is a fixed effect

#season emmean    SE   df lower.CL upper.CL
#a        5.73 0.535 1.07  -0.0567     11.5
#b        8.61 0.535 1.07   2.8254     14.4
#
#Degrees-of-freedom method: kenward-roger 
#Confidence level used: 0.95 

更多关于随机与固定的讨论可以在Cross Validated找到。

【讨论】：

看来我可能应该使用 Turtle 作为固定效果。我试过了，但 R 仍然需要很长时间来处理代码，因为有 20,000 个观察值。有什么建议吗？

哪一部分需要“永远”？如果是lmer，那你应该知道，这确实很辛苦。加速的一种方法是为效果定义良好的起始值，另一种方法是尝试不同的优化器。例如，可以从随机子样本中猜测起始值。

【参考方案2】：

你不一定让emmeans 直接做你想做的事，但某种合理的计算是可能的。

最简单的方法是使用跨季节的平均值对每只海龟进行平均预测：

ref_grid <- with(turtledata, 
   expand.grid(season=unique(season), turtle=unique(Turtle)))
pp <- predict(model1, newdata = ref_grid)
aggregate(pp, by=ref_grid$turtle, FUN=mean)

置信区间更难......