Sparowany test t jako specjalny przypadek liniowego modelowania z mieszanym efektem

20

Wiemy, że sparowany test t jest tylko specjalnym przypadkiem jednokierunkowej ANOVA z powtarzanymi pomiarami (lub wewnątrz podmiotu), a także liniowym modelem mieszanego efektu, który można zademonstrować za pomocą funkcji lme () pakiet nlme w języku R jak pokazano niżej.

#response data from 10 subjects under two conditions
x1<-rnorm(10)
x2<-1+rnorm(10)

# Now create a dataframe for lme
myDat <- data.frame(c(x1,x2), c(rep("x1", 10), rep("x2", 10)), rep(paste("S", seq(1,10), sep=""), 2))
names(myDat) <- c("y", "x", "subj")

Kiedy uruchamiam następujący sparowany test t:

t.test(x1, x2, paired = TRUE)

Mam ten wynik (otrzymasz inny wynik z powodu losowego generatora):

t = -2.3056, df = 9, p-value = 0.04657

Dzięki podejściu ANOVA możemy uzyskać ten sam wynik:

summary(aov(y ~ x + Error(subj/x), myDat))

# the F-value below is just the square of the t-value from paired t-test:
          Df  F value Pr(>F)
x          1  5.3158  0.04657

Teraz mogę uzyskać ten sam wynik w lme przy użyciu następującego modelu, zakładając dodatnio określoną symetryczną macierz korelacji dla dwóch warunków:

summary(fm1 <- lme(y ~ x, random=list(subj=pdSymm(form=~x-1)), data=myDat))

# the 2nd row in the following agrees with the paired t-test
# (Intercept) -0.2488202 0.3142115  9 -0.7918878  0.4488
# xx2          1.3325786 0.5779727  9  2.3056084  0.0466

Lub inny model, zakładając złożoną symetrię dla macierzy korelacji dwóch warunków:

summary(fm2 <- lme(y ~ x, random=list(subj=pdCompSymm(form=~x-1)), data=myDat))

# the 2nd row in the following agrees with the paired t-test
# (Intercept) -0.2488202 0.4023431  9 -0.618428  0.5516
# xx2          1.3325786 0.5779727  9  2.305608  0.0466

Za pomocą sparowanego testu t i jednokierunkowej ANOVA z powtarzanymi pomiarami mogę zapisać tradycyjny model średnich komórek jako

Yij = μ + αi + βj + εij, i = 1, 2; j = 1, ..., 10

gdzie i indeksuje warunek, j indeksuje podmiot, Y ij jest zmienną odpowiedzi, μ jest stały dla ustalonego efektu dla średniej ogólnej, α i jest ustalonym efektem dla warunku, β j jest efektem losowym dla badanego po N (0, σ p 2 ) (σ p 2 to wariancja populacyjna), a ε ij jest resztkowe po N (0, σ 2 ) (σ 2 jest wariancją wewnątrzosobniczą).

Myślałem, że powyższy model średniej komórkowej nie byłby odpowiedni dla modeli lme, ale problem polega na tym, że nie mogę wymyślić rozsądnego modelu dla dwóch podejść lme () przy założeniu struktury korelacji. Powodem jest to, że model lme wydaje się mieć więcej parametrów dla składników losowych niż model średniej komórki powyżej oferuje. Przynajmniej model lme zapewnia dokładnie tę samą wartość F, stopnie swobody, a także wartość p, czego gls nie może. Mówiąc dokładniej, gls podaje nieprawidłowe DF ze względu na fakt, że nie uwzględnia faktu, że każdy badany ma dwie obserwacje, co prowadzi do znacznie zawyżonego DF. Model lme najprawdopodobniej jest sparametryzowany w określaniu losowych efektów, ale nie wiem, co to jest model i jakie są parametry. Tak więc problem nadal jest dla mnie nierozwiązany.

