L’analyse de variance (ANOVA) est une méthode statistique utilisée pour comparer les moyennes de trois groupes ou plus. Elle permet de déterminer si les moyennes sont globalement différentes, ou plus précisément si au moins deux moyennes sont significativement différentes.
Le principe de cette approche est de décomposer la dispersion totale des données en une part factorielle (c’est-à-dire imputable aux groupes – on l’appelle aussi la dispersion inter-groupe) et en une part résiduelle (le reste – appelée aussi dispersion intra-groupes). Puis de comparer la variance inter-groupe à la variance intra-groupe par l’intermédiaire d’un test F de Fisher Snedecor.
Si la variance inter-groupe est suffisamment supérieure à la variance intra-groupe, alors cela suggère que les moyennes des groupes sont globalement différentes. Autrement formulé, cela suggère qu’il existe un effet du groupe sur la variable dépendante étudiée.
La puissance du test statistique de F de Fisher Snedecor (qui permet de comparer la variance inter-groupe à la variance intra-groupe) est sa capacité à mettre en évidence une différence significative entre les moyennes, si cette différence existe réellement. Comme pour tous les tests statistiques, pour un niveau de risque alpha fixé, la puissance statistique du test F augmente avec l’augmentation du nombre d’unités expérimentales observées au sein de chacun des k groupes.
C’est une bonne nouvelle ! La moins bonne, c’est que généralement, les coûts humains et financiers augmentent aussi avec l’augmentation du nombre d’unités expérimentales. Il est donc nécessaire de faire un compromis.
C’est dans ce contexte que savoir calculer le nombre minimum d’unités expérimentales nécessaires pour détecter une différence significative entre les groupes étudiés (avec un niveau de confiance et une puissance de test spécifiés) est essentiel.
L’objectif de cet article est de vous montrer comment réaliser ce calcul du nombre minimum d’unités expérimentales nécessaires, dans le cadre d’une comparaison d’au moins 3 moyennes, en utilisant les résultats d’une expérimentation pilote.
Pour plus d’informations sur l’anova à un facteur et la puissance statistique, vous pouvez consulter les articles suivants :
Le calcul du nombre d’unités expérimentales nécessaires pour une ANOVA à un facteur implique plusieurs étapes importantes, notamment le calcul de l’effect size f.
L’effect size est une mesure standardisée qui permet de quantifier la différence entre les moyennes des groupes étudiés.
Il existe plusieurs effect size dans le cadre de l’ANOVA à un facteur, nommés eta carré et f. Leur calcul se fait facilement à partir de la table de résultats de l’ANOVA, ou en employant des fonctions spécifiques qui les retournent directement
Voici les différentes étapes à suivre :
1. Réaliser l’ANOVA à un facteur à partir des données de l’expérimentation pilote
2. Calculer l’effect size f
3. Réaliser le calcul du nombre minimum d’unités expérimentales nécessaires en employant l’effect size f, et en fixant le risque alpha et la puissance souhaités
Pour illustrer ce calcul, j’ai utilisé le jeu de données iris tronqué aux 5 premières lignes de chaque espèce.
Nous allons chercher à déterminer le nombre minimum de fleurs à observer à l’intérieur de chaque espèce pour avoir une probabilité de 80% de mettre en évidence une différence entre les moyennes de ces 3 espèces.
# Installer et charger les packages nécessaires
#install.packages("tidyverse")
library(tidyverse)
# Charger les données iris
data(iris)
# Utiliser la fonction slice_head() du package dplyr pour sélectionner les 5 premières lignes de chaque espèce
iris_short <- iris %>%
group_by(species) %>%
slice_head(n=5)
# affichage des données
iris_short
## Sepal.Length Sepal.Width Petal.Length Petal.Width Species
## 1 5.1 3.5 1.4 0.2 setosa
## 2 4.9 3.0 1.4 0.2 setosa
## 3 4.7 3.2 1.3 0.2 setosa
## 4 4.6 3.1 1.5 0.2 setosa
## 5 5.0 3.6 1.4 0.2 setosa
## 6 7.0 3.2 4.7 1.4 versicolor
## 7 6.4 3.2 4.5 1.5 versicolor
## 8 6.9 3.1 4.9 1.5 versicolor
## 9 5.5 2.3 4.0 1.3 versicolor
## 10 6.5 2.8 4.6 1.5 versicolor
## 11 6.3 3.3 6.0 2.5 virginica
## 12 5.8 2.7 5.1 1.9 virginica
## 13 7.1 3.0 5.9 2.1 virginica
## 14 6.3 2.9 5.6 1.8 virginica
## 15 6.5 3.0 5.8 2.2 virginica
Nous nous intéressons ici à l’effet de l’espèce sur la largeur des sépales
model <- aov(Sepal.Width~Species, data=iris_short)
summary(model)
## Df Sum Sq Mean Sq F value Pr(>F)
## Species 2 0.372 0.186 2.138 0.161
## Residuals 12 1.044 0.087
On peut facilement obtenir les effect size eta carré et f avec la fonction parameters()
du package parameters
:
library(parameters)
parameters(model,effectsize_type = c("eta", "f" ))
## Parameter | Sum_Squares | df | Mean_Square | F | p | Eta2 | Cohen's f
## ----------------------------------------------------------------------------
## Species | 0.37 | 2 | 0.19 | 2.14 | 0.161 | 0.26 | 0.60
## Residuals | 1.04 | 12 | 0.09 | | | |
##
## Anova Table (Type 1 tests)
On peut également les obtenir, à partir de la table de résultats, en employant les équations suivantes :
\[ \eta^2 = \frac{SCF}{SCT} = \frac{0.37}{1.04+0.37}=0.26 \]
\[ f=\sqrt{\eta^2 /\left(1-\eta^2\right)}=\sqrt{0.26 /\left(1-0.26\right)}=0.6\]
Pour cette étape, nous allons employer la fonction pwr.anova.test()
du package pwr
.
Nous allons fixer le risque alpha à 5% et la puissance désirée à 80%.
# Installer et charger les packages nécessaires
#install.packages("pwr")
library(pwr)
# Définir les paramètres de la puissance de l'analyse
alpha <- 0.05 # Niveau de signification
power <- 0.8 # Puissance souhaitée
f <- 0.6 # Effect size f
# Calculer le nombre d'unités expérimentales nécessaires
n <- pwr.anova.test(k = 3, f = f, sig.level = alpha, power = power)
n
##
## Balanced one-way analysis of variance power calculation
##
## k = 3
## n = 9.990612
## f = 0.6
## sig.level = 0.05
## power = 0.8
##
## NOTE: n is number in each group
Ici, les résultats nous indiquent que pour mettre en évidence un effect size f de 0.6 (observé sur l’expérimentation pilote), avec une probabilité (puissance) de 80%, sachant que le risque alpha est fixé à 5%, il est nécessaire de disposer d‘au moins 10 fleurs par espèce.
Afin de visualiser l’évolution de la puissance en fonction du nombre d’unités expérimentales, nous pouvons réaliser des calculs inverses : calculer la puissance théorique, pour un nombre de répliquats allant de 5 à 50 par pas de 5, par exemple :
# Définition d'une séquence de nombres allant de 5 à 50 par incréments de 1 et stockage dans l'objet "n"
n <- seq(5,50,by=1)
# Calcul de la puissance statistique pour une ANOVA à un facteur avec trois groupes en fonction de n
# Les valeurs de f et sig.level doivent être définies au préalable
# La fonction pwr.anova.test du package "pwr" est utilisée pour effectuer le calcul
puissance <- pwr.anova.test(k = 3, f = f, sig.level = alpha,n=n)$power
# Création d'un objet data frame "df" contenant les valeurs de n et de puissance
df <- data.frame(n, puissance)
# Création d'un graphique à partir de l'objet "df" en utilisant la bibliothèque ggplot2
# La fonction geom_point ajoute des points au graphique correspondant à chaque paire (n, puissance)
# La fonction geom_hline ajoute deux lignes horizontales au graphique à y = 0.8 (en bleu) et y = 0.9 (en rouge)
ggplot(df, aes(x=n, y=puissance))+
geom_point()+
geom_hline(yintercept=0.8, colour="blue", linetype=2)+
geom_hline(yintercept=0.9, colour="red",linetype=2)
Vous l’aurez compris, déterminer le nombre minimum d’unités expérimentales nécessaires est une étape cruciale dans la planification d’une étude utilisant une ANOVA à un facteur. Cette étape permet de s’assurer que l’étude sera suffisamment puissante pour détecter des différences significatives entre les moyennes des groupes étudiés.
Lorsque j’anime des formations en biostatistiques, c’est une question que l’on me pose très régulièrement, j’espère donc que cet article aidera le plus grand nombre !
Retrouver le planning et les programmes de mes formations ici 👇👇👇
Retrouver mes propositions de services ici 👇👇👇
C’est possible en faisant un don sur la page Tipeee du blog
Enregistrez vous pour recevoir gratuitement mes fiches “aide mémoire” (ou cheat sheets) qui vous permettront de réaliser facilement les principales analyses biostatistiques avec le logiciel R et pour être informés des mises à jour du site.
11 Responses
Merci Claire pour ce post, toujours très intéressant !
Pour la sélection des 5 premières valeurs de chaque espèce, nous pouvons utiliser la fonction slice_head du tidyverse :
iris%>%
group_by(Species)%>%
slice_head(n=5)
L’ensemble des fonctions slice et notamment slice_sample sont très pratiques !
Bonjour Lilian,
Vous venez de me prendre en flagrant délit de sur-utilisation des fonctions map ! Comme je ne sais pas bien me servir de purrr, et des fonctions map, j’essaye de les employer le plus possible !
Votre suggestion est tellement plus élégante et simple, que j’ai modifié le code de l’article.
Merci.
Bonjour Madame Claire.
Je suis très ravi d’apprendre ce calcul mais quelques questions me viennent à l’esprit.
Dans l’exemple donné ci-haut, les résultats obtenus viennent d’une simulation à partir des données existantes.
Alors comment obtenir le nombre minimum d’unités expérimentales si ces données là n’existaient pas.
Merci
Bonjour Benjamin,
Bonjour,
vous pouvez vous aider des seuils d’effect size classiquement employés (small : 0.1, medium : 0.25, large : 0.4) (https://imaging.mrc-cbu.cam.ac.uk/statswiki/FAQ/effectSize), et en utilisant une puissance de 0.8 (minimum classique) et 0.9.
J’espère que cela vous aide.
Bonne continuation.
Bonjour Claire,
J’ai beaucoup apprécié votre article et j’ai une question à vous poser. Je souhaite commencer une expérimentation mais je n’ai pas encore de données. Comment puis-je déterminer les paramètres de f et de puissance pour obtenir le nombre minimal d’unités expérimentales nécessaires ?
Je vous remercie d’avance pour votre aide et vous adresse mes salutations cordiales.
Bonjour,
vous pouvez vous aider des seuils d’effect size classiquement employés (small : 0.1, medium : 0.25, large : 0.4) (https://imaging.mrc-cbu.cam.ac.uk/statswiki/FAQ/effectSize), et en utilisant une puissance de 0.8 (minimum classique) et 0.9.
J’espère que cela vous aide.
Bonne continuation.
Bonjour et merci beaucoup pour cet article très intéressant, très clair, qui répond à une question cruciale qu’on se pose toujours. Je sais enfin comment faire pour y répondre et j’ai hâte de mettre cela en application ! Mille merci pour cela !
Je me posais une petite question par curiosité sur l’un des termes de la fonction pwr.anova.test où je vois qu’il faut renseigner que k = 3. A quoi correspond ce k et sa valeur est-elle constante ou bien dépend-t-elle de quelque chose ?
Merci beaucoup !
Salut, k corresponds au nombre de modalités/ traitements différents, ici on a 3 espèces
Merci bien
Merci bien à vous
Merci Claire