Comparez deux groupes avec Mann-Whitney

Une équipe produit lance un test A/B : la moitié des visiteurs voit l’ancienne page d’accueil (version A), l’autre moitié une nouvelle (version B). On mesure la durée de chaque session. Les deux groupes sont formés de personnes différentes : impossible d’apparier, comme au chapitre précédent. Et les durées de session sont typiquement asymétriques, avec quelques sessions très longues.

Le test de Mann-Whitney est l’outil standard de cette situation. Vous verrez comment il se calcule, ce qu’il teste exactement (ce n’est pas tout à fait une égalité de moyennes), et pourquoi il est, mot pour mot, l’aire sous la courbe ROC que vous utilisez pour évaluer un classifieur.

Le fil rouge : un test A/B sur la durée des sessions

Sur un petit pilote, on a relevé ces durées de session (en minutes) :

Version A ($n_1 = 6$)	3,2	5,1	2,4	7,8	4,0	3,6
Version B ($n_2 = 7$)	6,3	8,9	4,4	12,5	7,1	5,7	9,8

Les médianes valent 3,8 minutes pour A et 7,1 minutes pour B. La version B semble retenir les visiteurs plus longtemps. Mais avec 13 sessions, est-ce concluant ?

La somme des rangs

L’idée du test, proposée par Wilcoxon puis reformulée par Mann et Whitney, est la suivante : on mélange les deux échantillons, on classe les $N = n_1 + n_2$ valeurs, puis on regarde comment les rangs se répartissent entre les groupes. Si les deux versions sont équivalentes, les rangs de A devraient être éparpillés au hasard parmi les rangs 1 à 13. Si A est systématiquement plus court, ses rangs seront concentrés en bas.

Exemple

Classons les 13 valeurs :

Valeur	2,4	3,2	3,6	4,0	4,4	5,1	5,7	6,3	7,1	7,8	8,9	9,8	12,5
Groupe	A	A	A	A	B	A	B	B	B	A	B	B	B
Rang	1	2	3	4	5	6	7	8	9	10	11	12	13

Somme des rangs de A : $R_1 = 1 + 2 + 3 + 4 + 6 + 10 = 26$. Somme des rangs de B : $R_2 = 65$. On vérifie que $R_1 + R_2 = 13 \times 14 / 2 = 91$.

Les sessions A occupent surtout les premiers rangs (1 à 4 et 6), les sessions B les derniers : c'est ce déséquilibre des sommes de rangs que mesure le test.

La statistique U et sa loi

On transforme la somme des rangs en une statistique plus parlante :

$$U_1 = R_1 - \frac{n_1(n_1 + 1)}{2}$$

où $R_1$ est la somme des rangs du groupe 1 et $n_1(n_1+1)/2$ la plus petite valeur possible de $R_1$ (quand le groupe 1 occupe les rangs 1 à $n_1$). On définit $U_2$ de même, et $U_1 + U_2 = n_1 n_2$.

Ici, $U_1 = 26 - 21 = 5$ et $U_2 = 65 - 28 = 37$, avec bien $5 + 37 = 42 = 6 \times 7$.

Pourquoi cette transformation ? Parce que $U_1$ a une interprétation directe :

$$U_1 = \#\lbrace (i, j) : X_i > Y_j \rbrace$$

c’est le nombre de paires (une session A, une session B) dans lesquelles la session A est la plus longue (un ex aequo compte pour 1/2). Sur les $6 \times 7 = 42$ paires possibles, seules 5 donnent l’avantage à A : par exemple, la session A de 7,8 minutes bat celles de B qui durent 4,4 ; 5,7 ; 6,3 et 7,1 minutes, et celle de 5,1 minutes bat celle de 4,4.

Chacune des 42 paires (A, B) est un duel. A n'en gagne que 5 : U₁ = 5 et U₂ = 37.

Loi exacte

Sous $H_0$ (les deux échantillons viennent de la même loi continue), toutes les façons de choisir les $n_1$ rangs de A parmi les $N$ sont équiprobables. Il y en a $\binom{13}{6} = 1716$. En les énumérant, on trouve 19 configurations avec $U_1 \le 5$. La p-valeur bilatérale vaut donc :

