Rozumiem, że używamy modeli efektów losowych (lub efektów mieszanych), gdy uważamy, że niektóre parametry modelu zmieniają się losowo w zależności od czynnika grupującego. Chcę dopasować model, w którym odpowiedź została znormalizowana i wyśrodkowana (nie idealnie, ale całkiem blisko) w obrębie czynnika grupującego, ale zmienna niezależna x
nie została w żaden sposób dostosowana. Doprowadziło mnie to do następującego testu (przy użyciu sfabrykowanych danych), aby upewnić się, że znajdę efekt, którego szukałem, jeśli rzeczywiście tam jest. Uruchomiłem jeden model efektów mieszanych z przypadkowym przechwytywaniem (pomiędzy grupami zdefiniowanymi przez f
) i drugi model efektu stałego z czynnikiem f jako predyktorem efektu stałego. Użyłem pakietu R lmer
dla modelu efektu mieszanego i funkcji podstawowejlm()
dla modelu z efektem stałym. Poniżej znajdują się dane i wyniki.
Zauważ, że y
niezależnie od grupy, waha się w okolicach 0. I x
różni się w zależności od y
grupy, ale różni się znacznie w zależności od grupy niży
> data
y x f
1 -0.5 2 1
2 0.0 3 1
3 0.5 4 1
4 -0.6 -4 2
5 0.0 -3 2
6 0.6 -2 2
7 -0.2 13 3
8 0.1 14 3
9 0.4 15 3
10 -0.5 -15 4
11 -0.1 -14 4
12 0.4 -13 4
Jeśli jesteś zainteresowany pracą z danymi, oto dput()
wynik:
data<-structure(list(y = c(-0.5, 0, 0.5, -0.6, 0, 0.6, -0.2, 0.1, 0.4,
-0.5, -0.1, 0.4), x = c(2, 3, 4, -4, -3, -2, 13, 14, 15, -15,
-14, -13), f = structure(c(1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L,
4L, 4L, 4L), .Label = c("1", "2", "3", "4"), class = "factor")),
.Names = c("y","x","f"), row.names = c(NA, -12L), class = "data.frame")
Dopasowanie modelu efektów mieszanych:
> summary(lmer(y~ x + (1|f),data=data))
Linear mixed model fit by REML
Formula: y ~ x + (1 | f)
Data: data
AIC BIC logLik deviance REMLdev
28.59 30.53 -10.3 11 20.59
Random effects:
Groups Name Variance Std.Dev.
f (Intercept) 0.00000 0.00000
Residual 0.17567 0.41913
Number of obs: 12, groups: f, 4
Fixed effects:
Estimate Std. Error t value
(Intercept) 0.008333 0.120992 0.069
x 0.008643 0.011912 0.726
Correlation of Fixed Effects:
(Intr)
x 0.000
Zwracam uwagę, że komponent wariancji przechwytywania jest szacowany na 0 i, co dla mnie ważne, x
nie jest znaczącym predyktorem y
.
Następnie dopasowuję model z efektem stałym f
jako predyktor zamiast czynnika grupującego dla przypadkowego przechwytywania:
> summary(lm(y~ x + f,data=data))
Call:
lm(formula = y ~ x + f, data = data)
Residuals:
Min 1Q Median 3Q Max
-0.16250 -0.03438 0.00000 0.03125 0.16250
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -1.38750 0.14099 -9.841 2.38e-05 ***
x 0.46250 0.04128 11.205 1.01e-05 ***
f2 2.77500 0.26538 10.457 1.59e-05 ***
f3 -4.98750 0.46396 -10.750 1.33e-05 ***
f4 7.79583 0.70817 11.008 1.13e-05 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.1168 on 7 degrees of freedom
Multiple R-squared: 0.9484, Adjusted R-squared: 0.9189
F-statistic: 32.16 on 4 and 7 DF, p-value: 0.0001348
Teraz zauważam, że zgodnie z oczekiwaniami x
jest znaczącym predyktorem y
.
To, czego szukam, to intuicja dotycząca tej różnicy. W jaki sposób błędne jest moje myślenie? Dlaczego niepoprawnie oczekuję znalezienia znaczącego parametru dla x
obu tych modeli, ale faktycznie widzę go tylko w modelu ze stałym efektem?
źródło
x
zmienna nie jest znacząca. Podejrzewam, że jest to ten sam wynik (współczynniki i SE), który można uruchomićlm(y~x,data=data)
. Nie mam więcej czasu na diagnozę, ale chciałem to podkreślić.Odpowiedzi:
Tutaj dzieje się kilka rzeczy. To interesujące kwestie, ale wyjaśnienie tego zajmie sporo czasu / miejsca.
Po pierwsze, wszystko to staje się o wiele łatwiejsze do zrozumienia, jeśli wykreślimy dane . Oto wykres punktowy, w którym punkty danych są pokolorowane według grup. Dodatkowo dla każdej grupy mamy osobną linię regresji specyficzną dla grupy, a także prostą linię regresji (ignorowanie grup) pogrubioną czcionką:
Model z efektem stałym
lm()
Model mieszany
Oto współczynniki dla prostego modelu regresji (przerywana pogrubiona linia na wykresie):
Jak widać, współczynniki tutaj są identyczne z uzyskanymi w modelu mieszanym. Dokładnie tego się spodziewaliśmy, ponieważ, jak już zauważyłeś, mamy oszacowanie 0 wariancji dla losowych przechwyceń, co czyni wspomnianą wcześniej korelację stosunku / wewnątrz klasy 0. Więc szacunki modelu mieszanego w tym przypadku są tylko proste oszacowania regresji liniowej, a jak widać na wykresie, nachylenie tutaj jest znacznie mniej wyraźne niż nachylenia wewnątrz gromady.
To prowadzi nas do jednego końcowego zagadnienia koncepcyjnego ...
Dlaczego wariancja losowych przechwyceń jest szacowana na 0?
Odpowiedź na to pytanie może stać się trochę techniczna i trudna, ale postaram się, aby była tak prosta i nietechniczna, jak to tylko możliwe (ze względu na nas samych!). Ale może nadal będzie trochę rozwlekły.
Rozważany przez nas model mieszany nie korzysta z wewnątrzklasowej metody korelacji reprezentującej zależność w danych. Zamiast tego opisuje zależność w kategoriach składników wariancji . Wszystko jest w porządku, dopóki korelacja wewnątrz klasy jest dodatnia. W takich przypadkach korelację wewnątrz klasy można łatwo napisać w kategoriach składników wariancji, szczególnie jako wcześniej wspomniany stosunek losowej wariancji przechwytywania do całkowitej wariancji. (Zobacz stronę wiki dotyczącą korelacji wewnątrzklasowejwięcej informacji na ten temat.) Niestety modele składowe wariancji mają trudności z radzeniem sobie z sytuacjami, w których mamy ujemną korelację wewnątrz klasy. W końcu pisanie korelacji wewnątrzklasowej pod względem składników wariancji wymaga zapisania jej jako proporcji wariancji, a proporcje nie mogą być ujemne.
Więc co możemy zrobić?
Aby to zrobić, bierzemy naszex xb x xw x
Wreszcie, nadal mamy oszacowanie 0 dla wariancji przypadkowych przechwyceń, z powodów, które rozwinąłem w poprzedniej sekcji. Nie jestem do końca pewien, co możemy z tym zrobić, przynajmniej bez przełączania się na oprogramowanie inne niż
lmer()
, i nie jestem też pewien, w jakim stopniu będzie to miało negatywny wpływ na nasze szacunki w tym ostatecznym mieszanym modelu. Być może inny użytkownik może wtrącić się z przemyśleniami na temat tego problemu.Bibliografia
źródło
lme
domyślnie ograniczała się, wynosiła> = 0? Zobacz to pytanie i wybraną odpowiedź , tj. Dopasowanie korelacji złożonej symetrii poprzezgls
dopasowanie lub ustawieniecorrelation = corCompSymm(form = ~1|f)
wlme
Po dłuższej kontemplacji wydaje mi się, że odkryłem własną odpowiedź. Wierzę, że ekonometryczny zdefiniowałby moją zmienną niezależną jako endogenną, a zatem byłby skorelowany zarówno ze zmiennymi niezależnymi, jak i zależnymi. W takim przypadku zmienne te są pomijane lub nie są obserwowane . Obserwuję jednak grupy, między którymi zmienna pominięta powinna się zmieniać.
Wierzę, że ekonometryczny sugerowałby model z efektem stałym . Oznacza to, że model zawiera atrapę dla każdego poziomu grupowania (lub równoważną specyfikację, która warunkuje model tak, że wiele manekinów grupujących nie jest wymaganych) w tym przypadku. Dzięki modelowi o stałym efekcie można mieć nadzieję, że wszystkie nieobserwowane i niezmienne w czasie zmienne mogą być kontrolowane przez warunkowanie w obrębie grupy (lub pomiędzy poszczególnymi) odmianami. Rzeczywiście, drugi model w moim pytaniu jest modelem o ustalonym efekcie i jako taki daje oszacowanie, którego się spodziewam.
Z zadowoleniem przyjmuję komentarze, które dodatkowo wyjaśnią tę okoliczność.
źródło