|
| 1 | +--- |
| 2 | +title: "Module 10 : Analyse de la variance (ANOVA) et Kruskal-Wallis" |
| 3 | +subtitle: "Traitement des données I" |
| 4 | +author: "Guyliann Engels & Philippe Grosjean" |
| 5 | +output: |
| 6 | + learnr::tutorial |
| 7 | +tutorial: |
| 8 | + id: "sdd1.10b" |
| 9 | + version: 1.0.0 |
| 10 | +runtime: shiny_prerendered |
| 11 | +--- |
| 12 | + |
| 13 | +```{r setup, include=FALSE} |
| 14 | +library(learnr) |
| 15 | +library(knitr) |
| 16 | +SciViews::R() |
| 17 | +
|
| 18 | +options(tutorial.event_recorder = BioDataScience::record_sdd) |
| 19 | +tutorial_options(exercise.checker = BioDataScience::checker_sdd) |
| 20 | +tutorial_options(exercise.timelimit = 60) |
| 21 | +tutorial_options(exercise.cap = "Code R") |
| 22 | +knitr::opts_chunk$set(echo = FALSE, comment = NA) |
| 23 | +
|
| 24 | +library(BioDataScience) |
| 25 | +``` |
| 26 | + |
| 27 | +```{r, echo=FALSE} |
| 28 | +fixedRow( |
| 29 | + column(9, div( |
| 30 | + img(src = 'images/BioDataScience-128.png', align = "left"), |
| 31 | + h1("Science des données biologiques"), |
| 32 | + "Réalisé par le service d'Écologie numérique des Milieux aquatiques, Université de Mons (Belgique)" |
| 33 | + )), |
| 34 | + column(3, div( |
| 35 | + textInput("user", "Utilisateur :", ""), |
| 36 | + textInput("email", "Email :", "") |
| 37 | + )) |
| 38 | +) |
| 39 | +textOutput("user") # This is newer shown, but required to trigger an event! |
| 40 | +textOutput("email") # Idem! |
| 41 | +``` |
| 42 | + |
| 43 | +```{r, context="server"} |
| 44 | +output$user <- renderText({BioDataScience::user_name(input$user);""}) |
| 45 | +output$email <- renderText({BioDataScience::user_email(input$email);""}) |
| 46 | +updateTextInput(session, "user", value = BioDataScience::user_name()) |
| 47 | +updateTextInput(session, "email", value = BioDataScience::user_email()) |
| 48 | +``` |
| 49 | + |
| 50 | +## Préambule |
| 51 | + |
| 52 | +Si vous n'avez jamais utilisé de tutoriel "learnr", familiarisez-vous d'abord avec son interface [ici](http://biodatascience-course.sciviews.org/sdd-umons/learnr.html). |
| 53 | + |
| 54 | + |
| 55 | + |
| 56 | +> Conformément au RGPD ([Règlement Général sur la Protection des Données](https://ec.europa.eu/info/law/law-topic/data-protection/reform/rules-business-and-organisations/principles-gdpr_fr)), nous sommes tenus de vous informer de ce que vos résultats seront collecté afin de suivre votre progression. **Les données seront enregistrées au nom de l'utilisateur apparaissant en haut de cette page. Corrigez si nécessaire !** En utilisant ce tutoriel, vous marquez expressément votre accord pour que ces données puissent être collectées par vos enseignants et utilisées pour vous aider et vous évaluer. Après avoir été anonymisées, ces données pourront également servir à des études globales dans un cadre scientifique et/ou éducatif uniquement. |
| 57 | +
|
| 58 | +## Objectifs |
| 59 | + |
| 60 | +- Pouvoir comparer plus de deux populations simultanément en utilisant des techniques de décomposition de la variance |
| 61 | + |
| 62 | +- Savoir effectuer une analyse de variance (ANOVA) |
| 63 | + |
| 64 | +- Savoir effectuer des tests de comparaison multiples |
| 65 | + |
| 66 | +- Connaitre l’équivalent non paramétrique à un facteur (test de Kruskal-Wallis) |
| 67 | + |
| 68 | +## Croissance de plantes |
| 69 | + |
| 70 | +Le jeu de données portant sur la croissance de plantes est mis à votre disposition afin de comparer l'analyse de variance du test de Kruskal-Wallis. |
| 71 | + |
| 72 | +```{r, echo = TRUE} |
| 73 | +# Import packages |
| 74 | +SciViews::R |
| 75 | +# Import dataset |
| 76 | +(plant <- read("PlantGrowth", package = "datasets")) |
| 77 | +``` |
| 78 | + |
| 79 | +Réalisez un résumé des données montrant la moyenne de la masse sèche de plantes par traitement. Indiquez également l'écart-type et le nombre d'observations par groupe. |
| 80 | + |
| 81 | +```{r prepare_plant} |
| 82 | +plant <- read("PlantGrowth", package = "datasets") |
| 83 | +``` |
| 84 | +Le snippet correspondant à cette instructions est **.hmanova1desc** |
| 85 | + |
| 86 | +```{r, echo = TRUE, eval = FALSE} |
| 87 | +DF %>.% |
| 88 | + group_by(., XFACTOR) %>.% |
| 89 | + summarise(., mean = mean(YNUM), sd = sd(YNUM), count = sum(!is.na(YNUM))) |
| 90 | +``` |
| 91 | + |
| 92 | +```{r plant1, exercise = TRUE, exercise.setup = "prepare_plant"} |
| 93 | +
|
| 94 | +``` |
| 95 | + |
| 96 | +```{r plant1-solution} |
| 97 | +plant %>.% |
| 98 | + group_by(., group) %>.% |
| 99 | + summarise(., mean = mean(weight), sd = sd(weight), count = sum(!is.na(weight))) |
| 100 | +``` |
| 101 | + |
| 102 | +```{r plant1-check} |
| 103 | +#TODO |
| 104 | +``` |
| 105 | + |
| 106 | +Réalisez à présent un graphique de type nuage de point . Ajoutez également la moyenne et l'écart type. |
| 107 | + |
| 108 | +Il n'existe pas directement de snippet dédié à ce graphique cependant, vous êtes capable d'employer les instructions ci-dessous pour réaliser le graphique demandé. |
| 109 | + |
| 110 | +```{r, echo = TRUE, eval=FALSE} |
| 111 | +chart(data = DF, YNUM ~ XFACTOR) + |
| 112 | + geom_jitter(alpha = 0.4, width = 0.2) + |
| 113 | + stat_summary(geom = "point", fun.y = "mean") + |
| 114 | + stat_summary(geom = "errorbar", width = 0.1, |
| 115 | + fun.data = "mean_sdl", fun.args = list(mult = 1)) |
| 116 | +``` |
| 117 | + |
| 118 | +```{r plant2, exercise = TRUE, exercise.setup ="prepare_plant"} |
| 119 | +
|
| 120 | +``` |
| 121 | + |
| 122 | +```{r plant2-solution} |
| 123 | +chart(data = plant, weight ~ group) + |
| 124 | + geom_jitter(alpha = 0.4, width = 0.2) + |
| 125 | + stat_summary(geom = "point", fun.y = "mean") + |
| 126 | + stat_summary(geom = "errorbar", width = 0.1, |
| 127 | + fun.data = "mean_sdl", fun.args = list(mult = 1)) |
| 128 | +``` |
| 129 | + |
| 130 | +```{r plant2-check} |
| 131 | +#TODO |
| 132 | +``` |
| 133 | + |
| 134 | +Réalisez maintenant une boite de dispersion montre la variation de masse sèche en fonction du traiement. |
| 135 | + |
| 136 | +```{r, echo = TRUE, eval=FALSE} |
| 137 | +chart(data = DF, YNUM ~ XFACTOR) + |
| 138 | + geom_boxplot() |
| 139 | +``` |
| 140 | + |
| 141 | +```{r plant3, exercise = TRUE, exercise.setup ="prepare_plant"} |
| 142 | +
|
| 143 | +``` |
| 144 | + |
| 145 | +```{r plant3-solution} |
| 146 | +chart(data = plant, weight ~ group) + |
| 147 | + geom_boxplot() |
| 148 | +``` |
| 149 | + |
| 150 | +```{r plant3-check} |
| 151 | +#TODO |
| 152 | +``` |
| 153 | + |
| 154 | +Posez vous maintenant cette question avant de continuer l'analyse. Sur base du tableau de données et des trois graphiques précédents pensez vous qu'il y a une différence entre les traitements ? L'intuition n'est cependant pas un test statistique, réalisez donc une Anova et/ou kruskal-Wallis |
| 155 | + |
| 156 | +### Anova |
| 157 | + |
| 158 | +Vous devez tout d'abord vérifier l'homoscédasticité. |
| 159 | + |
| 160 | +Le snippet correspondant à cette instructions est **.hvbartlett** |
| 161 | + |
| 162 | +```{r, echo = TRUE, eval = FALSE} |
| 163 | +bartlett.test(data = DF, YNUM ~ XFACTOR) |
| 164 | +``` |
| 165 | + |
| 166 | +```{r plant4, exercise = TRUE, exercise.setup = "prepare_plant"} |
| 167 | +
|
| 168 | +
|
| 169 | +``` |
| 170 | + |
| 171 | +```{r plant4-solution} |
| 172 | +bartlett.test(data = plant, weight ~ group) |
| 173 | +``` |
| 174 | + |
| 175 | +```{r plant4-check} |
| 176 | +#TODO |
| 177 | +``` |
| 178 | + |
| 179 | +```{r quiz1} |
| 180 | +question("Y a t'il homoscédasticité des variances ?", |
| 181 | + answer("oui", correct = TRUE), |
| 182 | + answer("non")) |
| 183 | +``` |
| 184 | + |
| 185 | +Réalisez votre anova. Le snippet correspondant à cette instructions est **.hmanova1** |
| 186 | + |
| 187 | +```{r, eval= FALSE, echo= TRUE} |
| 188 | +anova(anova. <- lm(data = DF, YNUM ~ XFACTOR)) |
| 189 | +``` |
| 190 | + |
| 191 | +```{r plant5, exercise = TRUE, exercise.setup = "prepare_plant"} |
| 192 | +
|
| 193 | +``` |
| 194 | + |
| 195 | +```{r plant5-solution} |
| 196 | +anova(anova. <- lm(data = plant, weight ~ group)) |
| 197 | +``` |
| 198 | + |
| 199 | +```{r plant5-check} |
| 200 | +#TODO |
| 201 | +``` |
| 202 | + |
| 203 | +Vérifiez la normalité de vos résidus. Le snippet correspondant à cette instruction est **.hmanova1** |
| 204 | + |
| 205 | +```{r echo = TRUE, eval = FALSE} |
| 206 | +#plot(anova., which = 2) |
| 207 | +anova. %>.% |
| 208 | + chart(broom::augment(.), aes(sample = .std.resid)) + |
| 209 | + geom_qq() + |
| 210 | + #geom_qq_line(colour = "darkgray") + |
| 211 | + labs(x = "Theoretical quantiles", y = "Standardized residuals") + |
| 212 | + ggtitle("Normal Q-Q") |
| 213 | +``` |
| 214 | + |
| 215 | +```{r prepare_plant1} |
| 216 | +plant <- read("PlantGrowth", package = "datasets") |
| 217 | +anova. <- lm(data = plant, weight ~ group) |
| 218 | +``` |
| 219 | + |
| 220 | +```{r plant6, exercise = TRUE, exercise.setup = "prepare_plant1"} |
| 221 | +
|
| 222 | +``` |
| 223 | + |
| 224 | +```{r plant6-solution} |
| 225 | +#plot(anova., which = 2) |
| 226 | +anova. %>.% |
| 227 | + chart(broom::augment(.), aes(sample = .std.resid)) + |
| 228 | + geom_qq() + |
| 229 | + #geom_qq_line(colour = "darkgray") + |
| 230 | + labs(x = "Theoretical quantiles", y = "Standardized residuals") + |
| 231 | + ggtitle("Normal Q-Q") |
| 232 | +``` |
| 233 | + |
| 234 | +```{r plant6-check} |
| 235 | +#TODO |
| 236 | +``` |
| 237 | + |
| 238 | +Complétez votre anova par une analyse post hoc. Le snippet correspondant à cette instruction est **.hmanovamult** |
| 239 | + |
| 240 | +```{r, echo = TRUE, eval = FALSE} |
| 241 | +summary(anovaComp. <- confint(multcomp::glht(anova., |
| 242 | + linfct = multcomp::mcp(XFACTOR = "Tukey")))) # Add a second factor if you want |
| 243 | +.oma <- par(oma = c(0, 5.1, 0, 0)); plot(anovaComp.); par(.oma); rm(.oma) |
| 244 | +``` |
| 245 | + |
| 246 | +```{r plant7, exercise = TRUE, exercise.setup = "prepare_plant1"} |
| 247 | +
|
| 248 | +``` |
| 249 | + |
| 250 | +```{r plant7-solution} |
| 251 | +summary(anovaComp. <- confint(multcomp::glht(anova., |
| 252 | + linfct = multcomp::mcp(group = "Tukey")))) # Add a second factor if you want |
| 253 | +.oma <- par(oma = c(0, 5.1, 0, 0)); plot(anovaComp.); par(.oma); rm(.oma) |
| 254 | +``` |
| 255 | + |
| 256 | +```{r plant7-check} |
| 257 | +#TODO |
| 258 | +``` |
| 259 | + |
| 260 | +## Kruskal-Wallis |
| 261 | + |
| 262 | +Réalisez à présent un test de Kruskal-Wallis afin de comparer vos résultats. Le snippet correspondant à cette instruction est **.hnwilkindep** |
| 263 | + |
| 264 | +```{r, echo = TRUE, eval = FALSE} |
| 265 | +wilcox.test(data = DF, YNUM ~ XFACTOR, |
| 266 | + alternative = "two.sided", conf.level = 0.95) |
| 267 | +``` |
| 268 | + |
| 269 | +```{r plant8, exercise = TRUE, exercise.setup ="prepare_plant"} |
| 270 | +
|
| 271 | +``` |
| 272 | + |
| 273 | +```{r plant8-solution} |
| 274 | +kruskal.test(data = plant, weight ~ group) |
| 275 | +``` |
| 276 | + |
| 277 | +```{r plant8-check} |
| 278 | +#TODO |
| 279 | +``` |
| 280 | + |
| 281 | +Suite à ce test vous pouvez réaliser une analyse complémentaire tout comme avec une anova. Le snippet correspondant à cette instruction est **.hnkrusmult** |
| 282 | + |
| 283 | +```{r, echo = TRUE, eval = FALSE} |
| 284 | +summary(kw_comp. <- nparcomp::nparcomp(data = DF, YNUM ~ XFACTOR)) |
| 285 | +plot(kw_comp.) |
| 286 | +``` |
| 287 | + |
| 288 | +```{r plant9, exercise = TRUE, exercise.setup = "prepare_plant"} |
| 289 | +
|
| 290 | +``` |
| 291 | + |
| 292 | +```{r plant9-solution} |
| 293 | +summary(kw_comp. <- nparcomp::nparcomp(data = plant, weight ~ group)) |
| 294 | +plot(kw_comp.) |
| 295 | +``` |
| 296 | + |
| 297 | +```{r plant9-check} |
| 298 | +#TODO |
| 299 | +``` |
| 300 | + |
| 301 | +Posez vous les questions suivantes : |
| 302 | + |
| 303 | +- est ce que les résultats obtenus avec l'Anova et le Kruskal-Wallis coincident avec votre intuition ? |
| 304 | +- Les résultats obtenus avec l'Anova et le Kruskal-Wallis dans les deux cas sont similaire, est ce que cela était attendu ? |
| 305 | +- Les résultats obtenus avec les comparaisons multiples sont ils similaires ? Lors de l'utilisation du test de comparaison multiple non paramétrique vous ne retrouvez pas de différence significative entre vos groupes, cela ne coincide pas avec ce que vous renvoie le Kruskal-Wallis. Que pouvez vous en déduire ? L'attitude d'un scientifique rigoureux doit être de recommencer l'expérience afin d'avoir un plus grand nombre de données pour réussir à confirmer ou infermer ces hypothèses de départ. La différence détectée par le kruskal wallis ne l'est plus pour un test de comparaison multiple, ce qui laisse planer le doute sur l'aspect significatif de notre test. Soyez vigilant |
| 306 | + |
| 307 | +## Conclusion |
| 308 | + |
| 309 | +Bravo! Vous venez de terminer votre séance d'exercices dans un tutoriel "learnr". |
| 310 | + |
| 311 | +Laissez nous vos impressions sur cet outil pédagogique ou expérimentez encore dans la zone ci-dessous. Rappelez-vous que pour placer un commentaire dans une zone de code R, vous devez utilisez un dièse (`#`) devant vos phrases. |
| 312 | + |
| 313 | +```{r comm, exercise=TRUE, exercise.lines = 8} |
| 314 | +# Ajout de commentaires |
| 315 | +# ... |
| 316 | +``` |
| 317 | + |
| 318 | +```{r comm-check} |
| 319 | +# Not yet... |
| 320 | +``` |
0 commit comments