Anova 2-czynnikowa

W Anovie 2-czynnikowej badamy 2 hipotezy. Pierwsza dotyczy ewentualnej różnicy między wynikami (score) ankiety nt. satysfakcji z pracy osób o różnej płci. Druga dotyczy ewentualnych różnic dla osob z różnym wykształceniem.

Badać możemy również trzecią - dodatkową hipotezę zerową: interakcje między czynnikami (płcią i wykształceniem).

Statystyki opisowe

Zobaczmy jak wygląda rozkład wyników - satysfakcji z pracy pracowników o różnej płci i wykształceniu.

ggplot(jobsatisfaction, aes(x=score)) +
  geom_histogram(bins=10) +
  facet_grid(gender~education_level)

pudelkowy<- ggboxplot(jobsatisfaction, x="gender",
          y="score", color="education_level")
pudelkowy

Zerknijmy, jak wyglądają poszczególne statystyki opisowe wg grup płci i wykształcenia dla ocen z satysfakcji z pracy:

Statystyki opisowe
Mężczyzna.Średnie (N = 9) Kobieta.Średnie (N = 10) Mężczyzna.Licencjat (N = 9) Kobieta.Licencjat (N = 10) Mężczyzna.Mgr i wyższe (N = 10) Kobieta.Mgr i wyższe (N = 10)
Min 4.78 4.93 5.58 5.65 8.70 6.52
Max 5.94 6.38 6.74 7.10 10.00 9.57
Q1 5.22 5.51 6.01 6.23 9.03 8.04
Mediana 5.51 5.72 6.30 6.45 9.20 8.48
Q3 5.65 6.05 6.45 6.78 9.50 9.06
Średnia 5.43 5.74 6.22 6.46 9.29 8.41
Odch. std. 0.36 0.47 0.34 0.47 0.44 0.94
Skośność -0.29 -0.12 -0.35 -0.13 0.45 -0.58
Kurtoza -1.20 -1.30 -0.89 -1.30 -1.25 -0.81

Założenia

Sprawdźmy, czy mamy prawo używać Anovy parametrycznej.

Obserwacje odstające

jobsatisfaction %>% 
  group_by(gender, education_level) %>%
  identify_outliers(score)

Brak obserwacji odstających.

Normalność

jobsatisfaction %>% 
  group_by(gender, education_level) %>%
  shapiro_test(score)

Nie ma podstaw, aby odrzucić hipotezę o normalności rozkładów ocen z satysfakcji z pracy dla każdej badanej grupy. Wyniki mają rozkłady normalne (p>0.05).

Zobaczmy, jak wyglądają wykresy normalności wg podgrup:

ggqqplot(jobsatisfaction, "score") +
  facet_grid(gender~education_level)

Wszystkie punkty (odpowiedzi ankietowanych) leżą na linii referencyjnej normalności - dla każdej z podgrup, a więc możemy uznać, że rozkłady satysfakcji z pracy są normalne.

Jednorodna wariancja

jobsatisfaction %>%
  levene_test(score~gender*education_level)

p-value > 0.05, a więc stwierdzamy brak istotnej różnicy wariancji (nie możemy odrzucić hipotezy zerowej). Zatem wariancje są jednorodne.

Anova

wyniki <- jobsatisfaction %>% anova_test(score~gender*education_level)
wyniki

W powyższej tabeli ANOVY widzimy, że istotne różnice w ocenie satysfakcji z pracy występują dla różnego poziomu wykształcenia. Dodatkowo zauważamy, że istnieją istotne interakcje między płcią a wykształceniem.

ges - oznacza uogólniony (generalized) efekt standaryzowany - czyli ilość wariancji, zróżnicowania objaśnianego przez czynnik.

Testy post-hoc

W przypadku odrzucenia hipotez zerowych w Anovie zwykle przeprowadza się tzw. testy post-hoc aby sprawdzić, co sprawiło, iż pojawia się istotna różnica między średnimi - czyli aby zobaczyć, które dokładnie z nich różnią się istotnie.

Sprawdźmy, jak dokładnie wyniki ocen satysfakcji z pracy różnią się wg poziomu wykształcenia i płci:

wynik2<- jobsatisfaction %>%
  group_by(gender) %>%
  tukey_hsd(score~education_level)
wynik2

Jak widzimy powyżej, dla prawie wszystkich par średnich (dla kobiet i mężczyzn vs. wykształcenie) różnice w ocenach satysfakcji z pracy okazały się istotne (p<0.05). Dla jednej pary różnice nie są istotne (ns - non significant) - dla kobiet z wykształceniem college vs. school.

Podsumowanie

