Как сообщить о случайном эффекте в модели смешанных эффектов

Я использую набор данных iris в R в качестве примера ниже.

library(lme4)
mixed.fit.b <- lmer(Sepal.Width ~ Sepal.Length + (1+ Sepal.Length|Species), data = iris)
summary(mixed.fit.b)

ranef(mixed.fit.b)$Species
coef(mixed.fit.b)$Species
predict(mixed.fit.b)

Случайный перехват и наклон показаны ниже.

           (Intercept) Sepal.Length
setosa     -0.49549054   0.78331501
versicolor  1.19197858   0.26689317
virginica   1.17260303   0.27282273 

Это означает, что точка пересечения равна -0,49549054 (фиксированная + случайная точка пересечения), а наклон равен 0,78331501 (фиксированный + случайный наклон) для setosa, верно? Итак, есть три пары перехватов и наклонов. В общей линейной модели мы можем сказать, что y = точка пересечения + наклон, а y изменил наклон на x. А вот в смешанных моделях есть три три пары перехватов и наклонов. Как сообщить о них?

"Как сообщить о них?" Обычно нет. Если бы вы хотели сообщить о них, я бы спросил, почему вы считаете Species случайным эффектом. Обычно я сообщал о фиксированных эффектах и ​​дисперсиях случайных эффектов. (Случайный эффект только с тремя субъектами очень сомнителен. Этого недостаточно для надежной оценки дисперсии. Здесь следует использовать модель с фиксированными эффектами.)

Roland 26.01.2023 08:55

Я согласен с вашей точкой зрения о кластерах случайных эффектов, но я должен с уважением не согласиться с отчетами о случайных эффектах. Хотя это не является общепринятой практикой, это действительно должно быть так, и ниже я подробно расскажу, почему.

Shawn Hemelstrand 26.01.2023 12:16

@ShawnHemelstrand Я сказал, что следует сообщать об отклонениях (и да, также о корреляциях). Вопрос ОП более конкретен, чем вопрос о том, как сообщать о случайных эффектах. Я предпочитаю отображать прогнозы на предметном уровне вместо коэффициентов.

Roland 27.01.2023 07:11

Ах, я думаю, что я неправильно понял ваше намерение. Я изменил свой ответ.

Shawn Hemelstrand 27.01.2023 07:20

Уважаемый @Roland, извините за неуместный пример. Я использовал набор данных iris, потому что это встроенный простой набор данных. Я полностью согласен с вашим беспокойством. Я не пытаюсь проиллюстрировать конкретную проблему через это, извините за неуместное выражение.

li jiaqi 30.01.2023 05:29

Не беспокойся. Я думаю, это просто хорошо знать в будущем. Кстати, пакеты lme4 и lmerTest на самом деле содержат наборы данных, с которыми можно попрактиковаться.

Shawn Hemelstrand 30.01.2023 05:58
Стоит ли изучать PHP в 2023-2024 годах?
Стоит ли изучать PHP в 2023-2024 годах?
Привет всем, сегодня я хочу высказать свои соображения по поводу вопроса, который я уже много раз получал в своем сообществе: "Стоит ли изучать PHP в...
Поведение ключевого слова "this" в стрелочной функции в сравнении с нормальной функцией
Поведение ключевого слова "this" в стрелочной функции в сравнении с нормальной функцией
В JavaScript одним из самых запутанных понятий является поведение ключевого слова "this" в стрелочной и обычной функциях.
Приемы CSS-макетирования - floats и Flexbox
Приемы CSS-макетирования - floats и Flexbox
Здравствуйте, друзья-студенты! Готовы совершенствовать свои навыки веб-дизайна? Сегодня в нашем путешествии мы рассмотрим приемы CSS-верстки - в...
Тестирование функциональных ngrx-эффектов в Angular 16 с помощью Jest
В системе управления состояниями ngrx, совместимой с Angular 16, появились функциональные эффекты. Это здорово и делает код определенно легче для...
Концепция локализации и ее применение в приложениях React ⚡️
Концепция локализации и ее применение в приложениях React ⚡️
Локализация - это процесс адаптации приложения к различным языкам и культурным требованиям. Это позволяет пользователям получить опыт, соответствующий...
Пользовательский скаляр GraphQL
Пользовательский скаляр GraphQL
Листовые узлы системы типов GraphQL называются скалярами. Достигнув скалярного типа, невозможно спуститься дальше по иерархии типов. Скалярный тип...
1
6
65
1
Перейти к ответу Данный вопрос помечен как решенный

Ответы 1

Ответ принят как подходящий

Рекомендация

Я согласен с Роландом в том, что использование смешанной модели с тремя кластерами нецелесообразно (см. Gelman & Hill, 2007, ссылка ниже). Говоря о вашем основном вопросе, мы хотели бы знать, сколько случайных вариаций анализируется моделью, и на самом деле это действительно должно использоваться в качестве оправдания для использования модели в первую очередь (см. Meteyard & Davies, 2020). , где они предлагают подгонку моделей только со случайными эффектами в качестве предшественника вашей полной модели).

Имея это в виду, я бы рекомендовал сообщать об этих частях вашей модели для учета случайных эффектов (включая то, что вы упомянули). Хотя они не все необходимы, они, безусловно, информативны:

  • SD случайных пересечений/случайных наклонов
  • Корреляция любых случайных эффектов (и, если возможно, объяснение почему)
  • ICC вашей модели (это объяснит, насколько сильно происходит кластеризация)
  • Псевдо R2, который пытается объяснить, какая часть эффектов объясняется фиксированными эффектами, а какая — фиксированными и случайными эффектами.
  • Гусеничный сюжет случайных эффектов.

