Повторные измерения ANOVA и ссылки на модели смешанных эффектов в R

У меня возникла проблема при выполнении двустороннего среднеквадратичного дисперсионного анализа в R для следующих данных (ссылка: https://drive.google.com/open?id=1nIlFfijUm4Ib6TJoHUUNeEJnZnnNzO29):

subjectnbr — это идентификатор субъекта, а blockType и linesTTL — независимые переменные. RT2 — зависимая переменная

Сначала я выполнил rm ANOVA с помощью ezANOVA со следующим кодом:

ANOVA_RTS <- ezANOVA(
    data=castRTs
    , dv=RT2
    , wid=subjectnbr
    , within = .(blockType,linesTTL)
    , type = 2
    , detailed = TRUE
    , return_aov = FALSE
)
ANOVA_RTS

Результат правильный (я перепроверил с помощью статистики).

Однако, когда я выполняю rm ANOVA с помощью функции lme, я получаю другой ответ, и я понятия не имею, почему.

Вот мой код:

lmeRTs <- lme(
      RT2 ~ blockType*linesTTL, 
      random = ~1|subjectnbr/blockType/linesTTL, 
      data=castRTs)

anova(lmeRTs)  

Вот результаты ezANOVA и lme.

Я надеюсь, что я был достаточно ясен и дал вам всю необходимую информацию.

Я с нетерпением жду вашей помощи, поскольку я пытаюсь понять это как минимум 4 часа!

Заранее спасибо.


person user9935785    schedule 13.06.2018    source источник
comment
Ваше описание проблемы почти идеально. Не хватает только данных. Первых 15 строк недостаточно, чтобы воспроизвести вашу проблему. Не могли бы вы попробовать разработать минимальный экземпляр вашей проблемы с небольшим набором данных, который вы бы включили в свой вопрос? Это может даже помочь вам понять, что не так. См. stackoverflow.com/help/mcve, чтобы узнать, как создать минимальный, полный и проверяемый пример. Надеюсь это поможет ;-)   -  person Oli    schedule 13.06.2018
comment
Да, включите пример данных и ожидаемый результат от ezANOVA. Я думаю, что проблема связана с вашей спецификацией случайных эффектов в lme, где вам нужно учитывать, что случайные эффекты blockType:subjectnbr и linesTTL:subjectnbr пересекаются. Эти двухфакторные модели случайных эффектов могут быть довольно сложными для указания в lme. Чтобы дать более конкретную помощь, нам нужен воспроизводимый пример с данными.   -  person Maurits Evers    schedule 13.06.2018
comment
Большое спасибо за ваши комментарии. Вот ссылка с моими данными: drive.google.com/open?id=1nIlFfijUm4Ib6TJoHUUNeEJnZnnNzO29 Ура :)   -  person user9935785    schedule 13.06.2018
comment
@ user9935785 Я добавил ответ, показывающий, как воспроизвести результат ezANOVA с помощью nlme::lme. Я также немного изменил заголовок вашего вопроса, чтобы лучше отразить постановку проблемы в основной части вашего сообщения.   -  person Maurits Evers    schedule 14.06.2018


Ответы (1)


Вот пошаговый пример того, как воспроизвести результаты ezANOVA с помощью nlme::lme.

Данные

Мы читаем данные и убеждаемся, что все категориальные переменные равны factors.

# Read in data
library(tidyverse);
df <- read.csv("castRTs.csv");
df <- df %>%
    mutate(
        blockType = factor(blockType),
        linesTTL = factor(linesTTL));

Результаты от ezANOVA

В качестве проверки воспроизводим результаты ez::ezANOVA.

## ANOVA using ez::ezANOVA
library(ez);
model1 <- ezANOVA(
    data = df,
    dv = RT2,
    wid = subjectnbr,
    within = .(blockType, linesTTL),
    type = 2,
    detailed = TRUE,
    return_aov = FALSE);
model1;
#        $ANOVA
#              Effect DFn DFd          SSn       SSd           F            p
#1        (Intercept)   1  13 2047405.6654 34886.767 762.9332235 6.260010e-13
#2          blockType   1  13     236.5412  5011.442   0.6136028 4.474711e-01
#3           linesTTL   1  13    6584.7222  7294.620  11.7348665 4.514589e-03
#4 blockType:linesTTL   1  13    1019.1854  2521.860   5.2538251 3.922784e-02
#  p<.05         ges
#1     * 0.976293831
#2       0.004735442
#3     * 0.116958989
#4     * 0.020088855

Результаты от nlme::lme

Теперь мы запускаем nlme::lme

## ANOVA using nlme::lme
library(nlme);
model2 <- anova(lme(
    RT2 ~ blockType * linesTTL,
    random = list(subjectnbr = pdBlocked(list(~1, pdIdent(~blockType - 1), pdIdent(~linesTTL - 1)))),
    data = df))
model2;
#                   numDF denDF  F-value p-value
#(Intercept)            1    39 762.9332  <.0001
#blockType              1    39   0.6136  0.4382
#linesTTL               1    39  11.7349  0.0015
#blockType:linesTTL     1    39   5.2538  0.0274

Результаты/заключение

Мы видим, что результаты F-теста обоих методов идентичны. Несколько сложная структура определения эффекта random в lme возникает из-за того, что у вас есть два перекрещенных случайных эффекта. Здесь пересечение означает, что для каждой комбинации blockType и linesTTL существует наблюдение для каждого subjectnbr.

Некоторые дополнительные (необязательные) детали

Чтобы понять роль pdBlocked и pdIdent, нам нужно взглянуть на соответствующую двухуровневую модель смешанного эффекта.

введите здесь описание изображения

Переменные-предикторы введите здесь описание изображения — это ваши категориальные переменные blockType и linesTTL, которые обычно кодируются с использованием фиктивных переменных.

Матрица дисперсии-ковариации для случайных эффектов введите здесь описание изображения может принимать различные формы в зависимости от базовой структуры корреляции ваших коэффициентов случайного эффекта. Чтобы соответствовать предположениям двухуровневого повторного измерения ANOVA, мы должны указать матрицу блочно-диагональной дисперсии-ковариации pdBlocked, где мы создаем диагональные блоки для смещения ~1 и для категориальных переменных-предикторов blockType pdIdent(~blockType - 1) и linesTTL pdIdent(~linesTTL - 1), соответственно. Обратите внимание, что нам нужно вычесть смещение из двух последних блоков (поскольку мы уже учли смещение).

Некоторые актуальные/интересные ресурсы

Пиньейро и Бейтс, Модели смешанных эффектов в S и S-PLUS , Спрингер (2000)

Потвин и Шутц, Статистическая мощность для двухфакторных повторных измерений ANOVA, Методы исследования поведения , Инструменты и компьютеры, 32, 347–356 (2000)

Деминг Ми, Как понимать и применять смешанные модели эффектов, кафедра биостатистики, университет Вандербильта

person Maurits Evers    schedule 13.06.2018