De test werd gelijktijdig ontwikkeld door Wilcoxon en door Mann en Whitney. Om deze reden wordt de test dikwijls de Wilcoxon-Mann-Whitney (WMW) test genoemd. Soms wordt de test ook de Wilcoxon rank sum test of de Mann-Whitney U test genoemd.
De test werd ontwikkeld voor het testen van de nulhypothese (4) tegenover het alternatief \(H_1: \mu_1\neq \mu_2\) (of de eenzijdige versies). Eerst wordt er een distributionele veronderstelling gemaakt: het locatie-shift model, later relaxeren we deze aanname.
Stel dat \(Y_1\) en \(Y_2\) uitkomsten zijn uit respectievelijk de eerste en tweede behandelingsgroep, met respectievelijke verdelingen \(f_1\) en \(f_2\). Het locatie-shift model geldt als er een \(\Delta\) bestaat waarvoor geldt
\[f_1(y)=f_2(y-\Delta) \;\;\;\text{ voor alle } y.\]Locatie-shift betekent dat \(f_1\) en \(f_2\) dezelfde vorm hebben, maar ze mogen over \(\Delta\) verschoven zijn. De \(\Delta\) uit de definitie heeft als interpretatie: \(\Delta = \mu_1-\mu_2\). Door locatie-shift aan te nemen, zal het verwerpen van \(H_0: f_1=f_2\) de conclusie \(\mu_1\neq \mu_2\) impliceren.
De klassieke two-sample \(t\)-teststatistiek is gebouwd rond het verschil in steekproefgemiddelden \(\bar{Y}_1-\bar{Y}_2\). We beschouwen nu ook het verschil in steekproefgemiddelden, maar niet op basis van de oorspronkelijke uitkomsten, maar op basis van de rank-getransformeerde uitkomsten. De ranks zijn toegekend op basis van de gepoolde observaties (i.e. na samenvoegen van de uitkomsten uit groep 1 en groep 2); dus \(R_{ij}=R(Y_{ij})\) is de rank van uitkomst \(Y_{ij}\) in de gepoolde steekproef. Beschouw de teststatistiek
\[T = \frac{1}{n_1}\sum_{i=1}^{n_1} R(Y_{i1}) - \frac{1}{n_2}\sum_{i=1}^{n_2} R(Y_{i2}) .\]De statistiek vergelijkt dus de gemiddelde rank in groep 1 met de gemiddelde rank in groep 2.
Dit is een zinvolle teststatistiek, want
als \(H_0\) waar is, dan verwachten we dat de gemiddelde rank in de eerste groep ongeveer gelijk is aan de gemiddelde rank in de tweede groep en dus verwachten we dat \(T\) dicht bij nul ligt.
als \(H_1\) waar is dan verwachten we dat de gemiddelde ranks zullen verschillen en dus dat \(T\) niet dicht bij nul zal liggen.
Er kan echter worden aangetoond dat het volstaat het om
\[S_1=\sum_{i=1}^{n_1} R(Y_{i1})\]als teststatistiek te beschouwen. \(S_1\) is de som van de ranks van de observaties uit de eerste behandelingsgroep; dit verklaart de naam rank sum test.
\(S_1\) en \(S_2\) bevatten immers dezelfde informatie en zijn gerelateerd via
\[S_1+S_2 = \text{som van alle ranks} = 1+2+\cdots + n=\frac{1}{2}n(n+1).\]Nu we weten dat \(S_1\) (en \(S_2\)) een goede teststatistiek is, kan de permutatietestmethode toegepast worden om de exacte permutatienuldistributie op te stellen en de test uit te voeren. Voor een gegeven steekproefgrootte \(n\), en veronderstellend dat er geen ties zijn, nemen de rank-getransformeerde uitkomsten altijd de waarden \(1, 2, \ldots, n\) aan. Voor gegeven groepsgroottes \(n_1\) en \(n_2\), zal de permutatienuldistributie dan ook steeds dezelfde zijn! In de vorige eeuw (tot ongeveer de jaren 1980) werd dit als een groot voordeel beschouwd omdat de nuldistributies voor gegeven \(n_1\) en \(n_2\) getabuleerd konden worden (belangrijke kwantielen werden als tabellen in boeken gepubliceerd zodat ze konden gebruikt worden voor het bepalen van kritische waarden en \(p\)-waarden), waardoor de gebruiker geen nood had aan zware rekencapaciteit. Vandaag de dag speelt dit argument niet meer mee, maar toch blijven de rank testen erg populair, maar dan wel om andere, heel belangrijke redenen.
Niettegenstaande \(S_1\) en \(S_2\) perfect als teststatistieken gebruikt kunnen worden, wordt dikwijls gewerkt met de gestandaardiseerde teststatistiek
\[T = \frac{S_1-\text{E}_{0}\left[S_1\right]}{\sqrt{\text{Var}_{0}\left[S_1\right]}},\]met \(\text{E}_{0}\left[S_1\right]\) en \(\text{Var}_{0}\left[S_1\right]\) de verwachtingswaarde en variantie van \(S_1\) onder \(H_0\). Dit zijn dus het gemiddelde en variantie van de permutatienuldistributie van \(S_1\).
Onder \(H_0\) geldt
\[\text{E}_{0}\left[S_1\right]= \frac{1}{2}n_1(n+1) \;\;\;\;\text{ en }\;\;\;\; \text{Var}_{0}\left[S_1\right]=\frac{1}{12}n_1n_2(n+1).\]Verder kan men onder \(H_0\) en als \(\min(n_1,n_2)\rightarrow \infty\) opgaat aantonen dat,
\[T = \frac{S_1-\text{E}_{0}\left[S_1\right]}{\sqrt{\text{Var}_{0}\left[S_1\right]}} \rightarrow N(0,1).\]Asymptotisch volgt de gestandaardiseerde teststatistiek dus een standaardnormaal verdeling.
We illustreren de WMW test aan de hand van de R functie wilcox.test.
wilcox.test(cholest ~ group, data = chol)
##
## Wilcoxon rank sum exact test
##
## data: cholest by group
## W = 24, p-value = 0.01587
## alternative hypothesis: true location shift is not equal to 0
We zien dat we op basis van de test de nulhypothese kunnen verwerpen op het 5% significantie-niveau.
De output geeft de teststatistiek \(W=\) 24. In volgende lijnen berekenen we \(S_1\) en \(S_2\) manueel voor de dataset.
S1 <- chol %>%
mutate(cholRank=rank(cholest)) %>%
filter(group==1) %>%
pull(cholRank) %>%
sum
S2 <- chol %>%
mutate(cholRank = rank(cholest)) %>%
filter(group == 2) %>%
pull(cholRank) %>%
sum
S1
## [1] 39
S2
## [1] 16
Waar komt \(W=\) 24 vandaan? Dit wordt zodadelijk toegelicht.
De teststatistieken \(S_1\) en \(S_2\) werden voorgesteld door Wilcoxon, maar tezelfdertijd werd een equivalente test voorgesteld door Mann en Whitney. Hun teststatistiek wordt gegeven door[51]
\[U_1 = \sum_{i=1}^{n_1}\sum_{k=1}^{n_2} \text{I}\left\{Y_{i1}\geq Y_{k2}\right\}.\]waarbij \(\text{I}\left\{.\right\}\) een indicator is die 1 is als de uitdrukking waar is en 0 als dit niet het geval is. Er wordt voor elke observatie uit de eerste groep geteld hoeveel keer zij groter of gelijk is aan een observatie uit de tweede groep. We berekenen de Mann-Whitney statistiek nu manueel in R.
y1 <- chol %>%
filter(group==1) %>%
pull("cholest")
y2 <- chol %>%
filter(group==2) %>%
pull("cholest")
u1Hlp <- sapply(
y1,
function(y1i,y2) {y1i>=y2},
y2=y2)
colnames(u1Hlp) <- y1
rownames(u1Hlp) <- y2
u1Hlp
## 244 206 242 278 236
## 188 TRUE TRUE TRUE TRUE TRUE
## 212 TRUE FALSE TRUE TRUE TRUE
## 186 TRUE TRUE TRUE TRUE TRUE
## 198 TRUE TRUE TRUE TRUE TRUE
## 160 TRUE TRUE TRUE TRUE TRUE
U1 <- sum(u1Hlp)
U1
## [1] 24
Er kan worden aangetoond dat
\[U_1 = S_1 - \frac{1}{2}n_1(n_1+1).\]S1-nGroups[1]*(nGroups[1]+1)/2
## 1
## 24
Hieruit concluderen we (1) dat \(U_1\) en \(S_1\) dezelfde informatie bevatten, (2) dat \(U_1\) ook een rankstatistiek is en dat exacte testen gebaseerd op \(U_1\) en \(S_1\) equivalent zijn.
De statistiek \(U_1\) heeft als voordeel dat het een informatieve interpretatie heeft. Stel \(Y_j\) een willekeurige uitkomst uit behandelingsgroep \(j\) (\(j=1,2\)). Dan geldt
\[\begin{eqnarray*} \frac{1}{n_1n_2}\text{E}\left[U_1\right] &=& \text{P}\left[Y_1 \geq Y_2\right]. \end{eqnarray*}\]Intuïtief voelen we dit aan: Op basis van de steekproef kunnen we die kans schatten door het gemiddelde te berekenen van alle indicator waarden \(\text{I}\left\{Y_{i1}\geq Y_{k2}\right\}\). We voerden inderdaad \(n_1 \times n_2\) vergelijkingen uit.
mean(u1Hlp)
## [1] 0.96
U1/(nGroups[1]*nGroups[2])
## 1
## 0.96
De kans \(\text{P}\left[Y_1 \geq Y_2\right]\) wordt een probabilistische index (Engels: probabilistic index) genoemd. Het is de kans dat een uitkomst uit de eerste groep groter of gelijk is dan een uitkomst uit de tweede groep. Als \(H_0\) waar is, dan is \(\text{P}\left[Y_1 \geq Y_2\right]=\frac{1}{2}\).
De gestandaardiseerde Mann-Whitney statistiek is
\[T = \frac{U_1 - \frac{n_1n_2}{2}}{\sqrt{\frac{1}{12}n_1n_2(n+1)}}.\]De R functie wilcox.test geeft niet de Wilcoxon rank sum statistiek, maar wel de Mann-Whitney statistiek \(U_1\).
We weten echter dat exacte permutatietesten gebaseerd op \(U_1\), \(U_2\), \(S_1\) of \(S_2\) dezelfde resultaten geven.
We bekijken nogmaals de output
wTest <- wilcox.test(cholest~group, data = chol)
wTest
##
## Wilcoxon rank sum exact test
##
## data: cholest by group
## W = 24, p-value = 0.01587
## alternative hypothesis: true location shift is not equal to 0
U1
## [1] 24
probInd <- wTest$statistic/prod(nGroups)
probInd
## W
## 0.96
Aangezien \(p=\) 0.0159 \(<0.05\) besluiten we op het \(5\%\) significantieniveau dat de gemiddelde cholestorolconcentratie groter is bij hartpatiënten kort na een hartaanval dan bij gezonde personen. We nemen aan dat locatie-shift opgaat.
Nu we weten hoe \(U_1\) berekend wordt, weten we ook meteen dat een cholestorolwaarde van hartpatiënten met een kans van \(U1/(n_1\times n_2)=\) 96% groter is die van gezonde personen. Aangezien we het locatie-shift model veronderstellen, besluiten we ook dat de gemiddelde uitkomst uit de behandelingsgroep groter is dan de gemiddelde uitkomst uit de placebogroep.
We zouden de veronderstelling van de locatie-shift moeten nagaan, maar met slechts 5 observaties in elke behandelingsgroep is dit zinloos. Zonder verder theorie hierover te geven, geven we nog mee dat zonder de locatie-shift veronderstelling de conclusie in termen van de probabilistische index correct blijft en de conclusie ook zo zou moeten worden geformuleerd.
Dus wanneer we geen locatie-shift veronderstellen en een tweezijdige test uitvoeren testen we eigenlijk
\[H_0: F_1=F_2 \text{ vs P}(Y_1 \geq Y_2) \neq 0.5.\]Er is een significant verschil in de distributie van de cholestorolconcentraties bij hartpatiënten 2 dagen na hun hartaanval en gezonde individuen (\(p=\) 0.0159). Het is meer waarschijnlijk om hogere cholestorolconcentraties te observeren bij hartpatiënten dan bij gezonde individuen. De puntschatting voor deze kans bedraagt 96%.