Les tests multiples

FIGURE68

Contexte:

Le grand principe de la théorie des tests est de calculer la probabilité de rejeter  l’hypothèse H_0 tout en minimisant les chances de se tromper. Ainsi, quand nous comparons la p-valeur, associée à une statistique de test, à la loi de probabilité qu’elle suit et que nous rejetons H_0 au risque de 5\%, dans la forme, nous prenons une décision mais, sur le fond, nous prenons également le risque de nous tromper avec une probabilité de 5\%.

Imaginons désormais que pour une longue série de variables nous appliquons le même test au risque \alpha de 5\%. Il paraît évident que nous augmentons les chances de conclure à tort au rejet de H_0 au moins pour un certain nombre de variables testées.

Justification:

En effet, avec la formule de Poincaré nous démontrons d’ailleurs toute l’exactitude du phénomène décrit dans le précèdent paragraphe:

P(\cup_{k = 1} ^K A_k) = \sum_k P(A_k) - \sum_{k_1 < k_2} P(A_{k_1} \cap A_{k_2}) + \\ \hspace*{25mm} \sum_{k_1 < k_2 < k_3} P(A_{k_1} \cap A_{k_2} \cap A_{k_3}) + \cdots + (-1) ^{K - 1} \cdot P(\cap_k A_k)

Où l’événement A_k = « le test N° k rejette à tort H_0 au seuil \alpha« . Nous avons alors P(A_k) = \alpha et grâce à l’indépendance des tests nous pouvons développer l’équation ci-dessus en:

P(\cup_{k = 1} ^K A_k) = K \times \alpha - \sum_{k_1 < k_2} \alpha ^2 + \cdots + (-1) ^{K - 1} \cdot \alpha ^K

Nous avons donc définit l’équation qui dépend du nombre K de variables à tester et du seuil \alpha de rejet de H_0. En faisant varier ces deux paramètres nous pouvons tracer le graphe qui nous permet de voir l’évolution de la probabilité de rejetter à tort H_0:

FIGURE89

Nous voyons ainsi qu’au risque \alpha de 5\% et pour 5 variables testées uniquement, nous avons un peu plus de 20\% de chance de rejeter H_0 à tort pour l’une d’elles! Face à un résultat aussi fort et aussi inquiétant, il est une solution à laquelle énormément de chercheurs ont contribué (dont les principaux sont Carlo Emilio Bonferroni, Yoav Benjamini et Yosef Höchberg): la théorie des tests multiples (nommée également comparaisons multiples).

Le grand principe des tests multiples est donc de corriger la liste des p-valeurs établies afin de rattraper les hypothèses H_0 rejetées à tort.

Définition du problème:

Définissons K le nombre d’hypothèses nulles H_1, \cdots, H_K que nous pourrions également assimiler au test de l’hypothèse H_0 pour K variables X ^1, \cdots, X ^K. Remarquons que l’écriture H_k régulièrement utilisées est un abus de langage voulant désigner en réalité l’hypothèse H_0 testée pour la k-ième fois et d’écriture exacte H_0 ^k.

Enfin, nous poserons R le nombre d’hypothèses nulles rejetées et m_0 le nombre d’hypothèses nulles vraies. 

add

Avec nos notations, précisons que K - m_0 représente le nombre d’hypothèses alternatives vraies.

La visualisation du problème au travers du tableau ci-dessus permet d’écrire les deux taux d’erreur de type I et ainsi de définir les deux grandes familles de tests multiples:

– (FamilyWise Error Rate) FWER = P(FP > 0) = 1 - P(FP = 0)

– (False Discovery Rate) FDR = E[\frac{FP}{R} | R > 0]

Et le taux d’erreur de type II: (Non Discovery Rate) NDR = E[\frac{FN}{K - R}]

Par définition, le FWER consiste à estimer la probabilité d’avoir au moins un faux positif, tandis que le FDR détermine le taux moyen d’avoir un faux positif. 

Chacun des deux taux d’erreur de type I peut se décliner sous différentes variantes que nous verrons dans la partie suivante.

Les procédures correctives de p-valeurs: 

Nous présentons ici la liste des procédures de tests multiples avec l’algorithme utilisé et l’application sur l’exemple qui suit. Imaginons la liste de p-valeurs suivantes:

p = [p_1 < 0.0001, p_2 = 0.0004, p_3 = 0.0007, p_4 = 0.6969, p_5 = 0.8129]

= [<0.0001, 0.0004, 0.0007, 0.6969, 0.8126]

p_k représente la p-valeur associée à la statistique de test calculée. Enfin, notons \alpha le seuil de significativité fixé au préalable.

 La procédure de Bonferroni (FWER):

Historique: La procédure de Bonferroni a été inventée par Carlo Emilio Bonferroni  en 1936. Elle porte également le nom de correction de Bonferroni.

Principe: Rejet de H_k si p_k \leq \frac{\alpha}{K} \Rightarrow \alpha \geq K \cdot p_k.

Algorithme associé:

– Ranger p par ordre croissant

\forall k \in [1,K], p_k * = min(K \times p_k, 1)

p * = [p_1 *, \cdots, p_K *]

Exemple:

Nous avons donc p = [<0.0001, 0.0004, 0.0007, 0.6969, 0.8129]. Comme nous avons 5 p-valeurs à ajuster, K = 5. Donc,

p * = [min(<0.0001 \times 5, 1), min(0.0004 \times 5, 1), min(0.0007 \times 5, 1), min(0.6969 \times 5, 1), min(0.8129 \times 5, 1)]

= [min(<0.0001, 1), min(0.002, 1), min (0.0035, 1), min (4.0645, 1), min(3.4845, 1)]

D’où,

p * = [<0.0001, 0.002, 0.0035, 1, 1]

Constat: Procédure jugée parfois trop conservatrice et peu puissante.

♦ La procédure step-down ajustée de Bonferroni (FWER):

Historique: La procédure step-down de Bonferroni a été inventée par Carlo Emilio Bonferroni en 1979.

Principe: Rejet de H_1, \cdots, H_k pour le plus grand k tel que p_k \leq \frac{\alpha}{K} \Rightarrow \alpha \geq \alpha \cdot K \cdot p_k.

Algorithme associé:

– Ranger p par ordre croissant

p_1 * = min(K \cdot p_1, 1)

\forall k \in [2, K], p_k * = min(max(p_{k - 1} *, (K - k + 1) \cdot p_k), 1)

p * = [p_1 *, \cdots, p_K *]

Exemple:

Nous avons donc p = [<0.0001, 0.0004, 0.0007, 0.6969, 0.8129]. Comme nous avons 5 p-valeurs à ajuster, K = 5. Donc,

– Etape 1, k = 1,

p_1 * = min(5 \times "<0.0001", 1) = "<0.0001"

– Etape 2, k = 2,

p_2 * = min(max(<0.0001, 4 \times 0.0004), 1)

 = min(max(<0.0001, 0.0016), 1)

 = 0.00016

– Etape 3, k = 3,

p_3 * = min(max(0.0016, 3 \times 0.0007), 1)

 = min(max(0.0016, 0.0021), 1)

 = 0.0021

– Etape 4, k = 4,

p_4 * = min(max(0.0021, 2 \times 0.6969), 1) = min(max(0.0021, 1.3938), 1) = 1

– Etape 5, k = 5,

p_5 * = min(max(1, 1 \times 0.8129), 1) = min(max(1, 08129), 1) = 1

D’où,

p * = [<0.0001, 0.0016, 0.0021, 1, 1]

Constat: Les méthodes step-down sont réputées comme étant plus puissantes que les méthodes d’origine.

 La procédure de Holm (FWER):

Historique: La procédure de Holm a été inventée par Sture Holm en 1978 et se base sur la correction de Bonferroni. Elle porte également le nom de procédure step-down de Bonferroni ou de procédure de Holm-Bonferroni.

Principe: Rejet de H_1, \cdots, H_{k - 1} pour le plus petit k tel que p_k > \frac{\alpha}{K + 1 - k} \Rightarrow \alpha < (K + 1 - k) \cdot p_k.

Algorithme associé:

– Ranger p par ordre croissant

p_1 * = min(K \times p_1, 1)

\forall k \in [2,K], p_k * = min(max(p_{k - 1} *, (K + 1 - k) \times p_k), 1)

p * = [p_1 *, \cdots, p_K *]

Exemple:

Nous avons donc p = [<0.0001, 0.0004, 0.0007, 0.6969, 0.8129]. Comme nous avons 5 p-valeurs à ajuster, K = 5. Donc,

Etape 1, k = 1,

p_1 * = min("<0.0001" \times 5, 1) = "<0.0001"

Etape 2, k = 2,

