Meervoudige regressie

multiple regression

Vooraf, in één zin: een tabel met vier voorspellers onder elkaar nodigt uit om naar het sterretje te zoeken, maar de eerste vraag is hoeveel voorspellers er op hoeveel mensen losgelaten zijn — hieronder wordt dezelfde voorspeller pas significant in het model dat níéts verklaart, en waarom dat zo werkt staat onderaan bij Wat het niet zegt.

Wat je ziet

In het kort. Geen enkel getal, maar een tabel: een rijtje voorspellers onder elkaar, elk met een b, een standaardfout en een p, en eronder één regel voor het model als geheel — een F en een R². Dat is een meervoudige regressie / multiple regression: één uitkomst, meerdere voorspellers tegelijk. Je komt hem tegen zodra een onderzoek wil weten wat iets voorspelt terwijl er van alles anders meespeelt: leeftijd, geslacht, hoe ernstig iemand begon.

Zo ziet het eruit in het wild:

Voorbeeldfragment, geschreven in de stijl van een resultatensectie — het artikel bestaat niet.

Change in total difficulties (SDQ, T0 minus T1; higher = improvement) was regressed on baseline severity (centered), age (centered) and gender in the 24 adolescents who completed the VR intervention. The model did not explain a significant proportion of the variance, F(3, 20) = 1.58, p = .22, R² = .19. Baseline severity was the only predictor approaching significance, b = 0.38, SE = 0.19, t = 1.94, p = .066, 95% CI [−0.03, 0.78]. Neither age (b = 0.27, SE = 0.44, p = .55) nor gender (b = −1.16, SE = 1.95, p = .56) contributed.

Wat er staat

In het kort. Vierentwintig jongeren, drie voorspellers, en een model dat het niet doet: de F van 1.58 met p = .22 zegt dat deze drie samen niet meer verklaren dan toeval kan maken. De R² van .19 zegt dat ze op papier 19% van de verschillen in verbetering dekken — een getal dat groot lijkt en het hier niet is, want met 24 mensen en drie voorspellers haal je die 19% ook zonder dat er iets aan de hand is. En de enige regel die in de buurt komt, baseline, is precies de regel waar je het minst aan hebt. Waarom, dat is de rest van dit blok.

Nu rustig van voren af aan. Wat een b betekent naast andere b’s. In de enkelvoudige regressie was b het verschil in de uitkomst dat hoort bij één punt verschil in de voorspeller. Dat blijft zo — er komt één woord bij: terwijl de andere voorspellers gelijk blijven. Hier betekent b = 0.38 dus: twee jongeren die één SDQ-punt in ernst verschilden bij binnenkomst, verschilden gemiddeld 0.38 punt in verbetering — bij gelijke leeftijd en gelijk geslacht.

Let op de vorm van die zin. Er staat verschilden, geen stegen. Er is hier niets gemanipuleerd: niemand heeft iemands beginscore hoger gezet om te kijken wat er gebeurt. Twee jongeren vergelijken is iets anders dan één jongere veranderen, en de taal moet dat verschil vasthouden. Dat is geen muggenzifterij; het is het enige wat een lezer tegen de verleiding beschermt om er een oorzaak van te maken.

Wat “gecorrigeerd voor” doet. Elke b in deze tabel is berekend met de andere twee erin. Je kunt het lezen als: eerst haalt het model uit de verbetering weg wat leeftijd en geslacht ervan kunnen verklaren, en pas over wat er dan nog ligt, legt het de lijn voor baseline. Vandaar de Nederlandse formulering die je in artikelen ziet — gecorrigeerd voor leeftijd en geslacht, of adjusted for. Wat je overhoudt, is het stuk van de samenhang dat níét al door de andere voorspellers verteld werd.