bluepole
źródło
2
Nie jestem pewien, o co pytasz. Model, który zapisałeś, jest dokładnie modelem modelu efektów losowych; struktura korelacji jest indukowana przez efekt losowy.
Aaron - Przywróć Monikę
@Aaron: losowy efekt βj w modelu średniej komórkowej powinien być zgodny z N (0, σp2). Moje zamieszanie polega na tym, w jaki sposób ten termin (z tylko jednym parametrem σp2) jest powiązany ze strukturą korelacji określoną przez symetrię złożoną lub prostą macierz symetryczną w modelu Lme?
bluepole
Kiedy obliczasz korelację między dwiema obserwacjami tego samego przedmiotu, korelacja jest sigma_p ^ 2 / (sigma_p ^ 2 + sigma ^ 2), ponieważ mają one tę samą beta_j. Patrz Pinheiro / Bates str.8. Ponadto model efektu losowego, który napisałeś, jest równoważny złożonej symetrii; inne struktury korelacji są bardziej złożone.
Aaron - Przywróć Monikę
@Aaron: Dzięki! Przeczytałem już o tym książkę Pinheiro / Bates i nadal nie mogę zrozumieć szczegółów losowych efektów. Bardziej trafne strony wydają się być przykładem z P.160-161. Również wyniki efektów losowych z lme () przy założeniu złożonej symetrii nie wydają się zgadzać z korelacją σp2 / (σp2 + σ2) w modelu średniej komórkowej. Wciąż mam wątpliwości co do struktury modelu.
bluepole
Cóż, prawie równoważne złożonej symetrii; w CS korelacja może być ujemna, ale nie z efektami losowymi. Być może na tym polega różnica. Szczegółowe informacje można znaleźć na stronie stats.stackexchange.com/a/14185/3601 .
Aaron - Przywróć Monikę

Odpowiedzi:

16

Równoważność modeli można zaobserwować, obliczając korelację między dwoma obserwacjami tej samej osoby, w następujący sposób:

Yjajot=μ+αja+βjot+ϵjajotβjotN.(0,σp2))ϵjajotN.(0,σ2))doov(yjak,yjotk)=doov(μ+αja+βk+ϵjak,μ+αjot+βk+ϵjotk)=doov(βk,βk)=σp2)V.zar(yjak)=V.zar(yjotk)=σp2)+σ2)σp2)/(σp2)+σ2))

Zauważ, że modele nie są jednak całkowicie równoważne, ponieważ model efektu losowego wymusza korelację na dodatnią. Model CS i model t-test / anova nie.

EDYCJA: Istnieją również dwie inne różnice. Po pierwsze, modele CS i losowe zakładają normalność efektu losowego, ale model t-test / anova nie. Po drugie, modele CS i losowe są dopasowane przy użyciu maksymalnego prawdopodobieństwa, podczas gdy anova jest dopasowana przy użyciu średnich kwadratów; kiedy wszystko jest zrównoważone, zgodzą się, ale niekoniecznie w bardziej złożonych sytuacjach. Wreszcie, będę ostrożny w stosowaniu wartości F / df / p z różnych dopasowań jako miary tego, jak bardzo modele się zgadzają; zobacz słynny stół Douga Batesa na df po więcej szczegółów. (EDYCJA KOŃCOWA)

Problem z twoim Rkodem polega na tym, że nie określasz poprawnie struktury korelacji. Musisz używać glsze corCompSymmstrukturą korelacji.

Wygeneruj dane, aby uzyskać efekt podmiotowy:

set.seed(5)
x <- rnorm(10)
x1<-x+rnorm(10)
x2<-x+1 + rnorm(10)
myDat <- data.frame(c(x1,x2), c(rep("x1", 10), rep("x2", 10)), 
                    rep(paste("S", seq(1,10), sep=""), 2))
names(myDat) <- c("y", "x", "subj")

Oto jak dopasujesz losowe efekty i złożone modele symetrii.

library(nlme)
fm1 <- lme(y ~ x, random=~1|subj, data=myDat)
fm2 <- gls(y ~ x, correlation=corCompSymm(form=~1|subj), data=myDat)

Standardowe błędy z modelu efektów losowych to:

m1.varp <- 0.5453527^2
m1.vare <- 1.084408^2

Korelacja i wariancja rezydualna z modelu CS jest następująca:

m2.rho <- 0.2018595
m2.var <- 1.213816^2

I są równe oczekiwanym:

> m1.varp/(m1.varp+m1.vare)
[1] 0.2018594
> sqrt(m1.varp + m1.vare)
[1] 1.213816

Inne struktury korelacji zwykle nie pasują do efektów losowych, ale po prostu przez określenie pożądanej struktury; jednym powszechnym wyjątkiem jest model efektu losowego AR (1) +, który ma efekt losowy i korelację AR (1) między obserwacjami tego samego efektu losowego.