p_2 * = min(max(<0.0001, 0.0004 \times (5 - 2 + 1)), 1)

 = min(max(<0.0001, 0.0016), 1)

 = 0.0016

Etape 3, k = 3,

p_3 * = min(max(0.0016, 0.0007 \times (5 - 3 +1)), 1)

 = min(max(0.0016, 0.0021), 1)

 = 0.0021

Etape 4, k = 4,

p_4 * = min(max(0.0021, 0.6969 \times (5 - 4 + 1)), 1)

 = min(max(0.0021, 1.34938), 1)

 = 1

Etape 5, k = 5,

 p_5 * = min(max(1, 0.6969 \times (5 - 5 + 1)), 1) = min(max(1, 0.8129), 1) = 1

D’où, 

p * = [<0.0001, 0.0016, 0.0021, 1, 1]

Constat: Plus puissant que la correction de Bonferroni mais moins que les procédures d’Hommel ou encore d’Höchberg. Néanmoins, la procédure présente l’avantage de ne requérir aucune hypothèse sur l’indépendance.

♦ La procédure adaptative de Holm (FWER):

Historique: La procédure adaptative de Holm a été inventée par Sture Holm en 1990. Elle porte également le nom d’ajustement de la procédure step-up de Bonferroni.

Principe: Rejet de H_1, \cdots, H_{k - 1} pour le plus petit k tel que p_k > \frac{\alpha}{min(K - k + 1, \hat{m_0})} \Rightarrow \alpha < min(K - k + 1, \hat{m_0}) \cdot p_k.

Algorithme associé:

– Ranger p par ordre croissant

– Estimer \hat{m_0} le nombre d’hypothèses nulles vraies

p_1 * = min(min(K, \hat{m_0}) \cdot p_1, 1)

\forall k \in [2, K], p_k = min(max(p_{k - 1} *, min(K - k + 1, \hat{m_0}) \cdot p_k), 1)

p * = [p_1 *, \cdots, p_K *]

Exemple:

Nous avons donc p = [<0.0001, 0.0004, 0.0007, 0.6969, 0.8129]. Comme nous avons 5 p-valeurs à ajuster, K = 5.

Afin d’estimer \hat{m_0} plusieurs méthodes existent. Ci-dessous, l’ensemble des méthodes extraitent du site du suppot SAS:

addNous utiliserons la méthode par défaut du logiciel sus-nommé: decrease slope (ou pente décroissante) qui consiste à étudier la courbe des p-valeurs tracée par le vecteur q = 1 - p puis de calculer la totalité des b_i = a(q_K, \cdots, q_{K - i})a représente le coefficient de la fonction affine de type y = a \cdot x + b (y prenant la valeur de l’indice de la p-valeur).

Il s’agit alors de retenir de manière itérative décroissante le premier coefficient b_i qui répond à l’inégalité b_i > b_{i + 1}. Nous avons alors,

\hat{m_0} = ceil(\frac{1}{b_{i+1}} - 1)

, où ceil(.) représente la fonction renvoyant le premier entier supérieur à « . ».

Si nous appliquons cet algorithme à notre exemple nous obtenons,

b_1 = a(1 - p_5, 1 - p_5) = non calculable

b_2 = a(1 - p_5, 1 - p_4) = - \infty

b_3 = a(1 - p_5, 1 - p_4, 1 - p_3) = 0.0785

b_4 = a(1 - p_5, 1 - p_4, 1 - p_3, 1 - p_2) = 0.8502

b_5 = a(1 - p_5, 1 - p_4, 1 - p_3, 1 - p_2, 1 - p_1) = 0.2563

Il en ressort, si nous parcourons les b_i de manière décroissante, que b_4 > b_5 et donc pour l’indice i = 4 nous avons le premier cas qui répond à l’inégalité recherchée. Par conséquent,

\hat{m_0} = ceil(\frac{1}{0.2563} - 1) = ceil(2.901688) = 3

Nous pouvons donc dérouler l’algorithme adaptatif de Holm. Ainsi,

– Etape 1, k = 1,

p_1 * = min(min(5, 3) \times "<0.0001", 1)

 = min(3 \times "<0.0001", 1)

 = "<0.0001"

– Etape 2, k = 2,

p_2 * = min(max(<0.0001, min(5 - 2 +1, 3) \times 0.0004), 1)

 = min(max(<0.0001, 0.0012), 1)

 = 0.0012)

– Etape 3, k = 3,

p_3 * = min(max(0.0012, min(5 - 3 +1, 3) \times 0.0007), 1)

 = min(max(0.0012, 0.0021), 1)

 = 0.0021

– Etape 4, k = 4,

p_4 * = min(max(0.0021, min(5 - 4 + 1, 3) \times 0.6969), 1)

 = min(max(0.0021, 1.3938), 1)

 = 1

– Etape 5, k = 5,

p_5 * = min(max(1, min(5 - 5 +1, 3) \times 0.8129), 1)

 = min(max(1, 0.8129), 1)

 = 1

D’où,

p * = [<0.0001, 0.0012, 0.0021, 1, 1]

Constat: Les méthodes adaptatives sont réputées comme étant plus puissantes que les méthodes step-down/step-up à condition que le paramètre \hat{m_0} soit correctement estimé.

♦ La procédure d’Höchberg (FWER):

Historique: La procédure d’Höchberg  a été mise au point par Yosef Höchberg en 1988. Elle porte également le nom de procédure step-up d’Höchberg.

Principe: Rejet de H_1, \cdots, H_k pour R le plus grand k tel que p_k \leq \frac{\alpha}{K - k + 1} \Rightarrow \alpha \geq (K - k + 1) \cdot p_k

Algorithme associé:

– Ranger p par ordre croissant

p_K * = p_K

\forall k \in [K - 1, 1], p_k * = min((K - k + 1) \times p_k, p_{k + 1} *)

p * = [p_1 *, \cdots, p_K *]

Exemple:

Nous avons donc p = [<0.0001, 0.0004, 0.0007, 0.6969, 0.8129]. Comme nous avons 5 p-valeurs à ajuster, K = 5. Donc,

Etape 1, k = 5,

 p_5 * = p_5 = 0.8129

Etape 2, k = 4,

p_4 * = min(0.6969 \times (5 - 4 +1), 0.8129) = min(0.8129, 1.3938) = 0.8129

Etape 3, k = 3,

p_3 * = min(0.0007 \times (5 - 3 +1), 0.8129) = min(0.0021, 0.8129) = 0.0021

Etape 4, k = 2,

p_2 * = min(0.0004 \times (5 - 2 +1), 0.0021) = min(0.0016, 0.0021) = 0.0016

Etape 5, k = 1,

p_1 * = min("<0.0001" \times (5 - 1 +1), 0.0016)

 = min(<0.0001, 0.0016)

 = "<0.0001"

D’où,

p * = [<0.0001, 0.0016, 0.0021, 0.8129, 0.8129]

Constat: Plus puissant que la procédure de Holm. Néanmoins, requiert l’hypothèse d’indépendance des tests.

♦ La procédure adaptative de Höchberg (FWER):

Historique: La procédure adaptative de Höchberg a été inventée par Yosef Höchberg en 1990. Elle porte également le nom d’ajustement de la procédure step-down de Bonferroni.

Principe: Rejet de H_1, \cdots, H_k pour R le plus grand k tel que p_k \leq \frac{\alpha}{min(K - k + 1, \hat{m_0})} \Rightarrow \alpha \geq min((K - k + 1), \hat{m_0}) \cdot p_k

Algorithme associé:

– Ranger p par ordre croissant

– Estimer \hat{m_0} le nombre d’hypothèses nulles vraies

– Pour p_K = min(1, \hat{m_0}) \cdot p_K

\forall k \in [K - 1, 1], p_k = min(p_{k + 1} *, min(K - k + 1, \hat{m_0}) \cdot p_k)

p * = [p_1 *, \cdots, p_K *]

Exemple:

Nous avons donc p = [<0.0001, 0.0004, 0.0007, 0.6969, 0.8129]. Comme nous avons 5 p-valeurs à ajuster, K = 5.

Afin d’estimer \hat{m_0} plusieurs méthodes existent. Ci-dessous, l’ensemble des méthodes extraitent du site du suppot SAS:

addNous utiliserons la méthode par défaut du logiciel sus-nommé: decrease slope (ou pente décroissante) qui consiste à étudier la courbe des p-valeurs tracée par le vecteur q = 1 - p puis de calculer la totalité des b_i = a(q_K, \cdots, q_{K - i})a représente le coefficient de la fonction affine de type y = a \cdot x + b (y prenant la valeur de l’indice de la p-valeur).

Il s’agit alors de retenir de manière itérative décroissante le premier coefficient b_i qui répond à l’inégalité b_i > b_{i + 1}. Nous avons alors,