Рабочий пример: подгонка модели

Я буду продолжать использовать данные, которые у вас есть, в качестве примера, но я изменю их на что-то сходящееся, поскольку ваша первая модель была единственной, которую нельзя использовать при составлении отчетов. Я загрузил четыре библиотеки, необходимые для используемых функций и соответствующих модели.

#### Load Libraries ####
library(lmerTest)
library(lattice)
library(broom.mixed)
library(performance)

#### Fit Model ####
fit <- lmer(Petal.Length
                    ~ Petal.Width
                    + (1 + Petal.Width |Species), 
                    data = iris)

Суммируя дисперсию случайного эффекта

Далее мы можем обобщить модель, но я предпочитаю использовать функцию tidy из пакета broom.mixed для быстрого получения оценок RE.

##### Summarise Model ####
tidy(fit)

Показано ниже:

# A tibble: 6 × 8
  effect   group    term                      estim…¹ std.e…² stati…³    df p.value
  <chr>    <chr>    <chr>                       <dbl>   <dbl>   <dbl> <dbl>   <dbl>
1 fixed    NA       (Intercept)                 2.43    0.898    2.71  1.97   0.115
2 fixed    NA       Petal.Width                 1.10    0.414    2.65  2.00   0.118
3 ran_pars Species  sd__(Intercept)             1.52   NA       NA    NA     NA    
4 ran_pars Species  cor__(Intercept).Petal.W…  -0.408  NA       NA    NA     NA    
5 ran_pars Species  sd__Petal.Width             0.646  NA       NA    NA     NA    
6 ran_pars Residual sd__Observation             0.362  NA       NA    NA     NA    
# … with abbreviated variable names ¹​estimate, ²​std.error, ³​statistic

Основываясь на этой информации, мы теперь знаем:

  • Виды различаются по условной средней длине лепестков примерно на 1,52 SD (после учета фиксированных и случайных эффектов).
  • Наклоны по видам различаются меньше (0,646 SD), что указывает на то, что соотношение между шириной лепестка и длиной лепестка не сильно меняется в зависимости от вида.
  • Корреляция между шириной лепестка и длиной лепестка составляет около -0,41, что указывает на то, что по мере увеличения средней длины лепестка наклоны обычно уменьшаются.

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

График случайных эффектов Caterpillar

Мы можем сообщить эту информацию в статье, если это необходимо, но мы также можем изобразить эти эффекты с помощью графика гусеницы:

#### Plot Dotplot of Random Effects ####
dotplot(ranef(fit))

Отчасти поэтому Роланд предупредил вас о большем количестве кластеров случайных эффектов. Здесь вы видите только три пересечения и наклона, что не отражает большого количества вариаций, которые вы хотели бы получить от смешанной модели. Однако здесь мы можем сделать два ключевых вывода:

  1. Растение virginica в среднем имеет гораздо большую длину лепестков (случайные пересечения), чем другие виды.
  2. Разноцветное растение имеет более высокий коэффициент наклона для ширины лепестков, чем другие.

Другие оценки случайных эффектов

Мы также можем получить оценки ICC и псевдо R2 с помощью этого кода:

#### Check Performance ####
performance(fit)

Показано ниже:

# Indices of model performance

AIC     |    AICc |     BIC | R2 (cond.) | R2 (marg.) |   ICC |  RMSE | Sigma
-----------------------------------------------------------------------------
154.478 | 155.066 | 172.542 |      0.957 |      0.231 | 0.944 | 0.355 | 0.362

Что следует отметить из этого резюме:

  • ICC указывает, что существует множество кластеров случайных эффектов ... около 94,4%, если быть точным, что указывает на то, что большая часть объясненной дисперсии исходит от случайных эффектов в модели.
  • Условный R2 — это дисперсия, объясняемая как фиксированными, так и случайными эффектами, которая оказывается равной 95,7%.
  • Предельный R2 указывает, какая часть дисперсии приходится на фиксированные эффекты, что составляет около 23,1%.

Цитаты

Ниже приведены цитаты, о которых я упоминал ранее. Gelman & Hill — канонический источник информации о смешанных моделях. Статья Meteyard & Davies представляет собой руководство по передовому опыту работы со смешанными моделями. Дайте мне знать, если вы нашли этот ответ полезным.

  • Гельман, А., и Хилл, Дж. (2007). Анализ данных с использованием регрессионных и многоуровневых/иерархических моделей. Издательство Кембриджского университета.
  • Метейард, Л., и Дэвис, Р.А.И. (2020). Руководство по передовой практике для линейных моделей смешанных эффектов в психологии. Журнал памяти и языка, 112, 104092. https://doi.org/10.1016/j.jml.2020.104092

Большое спасибо за подробное объяснение и цитаты. Вот что я хочу знать. Извините за неуместный пример. Я использовал набор данных iris, потому что это встроенный простой набор данных. Я полностью согласен с беспокойством вас и Роланда. Руководство по передовой практике великолепно, и я буду продолжать читать его.

li jiaqi 30.01.2023 05:23

Без проблем. Рад, что был полезен!

Shawn Hemelstrand 30.01.2023 05:56

Другие вопросы по теме