En dan b tegenover β. In veel tabellen — en in vrijwel elk artikel dat je in je literatuurstudie tegenkomt — staat naast de b een tweede kolom: β (bèta), soms met gestandaardiseerd of standardized erboven. Het verschil is de meetlat. De b staat in de eenheid van de wereld: hier 0.38 SDQ-punt verbetering per SDQ-punt ernst bij binnenkomst. De β heeft die eenheden weggedeeld, net zoals Cohens d en de correlatie r dat doen — hij zegt hoeveel standaardafwijkingen de uitkomst opschuift per standaardafwijking voorspeller. Daarom is de b het getal dat je aan een patiënt of een collega uitlegt, en de β het getal waarmee je voorspellers onderling vergelijkt: leeftijd in jaren en ernst in SDQ-punten zijn anders onvergelijkbaar, want de ene meetlat is nu eenmaal grover dan de andere. Twee dingen om vast te houden. De p en het betrouwbaarheidsinterval horen bij dezelfde toets — β verandert de eenheid, niet het verhaal. En staat er in een tabel alléén een β-kolom, dan is de vertaalslag terug naar de meetlat van de wereld jouw werk, en je hebt er de standaardafwijkingen uit de methodesectie voor nodig.

Wat de F en R² samen zeggen. De R² van .19 is het aandeel van de verschillen in verbetering dat het hele model dekt: negentien procent. De F stelt daar de nuchtere vraag bij — is dat meer dan je met drie willekeurige voorspellers op 24 mensen sowieso zou halen? Antwoord: nee, p = .22. Kijk daarom altijd eerst naar die regel. Een tabel vol voorspellers waarvan het model als geheel niet significant is, is een tabel waarin je geen enkele individuele regel als bevinding mag oppakken — en dat is precies wat er in dit fragment gebeurt.

En de baseline-regel dan, met p = .066? Twee dingen. Ten eerste is .066 groter dan .05, en “randsignificant” is geen aparte categorie — het betrouwbaarheidsinterval [−0.03, 0.78] loopt over de nul heen, dus een verband van nul past nog steeds bij deze data. Ten tweede, en belangrijker: zelfs als deze regel wél een sterretje had gehad, zou hij je weinig hebben verteld. Waarom niet, staat bij Wat het niet zegt.

Waar het vandaan komt

Dit blok leunt op de boek-eenheid Meer dan één verschil. Daar wordt de stap gezet die dit hele blok draagt: van één rechthoekje naar meerdere. Zolang je één voorspeller hebt, is regressie een lijn door een wolk. Komt er een tweede bij, dan wordt het wegrekenen — je haalt uit de uitkomst eerst weg wat de ene voorspeller verklaart, en kijkt wat de andere over de rest nog te zeggen heeft. Het boek tekent dat als een tweeluik: dezelfde wolk, één keer met en één keer zonder wat de eerste voorspeller al had opgeslokt.

Twee oude bekenden komen daarbij terug. R² ken je als r² uit het correlatie-blok — hetzelfde idee, nu met meerdere verklaarders tegelijk. En de F ken je van de ANOVA: nog steeds een breuk van verklaarde en onverklaarde spreiding, alleen heet “de groepen” hier “de voorspellers samen”.

Zelf doen in JASP

⬇ Download de data — verschilscore_mra.csv

Bij dit blok hoort verschilscore_mra.csv: 24 rijen, zeven kolommen — id, tds_t0 en tds_t1 (de SDQ-totaalscore vóór en na), verbetering (T0 min T1, dus positief = minder klachten), leeftijd, geslacht en diagnose. Het zijn dezelfde 24 jongeren als in het t-toets-blok; er is alleen bij gekomen wat je voor dit model nodig hebt. De rijen zijn gegenereerd — er zit geen echte jongere in — en ze reproduceren precies de getallen uit het fragment.

OpmerkingIn JASP · JASP 0.98.1

1. Open het bestand in JASP (menu linksboven → Open).
2. Klik in de balk bovenin op Regression en kies Linear Regression.
3. Sleep verbetering naar Dependent Variable, tds_t0 en leeftijd naar Covariates, en geslacht naar Factors.
4. Vink onder Statistics aan: Estimates, Model fit en Confidence intervals.
5. Wil je de twee modellen naast elkaar zien: open het tabblad Model, maak een tweede Block aan en zet diagnose daarin. JASP toont dan model 1 en model 2 onder elkaar, met de R² van beide.

In de tabel rechts staat nu bij tds_t0 een b van 0.377 met SE 0.194, t = 1.944, p = .066 en een 95%-interval van [−0.03, 0.78]; bij leeftijd 0.266 (p = .55) en bij geslacht −1.155 (p = .56). In de Model Summary: R² = .192, gecorrigeerde R² = .071, en in de ANOVA-regel F(3, 20) = 1.58, p = .22. Leg het fragment ernaast — dezelfde getallen.