\hat{m_0} = ceil(\frac{1}{b_{i+1}} - 1)

, où ceil(.) représente la fonction renvoyant le premier entier supérieur à « . ».

Si nous appliquons cet algorithme à notre exemple nous obtenons,

b_1 = a(1 - p_5, 1 - p_5) = non calculable

b_2 = a(1 - p_5, 1 - p_4) = - \infty

b_3 = a(1 - p_5, 1 - p_4, 1 - p_3) = 0.0785

b_4 = a(1 - p_5, 1 - p_4, 1 - p_3, 1 - p_2) = 0.8502

b_5 = a(1 - p_5, 1 - p_4, 1 - p_3, 1 - p_2, 1 - p_1) = 0.2563

Il en ressort, si nous parcourons les b_i de manière décroissante, que b_4 > b_5 et donc pour l’indice i = 4 nous avons le premier cas qui répond à l’inégalité recherchée. Par conséquent,

\hat{m_0} = ceil(\frac{1}{0.2563} - 1) = ceil(2.901688) = 3

Nous pouvons donc dérouler l’algorithme adaptatif de Höchberg. Ainsi,

– Etape 1, k = 5,

p_5 * = min(1,3) \times 0.8129 = 1 \times 0.8129 = 0.8129

– Etape 2, k = 4,

p_4 * = min(0.8129, min(5 - 4 + 1, 3) \times 0.6969) = min (0.8129, 1.3938) = 0.8129

– Etape 3, k = 3,

p_3 * = min(0.8129, min(5 - 3 +1, 3) \times 0.0007) = min(0.8129, 0.0021) = 0.0021

– Etape 4, k = 2,

p_2 * = min(0.0021, min(5 - 2 + 1, 3) \times 0.0004) = min(0.0021, 0.0012) = 0.0012

– Etape 5, k = 1,

p_1 * = min(0.0021, min(5 - 1 + 1, 3) \times "<0.0001") = "<0.0001"

D’où,

p * = [<0.0001, 0.0012, 0.0021, 0.8129, 0.8129]

Constat: Les méthodes adaptatives sont réputées comme étant plus puissantes que les méthodes step-down/step-up à condition que le paramètre \hat{m_0} soit correctement estimé.

♦ La procédure de Sidak (FWER):

Historique: La procédure de Sidak a été mise au point par Sidak en 1967. Elle porte également le nom de procédure de Dunn-Sidak.

Principe:  Rejet de H_k si p_k \leq 1 - (1 - \alpha) ^{\frac{1}{K}} \Rightarrow \alpha \geq 1 - (1 - p_k) ^K.

Algorithme associé:

– Ranger p par ordre croissant

\forall k \in [1, K], p_k * = 1 - (1 - p_k) ^K

p * = [p_1 *, \cdots, p_K *]

Exemple:

Nous avons donc p = [<0.0001, 0.0004, 0.0007, 0.6969, 0.8129]. Comme nous avons 5 p-valeurs à ajuster, K = 5. Donc,

p_1 * = 1 - (1 - "<0.0001") ^5 = "<0.0001"

p_2 * = 1 - (1 - 0.0004) ^5 = 1 - 0.9996 ^5 = 0.001998401

p_3 * = 1 - (1 - 0.0007) ^5 = 1 - 0.9993 ^5 = 0.003495103

p_4 * = 1 - (1 - 0.6969) ^5 = 1 - 0.3031 ^5 = 0.9974418

p_5 * = 1 - (1 - 0.8129) ^5 = 1 - 0.1871 = 0.9997707

D’où,

p* = [<0.0001, 0.001998401, 0.003495103, 0.9974418, 0.9997707]

Constat: Légèrement plus puissant que la procédure de Bonferroni.

 La procédure step-down de Sidak (FWER):

Historique: La procédure step-dow de Sidak a été inventée par Sidak en 1979.

Principe: Rejet des hypothèses H_{k + 1}, \cdots, H_K pour k le plus petit tel que p_k > 1 - (1 - \alpha) ^{\frac{1}{K}} \Rightarrow \alpha < 1 - (1 - p_k) ^K.

Algorithme associé:

– Ranger p par ordre croissant

p_1 * = 1 - (1 - p_1) ^K

\forall i \in [2, K], p_i * = max(p_{i - 1} *, 1 - (1 - p_i) ^{K - i + 1})

p * = [p_1 *, \cdots, p_K *]

Exemple:

Nous avons donc p = [<0.0001, 0.0004, 0.0007, 0.6969, 0.8129]. Comme nous avons 5 p-valeurs à ajuster, K = 5. Donc,

Etape 1, k = 1,

p_1 * = 1 - (1 - "<0.0001") ^5 = "<0.0001"

Etape 2, k = 2,

p_2 * = max(<0.0001, 1 - (1 - 0.0004) ^{5 - 2 +1})

 = max(<0.0001, 1 - 0.9996 ^4)

 = max(<0.0001, 0.00159904)

 = 0.001559904

Etape 3, k = 3,

p_3 * = max(0.00159904, 1 - (1 - 0.0007) ^{5 - 3 +1})

 = max(0.00159904, 0.00209853)

 = 0.00209853

Etape 4, k = 4,

p_4 * = max(0.00209853, 1 - (1 - 0.6969) ^{5 - 4 + 1})

 = max(0.00209853, 0.9081304)

 = 0.9081304

Etape 5, k = 5,

p_5 * = max(0.9081304, 1 - (1 - 0.8129) ^{5 - 5 +1})

 = max(0.9081304, 0.8129)

 = 0.9081304

D’où,

p * = [<0.0001, 0.00159904, 0.00209853, 0.9081304, 0.9081304]

Constat: Les méthodes step-down sont réputées comme étant plus puissantes que les méthodes d’origine.

♦ La procédure par combinaison de Fisher (FWER):

Historique: La procédure est basée sur le test des combinaisons de Fisher lui-même inventé par Ronald Aymer Fisher.

Principe: Rejet de H_k si p_k \leq (- 2 \sum log(p)) qui suit une loi du \chi ^2 à 2 \cdot K degrés de liberté.

Algorithme associé:

– Ranger p par ordre croissant

– Pour toutes les combinaisons C_j ^K d’éléments de p, calculer les p ^{k_1, \cdots, k, \cdots, k_{C_j ^K}} = - 2 \sum_{i \in \lbrace \mbox{ensembles des combinaisons de } C_j ^K \rbrace} log(p_i) qui suit une loi du \chi ^2 de paramètre 2 \cdot j

p_k * = max(p_k, p ^{k_1, \cdots, k, \cdots, k_{C_j ^K}})

p * = [p_1 *, \cdots, p_K *]

Exemple:

Nous avons donc p = [<0.0001, 0.0004, 0.0007, 0.6969, 0.8129]. Comme nous avons 5 p-valeurs à ajuster, K = 5.

Calculons pour chaque combinatoire C_1 ^5, C_2 ^5, C_3 ^5, C_5 ^5 l’ensemble des p-valeurs associées aux statistiques de test de formule - 2 \cdot \sum log(p). Les 5 tableaux ci-dessous présentent ces deux valeurs pour les cinq listes de combinaisons possibles:

add

add

add

add

add

Appliquons l’algorithme décisionnel pour calculer le vecteur p * à partir des 5 tableaux ci-dessus.

p_1 = "<0.0001"

D’après les cinq tableaux nous avons:

p ^1 = p ^{1, 2} = p ^{1, 3} = p ^{1, 4} = p ^{1, 5} = p ^{1, 2, 3} = p ^{1 , 2, 4} = p ^{1, 2, 5} = p ^{1, 3, 4} = p ^{1, 3, 5} = p ^{1, 4, 5} = p ^{1, 2, 3, 4} = p ^{1, 2, 3, 5} = p ^{1, 2, 4, 5} = p ^{1, 2, 4, 5} = p ^{1, 2, 3, 4, 5} = "<0.0001"

Par conséquent,

p_1 * = max(p ^1, p ^{1, 2}, p ^{1, 3}, p ^{1, 4}, p ^{1, 5}, p ^{1, 2, 3}, p ^{1 , 2, 4}, p ^{1, 2, 5}, p ^{1, 3, 4}, p ^{1, 3, 5}, p ^{1, 4, 5}, p ^{1, 2, 3, 4}, p ^{1, 2, 3, 5}, p ^{1, 2, 4, 5}, p ^{1, 2, 4, 5}, p ^{1, 2, 3, 4, 5})

 = "<0.0001"

\Rightarrow p_1 * = "<0.00001"

p_2 = 0.0004

D’après les cinq tableaux nous avons:

