ANOVA v R: Jednocestný a obousměrný test s příklady

Co je ANOVA?

Analýza rozptylu (ANOVA) je statistická technika používaná k porovnání průměrů dvou nebo více skupin. Test funguje tak, že rozdělí celkovou variabilitu v měření na část vysvětlenou příslušností ke skupině a část, která zbývá jako náhodný šum. ANOVA v R vám tedy říká, zda se alespoň jeden průměr skupiny liší od ostatních, nikoli který. Je to přímé rozšíření t-test v situacích, kdy má faktorová proměnná více než dvě úrovně.

Před spuštěním testu je užitečné vědět, který člen rodiny ANOVA odpovídá vašemu návrhu.

Typy ANOVA testů v R

„ANOVA“ je spíše skupina testů než jeden postup. Výběr správného člena závisí na tom, kolik faktorů máte a jak byla data shromážděna.

Tento tutoriál se zabývá prvními třemi. Zbývající varianty používají stejné rozhraní aov(), takže jakmile načtete jednu výstupní tabulku, můžete načíst všechny.

ANOVA vs. T-test v R: Klíčové rozdíly

Oba testy porovnávají průměry, takže je vhodné být přesný v tom, kde jeden nahrazuje druhý.

U tří skupin je pokušení provést tři samostatné t-testy. Odolejte mu. Každý test má svou vlastní 5% chybovost, takže tři srovnání zvyšují pravděpodobnost falešně pozitivního výsledku na zhruba 14 procent. ANOVA odpovídá na stejnou otázku jediným testem a Tukey HSD poté zpracovává párové detaily s kontrolovanou chybovostí. Pro případ dvou skupin viz tutoriál k t-testu.

Jednosměrná ANOVA

Existuje mnoho situací, kdy potřebujete porovnat průměr mezi více skupinami. Marketingové oddělení chce například vědět, zda mají tři týmy stejný prodejní výkon.

• Tým: Faktor 3 úrovní: A, B a C
• Prodej: Měřítkem výkonu

Test ANOVA může zjistit, zda mají tyto tři skupiny podobné výkony.

Chcete-li objasnit, zda data pocházejí ze stejné populace, můžete provést a jednosměrná analýza rozptylu (dále jen jednofaktorová ANOVA). Stejně jako jakýkoli jiný statistický test poskytuje důkaz o tom, zda lze hypotézu H0 zamítnout. Všimněte si, že nezamítnutí H0 neznamená prokázání její pravdivosti.

Hypotéza v jednocestném testu ANOVA

• H0: Průměry mezi skupinami jsou stejné
• H1: Alespoň průměr jedné skupiny je jiný

Jinými slovy, pokud H0 nezamítneme, znamená to, že neexistuje dostatek důkazů k závěru, že se jakýkoli skupinový průměr liší od ostatních.

Tento test je podobný t-testu, ale ANOVA je správnou volbou, pokud existuje více než dvě skupiny. Pokud jsou skupiny přesně dvě, jsou oba testy ekvivalentní a F-statistika se rovná druhé mocnině t-statistiky.

Předpoklady

Jednocestná ANOVA spočívá na třech podmínkách: pozorování jsou náhodně vzorkována a na sobě nezávislá, rezidua v každé skupině jsou přibližně normálně rozdělena a rozptyl je v každé skupině stejný (homogenita rozptylu). Níže uvedená část o kontrole předpokladů ukazuje, jak každý z nich otestovat v R.

Interpretujte test ANOVA

F-statistika se používá k testování, zda data pocházejí z významně odlišných populací, tj. různých průměrů vzorků.

Chcete-li vypočítat F-statistiku, musíte vydělit meziskupinová variabilita přes variabilita v rámci skupiny.

Jedno meziskupina Variabilita vyjadřuje, jak daleko se průměr každé skupiny liší od celkového průměru. Porovnejte dva grafy níže, abyste si to lépe představili.

