Anova 3-czynnikowa
3 - czynnikowa ANOVA
Trójczynnikowa ANOVA jest rozszerzeniem dwukierunkowej ANOVA do oceny, czy istnieje efekt interakcji pomiędzy trzema niezależnymi zmiennymi kategorycznymi (jakościowymi) na ciągłą zmienną wynikową (ilościową).
WykorzystaMY zbiór danych headache [z pakietu datarium], który zawiera ankietę dot. miary bólu w epizodzie migrenowego bólu głowy u 72 uczestników leczonych trzema różnymi metodami. Wśród uczestników jest 36 mężczyzn i 36 kobiet. Mężczyźni i kobiety zostali dalej podzieleni na grupy o niskim lub wysokim ryzyku migreny.
Chcemy zrozumieć, jak każda niezależna zmienna (rodzaj leczenia, ryzyko migreny i płeć) współdziałają w przewidywaniu wyniku bólu.
Hipoteza zerową jest stwierdzenie, że rodzaj leczenia, ryzyko migreny i płeć nie wpływają na wyniki bólu.
Statystyki opisowe
Jak wyglądaja dane?
Statystyki opisowe
headache%>%
group_by(gender, risk, treatment) %>%
get_summary_stats(pain_score, type="mean_sd")## # A tibble: 12 x 7
## gender risk treatment variable n mean sd
## <fct> <fct> <fct> <chr> <dbl> <dbl> <dbl>
## 1 Mężczyzna wysokie lek X pain_score 6 92.7 5.12
## 2 Mężczyzna wysokie lek Y pain_score 6 82.3 5.00
## 3 Mężczyzna wysokie lek Z pain_score 6 79.7 4.05
## 4 Mężczyzna nieskie lek X pain_score 6 76.1 3.86
## 5 Mężczyzna nieskie lek Y pain_score 6 73.1 4.76
## 6 Mężczyzna nieskie lek Z pain_score 6 74.5 4.89
## 7 Kobieta wysokie lek X pain_score 6 78.9 5.32
## 8 Kobieta wysokie lek Y pain_score 6 81.2 4.62
## 9 Kobieta wysokie lek Z pain_score 6 81.0 3.98
## 10 Kobieta nieskie lek X pain_score 6 74.2 3.69
## 11 Kobieta nieskie lek Y pain_score 6 68.4 4.08
## 12 Kobieta nieskie lek Z pain_score 6 69.8 2.72Histogram
ggplot(headache, aes(x=pain_score)) +
geom_histogram(bins=10) +
facet_grid(gender~risk~treatment)
Założenia
Obserwacje odstające
Sprawdźmy czy istnieją obserwacje odstające
headache %>%
group_by(treatment, risk, gender) %>%
identify_outliers(pain_score)## # A tibble: 4 x 7
## gender risk treatment id pain_score is.outlier is.extreme
## <fct> <fct> <fct> <int> <dbl> <lgl> <lgl>
## 1 Kobieta wysokie lek X 57 68.4 TRUE TRUE
## 2 Kobieta wysokie lek Y 62 73.1 TRUE FALSE
## 3 Kobieta wysokie lek Z 67 75.0 TRUE FALSE
## 4 Kobieta wysokie lek Z 71 87.1 TRUE FALSEBadana próba zawiera 4 obserwacje odstajace, w tym jedna to obserwacja ekstremalnie odstająca.
Normalność
Sprawdźmy założenie normalności
headache %>%
group_by(treatment, risk, gender) %>%
shapiro_test(pain_score)## # A tibble: 12 x 6
## gender risk treatment variable statistic p
## <fct> <fct> <fct> <chr> <dbl> <dbl>
## 1 Mężczyzna wysokie lek X pain_score 0.958 0.808
## 2 Kobieta wysokie lek X pain_score 0.714 0.00869
## 3 Mężczyzna nieskie lek X pain_score 0.982 0.962
## 4 Kobieta nieskie lek X pain_score 0.933 0.600
## 5 Mężczyzna wysokie lek Y pain_score 0.902 0.384
## 6 Kobieta wysokie lek Y pain_score 0.939 0.654
## 7 Mężczyzna nieskie lek Y pain_score 0.920 0.507
## 8 Kobieta nieskie lek Y pain_score 0.927 0.555
## 9 Mężczyzna wysokie lek Z pain_score 0.955 0.784
## 10 Kobieta wysokie lek Z pain_score 0.971 0.901
## 11 Mężczyzna nieskie lek Z pain_score 0.924 0.535
## 12 Kobieta nieskie lek Z pain_score 0.958 0.801Hipoteza zerowa odrzucona jest jedynie dla jednej zmiennej - kobiety o wysokim ryzyku migreny oraz leczeniem lekiem X. Wszytskie pozostałe przypadki posiadają p>0.05 zatem nie ma dla nich podstaw do odrzucenia hipotez zerowych.
Wykres normalności wg podgrup
ggqqplot(headache, x = ("pain_score"), color = "blue", title = "Wykres norlamności według podgrup" )+
facet_grid(gender~risk~treatment)
Jednorodność wariancji
Sprawdźmy czy wariancja jest jednorodna
headache %>% levene_test(pain_score ~ gender*risk*treatment)## # A tibble: 1 x 4
## df1 df2 statistic p
## <int> <int> <dbl> <dbl>
## 1 11 60 0.179 0.998p-value > 0.05, Zatem wariancje są jednorodne.
Anova
Przejdźmy do głównej analizy wariancji
headache %>%
anova_test(pain_score~gender*risk*treatment)## Coefficient covariances computed by hccm()## ANOVA Table (type II tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 gender 1 60 16.196 0.000163000000000 * 0.213
## 2 risk 1 60 92.699 0.000000000000088 * 0.607
## 3 treatment 2 60 7.318 0.001000000000000 * 0.196
## 4 gender:risk 1 60 0.141 0.708000000000000 0.002
## 5 gender:treatment 2 60 3.338 0.042000000000000 * 0.100
## 6 risk:treatment 2 60 0.713 0.494000000000000 0.023
## 7 gender:risk:treatment 2 60 7.406 0.001000000000000 * 0.198Wniosek: Analizując tabele ANOVY widzimy, że istotne różnice w poziomie bólu głowy występują dla wszytskich trzech czynników: płci, ryzyka oraz podawanego leku. Dodatkowo zauważamy, że istnieją istotne interakcje między płcią a otrzymywanym lekiem oraz trójczynnikowa interakcja płci, leku i ryzyka.
** WYkres ANOWY**
Testy post-hoc
Sprawdźmy co sprawiło, że pojawiły się te istotne różnice między zmiennymi
lm(pain_score~gender*risk*treatment,data=headache)##
## Call:
## lm(formula = pain_score ~ gender * risk * treatment, data = headache)
##
## Coefficients:
## (Intercept)
## 92.74
## genderKobieta
## -13.87
## risknieskie
## -16.69
## treatmentlek Y
## -10.40
## treatmentlek Z
## -13.06
## genderKobieta:risknieskie
## 11.98
## genderKobieta:treatmentlek Y
## 12.71
## genderKobieta:treatmentlek Z
## 15.23
## risknieskie:treatmentlek Y
## 7.48
## risknieskie:treatmentlek Z
## 11.46
## genderKobieta:risknieskie:treatmentlek Y
## -15.59
## genderKobieta:risknieskie:treatmentlek Z
## -18.01headache %>%
group_by(gender) %>%
anova_test(pain_score~risk*treatment)## Coefficient covariances computed by hccm()
## Coefficient covariances computed by hccm()## # A tibble: 6 x 8
## gender Effect DFn DFd F p `p<.05` ges
## * <fct> <chr> <dbl> <dbl> <dbl> <dbl> <chr> <dbl>
## 1 Mężczyzna risk 1 30 45.0 0.000000195 "*" 0.6
## 2 Mężczyzna treatment 2 30 9.15 0.000788 "*" 0.379
## 3 Mężczyzna risk:treatment 2 30 4.72 0.016 "*" 0.24
## 4 Kobieta risk 1 30 48.2 0.000000104 "*" 0.616
## 5 Kobieta treatment 2 30 0.542 0.587 "" 0.035
## 6 Kobieta risk:treatment 2 30 3.23 0.054 "" 0.177Różnice w poziomie bólu głowy są instotne dla 4 zmiennych: mężczyzna vs risk, mężczyzna vs treatment, mężczyna vs risk;treatment, konieta vs risk.