Eén regel van die tabel bevat vijf getallen die je in die volgorde leest, en ze hangen aan elkaar: de b is de schatting, de SE is hoe hard hij meewiebelt, de t is niets anders dan die twee op elkaar gedeeld (0.377 gedeeld door 0.194 is 1.94 — op de afronding na precies de t die JASP toont, en de betekenis is: hoeveel standaardfouten ligt de schatting van nul vandaan), de p hoort bij die t, en het 95%-interval is dezelfde informatie in de eenheid van de uitkomst. Zie je in een artikel de vorm t(20) = 1.94 staan, dan is dat getal tussen haakjes het aantal vrijheidsgraden: hier 24 deelnemers min 3 voorspellers min 1 voor de intercept, dezelfde 20 als in de F(3, 20). Eén regel, vijf getallen, één verhaal.

Twee vakken, twee woorden. JASP noemt die vakken Covariates en Factors, en dat zijn precies de woorden die je in artikelen en in je bijeenkomsten tegenkomt. Een covariabele (of covariaat) is een voorspeller die op een schaal gemeten is — leeftijd, beginscore — en die meestal wordt meegenomen om ervoor te corrigeren; dat is hetzelfde ding als wat wij hier “voorspeller” noemen, en het is de andere kant van de uitdrukking gecorrigeerd voor leeftijd. Zie je in een artikel een lijstje onder het kopje Covariabelen, dan lees je daar dus: dit is wat er is weggerekend, en alles wat er niet staat, is niet weggerekend. Een factor is een voorspeller die uit categorieën bestaat — geslacht, conditie, diagnose. Daar maakt JASP zelf dummy’s van, en één categorie wordt dan de maatstaf waar de andere tegen worden afgezet. Zeg je “factor”, dan noem je het vak in JASP; zeg je “dummy”, dan noem je wat er onder de motorkap gebeurt. Hetzelfde ding, twee lagen.

Eén ding valt op als je goed kijkt: het fragment noemt baseline centered en jij sleepte gewoon tds_t0 naar binnen. Toch komen de b, de p en het interval exact overeen. Dat klopt: centreren verschuift alleen het nulpunt, en daarmee de intercept — de helling verandert niet. Wat centreren wél doet, is de intercept betekenis geven (de verwachte verbetering bij een gemiddelde beginscore in plaats van bij een SDQ van 0, wat niemand heeft). Het is dus een leesbaarheidstruc, geen ander model.

Zelf opschrijven

OpmerkingSchrijf het op — de meervoudige regressie

Voor je eigen Nederlandstalige literatuurstudie:

De verandering in totale probleemscore werd voorspeld uit de beginscore, leeftijd en geslacht. Het model verklaarde geen significant deel van de verschillen in verandering, F(3, 20) = 1.58, p = .22, R² = .19. Jongeren die één punt hoger begonnen, verbeterden gemiddeld 0.38 punt meer, maar dit verband was niet significant, b = 0.38, SE = 0.19, p = .066, 95% BI [−0.03, 0.78].

Drie werkafspraken daarbij:

• Rapporteer het model vóór de regels. F, R² en het aantal voorspellers eerst; pas daarna de losse coëfficiënten. Wie met een significante regel begint, verzwijgt de vraag of het model überhaupt iets deed.
• Noem N én het aantal voorspellers. Vier voorspellers op 24 mensen is een ander onderzoek dan vier op 400, ook als de tabellen identiek ogen. Zonder die twee getallen kan je lezer de tabel niet wegen.
• Verschil-taal, geen verandering-taal. “Jongeren die hoger begonnen, verbeterden meer” — niet “hoger beginnen zorgt voor meer verbetering”. Er is niets gemanipuleerd; er zijn alleen mensen met elkaar vergeleken.

Wat het niet zegt

De verleidelijkste redenering die je over dit fragment kunt opschrijven is: “baseline voorspelt verbetering (p = .066, bijna significant), dus wie zieker binnenkomt, heeft meer baat bij de VR-training.” Eén zin, drie fouten.

“Bijna significant.” Die categorie bestaat niet. Het interval loopt van −0.03 tot 0.78: nul zit er nog in, dus “geen verband” is een uitkomst die bij deze data past. En het model als geheel deed het niet (p = .22). Wat je mag opschrijven is dat dit onderzoek geen significante voorspellers van verandering vond — niet dat het er bijna één vond.