EDYCJA 2: Kiedy dopasuję trzy opcje, otrzymuję dokładnie takie same wyniki, z tym wyjątkiem, że gls nie próbuje zgadnąć df dla interesującego nas terminu.

> summary(fm1)
...
Fixed effects: y ~ x 
                 Value Std.Error DF   t-value p-value
(Intercept) -0.5611156 0.3838423  9 -1.461839  0.1778
xx2          2.0772757 0.4849618  9  4.283380  0.0020

> summary(fm2)
...
                 Value Std.Error   t-value p-value
(Intercept) -0.5611156 0.3838423 -1.461839  0.1610
xx2          2.0772757 0.4849618  4.283380  0.0004

> m1 <- lm(y~ x + subj, data=myDat)
> summary(m1)
...
            Estimate Std. Error t value Pr(>|t|)   
(Intercept)  -0.3154     0.8042  -0.392  0.70403   
xx2           2.0773     0.4850   4.283  0.00204 **

(Punkt przecięcia jest tutaj inny, ponieważ przy domyślnym kodowaniu nie jest to średnia z wszystkich przedmiotów, ale średnia z pierwszego przedmiotu).

Warto również zauważyć, że nowszy lme4pakiet daje takie same wyniki, ale nawet nie próbuje obliczyć wartości p.

> mm1 <- lmer(y ~ x + (1|subj), data=myDat)
> summary(mm1)
...
            Estimate Std. Error t value
(Intercept)  -0.5611     0.3838  -1.462
xx2           2.0773     0.4850   4.283
Aaron - Przywróć Monikę
źródło
Jeszcze raz dziękuję za pomoc! Znam tę część z perspektywy modelu średniej komórki. Jednak z następującym wynikiem lme () ze złożoną symetrią: Efekty losowe: Wzór: ~ x - 1 | Subj Struktura: Compound Symmetry StdDev xx1 1,1913363 xx2 1,1913363 Corr: -0,036 Resztkowe 0,4466733. Nadal nie mogę pogodzić tych liczb z modelem średniej komórkowej. Może pomożesz mi uporządkować te liczby?
bluepole
A także jakieś przemyślenia na temat formułowania modelu z innymi strukturami korelacji, takimi jak prosta macierz symetryczna?
bluepole
Widzę! Powinienem był uważniej przeczytać twoją odpowiedź w drugim wątku. Myślałem wcześniej o użyciu gls (), ale nie udało mi się ustalić specyfikacji korelacji. Interesujące jest to, że lme () ze złożoną strukturą symetrii dla efektu losowego nadal renderuje tę samą wartość t, ale wydaje się, że wariancje efektów losowych nie są bezpośrednio interpretowalne. Doceniam twoją pomoc!
bluepole
Po drugiej przemyśleniu czuję, że moje pierwotne zamieszanie jest nadal nierozwiązane. Tak, gls można wykorzystać do zademonstrowania struktury korelacji i średnich kwadratów rumów, ale model pod nią nie jest dokładnie taki sam jak test sparowany t (lub ogólnie ANOVA z powtarzanymi pomiarami w jednym kierunku) i taka ocena jest dalej wspierane przez nieprawidłowe DF i wartość p z gls. W przeciwieństwie do tego, moje polecenie lme z symetrią złożoną zapewnia te same wartości F, DF i p. Jedyne, co mnie zastanawia, to sposób parametryzacji modelu lme, jak podano w moim oryginalnym poście. Czy jest tam jakaś pomoc?
bluepole
Nie jestem pewien, jak ci pomóc. Czy możesz napisać, co według ciebie są dwoma różnymi modelami? Coś jest nie tak z tym, jak myślisz o jednym z nich.
Aaron - Przywróć Monikę
3

Możesz również rozważyć użycie funkcji mixedw pakiecie afexdo zwrócenia wartości p z aproksymacją df Kenwarda-Rogera, która zwraca identyczne wartości p jako sparowany test t:

library(afex)
mixed(y ~ x + (1|subj), type=3,method="KR",data=myDat) 

Lub

library(lmerTest)
options(contrasts=c('contr.sum', 'contr.poly'))
anova(lmer(y ~ x + (1|subj),data=myDat),ddf="Kenward-Roger")
Tom Wenseleers
źródło