p ^{1, 2} = p ^{1, 2, 3} = p ^{1, 2, 4} = p ^{1, 2, 5} = p ^{1, 2, 3, 4} = p ^{1, 2, 3, 4, 5} = p ^{1, 2, 4, 5} = p ^{1, 2, 3, 4, 5} = p ^{2, 3} = p ^{2, 3, 4} = p ^{2, 3, 5} = "<0.0001" < p_2

Ainsi que, p ^{2, 4} = 0.002560453, p ^{2, 5} = 0.002936583, p ^{2, 4, 5} = 0.01010829 et p ^{2, 3, 4, 5} = 0.0001235494.

Par conséquent, 

p_2 * = max(p_2, p ^{2, 4}, p ^{2, 5}, p ^{2, 4, 5}, p ^{2, 3, 4, 5}) = 0.0101829

p_3 = 0.0007

D’après les cinq tableaux nous avons:

p ^{1, 3} = p ^{1, 2, 3} = p ^{1, 3, 4} = p ^{1, 3, 5} = p ^{1, 2, 3, 4} = p ^{1, 2, 3, 4, 5} = p ^{1, 3, 4, 5} = p ^{1, 2, 3, 4, 5} = p ^{2, 3} = p ^{2, 3, 4} = p ^{2, 3, 5} = "<0.0001" < p_3

Ainsi que, p ^{3, 4} = 0.004207799, p ^{3, 5} = 0.004820582, p ^{3, 4, 5} = 0.01566726 et p ^{2, 3, 4, 5} = 0.0001235494.

Par conséquent, 

p_3 * = max(p_3, p ^{3, 4}, p ^{3, 5}, p ^{3, 4, 5}, p ^{2, 3, 4, 5}) = 0.01566726

p_4 = 0.6969

D’après les cinq tableaux nous avons:

 p ^{1, 4} = p ^{1, 2, 4} = p ^{1, 3, 4} = p ^{1, 4, 5} = p ^{2, 3, 4} = p ^{1, 2, 3, 4} = p ^{1, 2, 4, 5} = p ^{1, 3, 4, 5} = p ^{1, 2, 3, 4, 5} = "<0.0001" < p_4

Ainsi que, p ^{2, 4} = 0.002560453, p ^{3, 4} = 0.004207799, p ^{4, 5} = 0.8884353, p ^{2, 4, 5} = 0.01010829, p ^{2, 3, 5} = 0.01566726 et p ^{2, 3, 4, 5} = 0.0001235494.

Par conséquent, 

p_4 * = max(p_4, p ^{2, 4}, p ^{3, 4}, p ^{4, 5}, p ^{2, 4, 5}, p ^{3, 4, 5}, p ^{2, 3, 4, 5}) = 0.8884353

p_5 = 0.8129

D’après les cinq tableaux nous avons:

 p ^{1, 5} = p ^{1, 2, 5} = p ^{1, 3, 5} = p ^{1, 4, 5} = p ^{2, 4, 5} = p ^{2, 3, 5} = p ^{1, 2, 3, 5} = p ^{1, 2, 4, 5} = p ^{1, 2, 3, 4, 5} = "<0.0001" < p_5

Ainsi que, p ^{2, 5} = 0.002936583, p ^{3, 5} = 0.004820582, p ^{4, 5} = 0.8884353, p ^{2, 4, 5} = 0.01010829, p ^{3, 4, 5} = 0.01566726 et p ^{2, 3, 4, 5} = 0.0001235494.

Par conséquent, 

p_5 * = max(p_5, p ^{2, 5}, p ^{3, 5}, p ^{4, 5}, p ^{2, 4, 5}, p ^{3, 4, 5}, p ^{2, 3, 4, 5}) = 0.8884353

D’où,

p * = [<0.0001, 0.0101829, 0.1566726, 0.8884353, 0.8884353]

Constat: Aucun.

♦ La procédure de Hommel (FWER):

Historique: la procédure d’Hommel a été inventée par G. Hommel en 1988. Elle se base sur les tests de Simes.

Principe: Pour chaque ensemble de K hypothèses H_0 auquel il faut associer les p-valeurs p_1 \leq \cdots \leq p_K, le test de Simes revient à calculer la p-valeur le min(\frac{K}{1} \cdot p_1, \cdots, \frac{K}{K} \cdot p_K). L’ajustement d’Hommel pour le k-ème test est le maximum de l’ensemble des p-valeurs de Simes basées sur les tests conjoints incluant le k-ième.

Algorithme associé:

– Ranger p par ordre croissant

– Pour k = K, pa ^K = q ^K = min(\frac{K \cdot p_1}{1}, \cdots, \frac{K \cdot p_K}{K}), avec q q-valeur et pa vecteur des p-valeurs ajustées.

\forall k \in [K - 1, 2], poser i_1 ^k = K - k + 1 et  i_2 ^k = K - k + 2 et calculer q ^k = [A_{1 \rightarrow i_1 ^k} | B_{i_2 ^k \rightarrow K}] où,

A_{1 \rightarrow i_1 ^k} = [min(min(\frac{k \cdot p_{i_2 ^k}}{2}, \cdots \frac{k \cdot p_K}{k}), k \times p_1), \cdots, min(min(\frac{k \cdot p_{i_2 ^k}}{2}, \cdots \frac{k \cdot p_K}{k}), k \times p_{i_1 ^k})]

B_{i_2 ^k \rightarrow K} = [A_{1 \rightarrow i_1 ^k} ^{K - k + 1}, \cdots, A_{1 \rightarrow i_1 ^k} ^{K - k + 1}]

pa ^k = [max(pa_1 ^{k + 1}, q_1 ^k), \cdots, max(pa_K ^{k + 1}, q_K)]

Exemple:

Nous avons donc p = [<0.0001, 0.0004, 0.0007, 0.6969, 0.8129]. Comme nous avons 5 p-valeurs à ajustées, K = 5. Donc,

– Etape 1, q = min(\frac{5 \times "<0.0001"}{1}, \frac{5 \times 0.0004}{2}, \frac{5 \times 0.0007}{3}, \frac{5 \times 0.6969}{4}, \frac{5 \times 0.8129}{5}]

 = min(<0.0001, 0.001, 0.00166667, 0.871125, 0.8129)

 = "<0.0001"

\Rightarrow pa ^5 = q ^5 = [<0.0001, <0.0001, <0.0001, <0.0001, <0.0001]

– Etape 2, k = 4 \Rightarrow i_1 ^4 = 5 - 4 + 1 = 2 et i_2 ^4 = 5 - 4 + 2 = 3

Nous calculons q ^4 = [A_{1 \rightarrow 2} | B_{3 \rightarrow 5}] avec,

A_{1 \rightarrow 2} = [min(min(\frac{4 \times 0.0007}{2}, \frac{4 \times 0.6969}{3}, \frac{4 \times 0.8129}{4}), 4 \times "<0.0001"), min(min(\frac{4 \times 0.0007}{2}, \frac{4 \times 0.6969}{3}, \frac{4 \times 0.8129}{4}), 4 \times 0.0004)]

 = [min(min(0.0014, 0.9292, 0.8129), <0.0001), min(0.0014, 0.9292, 0.8129), 0.0016)]

 = [<0.0001, 0.0014]

B_{3 \rightarrow 5} = [A_{1 \rightarrow 2} ^{5 - 4 + 1}, A_{1 \rightarrow 2} ^{5 - 4 + 1}, A_{1 \rightarrow 2} ^{5 - 4 + 1}] = [0.0014, 0.0014, 0.0014]

Nous avons donc;

q ^4 = [<0.0001, 0.0014, 0.0014, 0.0014, 0.0014]

pa ^5 = [<0.0001, <0.0001, <0.0001, <0.0001, <0.0001]

Comme q_1 ^4 = pa_1 ^5, q_2 ^4 > pa_2 ^5, q_3 ^3 > pa_3 ^5, q_4 ^4 > pa_4 ^5, q_5 ^4 > pa_5 ^5 alors,

pa ^4 = [<0.0001, 0.0014, 0.0014, 0.0014, 0.0014]

– Etape 3, k = 3 \Rightarrow i_1 ^3 = 5 - 3 +1 = 3 et i_2 ^3 = 5 - 3 + 2 = 4

Nous calculons q ^3 = [A_{1 \rightarrow 3} | B_{4 \rightarrow 5}] avec,

A_{1 \rightarrow 3} = [min(min(\frac{3 \times 0.6969}{2}, \frac{3 \times 0.8129}{3}), 3 \times "<0.0001"), min(min(\frac{3 \times 0.6969}{2}, \frac{3 \times 0.8129}{3}), 3 \times 0.0004), , min(min(\frac{3 \times 0.6969}{2}, \frac{3 \times 0.8129}{3}), 3 \times 0.0007)]

 = [min(min(1.04535, 0.8129), <0.0001), min(1.04535, 0.8129), 0.0012), min(1.04535, 0.8129), 0.0021)]

 = [<0.0001, 0.0012, 0.0021]

