Wykorzystanie lmera do liniowego modelu mieszanego z powtarzanymi pomiarami

EDYCJA 2: Początkowo myślałem, że muszę uruchomić ANOVA dwuskładnikową z powtarzanymi pomiarami dla jednego czynnika, ale teraz myślę, że liniowy model mieszanego efektu będzie działał lepiej dla moich danych. Myślę, że prawie wiem, co musi się wydarzyć, ale wciąż jestem zdezorientowany kilkoma punktami.

Eksperymenty, które muszę przeanalizować, wyglądają tak:

• Osobników przydzielono do jednej z kilku grup terapeutycznych
• Pomiary każdego pacjenta wykonywano przez wiele dni
• Więc:
  • Obiekt jest zagnieżdżony w trakcie leczenia
  • Leczenie krzyżuje się z dniem

(każdy pacjent jest przypisany tylko do jednego zabiegu, a pomiary są wykonywane dla każdego pacjenta każdego dnia)

Mój zestaw danych zawiera następujące informacje:

• Temat = czynnik blokujący (czynnik losowy)
• Dzień = w obrębie czynnika lub współczynnika powtarzanych pomiarów (współczynnik stały)
• Leczenie = między czynnikiem podmiotowym (czynnik stały)
• Obs = zmienna zmierzona (zależna)

AKTUALIZACJA OK, więc poszedłem porozmawiać ze statystyką, ale on jest użytkownikiem SAS. Uważa, że ​​modelem powinien być:

Leczenie + dzień + pacjent (leczenie) + dzień * pacjent (leczenie)

Oczywiście jego notacja różni się od składni R, ale ten model powinien uwzględniać:

• Leczenie (naprawione)
• Dzień (ustalony)
• interakcja Leczenie * Dzień
• Obiekt zagnieżdżony w ramach leczenia (losowo)
• Dzień skrzyżowany z „Podmiotem w trakcie leczenia” (losowo)

Czy to jest poprawna składnia?

m4 <- lmer(Obs~Treatment*Day + (1+Treatment/Subject) + (1+Day*Treatment/Subject), mydata)

Jestem szczególnie zaniepokojony tym, czy Dzień skrzyżowany z częścią „Podmiot w trakcie leczenia” ma rację. Czy ktoś zna SAS lub jest przekonany, że rozumie, co dzieje się w jego modelu, jest w stanie skomentować, czy moja smutna próba dopasowania składni R jest zgodna?

Oto moje poprzednie próby zbudowania modelu i pisania składni (omówione w odpowiedziach i komentarzach):

m1 <- lmer(Obs ~ Treatment * Day + (1 | Subject), mydata)

Jak poradzić sobie z faktem, że pacjent jest zagnieżdżony w trakcie leczenia? Czym m1różni się od:

m2 <- lmer(Obs ~ Treatment * Day + (Treatment|Subject), mydata)
m3 <- lmer(Obs ~ Treatment * Day + (Treatment:Subject), mydata)

a są m2i m3ekwiwalent (a jeśli nie, to dlaczego)?

Czy muszę też używać nlme zamiast lme4, jeśli chcę określić strukturę korelacji (jak correlation = corAR1)? Zgodnie z Powtarzanymi pomiarami, dla analizy powtarzanych pomiarów z powtarzanymi pomiarami jednego czynnika, ważna jest struktura kowariancji (charakter korelacji między pomiarami tego samego pacjenta).

Kiedy próbowałem wykonać ANOVA z powtarzanymi pomiarami, zdecydowałem się użyć SS typu II; czy jest to nadal aktualne, a jeśli tak, to jak mam to określić?

Oto przykład tego, jak wyglądają dane:

mydata <- data.frame(
  Subject  = c(13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 29, 30, 31, 32, 33, 
               34, 35, 36, 37, 38, 39, 40, 62, 63, 64, 65, 13, 14, 15, 16, 17, 18, 
               19, 20, 21, 22, 23, 24, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 
               40, 62, 63, 64, 65, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 
               29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 62, 63, 64, 65), 
  Day       = c(rep(c("Day1", "Day3", "Day6"), each=28)), 
  Treatment = c(rep(c("B", "A", "C", "B", "C", "A", "A", "B", "A", "C", "B", "C", 
                      "A", "A", "B", "A", "C", "B", "C", "A", "A"), each = 4)), 
  Obs       = c(6.472687, 7.017110, 6.200715, 6.613928, 6.829968, 7.387583, 7.367293, 
                8.018853, 7.527408, 6.746739, 7.296910, 6.983360, 6.816621, 6.571689, 
                5.911261, 6.954988, 7.624122, 7.669865, 7.676225, 7.263593, 7.704737, 
                7.328716, 7.295610, 5.964180, 6.880814, 6.926342, 6.926342, 7.562293, 
                6.677607, 7.023526, 6.441864, 7.020875, 7.478931, 7.495336, 7.427709, 
                7.633020, 7.382091, 7.359731, 7.285889, 7.496863, 6.632403, 6.171196, 
                6.306012, 7.253833, 7.594852, 6.915225, 7.220147, 7.298227, 7.573612, 
                7.366550, 7.560513, 7.289078, 7.287802, 7.155336, 7.394452, 7.465383, 
                6.976048, 7.222966, 6.584153, 7.013223, 7.569905, 7.459185, 7.504068, 
                7.801867, 7.598728, 7.475841, 7.511873, 7.518384, 6.618589, 5.854754, 
                6.125749, 6.962720, 7.540600, 7.379861, 7.344189, 7.362815, 7.805802, 
                7.764172, 7.789844, 7.616437, NA, NA, NA, NA))

