Функция Emmeans — в справочной сетке нет переменной

#r #lme4 #emmeans

Вопрос:

Я пытаюсь запустить функцию emmeans на наборе данных lmer, но она не работает. Вот мои данные:

 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)

код эммануэля:

 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
 

Я не уверен, почему в нем говорится, что Черепахи нет, так как это символьная переменная в моем наборе данных.

В принципе, я просто хочу, чтобы эммеаны запускались, но я также боюсь, что этого не произойдет, потому что полный набор данных составляет 20 000 строк.

Комментарии:

1. Черепаха-это случайный эффект, чего вы ожидали от эмминов?

2. @GeorgeSavva Я хотел бы объяснить тот факт, что Черепаха-это повторная мера. Есть ли лучший способ сделать это?

3. Очевидно, что вы не получили сообщения, показанные с помощью модели, соответствующей показанному набору данных, потому что сообщение подразумевает, что у вас было более 20 000 наблюдений. Пожалуйста, приведите воспроизводимые примеры.

4. @RussLenth показанный набор данных состоял только из первых 6 строк. Я сделал dput(голова(df))

5. Тогда вы должны были сказать это в операции.

Ответ №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 
 

Более подробную информацию о случайном и фиксированном можно найти в разделе Перекрестная проверка.

Комментарии:

1. Похоже, мне, вероятно, следует использовать Черепаху в качестве фиксированного эффекта. Я попробовал это, и R все еще требует вечности для обработки кода, так как есть 20 000 наблюдений. какие-нибудь советы по этому поводу?

2. Какая часть занимает «вечность»? Если это так 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)
 

Доверительные интервалы сложнее …