Levý graf ukazuje velmi malou variabilitu mezi třemi skupinami, takže průměry všech tří skupin se blíží celkový znamenat.

Pravý graf zobrazuje tři rozdělení daleko od sebe bez překrytí, takže rozdíl mezi celkovým průměrem a průměrem každé skupiny je velký.

Jedno v rámci skupiny Variabilita měří, jak daleko se jednotlivá pozorování odchylují od průměru vlastní skupiny. Některé body se nacházejí daleko od průměru skupiny a vnitroskupinový člen zachycuje přesně toto rozptylování, což je výběrová chyba.

Chcete-li vizuálně porozumět konceptu variability v rámci skupiny, podívejte se na níže uvedený graf.

Levá část znázorňuje rozdělení tří různých skupin. Zvětšili jste rozptyl každého vzorku a je zřejmé, že individuální rozptyl je velký. F-statistika klesá, takže byste nulovou hypotézu nezamítli.

Pravá část ukazuje vzorky se stejnými průměry, ale s mnohem menším rozptylem. To zvyšuje F-statistiku a svědčí ve prospěch alternativní hypotézy.

K vytvoření F-statistiky můžete použít obě měření. Je velmi intuitivní pochopit F-statistiku. Pokud se čitatel zvýší, znamená to, že variabilita mezi skupinami je vysoká a je pravděpodobné, že skupiny ve vzorku pocházejí ze zcela odlišných distribucí.

Jinými slovy, nízká F-statistika naznačuje malý nebo žádný významný rozdíl mezi průměry skupin.

Příklad Jednosměrný test ANOVA

K implementaci jednosměrného testu ANOVA použijete sadu dat jedu. Soubor dat obsahuje 48 řádků a 3 proměnné:

• Čas: Doba přežití zvířete
• jed: Typ použitého jedu: úroveň faktoru: 1,2 a 3
• léčba: Typ použité léčby: úroveň faktoru: 1,2 a 3

Než začnete počítat test ANOVA, musíte připravit data následovně:

• Krok 1: Importujte data
• Krok 2: Odstraňte nepotřebnou proměnnou
• Krok 3: Převeďte variabilní jed na uspořádanou úroveň

library(dplyr)
PATH <- "https://raw.githubusercontent.com/guru99-edu/R-Programming/master/poisons.csv"
df <- read.csv(PATH) %>%
select(-X) %>% 
mutate(poison = factor(poison, ordered = TRUE))
glimpse(df)

Výstup:

## Observations: 48
## Variables: 3
## $ time   <dbl> 0.31, 0.45, 0.46, 0.43, 0.36, 0.29, 0.40, 0.23, 0.22, 0...
## $ poison <ord> 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3, 1, 1, 1, 1, 2, 2, 2...
## $ treat  <fctr> A, A, A, A, A, A, A, A, A, A, A, A, B, B, B, B, B, B, ...		

Naším cílem je otestovat následující předpoklad:

• H0: Mezi skupinami není žádný rozdíl v průměrné době přežití
• H1: Průměrná doba přežití se alespoň u jedné skupiny liší.

Jinými slovy, chcete vědět, zda existuje statistický rozdíl mezi průměrnou dobou přežití podle typu jedu podaného morčeti.

Budete postupovat následovně:

• Krok 1: Zkontrolujte formát proměnné jed
• Krok 2: Vytiskněte souhrnnou statistiku: počet, průměr a standardní odchylku
• Krok 3: Nakreslete krabicový graf
• Krok 4: Vypočítejte jednosměrný test ANOVA
• Krok 5: Spusťte párové porovnání s Tukey HSD

Krok 1) Zkontrolujte úrovně jedu pomocí níže uvedeného kódu. Měli byste vidět třímístné hodnoty, protože příkaz mutate převedl sloupec na seřazený faktor.

levels(df$poison)

Výstup:

## [1] "1" "2" "3"

Krok 2) Vypočítáte průměr a směrodatnou odchylku.