$$p = 2 \times \frac{19}{1716} \approx 0{,}022$$

On rejette $H_0$ au seuil de 5 % : les sessions B sont plus longues.

Approximation normale

Pour de plus grands échantillons, on utilise, sous $H_0$ :

$$E(U) = \frac{n_1 n_2}{2}, \qquad \mathrm{Var}(U) = \frac{n_1 n_2 (n_1 + n_2 + 1)}{12}$$

Ici, $E(U) = 21$, $\mathrm{Var}(U) = 49$, donc $Z = (5 - 21)/7 \approx -2{,}29$ et $p \approx 0{,}022$. En présence d’ex aequo, on corrige la variance ; les logiciels ajoutent aussi souvent une correction de continuité de 1/2.

Sur les 1 716 façons de choisir les rangs de A, 19 donnent U₁ ≤ 5 et 19 donnent U₁ ≥ 37 : la valeur observée est dans une queue, d'où p ≈ 0,022.

À retenir

Le test de Mann-Whitney classe les deux échantillons ensemble et compare leurs sommes de rangs. Sa statistique $U_1$ compte les paires où le groupe 1 l’emporte. Sous $H_0$, $U$ est centré sur $n_1 n_2 / 2$.

Ce que teste vraiment Mann-Whitney

On lit souvent que Mann-Whitney « compare les médianes ». C’est vrai seulement sous une hypothèse supplémentaire. Le test est sensible à la quantité :

$$\theta = P(X > Y) + \tfrac{1}{2} P(X = Y)$$

où $X$ est une observation tirée du groupe 1 et $Y$ une observation tirée du groupe 2. L’estimateur naturel de $\theta$ est $U_1/(n_1 n_2)$. Sous $H_0$, $\theta = 1/2$ : une observation de A a autant de chances de battre une observation de B que l’inverse.

• Si les deux lois ont la même forme et ne diffèrent que par un décalage (modèle de translation, location shift), alors $\theta \ne 1/2$ équivaut à une différence de médianes, et le test compare bien les médianes.
• Si les formes diffèrent (par exemple, des variances très différentes), le test peut rejeter alors que les médianes sont égales, ou l’inverse.

À gauche, un pur décalage : θ ≠ 1/2 équivaut à des médianes différentes. À droite, les médianes sont égales et pourtant P(B > A) ≈ 0,60 : Mann-Whitney peut rejeter sans qu'aucune médiane ne diffère.

Attention

Mann-Whitney n’est pas un remède universel aux variances inégales. Si les deux groupes ont des dispersions très différentes et que vous voulez tester $\theta = 1/2$ de façon fiable, le test de Brunner-Munzel (scipy.stats.brunnermunzel) est plus adapté.

Estimer le décalage : Hodges-Lehmann

Une p-valeur ne dit pas de combien B dépasse A. Sous le modèle de translation, l’estimateur de Hodges-Lehmann du décalage est la médiane des $n_1 n_2$ différences $Y_j - X_i$. Ici, il vaut 3,2 minutes : la version B allonge les sessions d’environ 3 minutes. C’est l’estimation que fournit R avec conf.int = TRUE, accompagnée d’un intervalle de confiance.

En pratique : Python et R

import numpy as np
from scipy import stats

A = np.array([3.2, 5.1, 2.4, 7.8, 4.0, 3.6])
B = np.array([6.3, 8.9, 4.4, 12.5, 7.1, 5.7, 9.8])

res = stats.mannwhitneyu(A, B)          # bilatéral ; loi exacte pour les petits échantillons sans ex aequo
print(res.statistic, res.pvalue)        # 5.0 et environ 0.022

# Estimateur de Hodges-Lehmann du décalage B - A
print(np.median(B[:, None] - A[None, :]))   # 3.2

# Test t de Welch, pour comparaison
print(stats.ttest_ind(B, A, equal_var=False).pvalue)   # environ 0.023

scipy renvoie la statistique $U$ du premier échantillon passé en argument. En R, la même fonction que pour les données appariées sert, sans paired :

