R zonder P

Ik zou willen dat bepaalde standaardfunctionaliteiten wereldwijd uit statistische software gesloopt worden.

In vrijwel alle toegepaste wetenschappen spelen statistische analyses een centrale rol. En in vrijwel al die analyses wordt ‘mate van evidentie’ aangegeven door middel van de p-waarde. Die p-waarde ligt al lange tijd zwaar onder vuur, om een heleboel redenen: ten eerste, de p-waarde wordt veelal gebruikt in situaties waar dat gewoon wiskundig niet correct is, of het niet duidelijk is wat zij betekent. Ten tweede: zelfs indien correct gebruikt, is p ontzettend moeilijk te interpreteren: uit meerdere studies blijkt dat de meeste artsen, die toch allemaal regelmatig vakliteratuur met p-waarden lezen, er niets van begrijpen. En helaas ben ik er in de loop van mijn carrière achtergekomen dat ook statistici zélf er vaak het fijne niet van begrijpen. Ten derde: p is ‘krakkemikkig’ – als er ook maar iets in de analyse anders ging dan vooraf gepland, dan verliest p zijn betekenis. Ten vierde – als gevolg hiervan – is het buitengewoon eenvoudig om aan p-hacking te doen: de resultaten mooier doen voorkomen dan ze zijn – iets dat vaak onbewust, maar toch, gebeurt. Er zijn de afgelopen jaren een aantal invloedrijke papers verschenen waarin het gebruik van p-waarden simpelweg sterk wordt afgeraden. Bovendien zijn er tegenwoordig goede alternatieven, waarover zo meteen meer.

Maar toch blijkt het in de praktijk vrijwel onmogelijk om onderzoekers van hun p-verslaving af te houden: er zijn miljoenen onderzoekers die p gebruiken, een alternatieve aanpak leren beheersen kost veel tijd, men moet referenten overtuigen dat de alternatieve aanpak valide is, enzovoort. Er is nog een reden: SPSS en R, ’s werelds meest gebruikte statistische softwareprogramma’s,  leveren ongevraagd p-waardes wanneer je bijvoorbeeld een regressie-analyse of een andere statistische toets doet (vaak krijg je p onder een andere naam te zien, bij SPSS staat er bijvoorbeeld ‘significance’). Maar dan maak je het onderzoekers natuurlijk wel héél makkelijk om gewoon te blijven p-(hack-)en. En wordt het natuurlijk wel héél moeilijk om de gangbare praktijk ooit te veranderen. Zeker voor R, de favoriete software van statistici zelf, (software met het predicaat ‘wetenschappelijk’ en ‘open-source’), is dit kwalijk.

Hier is een voorbeeld van een regressieanalyse van de Boston housing data. Dit is data over 506 huizen, waarbij voor elk huis een aantal gegevens zoals verkoopprijs, stikstofpercentages, misdaadcijfers en industriele activiteit in de directe omgeving worden aangegeven. Ik kijk of die laatste drie gerelateerd kunnen worden aan 506 totaal willekeurig gegenereerde (normaal verdeelde) punten. Ik vond inderdaad bij mijn allereerste poging een significant (‘p < 0.05’) resultaat (de nabijheid van industrie heeft voorspellende waarde voor mijn willekeurige data), wat nog maar eens laat zien hoe makkelijk het mis kan gaan met p-waardes. Maar dat is niet mijn hoofdpunt: het hoofdpunt is dat in de output van R maarliefst vier p-waarden getoond worden (hieronder in rood). Dat hoeft niet!

Call:
lm(formula = willekeur ~ stikstof + misdaad + industrie)

Residuals:
    Min      1Q  Median      3Q     Max
-3.7082 -0.6176  0.0087  0.6550  2.7724

             Estimate Std. Error t value Pr(>|t|) 
(Intercept) -0.386783   0.268072  -1.443   0.1497 
stikstof     1.138110   0.611801   1.860   0.0634 .
misdaad     -0.003860   0.005824  -0.663   0.5078 
industrie   -0.024558   0.010260  -2.394   0.0170 *

---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.01 on 502 degrees of freedom
Multiple R-squared:  0.01342,       Adjusted R-squared:  0.007527
F-statistic: 2.277 on 3 and 502 DF, p-value: 0.07888