Myślę, że twoje podejście jest prawidłowe. Model m1określa osobny punkt przecięcia dla każdego przedmiotu. Model m2dodaje osobne nachylenie dla każdego obiektu. Twoje nachylenie trwa przez kilka dni, ponieważ badani uczestniczą tylko w jednej grupie leczenia. Jeśli napiszesz model m2w następujący sposób, bardziej oczywiste jest, że modelujesz osobny punkt przecięcia i nachylenie dla każdego przedmiotu

m2 <- lmer(Obs ~ Treatment * Day + (1+Day|Subject), mydata)

Jest to równoważne z:

m2 <- lmer(Obs ~ Treatment + Day + Treatment:Day + (1+Day|Subject), mydata)

To znaczy główne efekty leczenia, dzień i interakcja między nimi.

Myślę, że nie musisz się martwić o zagnieżdżanie, dopóki nie powtórzysz identyfikatorów uczestników w grupach leczenia. Który model jest poprawny, tak naprawdę zależy od twojego pytania badawczego. Czy istnieje powód, by sądzić, że nachylenie badanych różni się oprócz efektu leczenia? Możesz uruchomić oba modele i porównać je, anova(m1,m2)aby sprawdzić, czy dane obsługują jeden z nich.

Nie jestem pewien, co chcesz wyrazić za pomocą modelu m3? Składnia zagnieżdżania używa /np (1|group/subgroup).

Nie sądzę, że musisz martwić się autokorelacją z tak małą liczbą punktów czasowych.

Nie czuję się wystarczająco komfortowo, aby komentować problem z autokorelowanymi błędami (ani różne implementacje w lme4 vs. nlme), ale mogę porozmawiać z resztą.

Twój model m1jest modelem z losowym przechwytywaniem, w którym uwzględniono interakcję między poziomami między leczeniem a dniem (wpływ dnia może się różnić w zależności od grupy leczenia). Aby umożliwić różnicę zmian w czasie między uczestnikami (tj. Wyraźne modelowanie indywidualnych różnic w zmianach w czasie), należy również pozwolić, aby efekt Dnia był losowy . Aby to zrobić, określ:

m2 <- lmer(Obs ~ Day + Treatment + Day:Treatment + (Day | Subject), mydata)

W tym modelu:

• Punkt przecięcia, jeśli przewidywany wynik dla kategorii odniesienia leczenia w Dniu = 0
• Współczynnik dla dnia jest przewidywaną zmianą w czasie dla każdego 1-jednostkowego wzrostu liczby dni dla kategorii odniesienia leczenia
• Współczynniki dla dwóch kodów fikcyjnych dla grup terapeutycznych (automatycznie tworzone przez R) są przewidywaną różnicą między każdą pozostałą grupą leczoną a kategorią odniesienia w dniu = 0
• Współczynniki dla dwóch terminów interakcji są różnicą wpływu czasu (dzień) na przewidywane wyniki między kategorią odniesienia a pozostałymi grupami leczenia.

Zarówno przechwyty, jak i wpływ Dnia na wynik są losowe (każdy osobnik może mieć inny przewidywany wynik w Dniu = 0 i inną liniową zmianę w czasie). Modelowana jest również kowariancja między przechwytywaniem a nachyleniem (dozwolone jest tworzenie kowariancji).

Jak widać, interpretacja współczynników dla dwóch zmiennych obojętnych jest uzależniona od Dnia = 0. Powiedzą ci, czy przewidywany wynik w dniu = 0 dla kategorii referencyjnej różni się znacząco od dwóch pozostałych grup leczenia. Dlatego ważne jest, gdzie zdecydujesz się wyśrodkować zmienną Dzień. Jeśli skoncentrujesz się na Dniu 1, wówczas współczynniki powiedzą ci, czy przewidywany wynik dla kategorii referencyjnej w Dniu 1 różni się znacząco od przewidywanego wyniku dwóch pozostałych grup. W ten sposób możesz sprawdzić, czy istnieją wcześniej istniejące różnice między grupami . Jeśli skupisz się na Dniu 3, wówczas współczynniki powiedzą ci, czy przewidywany wynik dla kategorii odniesienia w Dniu 3różni się znacząco od przewidywanego wyniku dwóch pozostałych grup. W ten sposób można sprawdzić, czy pod koniec interwencji występują różnice między grupami .

Na koniec zauważ, że pacjenci nie są zagnieżdżeni w ramach leczenia. Twoje trzy zabiegi nie są przypadkowymi poziomami populacji poziomów, do których chcesz uogólniać swoje wyniki - raczej, jak wspomniałeś, twoje poziomy są stałe i chcesz uogólniać swoje wyniki tylko na te poziomy. (Nie wspominając, nie powinieneś używać modelowania wielopoziomowego, jeśli masz tylko 3 jednostki wyższego poziomu; patrz Maas i Hox, 2005.) Zamiast tego leczenie jest predyktorem poziomu 2, tj. Predyktorem, który przyjmuje pojedynczą wartość w dniach (jednostki poziomu 1) dla każdego przedmiotu. Dlatego jest on tylko uwzględniony jako predyktor w twoim modelu.

Odniesienie:
Maas, CJM i Hox, JJ (2005). Wystarczające rozmiary próbek do modelowania wielopoziomowego. Metodologia: European Journal of Research Methods for the Behavioural and Social Sciences , 1 , 86-92.