Jakiej metody wielokrotnego porównania użyć w modelu Lmer: lsmeans czy glht?

16

Analizuję zestaw danych przy użyciu modelu efektów mieszanych z jednym ustalonym efektem (warunkiem) i dwoma efektami losowymi (uczestnik ze względu na projekt i parę wewnątrz przedmiotu). Model ten został wygenerowany z lme4pakietu: exp.model<-lmer(outcome~condition+(1|participant)+(1|pair),data=exp).

Następnie wykonałem test współczynnika wiarygodności tego modelu względem modelu bez ustalonego efektu (warunku) i mam znaczącą różnicę. W moim zestawie danych są 3 warunki, więc chcę wykonać wielokrotne porównanie, ale nie jestem pewien, której metody użyć . Znalazłem wiele podobnych pytań na CrossValidated i innych forach, ale nadal jestem dość zdezorientowany.

Z tego, co widziałem, ludzie sugerowali używanie

1.lsmeans pakiet - lsmeans(exp.model,pairwise~condition)który daje mi wyjście następujące:

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 opakowanie na dwa różne sposoby - z użyciem mcp glht(exp.model,mcp(condition="Tukey"))powodując

     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)

i używanie lsm glht(exp.model,lsm(pairwise~condition))w wyniku

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)

Jak widać, metody dają różne wyniki. Po raz pierwszy pracuję z R i statystykami, więc coś może pójść nie tak, ale nie wiedziałbym gdzie. Moje pytania to:

Jakie są różnice między prezentowanymi metodami? W odpowiedzi na powiązane pytania przeczytałem, że chodzi o stopnie swobody (w lsmeansporównaniu glht). Czy istnieją jakieś zasady lub zalecenia, kiedy użyć której, tj. Metoda 1 jest odpowiednia dla tego typu zestawu danych / modelu itp.? Który wynik powinienem zgłosić? Nie wiedząc lepiej, prawdopodobnie po prostu podałbym najwyższą wartość p, którą mogłem grać bezpiecznie, ale byłoby miło mieć lepszy powód. Dzięki

schvaba986
źródło

Odpowiedzi:

17

Nie jest to kompletna odpowiedź ...

Różnica między glht(myfit, mcp(myfactor="Tukey"))dwiema innymi metodami polega na tym, że w ten sposób wykorzystuje się statystykę „z” (rozkład normalny), podczas gdy inne wykorzystują statystykę „t” (rozkład Studenta). Statystyka „z” jest taka sama jak statystyka „t” z nieskończonym stopniem swobody. Ta metoda jest asymptotyczna i zapewnia mniejsze wartości p i krótsze przedziały ufności niż inne. Wartości p mogą być zbyt małe, a przedziały ufności mogą być zbyt krótkie, jeśli zestaw danych jest mały.

Po uruchomieniu lsmeans(myfit, pairwise~myfactor)pojawia się następujący komunikat:

Loading required namespace: pbkrtest

Oznacza to, że lsmeans(dla lmermodelu) używa pbkrtestpakietu, który implementuje metodę Kenwarda i Rogersa dla stopni swobody statystyki „t”. Ta metoda ma zapewnić lepsze wartości p i przedziały ufności niż asymptotyczne (nie ma różnicy, gdy stopień swobody jest duży).

Teraz, o różnicy między lsmeans(myfit, pairwise~myfactor)$contrastsi glht(myfit, lsm(pairwise~factor), mam tylko zrobić kilka testów i moje spostrzeżenia są następujące:

  • lsmjest interfejsem między lsmeanspakietem a multcomppakietem (patrz ?lsm)

  • dla zrównoważonego projektu nie ma różnicy między wynikami

  • w przypadku niezrównoważonego projektu zaobserwowałem małe różnice między wynikami (błędy standardowe i stosunek t)

Niestety nie wiem, co jest przyczyną tych różnic. Wygląda jak lsmwywołania lsmeanstylko w celu uzyskania liniowej macierzy hipotez i stopni swobody, ale lsmeanswykorzystuje inny sposób obliczania standardowych błędów.

Stéphane Laurent
źródło
Dzięki za szczegółową odpowiedź! Całkowicie przeoczyłem różnicę w statystyce testowej ... Wspominasz, że wartości mogą być zbyt małe, a CI zbyt wąskie dla metody asymptotycznej. Mój zestaw danych składa się z ~ 30 uczestników, więc chyba trzymam się statystyki t. Kiedy mówisz, że metoda Kenwarda i Rogersa prowadzi do lepszych wartości p, masz na myśli bardziej dokładne czy mniejsze? Tak więc różnice wynikają z różnic w metodach obliczeniowych df i SE, a nie z powodu niewłaściwego użycia jednej z nich w moim modelu, jeśli dobrze cię zrozumiałem. Czy jest tutaj sposób na wybranie „najlepszej” metody?
schvaba986,
11
(Jestem deweloperem pakietu lsmeans ) lsmeansużywa pakietu pbkrtest, który zapewnia (1) obliczenia Kenward-Rogers df i (2) skorygowaną macierz kowariancji ze zmniejszonym odchyleniem w oszacowaniach. Jeśli ustawisz po raz pierwszy lsm.options(disable.pbkrtest=TRUE), lsmeanswywołanie z adjust="mvt"przyniesie takie same wyniki, jak glht, z wyjątkiem niewielkich różnic wynikających z algorytmu randomizowanego używanego przez oba pakiety do dystrybucji wielowymiarowej.
Russ Lenth
3
Sugeruję jednak korektę „mvt” bez wyłączania pbkrtest, ze względu na korektę polaryzacji i fakt, że bez df asymptotyczne (z) wartości zasadniczo przyjmują nieskończone df, tym samym uzyskując nierealistycznie niskie wartości P.
Russ Lenth
3
Nawiasem mówiąc, summarymetoda glhtpozwala na różne metody testowania obniżającego, oprócz domyślnej korekcji mnogości w jednym kroku (jednoczesne CI). W zupełnie innym punkcie, jeśli masz więcej niż jeden czynnik, lsmmożesz dość łatwo tworzyć zwykłe typy porównań, ale mcpwcale nie możesz tego zrobić.
Russ Lenth