B_{4 \rightarrow 5} = [A_{1 \rightarrow 3} ^{5 - 3 +1}, A_{1 \rightarrow 3} ^{5 - 3 +1}] = [0.0021, 0.0021]

Nous avons donc;

q ^3 = [<0.0001, 0.0012, 0.0021, 0.0021, 0.0021]

pa ^4 = [<0.0001, 0.0014, 0.0014, 0.0014, 0.0014]

Comme q_1 ^3 = pa_1 ^4, q_2 ^3 < pa_2 ^4, q_3 ^3 > pa_3 ^4, q_4 ^3 > pa_4 ^4, q_5 ^3 > pa_5 ^4 alors,

pa ^4 = [<0.0001, 0.0014, 0.0021, 0.0021, 0.0021]

– Etape 4, k = 2 \Rightarrow i_1 ^2 = 5 - 2 + 1 = 4 et i_2 ^2 = 5 - 2 + 2 = 5

Nous calculons q ^2 = [A_{1 \rightarrow 4} | B_{5 \rightarrow 5}] avec,

A_{1 \rightarrow 4} = [min(\frac{2 \times 0.8129}{2}, 2 \times "<0.0001"), min(\frac{2 \times 0.8129}{2}, 2 \times 0.0004), min(\frac{2 \times 0.8129}{2}, 2 \times 0.0007), min(\frac{2 \times 0.8129}{2}, 2 \times 0.6969)]

 = [min(0.8129, <0.0001), min(0.8129, 0.0008), min(0.8129, 0.0014), min(0.8129, 1.3938)]

 = [<0.0001, 0.0008, 0.0014, 0.8129]

B_{5 \rightarrow 5} = A_{1 \rightarrow 4} ^{5 - 4 + 1} = 0.8129

Nous avons donc;

q ^2 = [<0.0001, 0.0008, 0.0014, 0.8129, 0.8129]

pa ^3 = [<0.0001, 0.0014, 0.0021, 0.0021, 0.0021]

Comme q_1 ^2 = pa_1 ^3, q_2 ^1 < pa_2 ^3, q_3 ^2 w pa_3 ^3, q_4 ^2 > pa_4 ^3, q_5 ^2 > pa_5 ^3 alors,

pa ^3 = [<0.0001, 0.0014, 0.0021, 0.8129, 0.8129]

D’où,

p * = [<0.0001, 0.0014, 0.0021, 0.8129, 0.8129]

Constat: Aucun.

♦ La procédure par bootstrap (FWER):

Historique: La procédure par bootstrap a été décrite par P. J. Westfall et S. S. Young en 1993. Elle porte également les noms de procédure bootstrap minP ou encore Westfall-Young bootstrap minP.

Principe: Ajuster les p-valeurs par bootstrap.

Algorithme associé:

Il est souvent préférable, dans le cas où X est continue, de centrer-réduire les données afin d’approcher de manière optimal la distribution de l’hypothèse nulle. En effet, le rééchantillonnage par bootstrap, selon son paramètrage, peut transformer le format des données en multimodal.

– Déterminer le nombre B de rééchantillonnage par bootstrap à exécuter ainsi que le nombre d’observations à tirer. Par définition du bootstrap, le tirage se fait avec remise.

– Pour b \in [1, \cdots, B] tirer X_b

– Et calculer les différents p_1 ^b, \cdots, p_K ^b, \forall b pour chaque hypothèse

– Enfin, \forall k \in [1, \cdots, K], p_k * = \frac{\sharp \lbrace b; p_k ^b \leq p_k \rbrace }{B}

p * = [p_1 *, \cdots, p_K *]

Constat: La procédure bootstrap permet de conserver les corrélations entre les p-valeurs ainsi que les particularités de leur distribution.

♦ La procédure step-down minP par bootstrap (FWER):

Historique: La procédure step-down minP par bootstrap a été décrite par Yongchao Ge, Sandrine Dudoit, and Terence P. Speed en 2003.

Principe: Ajuster les p-valeurs par bootstrap.

Algorithme associé:

Il est souvent préférable, dans le cas où X est continue, de centrer-réduire les données afin d’approcher de manière optimal la distribution de l’hypothèse nulle. En effet, le rééchantillonnage par bootstrap, selon son paramètrage, peut transformer le format des données en multimodal.

– Déterminer le nombre B de rééchantillonnage par bootstrap à exécuter ainsi que le nombre d’observations à tirer. Par définition du bootstrap, le tirage se fait avec remise.

– Pour b \in [1, \cdots, B] tirer X_b

– Et calculer les différents p_1 ^b, \cdots, p_K ^b, \forall b pour chaque hypothèse

– Calculer les q-valeurs \forall b \in [1, B]:

_____ q_K ^b = p_K ^b

_____ \forall k \in [K - 1, \cdots, 1], q_k ^b = min(q_{k + 1} ^b, p_k ^b)

\forall k \in [1, \cdots, K], p_k * = \frac{\sharp \lbrace b; q_k ^b \leq p_k \rbrace }{B}

– Enfin, p_1 * = p_1 * et \forall k \in [2, K], p_k * = max (p_{k - 1} *, p_k *)

p * = [p_1 *, \cdots, p_K *]

Constat: Les méthodes step-down sont réputées comme étant plus puissantes que les méthodes d’origine.

♦ La procédure step-down maxP par bootstrap (FWER):

Historique: La procédure step-down maxP par bootstrap a été décrite par Yongchao Ge, Sandrine Dudoit, and Terence P. Speed en 2003.

Principe: Ajuster les p-valeurs par bootstrap.

Algorithme associé:

Il est souvent préférable, dans le cas où X est continue, de centrer-réduire les données afin d’approcher de manière optimal la distribution de l’hypothèse nulle. En effet, le rééchantillonnage par bootstrap, selon son paramètrage, peut transformer le format des données en multimodal.

– Déterminer le nombre B de rééchantillonnage par bootstrap à exécuter ainsi que le nombre d’observations à tirer. Par définition du bootstrap, le tirage se fait avec remise.

– Pour b \in [1, \cdots, B] tirer X_b

– Et calculer les différents p_1 ^b, \cdots, p_K ^b, \forall b pour chaque hypothèse

– Calculer les q-valeurs \forall b \in [1, B]:

_____ q_{K,b} = p_{K,t}

_____ \forall k \in [K - 1, \cdots, 1], q_{k,b} = max (q_{k + 1,b}, p_{k,b})

\forall k \in [1, \cdots, K], p_k * = \frac{\sharp \lbrace b; q_{k,b} \geq p_k \rbrace }{T}

– Enfin, p_1 * = p_1 * et \forall k \in [2, K], p_k * = max (p_{k - 1} *, p_k *)

p * = [p_1 *, \cdots, p_K *]

Constat: Les méthodes step-down sont réputées comme étant plus puissantes que les méthodes d’origine.

♦ La procédure par permutations (FWER):

Historique: La procédure par permutations a été décrite par C. C. Brown et  T. R. Fears en 1981.

Principe: Ajuster les p-valeurs par permutations.

Algorithme associé:

– Déterminer le nombre T de séries de permutations à exécuter sur X.

– Pour t \in [1, \cdots, T] permuter les N observations de X

– Calculer les p_1 ^t, \cdots, p_K ^t, \forall t pour chaque hypothèse

– Enfin, \forall k \in [1, \cdots, K], p_k * = \frac{\sharp \lbrace b; p_k ^t \leq p_k \rbrace }{T}

p * = [p_1 *, \cdots, p_K *]

Constat: La procédure par permutations permet de conserver les corrélations entre les p-valeurs ainsi que la structure multivariée des données.

♦ La procédure step-down minP par permutations (FWER):

Historique: La procédure step-down minP par permutations a été décrite par Yongchao Ge, Sandrine Dudoit, and Terence P. Speed en 2003.

Principe: Ajuster les p-valeurs par permutations.

Algorithme associé:

– Déterminer le nombre T séries de permutations à exécuter sur X.

– Pour t \in [1, \cdots, T] permuter les N observations de X

– Calculer les p_1 ^t, \cdots, p_K ^t, \forall t pour chaque hypothèse

\forall t \in [1, T], ranger p ^t par ordre croissant

– Calculer les q-valeurs \forall t \in [1, T]:

_____ q_K ^t = p_K ^t

_____ \forall k \in [K - 1, \cdots, 1], q_k ^t = min(q_{k + 1} ^t, p_k ^t)

