Porównanie modelu mieszanego (podmiot jako efekt losowy) z prostym modelem liniowym (przedmiot jako efekt stały)

Kończę analizę dużej grupy danych. Chciałbym wziąć model liniowy zastosowany w pierwszej części pracy i ponownie go dopasować za pomocą liniowego modelu mieszanego (LME). LME byłby bardzo podobny, z tym wyjątkiem, że jedna ze zmiennych zastosowanych w modelu byłaby zastosowana jako efekt losowy. Te dane pochodzą z wielu obserwacji (> 1000) w małej grupie podmiotów (~ 10) i wiem, że modelowanie efektu podmiotu lepiej jest wykonać jako efekt losowy (jest to zmienna, którą chcę przesunąć). Kod R wyglądałby następująco:

my_modelB <- lm(formula = A ~ B + C + D)    
lme_model <- lme(fixed=A ~ B + C, random=~1|D, data=my_data, method='REML')

Wszystko działa dobrze, a wyniki są bardzo podobne. Byłoby miło, gdybym mógł użyć czegoś takiego jak RLRsim lub AIC / BIC, aby porównać te dwa modele i zdecydować, który jest najbardziej odpowiedni. Moi koledzy nie chcą zgłaszać LME, ponieważ nie ma łatwo dostępnego sposobu wyboru, który jest „lepszy”, mimo że uważam, że LME jest bardziej odpowiednim modelem. Jakieś sugestie?

— MudPhud
źródło

Odpowiedzi:

Ma to na celu dodanie odpowiedzi @ ocram, ponieważ jest zbyt długi, aby opublikować komentarz. Traktowałbym A ~ B + Cjak twój model zerowy, abyś mógł ocenić istotność statystyczną Dprzechwytywania losowego na poziomie w konfiguracji modelu zagnieżdżonego. Jak wskazał ocram, warunki regularności są naruszane, gdy , a statystyka testu wskaźnika prawdopodobieństwa (LRT) niekoniecznie będzie asymptotycznie rozłożona . Rozwiązaniem, którego nauczyłem, było ładowanie LRT (którego rozkład prawdopodobnie nie będzie wynosił ) parametrycznie i obliczanie wartości p ładowania w następujący sposób: $H_0: \sigma^2 = 0$ $\chi^2$ $\chi^2$

library(lme4)
my_modelB <- lm(formula = A ~ B + C)
lme_model <- lmer(y ~ B + C + (1|D), data=my_data, REML=F)
lrt.observed <- as.numeric(2*(logLik(lme_model) - logLik(my_modelB)))
nsim <- 999
lrt.sim <- numeric(nsim)
for (i in 1:nsim) {
    y <- unlist(simulate(mymodlB))
    nullmod <- lm(y ~ B + C)
    altmod <- lmer(y ~ B + C + (1|D), data=my_data, REML=F)
    lrt.sim[i] <- as.numeric(2*(logLik(altmod) - logLik(nullmod)))
}
mean(lrt.sim > lrt.observed) #pvalue

Odsetek LRT z ładowaniem początkowym bardziej ekstremalny niż obserwowany LRT jest wartością p.

— zablokowane
źródło

Dziękuję za wypełnienie mojej odpowiedzi. Czasami ludzie używają mieszanki chi-kwadratów zamiast rozkładu chi-kwadrat do statystyk testowych.

— ocram

@ocram +1 za komentarz dotyczący decyzji, czy traktować zmienną jako losową, czy ustaloną oddzielnie od analizy. @MudPhud Jeśli twój PI nie rozumie problemu i nalega na wartość p, to może po prostu pokaż mu wynik testu losowego efektu (który i tak byś uwzględnił w piśmie).

— zablokowane

Dzięki za kod. Kiedy wpadłem to wynik jest żaden z bootstrapped LRTs są większe niż obserwowane, więc oznacza to, że mogę trzymać LM bez efektów losowych lub nawet zmiennej oryginalnej rzucone w.

— MudPhud

@MudPhud: Czy wystąpiły jakieś błędy? Spróbuj wpisać, lrt.simaby upewnić się, że nie wszystkie są zerami. W takim przypadku najbardziej prawdopodobnym winowajcą byłoby to, że nie masz lme4zainstalowanego pakietu .

— zablokowane

Nie są zerowe, tylko bardzo małe (~ 1e-6) w porównaniu do zaobserwowanych (63,95).

— MudPhud,

$0$ $H_0:variance = 0$ $H_1: variance > 0$ ...

EDYTOWAĆ

Aby uniknąć nieporozumień: Wspomniany powyżej test jest czasem wykorzystywany do podjęcia decyzji, czy efekt losowy jest znaczący ... ale nie do podjęcia decyzji, czy należy go przekształcić w efekt stały.

— ocram
źródło

Pytanie brzmi: czy istnieją testy pozwalające zdecydować, czy zmienna powinna być modelowana jako efekt mieszany czy efekt losowy? W przeciwnym razie możesz wykonać opisany test, a następnie przetestować go za pomocą dystansu chi-kwadrat (nie jestem pewien, jaki byłby odpowiedni test).

— MudPhud

@MudPhud: Modelowanie zmiennej jako efektu ustalonego lub losowego powinno faktycznie zostać ustalone przed analizą, kiedy badanie jest planowane. Zależy to w szczególności od zakresu twoich wniosków. Losowe efekty pozwalają na większą uogólnienie. Mogłoby to również uniknąć pewnych trudności technicznych. Na przykład asymptotyki mogą się załamać, gdy wzrośnie liczba parametrów, jak to jest w przypadku, gdy zmienna kategoryczna o wielu poziomach jest uważana za zmienną stałą.

— ocram

Zgadzam się, ale kiedy próbowałem wyjaśnić to mojemu PI, po prostu się odwrócił i poprosił o jakąś wartość p. Chcę zawrzeć tę analizę w manuskrypcie, ale nie włoży jej, jeśli nie będzie bardziej konkretnego uzasadnienia.

— MudPhud

@MudPhud: Zgodnie z moją najlepszą wiedzą nie ma wartości p dla takiej decyzji. Jeśli odsetki koncentrują się na wpływie wybranych wybranych poziomów, należy je uznać za ustalone. Jeśli dostępne poziomy czynników są postrzegane jako próba losowa z większej populacji i że wnioski są wymagane dla większej populacji, efekt powinien być losowy.

— ocram