“Wie zieker binnenkomt, heeft meer baat.” Hier zit de val die in pre-post-onderzoek zonder controlegroep altijd klaarligt: regressie naar het gemiddelde. Wie hoog scoort op een meting, zit daar deels hoog door een slechte dag of door meetruis; bij een tweede meting zakt dat stukje ruis vanzelf weg. Hoog beginnen voorspelt daarom bijna mechanisch een grotere daling — ook als de behandeling niets doet, en zelfs als er helemaal geen behandeling is. De baseline-regel in deze tabel is dus de regel waar je het minst uit mag concluderen, niet de meest interessante. Wil je weten of de VR-training werkt, dan heb je een groep nodig die hem niet kreeg.

“Gecorrigeerd voor leeftijd en geslacht, dus het ligt niet daaraan.” Corrigeren haalt eruit wat je hebt gemeten, en niets anders. Alles wat niet in de tabel staat — motivatie, thuissituatie, of iemand er tegelijk een andere behandeling bij kreeg — zit nog gewoon in dat verband. “Gecorrigeerd voor” is geen synoniem van “dus causaal”; het is een lijst van wat er wél is meegenomen, en die lijst lees je als lezer juist op wat er níét in staat.

En dan de belangrijkste, want die kun je in dit databestand zelf zien gebeuren. Voeg aan dit model de diagnose toe — zeven categorieën, dus zes dummy’s erbij — en je hebt negen voorspellers op 24 mensen. Wat er dan gebeurt:

	model 1 (3 voorspellers)	model 2 (9 voorspellers)
R²	.19	.45
gecorrigeerde R²	.07	.10
F	F(3, 20) = 1.58, p = .22	F(9, 14) = 1.27, p = .33
baseline	b = 0.38, p = .066	b = 0.56, p = .018

Lees die tabel van links naar rechts en er gebeurt iets ongemakkelijks. De R² springt van .19 naar .45 — het model lijkt twee keer zo goed. Het model als geheel wordt er níét beter van (p gaat van .22 naar .330). En de baseline-regel, die net nog niet significant was, krijgt een sterretje: p = .018.

Er is niets veranderd aan de wereld. Er is alleen ruis toegevoegd: dummy’s voor categorieën waarvan er drie precies één jongere bevatten. Zo’n dummy kan die ene persoon perfect natekenen, en dat kost niets aan uitleg maar wel aan vrijheidsgraden — van 20 naar 14. R² kán daarbij alleen maar omhoog: elke voorspeller die je toevoegt, verklaart altijd iets, al is het toeval. Dáárvoor bestaat gecorrigeerde R², die een verwachte toevalswinst aftrekt; die blijft hier rond de tien procent hangen en verraadt zo dat de sprong lucht was.

Dit is geen verzonnen schrikbeeld. In het echte onderzoek waar dit cohort op gebaseerd is, gebeurde precies dit: het eerste model had negen voorspellers op 24 deelnemers, en daarin leek een diagnose-dummy significant (p = .017). Afslanken naar drie voorspellers en er bleef niets van over. Bij ons springt het sterretje op een andere regel, maar het mechanisme is hetzelfde.

De leesregel die je hieraan overhoudt: tel bij elke regressietabel eerst het aantal voorspellers en zoek de N erbij. Een vuistregel die in de literatuur rondgaat is tien tot vijftien deelnemers per voorspeller; negen voorspellers vragen dus honderd tot honderdvijftig mensen, en dit onderzoek had er 24. Zie je een tabel met veel voorspellers, weinig deelnemers en één sterretje — dan heb je iets gevonden voor het kopje “beperkingen”, geen bevinding.

Oefening

OpmerkingOefening

Leesopgave. Hieronder een voorbeeldfragment (verzonnen). Vertel in gewone Nederlandse zinnen wat er staat. Wat weet je na deze zinnen wél, en wat niet? En welke twee getallen zou je als eerste opzoeken voor je er iets over opschrijft?

A multiple regression with eight predictors (age, gender, education, duration of illness, medication use, comorbidity, therapist experience and baseline severity) on 31 patients explained 52% of the variance in treatment outcome, R² = .52. Therapist experience emerged as a significant predictor, b = 2.14, p = .041. The authors recommend matching experienced therapists to complex cases.