\forall k \in [1, \cdots, K], p_k * = \frac{\sharp \lbrace t; q_k ^t \leq p_k \rbrace }{T}

– Enfin, p_1 * = p_1 * et \forall k \in [2, K], p_k * = max (p_{k - 1} *, p_k *)

p * = [p_1 *, \cdots, p_K *]

Constat: Les méthodes step-down sont réputées comme étant plus puissantes que les méthodes d’origine.

♦ La procédure step-down maxP par permutations (FWER):

Historique: La procédure step-down maxP par permutations a été décrite par Yongchao Ge, Sandrine Dudoit, and Terence P. Speed en 2003.

Principe: Ajuster les p-valeurs par permutation.

Algorithme associé:

– Déterminer le nombre T séries de permutations à exécuter sur X.

– Pour t \in [1, \cdots, T] permuter les N observations de X

– Calculer les p_1 ^t, \cdots, p_K ^t, \forall t pour chaque hypothèse

\forall t \in [1, T], ranger p ^t par ordre croissant

– Calculer les q-valeurs \forall t \in [1, T]:

_____ q_K ^t = p_K ^t

_____ \forall k \in [K - 1, \cdots, 1], q_k ^t = max (q_{k + 1} ^t, p_{k} ^t)

\forall k \in [1, \cdots, K], p_k * = \frac{\sharp \lbrace t; q_k ^t \geq p_k \rbrace }{T}

– Enfin, p_1 * = p_1 * et \forall k \in [2, K], p_k * = max (p_{k - 1} *, p_k *)

p * = [p_1 *, \cdots, p_K *]

Constat: Les méthodes step-down sont réputées comme étant plus puissantes que les méthodes d’origine.

♦ La procédure False Discovery Rate par bootstrap (FDR):

Historique: La procédure FDR par bootstrap a été décrite par D. Yekateuli et Yoav Benjamini en 1999.

Principe: Ajuster les p-valeurs par bootstrap.

Algorithme associé:

Il est souvent préférable, dans le cas où X est continue, de centrer-réduire les données afin d’approcher de manière optimal la distribution de l’hypothèse nulle. En effet, le rééchantillonnage par bootstrap, selon son paramètrage, peut transformer le format des données en multimodal.

– Déterminer le nombre B de rééchantillonnage par bootstrap à exécuter ainsi que le nombre d’observations à tirer. Par définition du bootstrap, le tirage se fait avec remise.

\forall b \in [1, B], ranger p ^b par ordre croissant

\forall b \in [1, B], calculer \forall k \in [1, K], R (p_k) ^b = \sharp \lbrace 1 \mbox{ a } b; p_k ^b \leq p_k \rbrace

\forall k \in [1, K], r_{\beta} (p_k) = 100(1 - \beta)-ème quantile de \lbrace R (p_k) ^1, \cdots, R (p_k) ^B \rbrace avec \beta = 1 - \alpha

\forall k \in [1, K], r (p_k) = \sharp \lbrace 1 \mbox{ a } B; p_k ^b \leq p_k \rbrace

\forall k \in [1, K], calculer l’estimateur local ou de la limite supérieure,