A <- c(3.2, 5.1, 2.4, 7.8, 4.0, 3.6)
B <- c(6.3, 8.9, 4.4, 12.5, 7.1, 5.7, 9.8)
wilcox.test(A, B)                   # W = 5, p = 0.02214
wilcox.test(B, A, conf.int = TRUE)  # estimation du décalage (Hodges-Lehmann) et IC

R appelle la statistique $W$, mais c’est bien $U_1$. D’où les deux noms du test : « Mann-Whitney U » et « test de la somme des rangs de Wilcoxon » (Wilcoxon rank-sum test).

Piège

Ne confondez pas les deux tests de Wilcoxon. Le test des rangs signés (signed-rank, chapitre 2) s’applique à des données appariées. Le test de la somme des rangs (rank-sum), c’est-à-dire Mann-Whitney, s’applique à deux groupes indépendants. En R, seul l’argument paired les distingue.

Mann-Whitney et l’AUC, c’est la même chose

Changeons de décor. Un classifieur attribue un score à des exemples positifs (des patients malades, des transactions frauduleuses) et à des exemples négatifs. Son AUC (aire sous la courbe ROC, présentée dans le cours d’Apprentissage automatique) admet une interprétation probabiliste : c’est la probabilité qu’un positif tiré au hasard reçoive un score plus élevé qu’un négatif tiré au hasard, les égalités comptant pour 1/2.

Vous reconnaissez $\theta$ ! L’AUC est exactement $U/(n_+ n_-)$, où $U$ est la statistique de Mann-Whitney des scores des positifs contre ceux des négatifs.

Exemple

Quatre positifs ont les scores 0,92 ; 0,81 ; 0,64 ; 0,47 et cinq négatifs les scores 0,70 ; 0,52 ; 0,33 ; 0,21 ; 0,15. Comptons, pour chaque positif, les négatifs qu’il dépasse :

• 0,92 et 0,81 dépassent les 5 négatifs : $5 + 5 = 10$ ;
• 0,64 dépasse 0,52 ; 0,33 ; 0,21 ; 0,15 : 4 ;
• 0,47 dépasse 0,33 ; 0,21 ; 0,15 : 3.

$U = 17$ sur $4 \times 5 = 20$ paires, donc AUC $= 17/20 = 0{,}85$. Par les rangs : les positifs occupent les rangs 9, 8, 6 et 4 parmi les 9 scores, $R = 27$, et $U = 27 - 4 \times 5 / 2 = 17$.

En parcourant les scores du plus grand au plus petit, chaque positif fait monter la courbe et chaque négatif la fait avancer. L'aire ombrée vaut 17/20 = 0,85 : exactement U/(n₊n₋).

Cette équivalence a trois conséquences pratiques :

• l’AUC ne dépend que du classement des scores, pas de leurs valeurs : toute transformation croissante des scores la laisse inchangée ;
• le calcul par les rangs coûte $O(N \log N)$, bien moins que l’énumération des $n_+ n_-$ paires ;
• tester « AUC = 0,5 » (le modèle ne fait pas mieux que le hasard), c’est faire un test de Mann-Whitney.

Pour notre test A/B, $U_2/(n_1 n_2) = 37/42 \approx 0{,}88$ : une session B tirée au hasard a 88 % de chances d’être plus longue qu’une session A. Cette « AUC » est une taille d’effet facile à expliquer, qu’on appelle aussi probabilité de supériorité.

En entretien

« Mann-Whitney et l’AUC sont la même quantité : $U$ divisé par le nombre de paires positif-négatif, c’est la probabilité qu’un positif soit mieux classé qu’un négatif. Cela explique pourquoi l’AUC ne dépend que de l’ordre des scores. »

En résumé