Wat zou je dan wel moeten laten zien? Mijn eigen onderzoeksgroep heeft de laatste jaren intensief gewerkt aan e-waardes, een veel robuuster alternatief voor de p-waarde. Zowel p-waardes als e-waardes meten evidentie dat er echt iets aan de hand is, dus bijvoorbeeld ‘er is een echt verband tussen industrie en willekeur, wat we zagen is geen toeval’), of, in de geneeskunde `het verschil tussen de groep die het geneesmiddel krijgt en de groep die het placebo krijgt is zo groot; dit nieuwe geneesmiddel doet echt iets’). De precieze definitie van e-waarde is te ingewikkeld om hier te geven (maar nog altijd simpeler dan die van de p-waarde), maar heel ruwweg gezegd is de e-waarde de hoeveelheid geld die je wint in een bepaalde reeks weddenschappen die zo opgezet zijn, dat je verwacht geen geld te winnen als er niets aan de hand is.

Hoe werkt de e-waarde?

Hier is een kleine illustratie. Stel, je gaat naar het casino met 100 euro. Je vermoedt dat het geen echt casino is, en dat rood een grotere kans dan zwart heeft. Je kunt proberen hiervan te profiteren door steeds iets meer dan de helft van je geld op rood en iets minder op zwart te zetten, en het geld dat je na elke ronde krijgt te op dezelfde manier te herinvesteren (het is niet slim om alles op rood te zetten, want dan is er een kans dat je al je geld verliest en kun je niet herinvesteren). Als je dit spelletje nu een aantal keren speelt en je hebt op een gegeven moment 2000 euro in handen, dan is je beginkapitaal dus 20 keer zo groot geworden. We zeggen dan E=20, de evidentie dat dit geen gewoon casino is en dat rood echt vaker valt dan zwart is 20. Je kunt laten zien dat, als het wel degelijk een echt casino is, dan is de kans dat E ooit groter wordt dan 20, kleiner  dan 1/20, de kans dat E ooit groter wordt dan 100 kleiner dan 1/100, enzovoort. Het wonderbaarlijke is dat je elk statistisch toetsprobleem kunt vertalen naar zo’n casino-setting: de nulhypothese (`geneesmiddel doet niets’, `industrie-nabijheid heeft geen invloed op ‘willekeur’, `muntje is eerlijk’) correspondeert met `het is een echt casino’. De e-waarde is de hoeveelheid geld die je wint in een spelletje waarbij je voor elk data-punt gokt wat het wordt en meer geld verdient naarmate je er dichter bij zit. Als de nulhypothese waar is, blijft e met grote kans voor altijd klein. Als e groot wordt heb je (virtueel) veel geld verdiend en heb je evidentie dat de nulhypothese misschien niet klopt. Een cruciaal verschil met p-waardes is dat e-waardes geldig zijn,  wanneer je ook maar stopt – net zoals de kans dat je in een echt casino rijk wordt heel klein is, wat ook je regel is voor wanneer je naar huis gaat. Als je zegt ‘ik speel lekker net zolang door totdat ik meer dan 2000 euro verdiend heb (E > 20) of ze me eruit gooien omdat het casino sluit’, dan helpt dat niet. Bij p-waardes werkt dit niet zo: als je zegt ‘ik ga door met data punten bekijken tot p < 0.05, dan weet je zeker dat dat op een bepaald moment gaat lukken, ook al is de nulhypothese waar en is er niets aan de hand.

In 2019 waren e-waardes nog nauwelijks bekend, nu zijn er 100en papers aan gewijd in top-tijdschriften zoals The Annals of Statistics, The Journal of the Royal Statistical Society, Proceedings of the National Academy of Sciences of the USA, zijn er meerdere prestigieuze prijzen uitgedeeld aan jonge onderzoekers die hier wereldwijd aan werken, zijn er workshops georganiseerd, enzovoort, enzovoort. Maar we `vechten’ natuurlijk wel tegen een paradigma dat bijna 100 jaar oud is. Een van de belangrijkste dingen die nu moeten gebeuren, is dat e-waardes een plek krijgen in standaard software zoals R en talen zoals Python. Daarvoor hebben we goede programmeurs nodig met veel verstand van statistiek. Maar die zijn moeilijk te krijgen, want met zo’n profiel kun je zo aan de slag bij één van de grote tech-bedrijven, en dan kun je heel veel geld verdienen…

Laat ik dit echter nuanceren: ik droom ervan (en deze droom is wellicht zelfs bereikbaar) dat berekening van e-waardes mogelijk wordt in standaard software. Maar zelfs e-waardes zijn niet zaligmakend: als ik heel eerlijk ben zou ik willen dat noch e-waardes, noch p-waardes standaard getoond worden: je zou elke p- of e-waarde, voor elke variabele, expliciet moeten opvragen. Dan stimuleer je onderzoekers om na te denken over hun statistische analyse – dat gebeurt nu veel te weinig. Wat moet de output zijn, en waarom? Misschien helpt het een beetje als ze die output niet zomaar over zich uitgestort krijgen, maar expliciet moeten aangeven wat ze specifiek zouden willen weten over hun data.