_____ Q (p_k) = \frac{1}{B} \sum_{b = 1} ^B \frac{R (p_k) ^b}{R (p_k) ^b + r (p_k) - p_k \cdot K} (dans le cas de sa version locale la condition est: r (p_k) - r_{\beta} (p_k) \geq p_k \cdot K, celle dans sa version de la limite supérieure: r (p_k) - r_{\beta} (p_k) \geq 0

_____ Q (p_k) = \frac{\sharp \lbrace R (p_k) ^b \geq 1 \rbrace}{B} sinon (quelque soit la version de l’estimateur)

\forall k \in [1, K], la p-valeur ajustée devient:

_____ p_k * = Q (p_k) pour k = K

_____ p_k * = min(p_{k + 1} *, Q (p_k)) \forall k \in \lbrace K - 1, \cdots, 1 \rbrace

p * = [p_1 *, \cdots, p_K *]

Constat: La procédure FDR bootstrap permet de conserver les corrélations entre les p-valeurs ainsi que les particularités de leur distribution.

♦ La procédure False Discovery Rate par permutations (FDR):

Historique: La procédure FDR par permutations a été décrite par D. Yekateuli et Yoav Benjamini en 1999.

Principe: Ajuster les p-valeurs par permutations.

Algorithme associé:

– Déterminer le nombre T de séries de permutations à exécuter sur X.

\forall t \in [1, T], ranger p ^t par ordre croissant

\forall t \in [1, T], calculer \forall k \in [1, K], R (p_k) ^t = \sharp \lbrace 1 \mbox{ a } t; p_k ^t \leq p_k \rbrace

\forall k \in [1, K], r_{\beta} (p_k) = 100(1 - \beta)-ème quantile de \lbrace R (p_k) ^1, \cdots, R (p_k) ^T \rbrace avec \beta = 1 - \alpha

\forall k \in [1, K], r (p_k) = \sharp \lbrace 1 \mbox{ a } T; p_k ^t \leq p_k \rbrace

\forall k \in [1, K], calculer l’estimateur local ou de la limite supérieure,

_____ Q (p_k) = \frac{1}{T} \sum_{t = 1} ^T \frac{R (p_k) ^t}{R (p_k) ^t + r (p_k) - p_k \cdot K} (dans le cas de sa version locale la condition est: r (p_k) - r_{\beta} (p_k) \geq p_k \cdot K, celle dans sa version de la limite supérieure: r (p_k) - r_{\beta} (p_k) \geq 0

_____ Q (p_k) = \frac{\sharp \lbrace R (p_k) ^t \geq 1 \rbrace}{T} sinon (quelque soit la version de l’estimateur)

\forall k \in [1, K], la p-valeur ajustée devient:

_____ p_k * = Q (p_k) pour k = K

_____ p_k * = min(p_{k + 1} *, Q (p_k)) \forall k \in \lbrace K - 1, \cdots, 1 \rbrace

p * = [p_1 *, \cdots, p_K *]

Constat: La procédure par permutations permet de conserver les corrélations entre les p-valeurs ainsi que la structure multivariée des données.

 La procédure de Benjamini-Höchberg (FDR):

Historique: La procédure de Benjamini-Höchberg a été inventée par Yoav Benjamini et Yosef Höchberg en 1995. Elle porte également le nom de méthode linéraire step-up de Benjamini-Höchberg.

Principe: Rejet de H_1, \cdots, H_k pour le plus grand k tel que p_k \leq \frac{k}{K} \cdot \alpha \Rightarrow \alpha \geq \frac{K}{k} \cdot p_k.

Algorithme associé:

– Ranger p par ordre croissant

– p_K * = p_K

\forall k \in [K - 1, \cdots, 1], p_k * = min(p_{k + 1}*, \frac{K}{k} \cdot p_k)

p * = [p_1 *, \cdots, p_K *]

Exemple:

Nous avons donc p = [<0.0001, 0.0004, 0.0007, 0.6969, 0.8129]. Comme nous avons 5 p-valeurs à ajuster, K = 5.

– Etape 1: k = 5,

p_5 * = p_5 = 0.8129

– Etape 2: k = 4,

p_4 * = min(0.8129, \frac{5}{4} \times 0.6969) = min(0.8129, 0.871125) = 0.8129

– Etape 3: k = 3,

p_3 * = min(0.8129, \frac{5}{3} \times 0.0007) = min(0.8129, 0.001166667) = 0.001166667

– Etape 4: k = 2,

p_2 * = min(0.001166667, \frac{5}{2} \times 0.0004) = min(0.001166667, 0.001) = 0.001

– Etape 5: k = 1,

p_1 * = min(0.001, \frac{5}{1} \times "<0.0001") = "<0.0001"

D’où,

p * = (<0.0001, 0.001, 0.001166667, 0.8129, 0.8129)

Constat: A l’instar de la procédure de Bonferroni pour la famille FWER, la procédure de Benjamini-Höchberg a inspiré bon nombre de procédures visant à améliorer la puissance pour la famille FDR.

♦ La procédure de Benjamini-Höchberg adaptative (Adaptive-FDR):

Historique: La procédure adaptative de Benjamini-Höchberg a été inventée par Yoav Benjamini et Yosef Höchberg en 2000.

Principe: Rejet de H_1, \cdots, H_k pour le plus grand k tel que p_k \leq \frac{k}{\hat{m_0}} \cdot \alpha \Rightarrow \alpha \geq \frac{\hat{m_0}}{k} \cdot p_k.

Algorithme associé:

– Ranger p par ordre croissant

– Estimer \hat{m_0}, le nombre d’hypothèse nulle vraie

\forall k \in [K, \cdots, 1],

_________ Si k > \hat{m_0}, p_k * = p_k

_________ Si p_k \leq \hat{m_0}, p_k * = min(p_{k + 1} *, \frac{\hat{m_0}}{k} \cdot p_k)

p * = [p_1 *, \cdots, p_K *]

Exemple:

Nous avons donc p = [<0.0001, 0.0004, 0.0007, 0.6969, 0.8129]. Comme nous avons 5 p-valeurs à ajuster, K = 5.

Afin d’estimer \hat{m_0} plusieurs méthodes existent. Ci-dessous, l’ensemble des méthodes extraitent du site du suppot SAS:

addNous utiliserons la méthode par défaut du logiciel sus-nommé: decrease slope (ou pente décroissante) qui consiste à étudier la courbe des p-valeurs tracée par le vecteur q = 1 - p puis de calculer la totalité des b_i = a(q_K, \cdots, q_{K - i})a représente le coefficient de la fonction affine de type y = a \cdot x + b (y prenant la valeur de l’indice de la p-valeur).

Il s’agit alors de retenir de manière itérative décroissante le premier coefficient b_i qui répond à l’inégalité b_i > b_{i + 1}. Nous avons alors,

\hat{m_0} = ceil(\frac{1}{b_{i+1}} - 1)

, où ceil(.) représente la fonction renvoiyant le premier entier supérieur à « . ».

Si nous appliquons cet algorithme à notre exemple nous obtenons,

b_1 = a(1 - p_5, 1 - p_5) = non calculable par formule

b_2 = a(1 - p_5, 1 - p_4) = - \infty

b_3 = a(1 - p_5, 1 - p_4, 1 - p_3) = 0.0785

b_4 = a(1 - p_5, 1 - p_4, 1 - p_3, 1 - p_2) = 0.8502

b_5 = a(1 - p_5, 1 - p_4, 1 - p_3, 1 - p_2, 1 - p_1) = 0.2563

Il en ressort, si nous parcourons les b_i de manière décroissante, que b_4 > b_5 est le premier cas qui répond à l’inégalité recherchée. Par conséquent,

\hat{m_0} = ceil(\frac{1}{0.2563} - 1) = ceil(2.901688) = 3

Nous pouvons donc dérouler l’algorithme adaptatif de Benjamin-Höchberg. Ainsi,

p_5 * = p_5 = 0.8129 car k > \hat{m_0} = 3

p_4 * = p_4 = 0.6969 car k > \hat{m_0} = 3

p_3 * = min(p_4 *, \frac{3}{3} \times 0.0007) = min(0.6969, 0.0007) = 0.0007

p_2 * = min(p_3 *, \frac{3}{2} \times 0.0004) = min(0.0007, 0.0006) = 0.0006

p_1 * = min(p_2 *, \frac{3}{1} \times "<0.0001")

 = min(<0.0001, <0.0001)

 = "<0.0001"

D’où,

p * = [<0.0001, 0.0006, 0.0007, 0.6969, 0.8129]

Constat: Plus puissant que la procédure de Benjamini-Höchberg.

 La procédure de Benjamini-Yekutieli (Dependant-FDR):

Historique: La procédure de Benjamini-Yekutieli a été inventée par Yoav Benjamini et D. Yekutieli en 2001. Elle porte également le nom de procédure de Benjamini-Höchberg-Yekutieli.

Principe:  Rejet de H_1, \cdots, H_k pour k le plus grand k tel que p_k \leq \frac{k}{K \cdot c(K)} \cdot \alpha \Rightarrow \alpha \geq \frac{K \cdot c(K)}{k} \cdot p_k avec c(K) = \sum_{k = 1} ^K \frac{1}{k}.

A noter que si nous ne sommes pas en présence de l’hypothèse de dépendance, alors:

– dans le cas de corrélation positive: c(K) = 1

– dans le cas de corrélation négative: c(K) = ln(K) + \gamma\gamma est la constante d’Euler-Mascheroni

Algorithme associé:

– Ranger p par ordre croissant

– Calculer c(K) en fonction de l’hypothèse de dépendance

\forall k \in [1, K], p_k * min(\frac{c(K) \cdot K}{k} \cdot p_k, \cdots, \frac{c(K) \cdot p_K}{K}, 1)

p * = [p_1 *, \cdots, p_K *]

Exemple:

Commençons pas calculer c(K) = c(5) = \sum_{k = 1} ^5 \frac{1}{k} = 2.283333333.

Nous avons donc:

– Etape 1, k = 5,

p_5 * = min(\frac{2.283333 \times 5}{5} \times 0.8129, 1) = min(1.856122, 1) = 1

– Etape 2, k = 4,

p_4 * = min(\frac{2.283333 \times 5}{4} \times 0.6969, 1, 1) = min(1.989069, 1, 1) = 1

– Etape 3, k = 3,

p_3 * = min(\frac{2.283333 \times 5}{3} \times 0.0007, 1, 1, 1) = min(0.002663889, 1, 1, 1) = 0.002663889

– Etape 4, k = 2,

p_2 * = min(\frac{2.283333 \times 5}{2} \times 0.0004, 0.002663889, 1, 1, 1)

 = min(0.00228333, 0.002663889, 1, 1, 1)

 = 0.002883333

– Etape 1, k = 1,

p_1 * = min(\frac{2.283333 \times 5}{1} \times "<0.0001", 0.00228333, 0.002663889, 1, 1, 1)

 = min(<0.0001, 0.00228333, 0.002663889, 1, 1, 1)

 = "<0.0001"

D’où,

p * = [<0.0001, 0.00228333, 0.002663889, 1, 1]

Constat: Requiert l’hypothèse d’indépendance.

 La procédure de Storey (Positive-FDR):

Historique: La procédure de Storey a été inventée par John D. Storey en 2002. Elle porte également le nom de positive FDR.

Principe: Détermination de la q-valeur.

Algorithme associé:

– Ranger p par ordre croissant

– Déterminer \lambda et poser f \in \lbrace 0, 1\rbrace en fonction de l’indépendance ou non des p-valeurs (f = 1 si indépendance, f = 0 sinon), à noter que si \lambda = 0 la procédure de Storey devient celle de Benjaminii-Höchberg (FDR)

– Calculer N( p_{seuil} ) = \sharp \lbrace p_k \leq p_{seuil} \rbrace pour p_{seuil} = \lambda

– Calculer \hat{\pi_0} = min(\frac{K - N(\lambda) - f}{K \cdot (1 - \lambda)}, 1), proportion d’hypothèses nulles vraies

– Fixer \forall p_k > \lambda, \hat{q}(p_k) = 1 en posant k' le plus grand indice tel que p_k \leq \lambda

\forall k \in [k', 1], \hat{q} (p_k) = min(\frac{\hat{\pi_0} \cdot K \cdot p_k}{k}, \hat{q} (p_{k + 1}))

– Dans le cas où aucune des p-valeurs p_k > \lambda, pour le cas k = K il faut calculer \hat{q}(p_K) = \hat{\pi_0} \cdot p_K, l’algorithme ci-dessus est ensuite appliqué pour \forall k \in [K - 1, 1]

p * = \hat{q} (p)

Exemple:

Nous avons donc p = [<0.0001, 0.0004, 0.0007, 0.6969, 0.8129]. Comme nous avons 5 p-valeurs à ajuster, K = 5. Nous cherchons à appliquer la procédure de Storey à p pour \lambda = 0.0005 dans le cas indépendant (\Rightarrow f = 1). Donc,

\lambda = 0.0005 \Rightarrow N(\lambda) = 2 \Rightarrow \hat{\pi_0} = \frac{5 - 2 + 1}{5 \times (1 - 0.0005)} = \frac{4}{4.9975} = 0.8004002

Comme,

p_3, p_4, p_5 > \lambda \Rightarrow \hat{q}(p_3) = \hat{q}(p_4) = \hat{q}(p_5) = 1

Enfin,

– Etape 1, k = 2

\hat{q}(p_2) = min(\frac{0.8004002 \times 5 \times 0.0004}{2}, 1) = min(0.008004002, 1) = 0.008004002

– Etape 2, k = 1,

\hat{q}(p_1) = min(\frac{0.8004002 \times 5 \times}{1}, 0.008004002)

 = min(<0.0001, 0.008004002)

 = "<0.0001"

D’où,

p * = \hat{q}(p) = (<0.0001, 0.008004002, 1, 1, 1)

Constat: Estimation facilement interprétable du nombre de faux positifs, recquiert la connaissance sur la distribution des p-valeurs.

Autres procédures:

Six autres procédures correctives existent mais restent beaucoup moins répandues, nous présentons ici les documents dans lesquels elle sont décrites:

– La procédure de Simes (aucune famille): Multiplicité des tests et calculs de taille d’echantillon en recherche clinique de Jérémie Riou

– La procédure step-up de Rom (FWER): Approximation du calcul de la taille échantillonnale pour les tests à hypothèses multiples lorsque r parmi m hypothèses doivent être significatives par Philippe Delorme

– La procédure de Dunnett (FWER): Approximation du calcul de la taille échantillonnale pour les tests à hypothèses multiples lorsque r parmi m hypothèses doivent être significatives par Philippe Delorme

– La procédure step-down de Dunnett (FWER): Approximation du calcul de la taille échantillonnale pour les tests à hypothèses multiples lorsque r parmi m hypothèses doivent être significatives par Philippe Delorme

– La procédure de Dunnett et Tamhane (FWER): Approximation du calcul de la taille échantillonnale pour les tests à hypothèses multiples lorsque r parmi m hypothèses doivent être significatives par Philippe Delorme

– La procédure step-wise de Tamhane-Liu-Dunnett (FWER): Approximation du calcul de la taille échantillonnale pour les tests à hypothèses multiples lorsque r parmi m hypothèses doivent être significatives par Philippe Delorme

Application informatique:

 Procédure SAS: http://support.sas.com/documentation/cdl/en/statug/63347/HTML/default/viewer.htm#statug_multtest_sect021.htm

 Fonction R:

p.multtest = function(p) {

    # La fonction applique les méthodes de Holm, Höchberg, Bonferroni, Benjamini-Höchberg, Benjamini-Yektulieli, step-down ajustée de Bonferroni,
    # adaptative de Holm, méthode adaptative de Höchberg, Sidak, Step-Down de Sidak, combinaison de Fisher, AFDR de Benjamini-Höchberg,

library(stats)
library(combinat)
nb.p = length(p)
p = sort(p)

# Calcul de l’estimateur d’hypothèse nulle vraie m0 pour les procédures concernées
coeffs = matrix(0,nb.p-1,1)
y = 1:nb.p
for (a in 2:nb.p) {coeffs[a-1] = line(1-sort(p,decreasing=TRUE)[1:a],y=y[1:a])$coefficients[1]}
a = 1
b = FALSE
while ((a<nb.p) && (b==FALSE)) {
if (coeffs[a] <= coeffs[a+1]) {a = a + 1}
if (coeffs[a] > coeffs[a+1]) {b = TRUE}
}
bi = coeffs[a+1]
m0 = (1/bi – 1)
if (m0 != floor(m0)) {m0 = floor(m0) + 1} # test fonction ceiling

# Création de la matrice de résultats
p.correct = matrix(0,nb.p,14)
p.correct[,1] = p

 # On applique la fonction p.adjust qui implémente déjà les méthodes de Holm, Höchberg, Hommel, Bonferroni, Benjamini-Höchberg et Benjamini-Yekutieli
exist = c(« holm », »hochberg », »hommel », »bonferroni », »BH », »BY »)
for (m in 1:6) {p.correct[,m+1] = p.adjust(p,method=exist[m])}

# On applique la méthode Step-Down ajustée de Bonferroni
pa = matrix(0,nb.p,1)
pa[1] = min(nb.p*p[1],1)
for (k in 2:nb.p) {pa[k] = min(max(pa[k-1],(nb.p-k+1)*p[k]),1)}
p.correct[,8] = pa

 # On applique la méthode adaptative de Holm
pa = matrix(0,nb.p,1)
pa[1] = min(min(nb.p,m0)*p[1],1)
for (k in 2:nb.p) {pa[k] = min(max(pa[k-1],min(nb.p-k+1,m0)*p[k]),1)}
p.correct[,9] = pa

    # On applique la méthode adaptative de Höchberg
pa = matrix(0,nb.p,1)
pa[nb.p] = min(1,m0)*p[nb.p]
for (k in (nb.p-1):1) {pa[k] = min(pa[k+1],min(nb.p-k+1,m0)*p[k])}
p.correct[,10] = pa

# On applique la méthode de Sidak
pa = matrix(0,nb.p,1)
for (k in 1:nb.p) {pa[k] = 1 – (1 – p[k]) ^nb.p}
p.correct[,11] = pa

# On applique la méthode Step-Down de Sidak
pa = matrix(0,nb.p,1)
pa[1] = 1 – (1 – p[1]) ^5
for (k in 1:nb.p) {pa[k] = max(pa[k-1], 1-(1-p[k])^(nb.p-k+1))}
p.correct[,12] = pa

# On applique la méthode par combinaison de Fisher
Combinaisons = list(nb.p)
for (k in 1:nb.p) {
combiK = t(combn(nb.p,k))
nbComb = dim(combiK)[1]
pCF = matrix(0,nbComb,1)
if (k ==1) {for (i in 1:nbComb) {pCF[i] = pchisq(-2*sum(log(p[combiK[i]])),2*k,lower.tail=FALSE)}}
if (k > 1) {for (i in 1:nbComb) {pCF[i] = pchisq(-2*sum(log(p[combiK[i,]])),2*k,lower.tail=FALSE)}}
Combinaisons[[k]] = cbind(combiK,pCF)
}
pa = matrix(0,nb.p,1)
for (k in 1:nb.p) {
PCFselect = p[k]
for (k.c in 1:nb.p) {
nbComb = dim(Combinaisons[[k.c]])[1]
for (i.c in 1:nbComb) {if (sum(Combinaisons[[k.c]][i.c,1:k.c] == k) > 0) {PCFselect = c(PCFselect,Combinaisons[[k.c]][i.c,k.c+1])}}
}
pa[k] = max(PCFselect)
}
p.correct[,13] = pa

 # On applique la méthode AFDR de Benjamini-Höchberg
pa = matrix(0,nb.p,1)
for (k in nb.p:1) {
if (k > m0) {pa[k] = p[k]}
if (k <= m0) {pa[k] = min(pa[k+1], p[k]*m0/k)}
}
p.correct[,14] = pa

p.correct = as.data.frame(p.correct)
names(p.correct) = c(« p Brut »,exist, »SD bonferroni », »Holm Adaptative », »Hochberg Adaptative », »Sidak », »SD Sidak », »FisherC », »AFDR »)
return(p.correct)

}

Bibliographie:

– Resampling-based multiple testing for microarray data analysis de Yongchao Ge, Sandrine Dudoit et Terence P. Speed

– Approximation du calcul de la taille échantillonnale pour les tests à hypothèses multiples lorsque r parmis m hypothèses doivent être significatives de Philippe Delorme

– Discussion of “Multiple Testing for Exploratory Research” de J. J. Goeman et A. Solari

– Une introduction au problème des tests multiples de E. Roquain

– Tests multiples pour données biologiques en grande dimension de C. Friguet

– Théorica statistica delle classi e clcolo delle probabilità de Carlo Emilio Bonferroni

– A simple sequentially rejective multiple test procédure de Sture Holm

– A sharper Bonferroni procedure for multiple tests of significance de Yosef Höchberg

– Controlling the false discovery rate: a pratical and powerfull approach to multiple testing de Yoav Benjamini et Yosef Höchberg

– The control of the false discovery rate in multiple testing under dependency de Yoav Benjamini et D. Yekutieli

– On the adaptative control of the false discovery rate in multiple testing with independant statistics de Yoav Benjamini et Yosef Höchberg

– Rectangular confidence regions for the means of multvariate normal distributions de Sidak

– More powerful procedures for multiple signiface testing de Yoav Benjamini et Yosef Höchberg

– Modified sequentially rejective multiple test procedures de J. P. Schaffer

– An improved sequentially rejective Bonferroni test procedure de B. S. Holland et M. D. Copenhaver

– Multiple comparison procedures de Yoav Höchberg et A. C. Tamhane

– A comparison of two modified Bonferroni procedures de G. Hommel

– Développement d’outils de géo-calcul haute performance pour l’identification de régions du génome potentiellement soumises à la sélection naturelle: analyse spatiale de la diversité de panels de polymorphismes nucléotidiques à haute densité (800k) chez Bos taurus et B. indicus en Ouganda de Sylvie Stucki

– A Direct Approach to False Discovery Rates de John D. Storey

– Strong Control, Conservative Point Estimation, and Simultaneous Conservative Consistency of False Discovery Rates: A Unified Approach de John D. Storey, J. E. Taylor et David Siegmund

– Resampling based multiple testing de P. J. Westfall et S. S. Young

– Exact Significance Levels for Multiple Binomial Testing with Application to Carcinogenicity Screens deC. C. Brown et T. R. Fears

– Adjusting for Multiplicity of Statistical Tests in the Analysis of Carcinogenicity Studies de J. Heyse et D. Rom

– Resampling based false discovery rate controlling multiple test procedure for correlated test statistics de D. Yekateuli et Y. Benjamini