• Mann-Whitney compare deux groupes indépendants en classant toutes les observations ensemble.
• $U_1 = R_1 - n_1(n_1+1)/2$ compte les paires où le groupe 1 l’emporte ; $U_1 + U_2 = n_1 n_2$ ; sous $H_0$, $E(U) = n_1 n_2/2$ et $\mathrm{Var}(U) = n_1 n_2(n_1 + n_2 + 1)/12$.
• Le test porte sur $\theta = P(X > Y) + \frac{1}{2} P(X = Y)$ ; il compare les médianes seulement si les deux lois ne diffèrent que par un décalage.
• L’estimateur de Hodges-Lehmann (médiane des différences croisées) chiffre le décalage.
• L’AUC d’un classifieur vaut $U/(n_+ n_-)$ : c’est une statistique de Mann-Whitney.
• scipy.stats.mannwhitneyu ; en R, wilcox.test(x, y) sans paired.

Exercices

Exercice · Calculer U à la main

facile

Deux méthodes de révision sont testées sur deux groupes d’étudiants. Notes du groupe X : 12, 15, 9, 20. Notes du groupe Y : 18, 25, 22, 30, 16.

1. Calculez $R_X$, puis $U_X$, de deux façons (par les rangs et par les paires).
2. Sous $H_0$, $P(U_X \le 2) = 4/126$. Donnez la p-valeur bilatérale et concluez au seuil de 5 %.

Voir la solution

1. Classement des 9 notes : 9 (X), 12 (X), 15 (X), 16 (Y), 18 (Y), 20 (X), 22 (Y), 25 (Y), 30 (Y). Les rangs de X sont 1, 2, 3 et 6, donc $R_X = 12$ et $U_X = 12 - 4 \times 5/2 = 2$. Par les paires : seule la note 20 dépasse des notes de Y (16 et 18), soit 2 paires. Les deux méthodes concordent.
2. $p = 2 \times 4/126 \approx 0{,}063$. On ne rejette pas $H_0$ au seuil de 5 %, même si le groupe Y semble meilleur. Avec 4 et 5 observations, il faudrait une séparation presque parfaite pour conclure.

En code

Python

import numpy as np
from scipy import stats

X = np.array([12, 15, 9, 20])
Y = np.array([18, 25, 22, 30, 16])

rangs = stats.rankdata(np.r_[X, Y])
R_X = rangs[:len(X)].sum()
print(R_X, R_X - len(X) * (len(X) + 1) / 2)        # R_X = 12.0, U_X = 2.0
print((X[:, None] > Y[None, :]).sum())             # par les paires : 2

res = stats.mannwhitneyu(X, Y)                     # loi exacte
print(res.statistic, res.pvalue)                   # 2.0 0.0635 (= 2 x 4/126)

R

X <- c(12, 15, 9, 20)
Y <- c(18, 25, 22, 30, 16)

rangs <- rank(c(X, Y))
R_X <- sum(rangs[seq_along(X)])
c(R_X, R_X - 4 * 5 / 2)             # R_X = 12, U_X = 2
sum(outer(X, Y, ">"))               # par les paires : 2

wilcox.test(X, Y)                   # W = 2, p-value = 0.06349 (= 2 x 4/126)

Exercice · De U à l'AUC

moyen

Un modèle de détection de fraude est évalué sur 50 transactions frauduleuses et 200 transactions normales. La somme des rangs des scores des fraudes, parmi les 250 scores, vaut $R_+ = 8\,475$ (aucun ex aequo).

1. Calculez $U$ et l’AUC.
2. Que vaudrait l’AUC si l’on remplaçait chaque score $s$ par $\log(s)$ ?
3. Sous l’hypothèse d’un modèle sans pouvoir discriminant, calculez $E(U)$ et l’écart-type de $U$. L’AUC observée est-elle compatible avec le hasard ?

Voir la solution

1. $U = R_+ - n_+(n_+ + 1)/2 = 8\,475 - 50 \times 51/2 = 8\,475 - 1\,275 = 7\,200$. AUC $= 7\,200/(50 \times 200) = 0{,}72$.
2. Toujours 0,72 : le logarithme est croissant, il conserve l’ordre des scores, donc les rangs et $U$.
3. $E(U) = 50 \times 200/2 = 5\,000$ et $\mathrm{Var}(U) = 50 \times 200 \times 251/12 \approx 209\,167$, soit un écart-type d’environ 457. $Z = (7\,200 - 5\,000)/457 \approx 4{,}8$ : on est très loin de ce que produirait le hasard ($p$ de l’ordre de $10^{-6}$). Le modèle discrimine réellement, même si une AUC de 0,72 reste modeste.