df %>%
	group_by(poison) %>%
	summarise(
		count_poison = n(),
		mean_time = mean(time, na.rm = TRUE),
		sd_time = sd(time, na.rm = TRUE)
	)

Výstup:

## 
# A tibble: 3 x 4
##   poison count_poison mean_time    sd_time
##    <ord>        <int>     <dbl>      <dbl>
## 1      1           16  0.617500 0.20942779
## 2      2           16  0.544375 0.28936641
## 3      3           16  0.276250 0.06227627

Krok 3) Ve třetím kroku můžete graficky zkontrolovat, zda existuje rozdíl mezi distribucí. Všimněte si, že zahrnete roztřesenou tečku.

ggplot(df, aes(x = poison, y = time, fill = poison)) +
    geom_boxplot() +
    geom_jitter(shape = 15,
        color = "steelblue",
        position = position_jitter(0.21)) +
    theme_classic()

Výstup:

Krok 4) Jednosměrný test ANOVA můžete spustit příkazem aov. Základní syntaxe testu ANOVA je:

aov(formula, data)
Arguments:			
- formula: The equation you want to estimate
- data: The dataset used	

Syntaxe vzorce je:

y ~ X1+ X2+...+Xn # X1 +  X2 +... refers to the independent variables
y ~ . # use all the remaining variables as independent variables

Nyní můžete odpovědět na otázku: Existuje nějaký rozdíl v době přežití mezi morčaty, vzhledem k typu podaného jedu?

Uložte model do objektu a předejte ho funkci summary(), abyste získali čitelný výpis výsledků.

anova_one_way <- aov(time~poison, data = df)
summary(anova_one_way)

Code Vysvětlení

• aov(čas ~ jed, data = df): Spusťte test ANOVA s následujícím vzorcem
• Summary(anova_one_way): Vytiskne shrnutí testu

Výstup:

##             Df Sum Sq Mean Sq F value   Pr(>F)
## poison       2  1.033  0.5165   11.79 7.66e-05 ***
## Residuals   45  1.972  0.0438                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Hodnota p je 7.66e-05, což je hluboko pod obvyklou hranicí 0.05 a tři hvězdičky označují kód s nejvyšší významností. Můžete odmítnout H0 a dojít k závěru, že alespoň jedna skupina jedů má odlišnou průměrnou dobu přežití.

Jak zkontrolovat předpoklady ANOVA v R

P-hodnota ANOVA je důvěryhodná pouze tehdy, když jsou splněny tři výše uvedené podmínky. Každá z nich má přímou kontrolu v R a všechny běží na fitovaném modelovém objektu.

1. Nezávislost pozorování. Toto je vlastnost designu studie, nikoli dat, takže žádný test ji nemůže zachránit. Každé morče musí být změřeno jednou a náhodně přiřazeno do své skupiny. Pokud se stejný subjekt objevuje v několika řádcích, potřebujete místo toho model opakovaných měření.

2. Normálnost reziduí. ANOVA předpokládá, že rezidua, nikoli nezpracovaná data, jsou přibližně normální. Prohlédněte si QQ graf a potvrďte Shapiro-Wilkovým testem:

par(mfrow = c(2, 2))
plot(anova_one_way)          # four diagnostic plots

shapiro.test(residuals(anova_one_way))

Body přiléhající k diagonále normálního QQ grafu označují normální rezidua. Shapiro-Wilkova p-hodnota nad 0.05 znamená, že nelze normalitu odmítnout.

3. Homogenita rozptylu. Každá skupina by měla vykazovat podobné rozptylování. Graf reziduí vs. fitovaných hodnot by měl vypadat spíše jako plochý pás než jako trychtýř. Potvrďte to Leveneovým testem, který je robustnější vůči nenormalitě než Bartlettova metoda:

library(car)
leveneTest(time ~ poison, data = df)

bartlett.test(time ~ poison, data = df)

Hodnota p nad 0.05 podporuje stejné rozptyly.

