Jednofaktorová ANOVA
Zadanie: Chcete porovnať priemerný úbytok hmotnosti pre tri rôzne diétne plány (diéta A, diéta B a diéta C). Zhromažďujete údaje od skupín účastníkov dodržiavajúcich každý diétny plán počas jedného mesiaca. Otestujte na hladine významnosti 0.01 pomocou testu implementovaného priamo v Rku. Rozhodnite sa na základe P hodnoty testu. Ak je to potrebné, tak urobte aj post testy. Testujeme nulovú hypotézu oproti alternatívnej:
𝐻0: 𝜇1=𝜇2=𝜇3 𝐻1: ∃ 𝑖, 𝑗: 𝜇𝑖 ≠ 𝜇𝑗
# Dáta o úbytku hmotnosti pre tri diéty
diet_A <- c(2.5, 3.0, 2.8, 3.2, 2.7)
diet_B <- c(3.5, 3.8, 3.2, 4.0, 3.6)
diet_C <- c(2.0, 2.5, 2.2, 2.3, 2.1)
# Vytvorenie dátového rámca
df_diet <- data.frame(
weight_loss = c(diet_A, diet_B, diet_C),
diet = factor(rep(c("A", "B", "C"), each = 5))
)
# Jednofaktorová ANOVA
df_diet$diet<-as.factor(df_diet$diet)
head(df_diet, 15)
## weight_loss diet
## 1 2.5 A
## 2 3.0 A
## 3 2.8 A
## 4 3.2 A
## 5 2.7 A
## 6 3.5 B
## 7 3.8 B
## 8 3.2 B
## 9 4.0 B
## 10 3.6 B
## 11 2.0 C
## 12 2.5 C
## 13 2.2 C
## 14 2.3 C
## 15 2.1 C
Grafické zobrazenie:
library(ggplot2)
ggplot(df_diet, aes(x= diet, y= weight_loss))+
geom_boxplot(fill=c("#4472c4", "#FF6666", "#00cc99"))+
labs(x= "Diétny plán", y= "Úbytok hmotnosti", title = "Porovnanie úbytky hmotnosti v závislosti od diétneho plánu")+
theme(plot.title = element_text(hjust=0.5))
Z
grafu vidíme, že dieta C má nižši úbytok hmotnosti ako ostatné diety a
dieta A ma najväčší úbytok hmotnosti. Použitím analýzy rozptylu zistíme,
či je tento rozdiel štatisticky významný.
Normalita- použijeme Shapiro Wilkov test. 𝐻0: 𝑋𝑖 -ty výber je z normálneho rozdelenia, 𝑖=1,2,…,𝑘 , kde 𝑘 - počet úrovní faktora.
𝐻1: ∃ 𝑋𝑖 , ktorý nepochádza z normálneho rozdelenia.
Keďže pre všetky 3 úrovne faktora je p-hodnota väčšia ako 0.01, 𝐻0 nemôžeme zamietnuť.
## $A
##
## Shapiro-Wilk normality test
##
## data: X[[i]]
## W = 0.98998, p-value = 0.9796
##
##
## $B
##
## Shapiro-Wilk normality test
##
## data: X[[i]]
## W = 0.99174, p-value = 0.9854
##
##
## $C
##
## Shapiro-Wilk normality test
##
## data: X[[i]]
## W = 0.97872, p-value = 0.9276
Pre všetky 3 úrovne faktora má úbytok hmotnosti normálne rozdelenie. Prvý predpoklad použitia parametrickej ANOVY je splnený.
Rovnosť rozptylov. Použijeme Bartlettov test. 𝐻𝑂: 𝜎1=𝜎2=…=𝜎𝑘
𝐻1: ∃ 𝑖,𝑗: 𝜎𝑖≠𝜎𝑗
##
## Bartlett test of homogeneity of variances
##
## data: df_diet$weight_loss and df_diet$diet
## Bartlett's K-squared = 0.74088, df = 2, p-value = 0.6904
Na základe p-hodnoty 𝐻0 nemôžeme zmietnuť na hladine významnosti 𝛼=0.01. Teda aj druhý predpoklad je dodržaný, rozptyly nie sú štatisticky významne odlišné pre rôzne úrovne faktora.
Samotná ANOVA
## Df Sum Sq Mean Sq F value Pr(>F)
## df_diet$diet 2 4.921 2.4607 36.55 7.87e-06 ***
## Residuals 12 0.808 0.0673
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = df_diet$weight_loss ~ df_diet$diet)
##
## $`df_diet$diet`
## diff lwr upr p adj
## B-A 0.78 0.3421668 1.2178332 0.0012611
## C-A -0.62 -1.0578332 -0.1821668 0.0068724
## C-B -1.40 -1.8378332 -0.9621668 0.0000054
## Warning: пакет 'agricolae' был собран под R версии 4.4.2
## $statistics
## MSerror Df F Mean CV Scheffe CriticalDifference
## 0.06733333 12 6.926608 2.893333 8.968433 3.721991 0.61083
##
## $parameters
## test name.t ntr alpha
## Scheffe diet 3 0.01
##
## $means
## weight_loss std r se Min Max Q25 Q50 Q75
## A 2.84 0.2701851 5 0.116046 2.5 3.2 2.7 2.8 3.0
## B 3.62 0.3033150 5 0.116046 3.2 4.0 3.5 3.6 3.8
## C 2.22 0.1923538 5 0.116046 2.0 2.5 2.1 2.2 2.3
##
## $comparison
## NULL
##
## $groups
## weight_loss groups
## B 3.62 a
## A 2.84 b
## C 2.22 c
##
## attr(,"class")
## [1] "group"
Na základe oboch post testov je štatisticky významný rozdiel medzi diétnymi plánmi A, B a C. Diétny plán B je účinnejší ako diétny plán A a C a diétny plán A je účinnejší ako diétny plán C.
Jednofaktorová ANOVA s F testom
Zadanie: Chcete otestovať, či tri rôzne hnojivá (hnojivo A, B a C) vedú k rôznym výnosom plodín. Zhromažďujete údaje o výnosoch z polí, na ktorých sa aplikovalo každé hnojivo, a zaujíma vás, či existuje významný rozdiel medzi priemernými výnosmi. Otestujte na hladine významnosti 0.1 pomocou testu implementovaného priamo v Rku. Rozhodnite sa na základe F štatistiky. Ak je to potrebné, tak urobte aj post testy.
fertilizer_A <- c(10.2, 10.5, 9.8, 10.3, 10.4)
fertilizer_B <- c(11.1, 11.5, 11.3, 10.8, 11.0)
fertilizer_C <- c(9.5, 9.8, 9.7, 9.6, 9.9)
df_diet2<-data.frame(yield= c(fertilizer_A, fertilizer_B, fertilizer_C),
fertilizer = rep(c("A", "B", "C"), c(length(fertilizer_A), length(fertilizer_B), length(fertilizer_C)))
)
df_diet2$fertilizer<-as.factor(df_diet2$fertilizer)
head(df_diet2, 15)
## yield fertilizer
## 1 10.2 A
## 2 10.5 A
## 3 9.8 A
## 4 10.3 A
## 5 10.4 A
## 6 11.1 B
## 7 11.5 B
## 8 11.3 B
## 9 10.8 B
## 10 11.0 B
## 11 9.5 C
## 12 9.8 C
## 13 9.7 C
## 14 9.6 C
## 15 9.9 C
Grafické zobrazenie:
ggplot(df_diet2, aes(x= fertilizer, y= yield))+
geom_boxplot(fill=c("#4472c4", "#FF6666", "#00cc99"))+
labs(x= "Hnojivo", y= "Výnos plodín", title = "Porovnanie výnosov plodín v závislosti od hnojív")+
theme(plot.title = element_text(hjust=0.5))
Z
grafu vidíme, že hnojivo C má nižši výnos plodin ako ostatné hnojíva a
hnojivo A ma najväčší výnos plodin. Použitím analýzy rozptylu zistíme,
či je tento rozdiel štatisticky významný.
Predpoklady
normalita- použijeme Shapiro Wilkov test. 𝐻0: 𝑋𝑖 -ty výber je z normálneho rozdelenia, 𝑖=1,2,…,𝑘 , kde 𝑘 - počet úrovní faktora.
𝐻1: ∃ 𝑋𝑖 , ktorý nepochádza z normálneho rozdelenia.
Keďže pre všetky 3 úrovne faktora je p-hodnota väčšia ako 0.1, 𝐻0 nemôžeme zamietnuť.
## $A
##
## Shapiro-Wilk normality test
##
## data: X[[i]]
## W = 0.90309, p-value = 0.4272
##
##
## $B
##
## Shapiro-Wilk normality test
##
## data: X[[i]]
## W = 0.98998, p-value = 0.9796
##
##
## $C
##
## Shapiro-Wilk normality test
##
## data: X[[i]]
## W = 0.98676, p-value = 0.9672
Pre všetky 3 úrovne faktora má úbytok hmotnosti normálne rozdelenie. Prvý predpoklad použitia parametrickej ANOVY je splnený.
Rovnosť rozptylov. Použijeme Bartlettov test. 𝐻𝑂: 𝜎1=𝜎2=…=𝜎𝑘
𝐻1: ∃ 𝑖,𝑗: 𝜎𝑖≠𝜎𝑗
##
## Bartlett test of homogeneity of variances
##
## data: df_diet2$yield and df_diet2$fertilizer
## Bartlett's K-squared = 1.1857, df = 2, p-value = 0.5528
Na základe p-hodnoty 𝐻0 nemôžeme zmietnuť na hladine významnosti 𝛼=0.01. Teda aj druhý predpoklad je dodržaný, rozptyly nie sú štatisticky významne odlišné pre rôzne úrovne faktora.
Samotná ANOVA
## Df Sum Sq Mean Sq F value Pr(>F)
## df_diet2$fertilizer 2 5.292 2.646 46.42 2.25e-06 ***
## Residuals 12 0.684 0.057
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Podľa výsledku testu vidíme, že F-hodnota (46.42) je rádovo väčšia ako jednota. Teda hypotézu o nulovosti faktora môžeme zamietnuť. Môžeme teda povedať, že rozdiely vo výnosoch plodín v závislosti od hnojív sú štatisticky významné na hladine významnosti 0.1.
Na určenie, ktorý plan diety sa štatisticky významne odlišujé, budeme testovať kontrasty (alebo inak post hoc test, viacnásobné porovnávanie priemerov). Použijeme Tukeyov a Scheffeho test.
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = df_diet2$yield ~ df_diet2$fertilizer)
##
## $`df_diet2$fertilizer`
## diff lwr upr p adj
## B-A 0.90 0.4971614 1.3028386 0.0001806
## C-A -0.54 -0.9428386 -0.1371614 0.0098537
## C-B -1.44 -1.8428386 -1.0371614 0.0000017
an21<-aov(yield ~ fertilizer, data=df_diet2)
result<- scheffe.test(an21,"fertilizer", alpha = 0.1)
result
## $statistics
## MSerror Df F Mean CV Scheffe CriticalDifference
## 0.057 12 2.806796 10.36 2.304505 2.369302 0.3577567
##
## $parameters
## test name.t ntr alpha
## Scheffe fertilizer 3 0.1
##
## $means
## yield std r se Min Max Q25 Q50 Q75
## A 10.24 0.2701851 5 0.1067708 9.8 10.5 10.2 10.3 10.4
## B 11.14 0.2701851 5 0.1067708 10.8 11.5 11.0 11.1 11.3
## C 9.70 0.1581139 5 0.1067708 9.5 9.9 9.6 9.7 9.8
##
## $comparison
## NULL
##
## $groups
## yield groups
## B 11.14 a
## A 10.24 b
## C 9.70 c
##
## attr(,"class")
## [1] "group"
Na základe oboch post testov je štatisticky významný rozdiel medzi hnojivami A, B a C. Hnojivo B je účinnejšie ako hnojivo A a C a hnojivo A je účinnejšie ako hnojivo C.
Dvojfaktorová ANOVA
Zadanie: Predpokladajme, že študujete účinky typu cvičenia a typu stravy na chudnutie. Pre každý faktor máte dve úrovne: Typ cvičenia: Kardio a sila Typ stravy: Nízkosacharidový a nízkotučný Máte podozrenie, že medzi typom cvičenia a typom stravy môže existovať interakcia na chudnutie, t. j. účinnosť cvičenia sa môže líšiť v závislosti od typu dodržiavanej diéty.
weight_loss <- c(
# Cardio + Low-Carb
rnorm(10, mean = 5, sd = 1),
# Cardio + Low-Fat
rnorm(10, mean = 5, sd = 1),
# Strength + Low-Carb
rnorm(10, mean = 5, sd = 1),
# Strength + Low-Fat
rnorm(10, mean = 7, sd = 1)
)
exercise <- factor(rep(c("Cardio", "Strength"), each = 20))
diet <- factor(rep(c("Low-Carb", "Low-Fat"), times = 10, each = 10))
df_diet3 <- data.frame(
weight_loss,
exercise,
diet
)
df_diet3$exercise<-as.factor(df_diet3$exercise)
df_diet3$diet<-as.factor(df_diet3$diet)
Urobme si prvotnu grafickú analýzu:
ggplot(df_diet3, aes(x= diet, y= weight_loss))+
geom_boxplot(aes(fill=diet))+
scale_fill_brewer(palette = "Greys", name="")+
labs(x= "Diéta", y= "Сhudnutie [kg]", title = "Porovnanie chudnutia vzhľadom na diéty")+
theme(plot.title = element_text(hjust=0.5))
ggplot(df_diet3, aes(x= exercise, y= weight_loss))+
geom_boxplot(aes(fill=exercise))+
scale_fill_brewer(palette = "Greys", name="")+
labs(x= "Сvičenie", y= "Сhudnutie [kg]", title = "Porovnanie chudnutia vzhľadom na cvičenia")+
theme(plot.title = element_text(hjust=0.5))
ggplot(df_diet3, aes(x = diet, y= weight_loss, fill = exercise))+
geom_boxplot()+
scale_fill_brewer(palette = "Greys", name="Сvičenie")+
labs(x= "Diéta", y= "Сhudnutie [kg]", title = "Porovnanie možných interakcií")+
theme(plot.title = element_text(hjust=0.5))
Z boxplotov sa zdá že kombinácia nízkotučnej diety a silových cvičení dáva lepšie chudnutie.
## Df Sum Sq Mean Sq F value Pr(>F)
## diet 1 103.56 103.56 110.31 < 2e-16 ***
## exercise 1 18.28 18.28 19.47 1.69e-05 ***
## diet:exercise 1 61.77 61.77 65.80 5.36e-14 ***
## Residuals 196 184.00 0.94
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Čím je F vyššie, tým silnejšie tento faktor ovplyvňuje konečný výsledok, preto najsilnejšie pôsobí na chudnutie kombinácia stravy a cvičenia.
Na určenie, ktorý plan diety sa štatisticky významne odlišujé, budeme testovať kontrasty (alebo inak post hoc test, viacnásobné porovnávanie priemerov). Použijeme Tukeyov test.
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = weight_loss ~ diet + exercise + diet * exercise, data = df_diet3)
##
## $diet
## diff lwr upr p adj
## Low-Fat-Low-Carb 1.439157 1.16893 1.709385 0
##
## $exercise
## diff lwr upr p adj
## Strength-Cardio 0.6045807 0.334353 0.8748084 1.69e-05
##
## $`diet:exercise`
## diff lwr upr p adj
## Low-Fat:Cardio-Low-Carb:Cardio 0.3276983 -0.1744235 0.829820053 0.3311285
## Low-Carb:Strength-Low-Carb:Cardio -0.5068783 -1.0090000 -0.004756484 0.0469135
## Low-Fat:Strength-Low-Carb:Cardio 2.0437380 1.5416162 2.545859732 0.0000000
## Low-Carb:Strength-Low-Fat:Cardio -0.8345765 -1.3366983 -0.332454757 0.0001528
## Low-Fat:Strength-Low-Fat:Cardio 1.7160397 1.2139179 2.218161459 0.0000000
## Low-Fat:Strength-Low-Carb:Strength 2.5506162 2.0484944 3.052737996 0.0000000
Vidíme, že najlepšie výsledky prináša nízkotučná dieta v kombinácii so silovými cvičeniami. (má najmenšie p-hodnoty a kladny rozdiel)