Podsumujmy testy Anova za pomocą wykresu ujmującego razem wszystkie poziomy badanych czynników z etykietami i podpisami.

wynik2 <- wynik2 %>% add_xy_position(x="gender")
pudelkowy +
  stat_pvalue_manual(wynik2) +
  labs(
    subtitle=get_test_label(wyniki,detailed=TRUE),
    caption=get_pwc_label(wynik2)
  )

Jak widać, interakcje między czynnikami są istotne (F=7.34, p=0.002).

attach(jobsatisfaction)
interaction.plot(education_level,gender,score,type="b",col=c(2:3),leg.bty="o ",leg.bg="beige",lwd=2,pch=c(18,24) ,xlab="Wykształcenie",ylab="Satysfakcja z pracy",main="Wykres interakcji")

Challenge

Wykorzystamy zestaw danych dotyczących bólu głowy [pakiet datarium], który zawiera miary wyniku bólu epizodu migrenowego bólu głowy u 72 uczestników leczonych trzema różnymi metodami.

Wśród uczestników znalazło się 36 mężczyzn i 36 kobiet. Mężczyźni i kobiety byli dalej podzieleni w podgrupy, zależnie od tego czy byli na niskim lub wysokim ryzyku migreny.

Chcemy zrozumieć, jak każda zmienna niezależna (rodzaj zabiegów, ryzyko migreny i płeć) współdziałają w przewidywaniu wyniku bólu.

Załaduj dane i skontroluj jeden losowy rząd według kombinacji grup:

set.seed(123)
data("headache", package = "datarium")
headache %>% sample_n_by(gender, risk, treatment, size = 1)

Wizualizacje

pudelkowy <- ggboxplot(
  headache, x = "treatment", y = "pain_score", 
  color = "risk", palette = "jco", facet.by = "gender"
  )
pudelkowy

Testy założeń

ggqqplot(headache, "pain_score", ggtheme = theme_bw()) +
  facet_grid(gender + risk ~ treatment, labeller = "label_both")

headache %>% levene_test(pain_score ~ gender*risk*treatment)

Jak widać powyżej, nie mamy problemu z brakiem normalności czy różnicą w wariancjach badanych podgrup danych.

Anova

wyniki <- headache %>% anova_test(pain_score~gender*risk*treatment)
wyniki

Testy post-hoc

install.packages("emmeans")
## Instalowanie pakietu w 'C:/Users/anisz/OneDrive/Documents/R/win-library/4.1'
## (ponieważ 'lib' nie jest określony)
## package 'emmeans' successfully unpacked and MD5 sums checked
## 
## The downloaded binary packages are in
##  C:\Users\anisz\AppData\Local\Temp\RtmpqMssPM\downloaded_packages
library(emmeans)
## Warning: pakiet 'emmeans' został zbudowany w wersji R 4.1.3
parami <- headache %>%
  group_by(gender, risk) %>%
  emmeans_test(pain_score ~ treatment, p.adjust.method = "bonferroni") %>%
  select(-df, -statistic, -p) 

parami %>% filter(gender == "male", risk == "high")

Przeprowadzono trójczynnikową ANOVA w celu określenia wpływu płci, ryzyka i leczenia na pain_score epizodu migrenowego bólu głowy.

Istniała statystycznie istotna trójstronna interakcja między płcią, ryzykiem i leczeniem, F(2, 60) = 7,41, p = 0,001.

Istniała statystycznie istotna prosta dwukierunkowa interakcja między ryzykiem a leczeniem dla mężczyzn, F(2, 60) = 5,2, p = 0,008, ale nie dla kobiet, F(2, 60) = 2,8, p = 0,065.

Istniał statystycznie istotny prosty efekt główny leczenia dla mężczyzn z wysokim ryzykiem migreny, F(2, 60) = 14,8, p < 0,0001), ale nie dla mężczyzn z niskim ryzykiem migreny, F(2, 60) = 0,66, p = 0,521.

Wystąpiła statystycznie istotna średnia różnica między leczeniem X a leczeniem Y. Jednak różnica między leczeniem Y a leczeniem Z, nie była statystycznie istotna.

Wizualizacja Anovy

parami <- parami %>% add_xy_position(x = "treatment")
filtrowane_dane <- parami %>% filter(gender == "male", risk == "high")
pudelkowy +
  stat_pvalue_manual(
    filtrowane_dane, color = "risk", linetype = "risk", hide.ns = TRUE,
    tip.length = 0, step.increase = 0.1, step.group.by = "gender"
  ) +
  labs(
    subtitle = get_test_label(wyniki, detailed = TRUE),
    caption = get_pwc_label(parami)
    )