Dlaczego wartości p zmieniają się w znaczeniu przy zmianie kolejności zmiennych towarzyszących w modelu aov?

10

Mam zestaw danych z 482 obserwacji.

data=Populationfull

Zamierzam przeprowadzić analizę asocjacji genotypu dla 3 SNP. Próbuję zbudować model do mojej analizy i używam aov (y ~ x, data = ...). Dla jednej cechy mam kilka ustalonych efektów i zmiennych towarzyszących, które zawarłem w modelu, na przykład:

Starts <- aov(Starts~Sex+DMRT3+Birthyear+Country+Earnings+Voltsec+Autosec, data=Populationfull)

summary(Starts)
                Df Sum Sq Mean Sq F value   Pr(>F)    
Sex              3  17.90    5.97  42.844  < 2e-16 ***
DMRT3            2   1.14    0.57   4.110    0.017 *  
Birthyear        9   5.59    0.62   4.461 1.26e-05 ***
Country          1  11.28   11.28  81.005  < 2e-16 ***
Earnings         1 109.01  109.01 782.838  < 2e-16 ***
Voltsec          1  12.27   12.27  88.086  < 2e-16 ***
Autosec          1   8.97    8.97  64.443 8.27e-15 ***
Residuals      463  64.48    0.14                     
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1   1

Odkryłem, że jeśli zmieniłem kolejność zmiennych w modelu, otrzymałem różne wartości p, patrz poniżej.

Starts2 <- aov(Starts~Voltsec+Autosec+Sex+DMRT3+Birthyear+Country+Earnings, data=Populationfull)

summary(Starts2)

                Df Sum Sq Mean Sq F value   Pr(>F)    
Voltsec   1   2.18    2.18  15.627 8.92e-05 ***
Autosec   1 100.60  100.60 722.443  < 2e-16 ***
Sex              3  10.43    3.48  24.962 5.50e-15 ***
DMRT3            2   0.82    0.41   2.957  0.05294 .  
Birthyear        9   3.25    0.36   2.591  0.00638 ** 
Country          1   2.25    2.25  16.183 6.72e-05 ***
Earnings      1  46.64   46.64 334.903  < 2e-16 ***
Residuals      463  64.48    0.14                     
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1   1

Dlaczego otrzymuję różne wartości p w zależności od tego, w jakiej kolejności kodowane są zmienne / czynniki / zmienne towarzyszące / ustalone efekty (?)? Czy istnieje sposób, aby to „poprawić”? Czy to możliwe, że używam niewłaściwego modelu? Nadal jestem całkiem nowy w R, więc jeśli możesz mi w tym pomóc, proszę, zastanów się naprawdę, abym mógł zrozumieć odpowiedź hehe ... Dziękuję, mam nadzieję, że ktoś może mi pomóc to zrozumieć!

Rbeginner
źródło
3
Podaj przykładowe dane, Populationfullaby problem mógł być powtarzalny . Nie dzieje się tak w przypadku przykładu ze strony aov()pomocy. summary(aov(yield ~ block + N + P + K, npk)); summary(aov(yield ~ K + P + block + N , npk))
MrFlick
Wartości p zmieniają się, ponieważ zmieniło się całe pole wartości. twój pierwszy bieg Earnings 1 109.01 109.01 782.838 < 2e-16 ***twój drugi bieg Earnings 1 46.64 46.64 334.903 < 2e-16 ***. Twoje wyniki nie są takie same. Zacznij od sprawdzenia, czy nie zrobiłeś nic więcej niż zmiana kolejności zmiennych.
1
RÓWNIEŻ: w drugim modelu używasz Zarabiaj, a nie Zarabiaj ... jeśli istnieją dwie zmienne o różnych nazwach, może to być problem, jeśli nie jest to literówka podczas kopiowania do przestrzeni pytań SO.
Tak, wartości się zmieniają, ale dlaczego? Użyłem dokładnie tych samych kolumn z dokładnie tej samej ramki danych w obu modelach (Zysk vs. Zarobek w drugim modelu polega na tym, że źle napisałem, poprawiłem go teraz).
Rbeginner
1
Dzieje się tak, ponieważ masz niezrównoważony projekt. Znajdziesz dużo pomocy na ten temat, jeśli przeszukujesz to forum lub po prostu Google „niezrównoważona ANOVA w R”. Polecam carzajrzeć do pakietu - implementuje on ANOVA typu II i typu III, które nie zależą od kolejności zmiennych, natomiast aovANOVA typu I.
Slow loris

Odpowiedzi:

15

Problem pochodzi ze sposobu, w jaki aov()wykonuje domyślne testowanie istotności. Wykorzystuje analizę ANOVA o nazwie „Typ I”, w której testowanie odbywa się w kolejności, w której zmienne są określone w modelu. Tak więc w pierwszym przykładzie określa, ile wariancji jest wyjaśnione sexi testuje jej znaczenie, a następnie jaka część pozostałej wariancji jest wyjaśniona DMRT3i testuje jej znaczenie pod względem tej pozostałej wariancji i tak dalej. W drugim przykładzie DMRT3jest oceniany tylko po Voltsec, Autoseci sexw tej kolejności, więc pozostało mniej wariancji DMRT3do wyjaśnienia.

Jeśli dwie zmienne predykcyjne są skorelowane, pierwsza wprowadzona do modelu uzyska pełny „kredyt”, pozostawiając mniejszą wariancję, którą można „wyjaśnić” drugą, co może wydawać się mniej „statystycznie znaczące” niż pierwsza, nawet jeśli jest nie funkcjonalnie. To pytanie i jego odpowiedź wyjaśniają różne typy analiz ANOVA.

Jednym ze sposobów, aby ominąć to, aby wyodrębnić się od zwężeń klasycznej analizy wariancji i używa prostego modelu liniowego, ze lm()w R, zamiast aov(). To skutecznie analizuje wszystkie predyktory równolegle, „korygując” wszystkie predyktory jednocześnie. W takim przypadku dwa skorelowane predyktory mogą mieć duże błędy standardowe szacowanych współczynników regresji, a ich współczynniki mogą różnić się między różnymi próbkami z populacji, ale przynajmniej kolejność wprowadzania zmiennych do specyfikacji modelu nie będzie miała znaczenia.

Jeśli twoja zmienna odpowiedzi jest pewnego rodzaju zmienną zliczającą, jak Startssugeruje jej nazwa , prawdopodobnie i tak nie powinieneś używać ANOVA, ponieważ jest mało prawdopodobne, aby reszty były normalnie rozłożone, jak wymaga interpretacja wartości p . Zmienne zliczania są lepiej obsługiwane przy uogólnionych modelach liniowych (np. glm()W R), które można traktować jako uogólnienie lm()dla innych rodzajów struktur błędów resztkowych.

EdM
źródło