En code

Python

import numpy as np
from scipy import stats

n_pos, n_neg, R_pos = 50, 200, 8475
U = R_pos - n_pos * (n_pos + 1) / 2
print(U, U / (n_pos * n_neg))                    # 7200.0 0.72

# Sous H0 : moyenne et écart-type de U, puis score z
E_U = n_pos * n_neg / 2
sd_U = np.sqrt(n_pos * n_neg * (n_pos + n_neg + 1) / 12)
z = (U - E_U) / sd_U
print(E_U, sd_U, z, 2 * stats.norm.sf(z))        # 5000.0 457.3 4.81 1.5e-06

R

n_pos <- 50; n_neg <- 200; R_pos <- 8475
U <- R_pos - n_pos * (n_pos + 1) / 2
c(U, U / (n_pos * n_neg))                        # 7200 et 0.72

# Sous H0 : moyenne et écart-type de U, puis score z
E_U <- n_pos * n_neg / 2
sd_U <- sqrt(n_pos * n_neg * (n_pos + n_neg + 1) / 12)
z <- (U - E_U) / sd_U
c(E_U, sd_U, z, 2 * pnorm(-z))                   # 5000 457.3 4.81 1.5e-06

Exercice · Pourquoi U compte les paires

difficile

Montrez que $U_1 = R_1 - n_1(n_1+1)/2$ est égal au nombre de paires $(X_i, Y_j)$ telles que $X_i > Y_j$ (on suppose qu’il n’y a pas d’ex aequo).

Indication : exprimez le rang de $X_{(k)}$, la $k$-ième plus petite valeur du groupe 1, dans l’échantillon mélangé.

Voir la solution

Dans l’échantillon mélangé, le rang de $X_{(k)}$ est égal au nombre de valeurs inférieures ou égales à $X_{(k)}$. Parmi elles, il y a exactement $k$ valeurs du groupe 1 ($X_{(1)}, \dots, X_{(k)}$), et un certain nombre $c_k$ de valeurs du groupe 2, celles qui sont plus petites que $X_{(k)}$. Donc $\text{rang}(X_{(k)}) = k + c_k$.

En sommant sur $k = 1, \dots, n_1$ :

$$R_1 = \sum_{k=1}^{n_1} (k + c_k) = \frac{n_1(n_1+1)}{2} + \sum_{k=1}^{n_1} c_k$$

Or $\sum_k c_k$ compte, pour chaque $X$, le nombre de $Y$ qu’il dépasse : c’est exactement le nombre de paires avec $X_i > Y_j$. D’où $U_1 = R_1 - n_1(n_1+1)/2 = \#\lbrace (i, j) : X_i > Y_j \rbrace$. Ce résultat justifie à la fois l’interprétation de $U$ et le calcul rapide de l’AUC par les rangs.

En code

Python

import numpy as np
from scipy import stats

# Vérification numérique sur des échantillons aléatoires sans ex aequo
rng = np.random.default_rng(0)
for n1, n2 in [(4, 5), (10, 7), (30, 40)]:
    X, Y = rng.normal(size=n1), rng.normal(0.5, 1, size=n2)
    R1 = stats.rankdata(np.r_[X, Y])[:n1].sum()
    U_rangs = R1 - n1 * (n1 + 1) / 2
    U_paires = (X[:, None] > Y[None, :]).sum()
    print(n1, n2, U_rangs, U_paires)     # les deux colonnes sont toujours égales

R

# Vérification numérique sur des échantillons aléatoires sans ex aequo
set.seed(0)
for (n in list(c(4, 5), c(10, 7), c(30, 40))) {
  X <- rnorm(n[1]); Y <- rnorm(n[2], 0.5)
  R1 <- sum(rank(c(X, Y))[1:n[1]])
  U_rangs <- R1 - n[1] * (n[1] + 1) / 2
  U_paires <- sum(outer(X, Y, ">"))
  cat(n, U_rangs, U_paires, "\n")      # les deux dernières valeurs sont toujours égales
}