Co dělat, když předpoklad selže. Pokud jsou rozptyly nerovnaké, spusťte oneway.test(time ~ poison, data = df, var.equal = FALSE), Welchovu korekci. Pokud jsou rezidua zjevně nenormální a vzorek je malý, přepněte na Kruskal-Wallisův pořadový test, krusalk.test(time ~ poison, data = df). U velkých vyvážených vzorků je ANOVA poměrně robustní až do mírných odchylek od normality, takže hraniční Shapiro-Wilkův výsledek je zřídka fatální.

Párové srovnání

Signifikantní F-test vám říká, že průměry skupin nejsou všechny stejné, ale ne který pár se liší. Tukeyho test signifikantních rozdílů na to odpovídá porovnáním každého páru a zároveň kontrolou míry chyb v rámci rodiny.

TukeyHSD(anova_one_way)

Výstup:

Přečtěte si výstup, jeden řádek pro každou dvojici. diff sloupec obsahuje rozdíl mezi průměry obou skupin, lwr a licence ohraničil 95% interval spolehlivosti pro tento rozdíl a p příd. l. je p-hodnota upravená pro vícenásobná srovnání. Dvojice se významně liší, když její interval vylučuje nulu, což odpovídá případu, kdy je p adj nižší než 0.05. V tomto souboru dat jsou významná srovnání zahrnující jed 3, což odpovídá krabicovému grafu: skupina 3 má zřetelně nižší průměrnou dobu přežití než skupiny 1 a 2, zatímco skupiny 1 a 2 jsou od sebe statisticky nerozeznatelné.

Obousměrná ANOVA

Dvoucestná ANOVA přidává do vzorce druhý faktor. Funguje přesně jako jednocestný test, mění se pouze vzorec:

y ~ x1 + x2

Zde y je kvantitativní proměnná odezvy, zatímco x1 a x2 jsou oba kategoriální faktory.

Hypotéza ve dvoucestném testu ANOVA

• H0: Průměry skupiny jsou pro obě faktorové proměnné stejné
• H1: Alespoň jeden skupinový průměr se liší u alespoň jednoho ze dvou faktorů.

Do modelu přidáte proměnnou „pamlsek“. Tato proměnná zaznamenává léčbu poskytnutou morčeti. Níže uvedený aditivní vzorec testuje, zda každý faktor ovlivňuje dobu přežití samostatně, po započtení ostatních faktorů.

Upravte kód přidáním proměnné treat vedle první nezávislé proměnné.

anova_two_way <- aov(time~poison + treat, data = df)
summary(anova_two_way)

Výstup:

##             Df Sum Sq Mean Sq F value  Pr(>F)    
## poison       2 1.0330  0.5165   20.64 5.7e-07 ***
## treat        3 0.9212  0.3071   12.27 6.7e-06 ***
## Residuals   42 1.0509  0.0250                    
## ---

Obě p-hodnoty (5.7e-07 pro jed a 6.7e-06 pro léčbu) jsou hluboko pod 0.05, takže H0 pro oba faktory zamítnete a dojdete k závěru, že změna buď jedu, nebo léčby ovlivňuje dobu přežití.

Přidání interakčního termínu

Výše uvedený aditivní model předpokládá, že účinek jedu je stejný bez ohledu na léčbu. Pro ověření tohoto předpokladu nahraďte znaménko plus hvězdičkou, která odpovídá jak hlavním účinkům, tak jejich interakci:

anova_interaction <- aov(time~poison * treat, data = df)
summary(anova_interaction)

Pokud řádek poison:treat není významný, je aditivní model lepší volbou, protože využívá méně stupňů volnosti.

ANOVA v R: Rychlý referenční test

V následující tabulce jsou uvedeny všechny výše použité testy, volání jazyka R, které je spouští, a hypotéza, kterou vyhodnocuje:

aov(y ~ X, data = df)

TukeyHSD(ANOVA summary)

aov(y ~ X1 + X2, data = df)

Nejčastější dotazy