Porównanie modeli efektów mieszanych i modeli efektów stałych (testowanie istotności efektów losowych)

Podane trzy zmienne, yi x, które są pozytywne i ciągły z, co jest kategoryczny, mam dwa modele kandydujące dana przez:

fit.me <- lmer( y ~ 1 + x + ( 1 + x | factor(z) ) )

fit.fe <- lm( y ~ 1 + x )

Mam nadzieję porównać te modele, aby ustalić, który model jest bardziej odpowiedni. Wydaje mi się, że w pewnym sensie fit.fejest zagnieżdżone fit.me. Zazwyczaj, gdy utrzymuje się ten ogólny scenariusz, można wykonać test chi-kwadrat. W Rmożemy wykonać ten test za pomocą następującego polecenia,

anova(fit.fe,fit.me)

Gdy oba modele zawierają efekty losowe (generowane przez lmerz lme4pakietu), anova()polecenie działa poprawnie. Ze względu na parametry brzegowe zwykle zaleca się przetestowanie uzyskanej statystyki Chi-Square za pomocą symulacji, niemniej jednak nadal możemy korzystać ze statystyki w procedurze symulacji.

Gdy oba modele zawierają tylko efekty stałe, to podejście --- i powiązane anova()polecenie --- działają dobrze.

Jednak gdy jeden model zawiera efekty losowe, a model zredukowany zawiera tylko efekty stałe, tak jak w powyższym scenariuszu, anova()polecenie nie działa.

Mówiąc dokładniej, pojawia się następujący błąd:

 > anova(fit.fe, fit.me)
 Error: $ operator not defined for this S4 class

Czy jest coś złego w podejściu Chi-Square z góry (z symulacją)? Czy jest to po prostu problem polegający na tym, że anova()nie umiemy radzić sobie z modelami liniowymi generowanymi przez różne funkcje?

Innymi słowy, czy właściwe byłoby ręczne wygenerowanie statystyki Chi-Square uzyskanej z modeli? Jeśli tak, jakie są odpowiednie stopnie swobody przy porównywaniu tych modeli? Według mojej oceny:

F = \frac{((S S E_{r e d u c e d} - S S E_{f u l l}) / (p - k))}{((S S E_{f u l l}) / (n - p - 1))} \sim F_{p - k, n - p - 1}

$F = \frac{\left((SSE_{reduced}-SSE_{full})/(p-k)\right)}{\left((SSE_{full})/(n-p-1)\right)} \sim F_{p-k,n-p-1}$

Szacujemy dwa parametry w modelu efektów stałych (nachylenie i przecięcie) i dwa kolejne parametry (parametry wariancji dla losowego nachylenia i losowego przechwytywania) w modelu efektów mieszanych. Zazwyczaj parametr przecięcia nie jest liczona w obliczeniach stopni swobody, tak, że oznacza, że , a ; powiedziawszy, że nie jestem pewien, czy parametry wariancji parametrów efektów losowych należy uwzględnić w obliczeniach stopni swobody; szacunki wariancji dla parametrów o stałym efekcie nie są brane pod uwagę , ale uważam, że dzieje się tak, ponieważ zakłada się, że oszacowania parametrów dla stałych efektów są nieznane, podczas gdy są one uważane za niepoznawalne zmienne losowe $k=1$ $p=k+2=3$ dla mieszanych efektów. Byłbym wdzięczny za pomoc w tej sprawie.

Wreszcie, czy ktoś ma bardziej odpowiednie ( Roparte) rozwiązanie do porównywania tych modeli?

r hypothesis-testing mixed-model lme4-nlme

— użytkownik9171
źródło

Jeśli zastąpi lm()się gls()z nlmeopakowania i lmer()z lme()(znowu od nlmepakietu), wszystko będzie działać dobrze. Pamiętaj jednak, że otrzymasz test zachowawczy (zbyt duże wartości p ), ponieważ parametry prostszego modelu znajdują się na granicy przestrzeni parametrów. I tak naprawdę wybór, czy uwzględnić efekty losowe, powinien opierać się na teorii (np. Planie próbkowania), a nie na teście statystycznym.

— Karl Ove Hufthammer

Co chcesz zrobić z modelami? Jeden model może być lepszy do niektórych celów, a drugi model lepszy do innych celów. Wszystkie modele są błędne, więc pytanie nie brzmi, który model jest odpowiedni, ale który jest bardziej przydatny w przypadku konkretnego problemu.

— Kodiolog,

@Kodiologist Zasadniczo chcę upewnić się, że oszacowania parametrów dla ustalonych efektów są wiarygodne. Standardowe błędy, które mogą być niewiarygodne, jeżeli zakłada się, że obserwacje są niezależne. Ponadto byłoby miło powiedzieć kilka słów o tym, jak zmienna jest efektem losowym, ale myślę, że nie jest to tak istotne.

— user9171

@ user9171 Dobrym sposobem sprawdzenia stabilności (niezawodności) w oszacowaniach parametrów modelu jest użycie ładowania początkowego. Wykresy rozkładu ładowania początkowego dla każdego parametru dzielonego przez dwa modele, z jednym wykresem na parametr i model. Węższe rozkłady oznaczają większą stabilność. Prawdopodobnie okaże się, że prostszy model daje bardziej stabilne oszacowania, ponieważ mniej parametrów pozwala na bardziej precyzyjne oszacowanie każdego parametru.

— Kodiolog

Technicznie możesz to uruchomić, zmieniając kolejność parametrów:

> anova(fit.me, fit.fe)

Będzie działać dobrze. Jeśli miniesz obiekt wygenerowany przez lmerpierwszy, anova.merModzostanie wywołany zamiast anova.lm(co nie wie, jak obsługiwać lmerobiekty). Widzieć:

?anova.merMod

Chociaż wybór modelu mieszanego lub modelu stałego jest wyborem modelowania, który musi uwzględniać projekt eksperymentalny, a nie problem wyboru modelu. Zobacz @ BenBolkera https://bbolker.github.io/mixedmodels-misc/glmmFAQ.html#testing-signiance-of-landom-effects po więcej szczegółów:

Rozważ nie testowanie znaczenia efektów losowych.

— witek
źródło

+1. Zezwoliłem na wstawienie linku do FAQ @ BenBolkera, który zawiera dalsze dyskusje i odniesienia.

— ameba