testy statystyczne
Anova 1-czynnikowa
niejednorodność wariancji
Przykad 1 LAB
Zmienną badaną jest mpg - liczba mil, którą można przejechać na galonie paliwa. Chcielibyśmy sprawdzić, czy rozkład mpg różni się istotnie wg różnych charakterystyk samochodów.
H0: Liczba mpg nie różni się istotnie dla różnej liczby cylindrów.
Ha: Co najmniej jedna para średnich różni się istotnie.
ggboxplot
ggboxplot(mtcars, x = "cyl", y = "mpg",
color = "cyl", palette = c("#00AFBB", "#E7B800", "#FC4E07"),
order = c("4", "6", "8"),
ylab = "Mil na galon", xlab = "Liczba cylindrów")1. Norlamność
mtcars %>%
group_by(cyl) %>%
shapiro_test(mpg)## # A tibble: 3 x 4
## cyl variable statistic p
## <dbl> <chr> <dbl> <dbl>
## 1 4 mpg 0.912 0.261
## 2 6 mpg 0.899 0.325
## 3 8 mpg 0.932 0.323
Wniosek: Dla każdego cyl rozkład jest normalny.
2. jednorodność wariancji
leveneTest(mpg~cyl)
Wniosek: p<<alfa, co świadczy o niejednorodnej wariancji w badanych próbach.
3. ANOVA 1-czynnikowa dla niejednorodnych wariancji
oneway.test(mpg ~ cyl,
data = mtcars,
var.equal = FALSE # wariancje niejednorodne
)##
## One-way analysis of means (not assuming equal variances)
##
## data: mpg and cyl
## F = 31.624, num df = 2.000, denom df = 18.032, p-value = 1.271e-06
Wniosek: Odrzucamy H0 i stwierdzamy, że co najmniej jedna para średnich mpg różni się istotnie wg liczby cyl.
4. Post-hoc*
model <- aov(mpg~cyl,data=mtcars)
post_hoc<-TukeyHSD(model,“cyl”,ordered=TRUE)
print(post_hoc)
Wniosek: Dla wszystkich par średnich poszczególne testy t odrzuciły H0 o różnicy mpg wg liczby cylindrów. p<0,05
ANOVA 2-czynnikowa
Próby niezależne
Przykład 1 LAB
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(gender). Druga dotyczy ewentualnych różnic dla osob z różnym wykształceniem(labels).
ggboxplot
attach(jobsatisfaction)
pudelkowy<- ggboxplot(jobsatisfaction, x="gender",
y="score", color="education_level")
pudelkowy1. Normalność
jobsatisfaction %>%
group_by(gender, education_level) %>%
shapiro_test(score)## # A tibble: 6 x 5
## gender education_level variable statistic p
## <fct> <fct> <chr> <dbl> <dbl>
## 1 male school score 0.980 0.966
## 2 male college score 0.958 0.779
## 3 male university score 0.916 0.323
## 4 female school score 0.963 0.819
## 5 female college score 0.963 0.819
## 6 female university score 0.950 0.674
Wykresy normalności wg grup
ggqqplot(jobsatisfaction, "score") +
facet_grid(gender~education_level)Wniosek: 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).
2. Obserwacje odstające
jobsatisfaction %>%
group_by(gender, education_level) %>%
identify_outliers(score)## [1] gender education_level id score
## [5] is.outlier is.extreme
## <0 wierszy> (lub 'row.names' o zerowej długości)
Wniosek: Brak obserwacji odstajacych
3. Jednorodne wariancje
jobsatisfaction %>%
levene_test(score~gender*education_level)## # A tibble: 1 x 4
## df1 df2 statistic p
## <int> <int> <dbl> <dbl>
## 1 5 52 2.20 0.0686
Wniosek: p-value > 0.05, a więc stwierdzamy brak istotnej różnicy wariancji (nie możemy odrzucić hipotezy zerowej). Zatem wariancje są jednorodne.
4. ANOVA 2-czynnikowa
wyniki <- jobsatisfaction %>% anova_test(score~gender*education_level)## Coefficient covariances computed by hccm()
wyniki## ANOVA Table (type II tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 gender 1 52 0.745 3.92e-01 0.014
## 2 education_level 2 52 187.892 1.60e-24 * 0.878
## 3 gender:education_level 2 52 7.338 2.00e-03 * 0.220
Wniosek: Istotne różnice w ocenie satysfakcji z pracy występują dla różnego poziomu wykształcenia. + Istotne interakcje między płcią a wykształceniem.
5. POst-hoc
wynik2<- jobsatisfaction %>%
group_by(gender) %>%
tukey_hsd(score~education_level)
wynik2## # A tibble: 6 x 10
## gender term group1 group2 null.value estimate conf.low conf.high p.adj
## * <fct> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 male educati~ school college 0 0.797 0.341 1.25 5.59e- 4
## 2 male educati~ school univer~ 0 3.87 3.42 4.31 3.28e-14
## 3 male educati~ colle~ univer~ 0 3.07 2.62 3.51 3.99e-14
## 4 female educati~ school college 0 0.722 -0.0163 1.46 5.62e- 2
## 5 female educati~ school univer~ 0 2.66 1.93 3.40 4.3 e- 9
## 6 female educati~ colle~ univer~ 0 1.94 1.20 2.68 1.58e- 6
## # ... with 1 more variable: p.adj.signif <chr>
p.value parami: ns - nie istotne * 0,05 ** 0,01 *** 0,001
Wniosek: 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.
Próby zależne
Przykład 1 LAB
W zbiorze “depression” zawarto wyniki poziom stanu depresjnego u 24 osób w czterech punktach czasowych. połowa badanych otrzymała lek, a pozostali są poddani kontroli. Badanie ma na celu orzec, czy na wynik depresji wpływa czas i leczenie.
data("depression", package = "datarium")
head(depression, 3)## # A tibble: 3 x 6
## id treatment t0 t1 t2 t3
## <fct> <fct> <dbl> <dbl> <dbl> <dbl>
## 1 1 ctr 296 175 187 242
## 2 2 ctr 376 329 236 126
## 3 3 ctr 309 238 150 173
Przekształcenie ramki w potzrebną do wykonania obliczeń ramkę danych na długi format:
depression <- depression %>%
gather(key = "time", value = "score", t0, t1, t2, t3) %>%
convert_as_factor(id, time, treatment)
head(depression, 3)## # A tibble: 3 x 4
## id treatment time score
## <fct> <fct> <fct> <dbl>
## 1 1 ctr t0 296
## 2 2 ctr t0 376
## 3 3 ctr t0 309
1. Normalność
depression%>%
group_by(time,treatment)%>%
shapiro_test(score)## # A tibble: 8 x 5
## treatment time variable statistic p
## <fct> <fct> <chr> <dbl> <dbl>
## 1 ctr t0 score 0.945 0.572
## 2 treated t0 score 0.836 0.0248
## 3 ctr t1 score 0.918 0.268
## 4 treated t1 score 0.953 0.680
## 5 ctr t2 score 0.951 0.652
## 6 treated t2 score 0.957 0.734
## 7 ctr t3 score 0.873 0.0706
## 8 treated t3 score 0.932 0.398
Wykres normalności:
ggqqplot(depression, "score") + facet_grid(time~treatment)Wniosek: Tylko dla jednej zmiennej (osoba leczona w czasie 0) hipoteza zerowa zostaje odrzucona. Pozostałe podlegają normalności.
2. Obserwacje odstające
depression%>%
group_by(time,treatment)%>%
identify_outliers(score)## # A tibble: 4 x 6
## treatment time id score is.outlier is.extreme
## <fct> <fct> <fct> <dbl> <lgl> <lgl>
## 1 ctr t0 5 150 TRUE FALSE
## 2 ctr t0 8 447 TRUE FALSE
## 3 treated t0 18 138 TRUE FALSE
## 4 treated t0 22 150 TRUE FALSE
Wniosek: Dane depression posiadają 4 zmienne odstają, które nie odstają na szczęście ekstremalnie. Ten sam wniosek, można wydedukować dzięki zamieszczonemu wyżej wykresowi ramka-wąs
3. Anova dla zmniennych zależnych
anova <- anova_test(data=depression, formula = score~treatment*time,wid = id, within = time)## Coefficient covariances computed by hccm()
get_anova_table(anova)## ANOVA Table (type II tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 treatment 1 88 43.423 3.13e-09 * 0.330
## 2 time 3 88 19.713 7.31e-10 * 0.402
## 3 treatment:time 3 88 3.755 1.40e-02 * 0.113
Wniosek: Zgodnie z testem ANOVA na wynik depresji istotnie wpływa czas oraz leczenie.
Testy post-hoc
Co sprawiło, że nastapiła widoczna różnica?
depression %>%
group_by(treatment) %>%
emmeans_test(score~time)## # A tibble: 12 x 10
## treatment term .y. group1 group2 df statistic p p.adj
## * <chr> <chr> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl>
## 1 ctr time score t0 t1 88 1.60 0.113 0.676
## 2 ctr time score t0 t2 88 2.98 0.00375 0.0225
## 3 ctr time score t0 t3 88 3.60 0.000520 0.00312
## 4 ctr time score t1 t2 88 1.38 0.172 1
## 5 ctr time score t1 t3 88 2.00 0.0485 0.291
## 6 ctr time score t2 t3 88 0.625 0.534 1
## 7 treated time score t0 t1 88 6.25 0.0000000142 0.0000000855
## 8 treated time score t0 t2 88 6.05 0.0000000338 0.000000203
## 9 treated time score t0 t3 88 5.84 0.0000000854 0.000000513
## 10 treated time score t1 t2 88 -0.196 0.845 1
## 11 treated time score t1 t3 88 -0.409 0.683 1
## 12 treated time score t2 t3 88 -0.213 0.832 1
## # ... with 1 more variable: p.adj.signif <chr>
Wniosek: Dla pięciu z dwunastu powiązanych zmiennych różniece w wynikach poziomu stanu depresyjnego okazują sie istotne (p_value<0.5). Są to: 1. osoba pod obserwacjąT0T2 2. osoba pod obserwacjąT0T3 3. osoba leczonaT0T1 4. osoba leczonaT0T2 5. osoba leczonaTOT3
ANOVA 3-czynnikowa
ANOVA do oceny, czy istnieje efekt interakcji pomiędzy trzema niezależnymi zmiennymi kategorycznymi (jakościowymi) na ciągłą zmienną wynikową (ilościową).
Próby niezależne
Przykład 1 LAB
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.
ggboxplot
ggboxplot(depression, x = "time", y = "score", color = "treatment", palette = "jco")Opisówka:
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 male high X pain_score 6 92.7 5.12
## 2 male high Y pain_score 6 82.3 5.00
## 3 male high Z pain_score 6 79.7 4.05
## 4 male low X pain_score 6 76.1 3.86
## 5 male low Y pain_score 6 73.1 4.76
## 6 male low Z pain_score 6 74.5 4.89
## 7 female high X pain_score 6 78.9 5.32
## 8 female high Y pain_score 6 81.2 4.62
## 9 female high Z pain_score 6 81.0 3.98
## 10 female low X pain_score 6 74.2 3.69
## 11 female low Y pain_score 6 68.4 4.08
## 12 female low Z pain_score 6 69.8 2.72
1. Normalność
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 male high X pain_score 0.958 0.808
## 2 female high X pain_score 0.714 0.00869
## 3 male low X pain_score 0.982 0.962
## 4 female low X pain_score 0.933 0.600
## 5 male high Y pain_score 0.902 0.384
## 6 female high Y pain_score 0.939 0.654
## 7 male low Y pain_score 0.920 0.507
## 8 female low Y pain_score 0.927 0.555
## 9 male high Z pain_score 0.955 0.784
## 10 female high Z pain_score 0.971 0.901
## 11 male low Z pain_score 0.924 0.535
## 12 female low Z pain_score 0.958 0.801
Wniosek: Hipoteza 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)2. 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 female high X 57 68.4 TRUE TRUE
## 2 female high Y 62 73.1 TRUE FALSE
## 3 female high Z 67 75.0 TRUE FALSE
## 4 female high Z 71 87.1 TRUE FALSE
Wniosek: Badana próba zawiera 4 obserwacje odstajace, w tym jedna to obserwacja ekstremalnie odstająca.
3. Jednrorodnośc wariancji
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.998
Wniosek: p-value > 0.05, Zatem wariancje są jednorodne.
4. ANOVA 3-czynnikowa
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 1.63e-04 * 0.213
## 2 risk 1 60 92.699 8.80e-14 * 0.607
## 3 treatment 2 60 7.318 1.00e-03 * 0.196
## 4 gender:risk 1 60 0.141 7.08e-01 0.002
## 5 gender:treatment 2 60 3.338 4.20e-02 * 0.100
## 6 risk:treatment 2 60 0.713 4.94e-01 0.023
## 7 gender:risk:treatment 2 60 7.406 1.00e-03 * 0.198
Wniosek:
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.
5. Post-hoc lm ()
lm(pain_score~gender*risk*treatment,data=headache)##
## Call:
## lm(formula = pain_score ~ gender * risk * treatment, data = headache)
##
## Coefficients:
## (Intercept) genderfemale
## 92.739 -13.874
## risklow treatmentY
## -16.687 -10.397
## treatmentZ genderfemale:risklow
## -13.058 11.978
## genderfemale:treatmentY genderfemale:treatmentZ
## 12.708 15.228
## risklow:treatmentY risklow:treatmentZ
## 7.484 11.462
## genderfemale:risklow:treatmentY genderfemale:risklow:treatmentZ
## -15.589 -18.009
headache %>%
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 male risk 1 30 45.0 0.000000195 "*" 0.6
## 2 male treatment 2 30 9.15 0.000788 "*" 0.379
## 3 male risk:treatment 2 30 4.72 0.016 "*" 0.24
## 4 female risk 1 30 48.2 0.000000104 "*" 0.616
## 5 female treatment 2 30 0.542 0.587 "" 0.035
## 6 female risk:treatment 2 30 3.23 0.054 "" 0.177
Wniosek: Róż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.
Próby zależne
Przykład 2
Czy osiąnięcia uczestnika egzaminu istotnie różnią się w zależności od czasu, płci oraz poziomu stresu?
library(datarium)
data("performance")Konwertacja:
performance <- performance %>%
gather(key = "time", value = "score", t1, t2) %>%
convert_as_factor(id, gender, time)1. Normalność
performance %>%
group_by(gender,stress,time) %>%
shapiro_test(score)## # A tibble: 12 x 6
## gender stress time variable statistic p
## <fct> <fct> <fct> <chr> <dbl> <dbl>
## 1 male low t1 score 0.942 0.574
## 2 male low t2 score 0.966 0.849
## 3 male moderate t1 score 0.848 0.0547
## 4 male moderate t2 score 0.958 0.761
## 5 male high t1 score 0.915 0.319
## 6 male high t2 score 0.925 0.403
## 7 female low t1 score 0.898 0.207
## 8 female low t2 score 0.886 0.154
## 9 female moderate t1 score 0.946 0.626
## 10 female moderate t2 score 0.865 0.0880
## 11 female high t1 score 0.989 0.996
## 12 female high t2 score 0.930 0.452
Wniosek: Dla każdej zmiennej p.value > 0.05 więc nie ma podstaw do odrzucenia hipotezy zerowej.
2. Wartości odstające
performance %>%
group_by(gender, stress, time)%>%
identify_outliers(score)## # A tibble: 1 x 7
## gender stress time id score is.outlier is.extreme
## <fct> <fct> <fct> <fct> <dbl> <lgl> <lgl>
## 1 female low t2 36 6.15 TRUE FALSE
Tylko jedna zmienna jest wartością odstającą - kobieta o niskim poziomie stresu i czasie t2
3. ANOVA dla prob zależnych
wyniki<-anova_test(data=performance, formula = score~gender*time*stress,wid = id, within = time)## Coefficient covariances computed by hccm()
get_anova_table(wyniki)## ANOVA Table (type II tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 gender 1 108 2.487 1.18e-01 0.023000
## 2 time 1 108 0.061 8.06e-01 0.000564
## 3 stress 2 108 21.878 1.05e-08 * 0.288000
## 4 gender:time 1 108 4.571 3.50e-02 * 0.041000
## 5 gender:stress 2 108 1.607 2.05e-01 0.029000
## 6 time:stress 2 108 1.759 1.77e-01 0.032000
## 7 gender:time:stress 2 108 5.896 4.00e-03 * 0.098000
Tabela anovy przedstawia widoczne różnice w wynikach na sprawdzianie, zależne od poziomu stresu. Istotna jest też interakcja płci i czasu. Ponadto osiąnięcia uczestnika egzaminu istotnie różnią się w zależności od czasu, płci oraz poziomu stresu.
post-hoc
performance %>%
group_by(gender, stress) %>%
emmeans_test(score~time)## # A tibble: 6 x 11
## gender stress term .y. group1 group2 df statistic p p.adj
## * <chr> <chr> <chr> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl>
## 1 female high time score t1 t2 108 2.13 0.0351 0.0351
## 2 female low time score t1 t2 108 -1.93 0.0561 0.0561
## 3 female moderate time score t1 t2 108 2.72 0.00767 0.00767
## 4 male high time score t1 t2 108 -1.95 0.0537 0.0537
## 5 male low time score t1 t2 108 0.266 0.791 0.791
## 6 male moderate time score t1 t2 108 -0.631 0.529 0.529
## # ... with 1 more variable: p.adj.signif <chr>
Jak widzimy powyżej, dla TYLKO DWÓCH par średnich (dla kobiet vs. poziom stresu) 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).