Какой метод множественного сравнения использовать для модели lmer: lsmeans или glht?

15

Я анализирую набор данных, используя модель смешанных эффектов с одним фиксированным эффектом (условием) и двумя случайными эффектами (участник из-за дизайна объекта и пары). Модель была сгенерирована с lme4пакетом: exp.model<-lmer(outcome~condition+(1|participant)+(1|pair),data=exp).

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

Из того, что я видел, люди предложили использовать

1.lsmeans пакет - lsmeans(exp.model,pairwise~condition)который дает мне выход следующий:

condition     lsmean         SE    df  lower.CL  upper.CL
 Condition1 0.6538060 0.03272705 47.98 0.5880030 0.7196089
 Condition2 0.7027413 0.03272705 47.98 0.6369384 0.7685443
 Condition3 0.7580522 0.03272705 47.98 0.6922493 0.8238552

Confidence level used: 0.95 

$contrasts
 contrast                   estimate         SE    df t.ratio p.value
 Condition1 - Condition2 -0.04893538 0.03813262 62.07  -1.283  0.4099
 Condition1 - Condition3 -0.10424628 0.03813262 62.07  -2.734  0.0219
 Condition2 - Condition3 -0.05531090 0.03813262 62.07  -1.450  0.3217

P value adjustment: tukey method for comparing a family of 3 estimates 

2.multcomp пакет двумя различными способами - с использованием в mcp glht(exp.model,mcp(condition="Tukey"))результате

     Simultaneous Tests for General Linear Hypotheses

Multiple Comparisons of Means: Tukey Contrasts


Fit: lmer(formula = outcome ~ condition + (1 | participant) + (1 | pair), 
    data = exp, REML = FALSE)

Linear Hypotheses:
                             Estimate Std. Error z value Pr(>|z|)  
Condition2 - Condition1 == 0  0.04894    0.03749   1.305    0.392  
Condition3 - Condition1 == 0  0.10425    0.03749   2.781    0.015 *
Condition3 - Condition2 == 0  0.05531    0.03749   1.475    0.303  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1   1
(Adjusted p values reported -- single-step method)

и используя в lsm glht(exp.model,lsm(pairwise~condition))результате

Note: df set to 62

     Simultaneous Tests for General Linear Hypotheses

Fit: lmer(formula = outcome ~ condition + (1 | participant) + (1 | pair), 
    data = exp, REML = FALSE)

Linear Hypotheses:
                             Estimate Std. Error t value Pr(>|t|)  
Condition1 - Condition2 == 0 -0.04894    0.03749  -1.305   0.3977  
Condition1 - Condition3 == 0 -0.10425    0.03749  -2.781   0.0195 *
Condition2 - Condition3 == 0 -0.05531    0.03749  -1.475   0.3098  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1   1
(Adjusted p values reported -- single-step method)

Как видите, методы дают разные результаты. Я впервые работаю с R и статистикой, поэтому что-то может пойти не так, но я не знаю где. Мои вопросы:

Каковы различия между представленными методами? Я прочитал в ответе на связанные вопросы, что это о степенях свободы ( lsmeansпротив glht). Существуют ли какие-либо правила или рекомендации, когда использовать какой-либо метод, т. Е. Метод 1, подходящий для этого типа набора данных / модели и т. Д.? О каком результате я должен сообщить? Не зная лучше, я бы, вероятно, просто пошел и сообщил о самом высоком p-значении, которое я получил, чтобы не рисковать, но было бы неплохо иметь более вескую причину. Благодарность

schvaba986
источник

Ответы:

17

Не полный ответ ...

Разница между glht(myfit, mcp(myfactor="Tukey"))этими двумя и двумя другими методами заключается в том, что в этом способе используется статистика "z" (нормальное распределение), в то время как другие используют статистику "t" (распределение Стьюдента). Статистика "z" это то же самое, что статистика "t" с бесконечной степенью свободы. Этот метод является асимптотическим и обеспечивает меньшие значения p и более короткие доверительные интервалы, чем другие. Значения p могут быть слишком малыми, а доверительные интервалы могут быть слишком короткими, если набор данных мал.

При запуске lsmeans(myfit, pairwise~myfactor)появляется следующее сообщение:

Loading required namespace: pbkrtest

Это означает, что lsmeans(для lmerмодели) используется pbkrtestпакет, который реализует метод Кенварда и Роджерса для степеней свободы статистики "t". Этот метод предназначен для обеспечения лучших p-значений и доверительных интервалов, чем асимптотический (нет разницы, когда степень свободы велика).

Теперь о разнице между lsmeans(myfit, pairwise~myfactor)$contrastsи glht(myfit, lsm(pairwise~factor), я только что провел несколько тестов, и мои наблюдения следующие:

  • lsmинтерфейс между lsmeansпакетом и multcompпакетом (см. ?lsm)

  • для сбалансированного дизайна нет разницы между результатами

  • для несбалансированного дизайна я заметил небольшие различия между результатами (стандартные ошибки и отношение т)

К сожалению, я не знаю, что является причиной этих различий. Это выглядит как lsmвызовы lsmeansтолько для получения матрицы линейных гипотез и степеней свободы, но lsmeansиспользует другой способ для вычисления стандартных ошибок.

Стефан Лоран
источник
Спасибо за подробный ответ! Я полностью пропустил разницу в статистике теста ... Вы упоминаете, что значения могут быть слишком маленькими, а CI слишком узкими для асимптотического метода. Мой набор данных состоит из ~ 30 участников, поэтому, я думаю, я буду придерживаться t-статистики. Когда вы говорите, что метод Кенварда и Роджерса приводит к лучшим p-значениям, вы имеете в виду более точный или меньший? Таким образом, различия связаны с различиями в методах расчета df и SE, а не с неправильным использованием одного из них в моей модели, если я вас правильно понял. Есть ли способ выбрать «лучший» метод здесь?
schvaba986
11
(Я - разработчик пакета lsmeans ) lsmeansиспользует пакет pbkrtest, который обеспечивает (1) вычисления Kenward-Rogers df и (2) скорректированную ковариационную матрицу с уменьшенным смещением в оценках. Если вы сначала установите lsm.options(disable.pbkrtest=TRUE), то lsmeansвызов с adjust="mvt"выдаст те же результаты, что glhtи за исключением небольших различий из-за рандомизированного алгоритма, используемого обоими пакетами для многомерного t-распределения.
Расс Лент
3
Тем не менее, я предлагаю настройку "mvt" без отключения pbkrtest из-за настройки смещения и того факта, что при отсутствии df асимптотические (z) значения по существу предполагают бесконечное df, что приводит к нереально низким значениям P.
Расс Лент
3
Кстати, этот summaryметод glhtпозволяет использовать различные методы пошагового тестирования, помимо настройки множественности одношаговых (одновременных КИ) по умолчанию. С другой стороны, если у вас есть более одного фактора, вы lsmможете легко создавать обычные типы сравнений, но mcpне можете делать это вообще.
Расс Лент