Doe-opgave. Draai op verschilscore_mra.csv eerst het model met drie voorspellers en daarna dat met negen (diagnose erbij, als tweede Block). Noteer van allebei R², gecorrigeerde R², F met vrijheidsgraden, en de regel van tds_t0. Beantwoord dan drie vragen. (1) Welk getal in de tweede tabel zou je op het verkeerde been zetten als je alleen díé tabel zag? (2) Waarom is gecorrigeerde R² nauwelijks meegestegen met R²? (3) De vrijheidsgraden van de F gingen van (3, 20) naar (9, 14) — reken na waar die 14 vandaan komt, en wat dat betekent voor hoe scherp dit onderzoek nog iets kán aantonen.

OpgeletUitwerkingen — klik open

Uitwerking leesopgave

Wat er staat: acht voorspellers op 31 patiënten verklaren samen 52% van de verschillen in behandeluitkomst, en één daarvan — ervaring van de behandelaar — heeft een p onder de .05. Wat je wél weet: in déze 31 patiënten hangt ervaring samen met uitkomst, in dit model. Wat je niet weet: zo’n beetje al het andere. Acht voorspellers op 31 mensen is ongeveer vier mensen per voorspeller, ver onder elke vuistregel; een R² van .52 is bij die verhouding precies wat je van toeval mag verwachten, want elke voorspeller verklaart altijd íéts.

De twee getallen die je als eerste opzoekt: de gecorrigeerde R² (die staat er niet, en zal fors lager liggen) en de F van het hele model met zijn vrijheidsgraden en p (die staat er ook niet — en zonder die regel weet je niet eens of het model als geheel iets deed). Ontbreken ze allebei, dan is de enige eerlijke zin in je literatuurstudie dat het onderzoek een verband rapporteert dat op deze steekproefgrootte niet te onderscheiden is van toeval. En de aanbeveling van de auteurs — ervaren behandelaars koppelen aan complexe casus — is een beleidsadvies op één p-waarde uit een overladen model, in een onderzoek waarin niemand behandelaars heeft toegewezen. Dat is twee stappen te ver: van samenhang naar oorzaak, en van één regel naar beleid.

Uitwerking doe-opgave

De twee modellen naast elkaar: model 1 geeft R² = .192, adjusted R² = .071, F(3, 20) = 1.58, p = .22, met tds_t0 op b = 0.377, SE = 0.194, p = .066. Model 2 geeft R² = .450, gecorrigeerde R² = .097, F(9, 14) = 1.27, p = .33, met tds_t0 op b = 0.563, SE = 0.211, p = .018.

1. Het getal dat je op het verkeerde been zet, is die p = .018 — of, als je alleen naar het model kijkt, de R² van .45. Beide zien er sterker uit dan model 1, terwijl het model als geheel juist minder significant werd. (2) Gecorrigeerde R² stijgt nauwelijks omdat hij precies corrigeert voor wat hier gebeurde: hij trekt van R² af wat je van puur toeval mocht verwachten gegeven het aantal voorspellers. Zes extra dummy’s leveren op papier veel op en na aftrek bijna niets — dat gaatje tussen .45 en .10 ís de overfit, zichtbaar gemaakt. (3) De 14: 24 deelnemers min 9 voorspellers min 1 voor de intercept. Elke voorspeller kost een vrijheidsgraad, en met 14 over is de toets bot geworden — dezelfde data, minder scherpte. Dat is de prijs van voorspellers toevoegen, en die betaal je ook als ze niets verklaren.

◠ Zelfde vorm in het boek

Hier: meerdere voorspellers tegelijk, elke b gelezen als “bij gelijke andere voorspellers” — en de vraag hoeveel voorspellers er op hoeveel mensen los mogen.

In het boek: W8 — Meer dan één verschil.

Daar gaat het van één rechthoekje naar meerdere: eerst wegrekenen wat de ene voorspeller al verklaart, dan kijken wat de volgende over de rest nog te zeggen heeft.

TODO: stem-pas — dit blok is nog niet langs Bens transcript-units gelegd (Bib-query op meervoudige regressie / confounding / overfit). Data zijn expliciet synthetisch (verschilscore_mra.csv, gegenereerd uit hetzelfde cohort als het t-toets-blok); het anker en de bewuste afwijking staan in de kop van genereer_verschilscore_mra.R.