Chapitre 4 Masse des oeufs
La démarche sera la même que celle du chapitre précédent sur la matière grasse. Il se peut qu’il y ait moins de commentaires.
Nous utiliserons le jeu de données oeuf.csv
qui contient différentes mesures dont les masses des oeufs mesurées en 5 séances. Mêmes traitements (régimes).
La question est de savoir si les différents régimes induisent des masses d’oeufs significativement différentes avec le temps. Mais chaque traitement n’ayant pas été appliqué sur chaque groupe d’oiseaux, l’ANOVA à mesures répétées ne pourrait pas être appliquée. Nous comparerons les effets des traitements séance par séance, puis à l’aide d’une figure on appréciera s’il y a une évolution en fonction du temps.
4.1 Les données
mo <- read_csv("data/oeuf.csv")
mo <- mo %>%
select(seance, regime, no_oeuf, masse_oeuf) %>%
mutate(id = rep(1:30, 5), .before = 1,
id = factor(id),
seance = factor(seance),
regime = factor(regime))
Le tableau a été préalablement structuré en format long en Excel. J’ai ajouté un identifiant (id
) pour les échantillons des séances.
## Rows: 150
## Columns: 5
## $ id <fct> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, ~
## $ seance <fct> seance 1, seance 1, seance 1, seance 1, seance 1, seance 1,~
## $ regime <fct> "Ba 0,25", "Ba 0,25", "Ba 0,25", "Ba 0,50", "Ba 0,50", "Ba ~
## $ no_oeuf <dbl> 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2,~
## $ masse_oeuf <dbl> 10.762, 10.807, 10.650, 10.268, 10.368, 10.620, 11.214, 11.~
4.2 Visualisation boxplots
bxp <- ggplot(mo, aes(x = regime, y = masse_oeuf, fill = regime)) +
geom_boxplot() +
facet_grid(seance ~ .) +
theme(axis.text.x = element_text(angle = 90, color = "black", vjust = 0.5, hjust = 1)) +
theme_bw()
bxp
Il y a des variations notables d’une séance à l’autre pour certains régimes.
4.3 Détection des valeurs aberrantes extrêmes
## [1] seance regime id no_oeuf masse_oeuf is.outlier is.extreme
## <0 rows> (or 0-length row.names)
=> Pas de valeurs aberrantes extrêmes pour toutes les séances.
4.4 Conditions de l’ANOVA
4.4.1 Normalité
Si les données sont normalement distribuées, la p-value de Shapiro-Wilk doit être supérieure à 0,05 pour chaque régime.
## # A tibble: 5 x 4
## seance variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 seance 1 masse_oeuf 0.945 0.126
## 2 seance 2 masse_oeuf 0.603 0.0000000787
## 3 seance 3 masse_oeuf 0.962 0.342
## 4 seance 4 masse_oeuf 0.948 0.146
## 5 seance 5 masse_oeuf 0.963 0.363
=> Normalité confirmée pour toutes les séances sauf la 2.
Créer des QQ-plots pour chaque point par séance
Tous les points se situent approximativement le long de la ligne de référence sauf pour la séance 2 où 1 point se démarque. Il n’a cependant pas été identifié comme aberrant. Mais on pourra l’exclure et apprécier au moment venu.
4.4.2 Homogénéité des variances
## # A tibble: 5 x 5
## seance df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 seance 1 9 20 1.19 0.353
## 2 seance 2 9 20 0.887 0.553
## 3 seance 3 9 20 0.805 0.617
## 4 seance 4 9 20 0.718 0.687
## 5 seance 5 9 20 0.383 0.930
=> Toutes les valeurs p sont > 0.05 => toutes les variances sont homogènes.
Les conditions de la validité d’une ANOVA étant remplies, les interprétations seront donc valides. Le cas de la séance 2 vera une vérification supplémentaire.
4.5 ANOVA à 1 facteur séance par séance
4.5.1 Séance 1
4.5.1.1 Le modèle
## Anova Table (Type II tests)
##
## Response: masse_oeuf
## Sum Sq Df F value Pr(>F)
## regime 20.9072 9 5.7365 0.0005679 ***
## Residuals 8.0991 20
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
La p-value < 0.01 => différence très significative entre les effets d’au moins 2 régimes sur la masse de l’oeuf à la séance 1.
4.5.1.2 Comparaisons par paires
Comparaisons des moyennes par paires, Student - Newman - Keuls.
cm1 <- (SNK.test(lm1, "regime", group = TRUE))$groups %>%
mutate(regime = rownames(.)) %>%
select(regime, masse_oeuf, groups) %>%
as_tibble()
cm1
## # A tibble: 10 x 3
## regime masse_oeuf groups
## <chr> <dbl> <chr>
## 1 Ba 10 12.2 a
## 2 Ba 0,75 11.5 ab
## 3 WC 11.5 ab
## 4 Ba 1 11.2 ab
## 5 YC 11.1 ab
## 6 Ba 2,5 10.9 ab
## 7 Ba 0,25 10.7 ab
## 8 Ba 0,50 10.4 b
## 9 Ba 5 10.0 bc
## 10 Ba 7,5 9.03 c
4.5.1.3 Visualisation des groupes, bareplots avec labels
Figure pas nécessaire.
ggplot(data = cm1, mapping = aes(x = regime, y = masse_oeuf)) +
geom_bar(stat = "identity", color = "blue", fill = "grey", width = 0.6) +
geom_text(aes(label = groups), vjust = -0.5, size = 4) +
ylim(0, 13) +
xlab("Régimes") + ylab("Masse de l'oeuf (g)") +
theme(axis.text.x = element_text(angle = 45, color = "black", vjust = 1, hjust = 1)) +
theme_bw()
4.5.2 Séance 2
4.5.2.1 Le modèle
## # A tibble: 1 x 7
## id seance regime no_oeuf masse_oeuf is.outlier is.extreme
## <fct> <fct> <fct> <dbl> <dbl> <lgl> <lgl>
## 1 19 seance 2 Ba 7,5 1 19.9 TRUE TRUE
Excluons cette observation pour apprécier.
## Anova Table (Type II tests)
##
## Response: masse_oeuf
## Sum Sq Df F value Pr(>F)
## regime 17.0023 9 5.7435 0.0006788 ***
## Residuals 6.2494 19
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
La p-value est < 0.01 => Différence très significative entre les effetes d’au moins 2 régimes.
4.5.2.2 Comparaisons par paires
cm2 <- (SNK.test(lm2, "regime", group = TRUE))$groups %>%
mutate(regime = rownames(.)) %>%
select(regime, masse_oeuf, groups) %>%
as_tibble()
cm2
## # A tibble: 10 x 3
## regime masse_oeuf groups
## <chr> <dbl> <chr>
## 1 Ba 0,75 12.1 a
## 2 Ba 7,5 11.3 ab
## 3 YC 11.2 abc
## 4 Ba 1 10.9 abc
## 5 Ba 0,25 10.5 bc
## 6 Ba 5 10.4 bc
## 7 Ba 10 10.0 bc
## 8 WC 9.83 bc
## 9 Ba 2,5 9.76 bc
## 10 Ba 0,50 9.60 c
4.5.2.3 Visualisation des groupes
ggplot(data = cm2, mapping = aes(x = regime, y = masse_oeuf)) +
geom_bar(stat = "identity", color = "blue", fill = "grey", width = 0.6) +
ylim(0, 13) +
geom_text(aes(label = groups), vjust = -0.5, size = 4) +
xlab("Régimes") + ylab("Masse de l'oeuf (g)") +
theme(axis.text.x = element_text(angle = 45, color = "black", vjust = 1, hjust = 1)) +
theme_bw()
4.5.3 Séance 3
4.5.3.1 Le modèle
## Anova Table (Type II tests)
##
## Response: masse_oeuf
## Sum Sq Df F value Pr(>F)
## regime 34.950 9 8.5844 3.592e-05 ***
## Residuals 9.047 20
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
La p-value < 0.05 => Différence significative entre les effetes d’au moins 2 régimes.
4.5.3.2 Comparaisons par paires
cm3 <- (SNK.test(lm3, "regime", group = TRUE))$groups %>%
mutate(regime = rownames(.)) %>%
select(regime, masse_oeuf, groups) %>%
as_tibble()
cm3
## # A tibble: 10 x 3
## regime masse_oeuf groups
## <chr> <dbl> <chr>
## 1 WC 12.9 a
## 2 Ba 1 11.9 ab
## 3 Ba 0,25 11.3 bc
## 4 Ba 0,75 11.3 bc
## 5 YC 11.2 bc
## 6 Ba 7,5 10.7 bc
## 7 Ba 10 10.5 bc
## 8 Ba 2,5 10.1 bc
## 9 Ba 0,50 9.96 c
## 10 Ba 5 8.74 d
4.5.3.3 Visualisation des groupes
ggplot(data = cm3, mapping = aes(x = regime, y = masse_oeuf)) +
geom_bar(stat = "identity", color = "blue", fill = "grey", width = 0.6) +
geom_text(aes(label = groups), vjust = -0.5, size = 4) +
ylim(0, 14) +
xlab("Régimes") + ylab("Masse de l'oeuf (g)") +
theme(axis.text.x = element_text(angle = 45, color = "black", vjust = 1, hjust = 1)) +
theme_bw()
4.5.4 Séance 4
4.5.4.1 Le modèle
## Anova Table (Type II tests)
##
## Response: masse_oeuf
## Sum Sq Df F value Pr(>F)
## regime 10.145 9 2.105 0.07953 .
## Residuals 10.710 20
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
La p-value > 0.05 => différence pas significative entre les effets des régimes.
4.5.4.2 Comparaisons par paires, séance 4
cm4 <- (SNK.test(lm4, "regime", group = TRUE))$groups %>%
mutate(regime = rownames(.)) %>%
select(regime, masse_oeuf, groups) %>%
as_tibble()
cm4
## # A tibble: 10 x 3
## regime masse_oeuf groups
## <chr> <dbl> <chr>
## 1 YC 11.4 a
## 2 Ba 1 11.2 ab
## 3 Ba 2,5 11.0 ab
## 4 WC 10.5 ab
## 5 Ba 10 10.5 ab
## 6 Ba 0,75 10.3 ab
## 7 Ba 5 10.3 ab
## 8 Ba 0,50 10.3 ab
## 9 Ba 0,25 10.2 ab
## 10 Ba 7,5 9.23 b
Bizare ! j’ai essayé avec les comparaisons de Tukey qui me sortent également ces différences, pendant que l’ANOVA estime le contraire.
4.5.4.3 Visualisation des groupes
ggplot(data = cm4, mapping = aes(x = regime, y = masse_oeuf)) +
geom_bar(stat = "identity", color = "blue", fill = "grey", width = 0.6) +
ylim(0, 12) +
geom_text(aes(label = groups), vjust = -0.5, size = 4) +
xlab("Régimes") + ylab("Masse de l'oeuf (g)") +
theme(axis.text.x = element_text(angle = 45, color = "black", vjust = 1, hjust = 1)) +
theme_bw()
4.5.5 Séance 5
4.5.5.1 Le modèle
## Anova Table (Type II tests)
##
## Response: masse_oeuf
## Sum Sq Df F value Pr(>F)
## regime 13.056 9 1.5069 0.2124
## Residuals 19.254 20
La p-value > 0.05 => différence pas significative entre les effets des régimes.
4.5.5.2 Comparaisons par paires, séance 5
cm5 <- (SNK.test(lm5, "regime", group = TRUE))$groups %>%
mutate(regime = rownames(.)) %>%
select(regime, masse_oeuf, groups) %>%
as_tibble()
cm5
## # A tibble: 10 x 3
## regime masse_oeuf groups
## <chr> <dbl> <chr>
## 1 Ba 10 12.1 a
## 2 Ba 1 11.7 a
## 3 YC 11.6 a
## 4 Ba 0,25 11.4 a
## 5 Ba 2,5 10.9 a
## 6 Ba 7,5 10.5 a
## 7 WC 10.5 a
## 8 Ba 0,75 10.4 a
## 9 Ba 0,50 10.3 a
## 10 Ba 5 10.1 a
4.5.5.3 Visualisation
ggplot(data = cm5, mapping = aes(x = regime, y = masse_oeuf)) +
geom_bar(stat = "identity", color = "blue", fill = "grey", width = 0.6) +
ylim(0, 13) +
geom_text(aes(label = groups), vjust = -0.5, size = 4) +
xlab("Régimes") + ylab("Masse de l'oeuf (g)") +
theme(axis.text.x = element_text(angle = 45, color = "black", vjust = 1, hjust = 1)) +
theme_bw()
4.6 Évolution de la masse de l’oeuf par régime au cours du temps
4.6.1 Sommaire
mo_ic <- summarySE(mo,
measurevar = "masse_oeuf",
groupvars = c("seance", "regime"),
na.rm = TRUE)
mo_ic
## seance regime N masse_oeuf sd se ci
## 1 seance 1 Ba 0,25 3 10.739667 0.08084759 0.04667738 0.20083656
## 2 seance 1 Ba 0,50 3 10.418667 0.18138725 0.10472398 0.45059090
## 3 seance 1 Ba 0,75 3 11.531667 0.31750013 0.18330879 0.78871405
## 4 seance 1 Ba 1 3 11.196000 1.17814600 0.68020291 2.92667691
## 5 seance 1 Ba 10 3 12.188000 0.00000000 0.00000000 0.00000000
## 6 seance 1 Ba 2,5 3 10.885667 0.10700156 0.06177738 0.26580660
## 7 seance 1 Ba 5 3 10.041333 1.00905963 0.58258085 2.50664308
## 8 seance 1 Ba 7,5 3 9.034000 0.17617889 0.10171693 0.43765262
## 9 seance 1 WC 3 11.451000 1.20384426 0.69503981 2.99051492
## 10 seance 1 YC 3 11.139000 0.10652230 0.06150068 0.26461606
## 11 seance 2 Ba 0,25 3 10.494000 0.29480163 0.17020380 0.73232784
## 12 seance 2 Ba 0,50 3 9.600333 0.40317531 0.23277338 1.00154300
## 13 seance 2 Ba 0,75 3 12.113667 0.64016899 0.36960174 1.59026793
## 14 seance 2 Ba 1 3 10.947333 0.11727034 0.06770606 0.29131568
## 15 seance 2 Ba 10 3 10.039333 0.18741487 0.10820403 0.46556435
## 16 seance 2 Ba 2,5 3 9.757000 0.28391372 0.16391766 0.70528078
## 17 seance 2 Ba 5 3 10.418333 1.48207602 0.85567699 3.68168094
## 18 seance 2 Ba 7,5 3 14.160333 4.96930542 2.86902982 12.34443899
## 19 seance 2 WC 3 9.828000 0.12134661 0.07005950 0.30144169
## 20 seance 2 YC 3 11.154000 0.32261897 0.18626415 0.80142995
## 21 seance 3 Ba 0,25 3 11.348000 1.15412608 0.66633500 2.86700811
## 22 seance 3 Ba 0,50 3 9.961333 0.59706142 0.34471357 1.48318278
## 23 seance 3 Ba 0,75 3 11.269000 0.85362697 0.49284176 2.12052694
## 24 seance 3 Ba 1 3 11.860000 0.69425716 0.40082956 1.72463039
## 25 seance 3 Ba 10 3 10.483000 0.00000000 0.00000000 0.00000000
## 26 seance 3 Ba 2,5 3 10.127667 0.36750011 0.21217629 0.91292089
## 27 seance 3 Ba 5 3 8.735667 0.25350016 0.14635839 0.62972932
## 28 seance 3 Ba 7,5 3 10.730000 1.06539617 0.61510677 2.64659080
## 29 seance 3 WC 3 12.892000 0.00000000 0.00000000 0.00000000
## 30 seance 3 YC 3 11.176000 0.53869101 0.31101340 1.33818264
## 31 seance 4 Ba 0,25 3 10.156333 0.52538874 0.30333333 1.30513799
## 32 seance 4 Ba 0,50 3 10.329667 0.64946234 0.37496726 1.61335390
## 33 seance 4 Ba 0,75 3 10.345667 1.11116710 0.64153263 2.76029211
## 34 seance 4 Ba 1 3 11.195667 0.03901709 0.02252653 0.09692383
## 35 seance 4 Ba 10 3 10.507000 0.49892384 0.28805381 1.23939553
## 36 seance 4 Ba 2,5 3 11.034667 1.01361548 0.58521117 2.51796043
## 37 seance 4 Ba 5 3 10.343000 0.00000000 0.00000000 0.00000000
## 38 seance 4 Ba 7,5 3 9.231333 1.35164949 0.78037519 3.35768346
## 39 seance 4 WC 3 10.525000 0.46992340 0.27131040 1.16735443
## 40 seance 4 YC 3 11.398000 0.31091317 0.17950580 0.77235113
## 41 seance 5 Ba 0,25 3 11.410333 1.64365031 0.94896195 4.08305371
## 42 seance 5 Ba 0,50 3 10.333333 0.44889011 0.25916683 1.11510486
## 43 seance 5 Ba 0,75 3 10.377333 1.10940044 0.64051264 2.75590347
## 44 seance 5 Ba 1 3 11.660000 0.56351664 0.32534648 1.39985293
## 45 seance 5 Ba 10 3 12.112000 0.00000000 0.00000000 0.00000000
## 46 seance 5 Ba 2,5 3 10.876333 0.59530860 0.34370158 1.47882855
## 47 seance 5 Ba 5 3 10.143000 1.10098683 0.63565504 2.73500290
## 48 seance 5 Ba 7,5 3 10.490000 1.24203019 0.71708647 3.08537404
## 49 seance 5 WC 3 10.471333 0.84060831 0.48532543 2.08818681
## 50 seance 5 YC 3 11.643000 1.16600000 0.67319041 2.89650457
4.6.2 Visualisation
ggplot(mo_ic, aes(x = seance, y = masse_oeuf, colour = regime, group = regime)) +
geom_line(size = 1) +
geom_point(size = 2) +
ylab("Masse de l'oeuf (g)") +
theme_bw()
Il ne semble pas se dégager une tendance concrète entre les différentes séances.
Nous savons par les analyses pour chaque séance plus haut, que
- séance 1 : il existe des différences d’effet entre les régimes
- séance 2 : il existe des différences d’effet entre les régimes (1 valeur extrême exclue)
- séance 3 : il existe des différences d’effet entre les régimes
- séance 4 : pas de différences signicatives, mais des groupes sont constitués, p-value = 0.07
- séance 5 : pas de différences signicatives
Puisque les données ne répondent pas aux conditions pour évaluer les effets des régimes au cours du temps, on négligera l’effet des régimes pour évaluer globalement l’effet du temps sur les masses d’oeuf.
On pourrait se demander si les masses d’oeuf mesurées sur l’ensemble des sujets sont significativement différentes d’une séance à l’autre (c’est-à-dire avec le temps).
4.6.3 Effet du temps
4.6.3.1 boxplots, facteur temps
4.6.3.2 Valeurs aberrantes, facteur temps
mo <- mo %>% mutate(id2 = 1:nrow(.), .before = 1)
mo_out <- mo %>%
group_by(seance) %>%
identify_outliers(masse_oeuf)
mo_out
## # A tibble: 3 x 8
## seance id2 id regime no_oeuf masse_oeuf is.outlier is.extreme
## <fct> <int> <fct> <fct> <dbl> <dbl> <lgl> <lgl>
## 1 seance 2 49 19 Ba 7,5 1 19.9 TRUE TRUE
## 2 seance 4 109 19 Ba 7,5 1 8.46 TRUE FALSE
## 3 seance 4 110 20 Ba 7,5 2 8.44 TRUE FALSE
Il y a une observation aberrante extrême pour la séance 2. On pourrait l’exclure. Mais je l’ai conservé.
4.6.3.3 Homogénéité des variances et ANOVA, facteur temps
Les autres conditions ont déjà été vérifiées. La fonction anova_test()
réalise également le test de sphéricité de Mauchly.
lm <- anova_test(data = mo,
dv = masse_oeuf, # dependant variable, num
wid = id, # identificateur de cas/échantillon (facteur)
within = seance) # facteur de groupement intra-sujets
get_anova_table(lm)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 seance 2.18 63.36 0.562 0.588 0.015
=> C’est la p-value qui nous intéresse et elle est > 0.05 => pas de différence d’une séance à l’autre.
4.6.3.4 Comparaisons par paires, facteur temps
tph <- mo %>%
pairwise_t_test(masse_oeuf ~ seance,
paired = TRUE,
p.adjust.method = "bonferroni")
tph %>%
select(group1, group2, p, p.adj, p.adj.signif)
## # A tibble: 10 x 5
## group1 group2 p p.adj p.adj.signif
## <chr> <chr> <dbl> <dbl> <chr>
## 1 seance 1 seance 2 0.98 1 ns
## 2 seance 1 seance 3 0.986 1 ns
## 3 seance 1 seance 4 0.075 0.753 ns
## 4 seance 1 seance 5 0.697 1 ns
## 5 seance 2 seance 3 0.987 1 ns
## 6 seance 2 seance 4 0.432 1 ns
## 7 seance 2 seance 5 0.817 1 ns
## 8 seance 3 seance 4 0.181 1 ns
## 9 seance 3 seance 5 0.692 1 ns
## 10 seance 4 seance 5 0.039 0.393 ns