Vierde beweging · Verschil onthechten
Oefening 9.1 · Dummy-codering
Een groep als 0 en 1 · dummy coding
Data bij deze oefening — penguins.csv 342 pinguïns · diamantjes.csv de twaalf steentjes uit het boek · wat elke kolom betekent
prikkel
Een voorspeller waar woorden in staan
groepslabel · voorspeller · regressielijn
Een regressielijn is een rekensom, en verder niets. Je stopt er een getal in en er komt een getal uit. Bij de twaalf steentjes stopte je karaat erin en kreeg je glans terug, en de helling vertelde wat één karaat erbij waard was in glans.
Kijk nu eens naar de kolom soort in het pinguïn-bestand.1 Daar staan geen getallen in. Daar staat Adelie, Chinstrap, Gentoo — drie woorden. Je kunt er niet mee vermenigvuldigen, je kunt er geen helling van maken, en “één soort erbij” betekent niks.
En toch. Stop de soort in een model en er rolt gewoon een tabel met getallen uit. Geen foutmelding, geen gemopper. Er zijn ergens getallen vandaan gekomen die in het bestand niet staan.
Waar komen die vandaan? Trek even mee.
speel het
Maak de dummy’s zelf, en laat R het daarna nog eens overdoen
dummy · k − 1 · lm()
d <- read.csv("diamantjes.csv")
# kolom D staat er al: 1 voor een sterretje, 0 voor een golfje
table(d$D, d$merk) # → de zes golfjes op 0, de zes sterretjes op 1
# en zo'n nul-of-een-kolom mag gewoon de regressie in
lm(glans ~ D, data = d) # → (Intercept) 40 D 20lm(), t.test() en aov() zijn in beide talen precies dezelfde functie — daar valt niets te vertalen. Het verschil tussen base R en dplyr zit in het klaarzetten van de data, niet in het model.
Twee groepen, één kolom met nullen en enen: dat is een dummy. En kijk nu wat de twee getallen zeggen. Het intercept is 40, en dat is precies de gemiddelde glans van de golfjes. De helling is 20, en dat is precies wat er bij de sterretjes bovenop komt: 40 + 20 = 60, hun gemiddelde. Reken het na met tapply(d$glans, d$merk, mean).
Die kolom D hoef je trouwens niet cadeau te krijgen. ifelse(d$merk == "sterretje", 1, 0) maakt hem in één regel na. En je hoeft hem zelfs helemaal niet te máken, want R doet het ook zelf als je hem de woordkolom geeft:
lm(glans ~ merk, data = d) # → (Intercept) 40 merksterretje 20Dezelfde twee getallen. R heeft in stilte precies dezelfde nul-of-een-kolom gebouwd, en hij zet zijn eigen naam erop: merksterretje, oftewel “de dummy die aanstaat bij sterretje”.
Eén ding om nu al te zien, want je komt het in de volgende oefening tegen. Vraag je diezelfde vergelijking als t-toets op — t.test(glans ~ merk, data = d, var.equal = TRUE) — dan staat er t(10) = −1.73. Vraag je ’m als dummy-regressie op, dan staat er t(10) = 1.73, zonder min. Zelfde data, zelfde toets, gespiegeld teken. De t-toets rekent golfje mín sterretje; de dummy rekent sterretje mín golfje. Welke kant het op gaat hangt af van wie je als vertrekpunt neemt — en dáár gaat deze hele oefening over.
En de vraag die je nu wilt stellen: wélke van de twee schrijf je dan op? Allebei mag, zolang je er maar bij zegt welke kant je op gerekend hebt. Een min-teken zonder die zin erbij is niet fout maar onleesbaar; met die zin erbij is het antwoord af.
Nu jij — dezelfde pinguïns die je kent, maar met een vraag die je nog niet aan ze gesteld hebt. Bij de steentjes waren er twee groepen en had je aan één dummy genoeg. Bij de pinguïns zijn het er drie. Bouw ze eerst met de hand, want dan zie je letterlijk wat er straks in de tabel komt te staan.
peng <- read.csv("penguins.csv")
# drie soorten, twee dummy's — de Gentoo's laat je weg
peng$adelie <- as.numeric(peng$soort == "Adelie")
peng$chinstrap <- as.numeric(peng$soort == "Chinstrap")
# kijk of ze staan waar je denkt dat ze staan
table(peng$adelie, peng$soort) # → 151 Adelies op 1, alle andere op 0
model.dummy <- lm(gewicht ~ adelie + chinstrap, data = peng)
round(coef(model.dummy), 2)
# (Intercept) adelie chinstrap
# 5076.02 -1375.35 -1342.93Die table()-regel is geen versiering. Een dummy die per ongeluk overal 0 is geeft geen foutmelding — hij geeft een model dat rustig doorrekent en nergens meer over gaat. Zet er dus altijd de kruistabel naast. Meer over factoren, niveaus en welke groep de referentie wordt staat bij nieuwe variabelen maken en hercoderen.
Nu de vragen. Ze staan los onder elkaar zodat je ze later kunt terugvinden; probeer ze eerst zelf, en klap daarna het antwoord open. Tussen haakjes staat telkens wat er in een compleet antwoord hoort — dan kun je je eigen antwoord nakijken zonder het onze te openen. De laatste heet met opzet Z) en niet l): dat is geen volgende stap maar de vraag die je bij élk model opnieuw stelt, en daarom staat hij achteraan. Z is de laatste letter van het alfabet, dus je ziet in één oogopslag dat er niets meer achter komt — en het blíjft Z, of het rijtje ervoor nu bij d) ophoudt of bij k).
a) Hoeveel dummy’s heb je nodig voor drie soorten, welke soort heb je weggelaten, en hoe heet die groep? (het aantal, de soort, en de term)
b) Wat zijn de nulhypothese en de alternatieve hypothese voor het model als geheel, in woorden én in symbolen — en gaat je uitspraak over deze 342 pinguïns, of over alle pinguïns waar zij een steekproef uit zijn? (twee hele zinnen met populatie of steekproef erin, en die twee ook in symbolen)
c) Kan die nulhypothese verworpen worden? (rapporteer F, df en p, en zeg wát je verwerpt)
d) Hoeveel procent van het verschil in gewicht tussen deze pinguïns verklaart de soort? (rapporteer R², en zeg wat de rest is)
e) Wat zijn de nulhypothese en de alternatieve hypothese voor de twee losse regressiegewichten, weer in woorden én in symbolen — en waarin verschillen die van b)? (twee verschillen, en de hypothesen in symbolen)
f) Zijn de twee regressiegewichten, −1375.35 en −1342.93, significant? (rapporteer per gewicht b, t, df en p)
g) Schrijf de geschatte regressievergelijking op. (met de namen van je eigen dummy’s erin)
h) Welk getal in de uitvoer is het gemiddelde gewicht van de Gentoo’s — en welke schatter gebruik je daarvoor? (het getal, en hoe die schatter heet)
i) Reken de andere twee groepsgemiddelden uit de vergelijking uit. (rekenmachine erbij, het zijn twee sommen)
j) Drie soorten zijn drie paren om te vergelijken, maar je model geeft er maar twee. Welk paar mis je, hoe kom je aan dat verschil, en waarom zit er nog geen toets bij? (een getal in grammen, en één zin over die toets)
k) Geef een inhoudelijke conclusie in niet meer dan vijf zinnen. (over pinguïns, niet over getallen)
Z) Mag je alles hierboven geloven? Loop de aannames na:
- hoe ver lopen de spreidingen per soort uiteen?
- telt elke pinguïn maar één keer mee?
- liggen de residuen ongeveer symmetrisch?
En dan de vraag die er altijd achteraan hoort: gaat er iets mis, welke antwoorden hierboven doe je dan over en welke mogen blijven staan?
a) Twee: altijd eentje minder dan het aantal groepen, k − 1. De Gentoo’s heb je weggelaten, dus Gentoo is de referentiecategorie.
b) Over de populatie, niet over deze 342. Die heb je gewogen; hun gemiddelden kun je gewoon aflezen, dus daar valt niets aan te toetsen.
De nulhypothese: in de populatie pinguïns waar deze 342 een steekproef uit zijn, hangt de soort niet samen met het gewicht — beide regressiegewichten zijn daar nul.
De alternatieve hypothese: minstens één van die twee is niet nul, dus érgens hangt de soort wél met het gewicht samen.
Diezelfde twee zinnen in symbolen:
\[H_0:\ \beta_1 = \beta_2 = 0 \qquad H_1:\ \text{minstens één } \beta_j \neq 0\]
Grieks (β), want het gaat over de populatie; de b’s in je uitvoer zijn de schattingen daarvan uit déze 342. Je mag hem ook over de verklaarde variantie schrijven — \(H_0:\ \rho^2 = 0\) tegen \(H_1:\ \rho^2 > 0\) — want dat zegt hetzelfde. Eén van de twee is genoeg.
c) Ja: F(2, 339) = 343.63, p < .001. Wat je verwerpt is de nulhypothese dat de soort in de populatie niets met het gewicht te maken heeft. Wat je er níét mee weet, is wélke soorten van elkaar verschillen — daarvoor moet je naar de losse gewichten kijken, bij f) en j).
d) R² = .67, dus 67%. De andere 33% is verschil in gewicht dat de soort níét verklaart; dat zit bínnen de soorten, tussen pinguïns van dezelfde soort. Je leest hem af onderaan summary(model.dummy), bij Multiple R-squared.
e) Bij b) ging het over de twee gewichten samen, hier gaat het per gewicht apart. De nulhypothese voor de Adélie-dummy: in de populatie is het verschil in gewicht tussen Adélies en Gentoo’s nul; de alternatieve: dat verschil is niet nul. En hetzelfde nog eens voor de Chinstrap-dummy. Twee verschillen met b) dus: b) toetst ze in één keer en zegt alleen of er ergens iets zit, deze twee toetsen elk één verschil en zeggen ook welk.
In symbolen, nu per gewicht apart:
\[H_0:\ \beta_j = 0 \qquad H_1:\ \beta_j \neq 0\]
En een derde verschil krijg je er cadeau, want je ziet het aan de notatie: deze staat twee kanten op. Een los verschil kan ook de andere kant uit vallen, terwijl de F van b) alleen maar groter dan nul kan uitpakken.
f) Ja, allebei. Adélies tegenover Gentoo’s: b = −1375.35 g, t(339) = −24.50, p < .001. Chinstraps tegenover Gentoo’s: b = −1342.93 g, t(339) = −19.22, p < .001. Je vindt ze met round(summary(model.dummy)$coefficients, 2).
g) gewicht = 5076.02 − 1375.35 × adelie − 1342.93 × chinstrap. De twee minnen zitten in de gewichten zelf, dus je mag ze ook als plus opschrijven met een negatief getal erachter — als de uitkomst maar hetzelfde is.
h) 5076.02, en de schatter die je daarvoor gebruikt is het intercept. Vul in de vergelijking bij adelie en chinstrap allebei een 0 in en er blijft niets over dan dat intercept — en een pinguïn die op allebei de dummy’s nul scoort ís een Gentoo.
i) Adélie = 5076.02 − 1375.35 = 3700.67; Chinstrap = 5076.02 − 1342.93 = 3733.09. Kijk ze na:
round(tapply(peng$gewicht, peng$soort, mean), 2)
# Adelie Chinstrap Gentoo
# 3700.66 3733.09 5076.02Jouw Adélie-som landt dus één honderdste te hoog. Dat is geen rekenfout van jou: de exacte waarden zijn 5076.016 en 1375.354, en R rekent met álle decimalen door. Rond af aan het eind, nooit onderweg.
j) Je mist Chinstrap tegenover Adélie. Je model meet álles af aan de Gentoo’s, dus die twee komen elkaar nergens tegen.
Het verschil krijg je door de gewichten van elkaar af te trekken: −1342.93 − (−1375.35) = 32.42 gram, de Chinstraps net iets zwaarder.
Maar een verschil is nog geen toets. Daarvoor heb je de standaardfout van dít verschil nodig, en die staat nergens in deze tabel — je krijgt hem pas als je de Adélies tot referentie maakt. Wil je alle drie de paren tegelijk, dan is daar een eigen gereedschap voor: de post-hoc-toets, in de oefening over de eenweg ANOVA.
k) Bijvoorbeeld: welke soort een pinguïn is, hangt sterk samen met wat hij weegt — de soort verklaart 67% van het verschil in gewicht tussen deze 342 pinguïns, F(2, 339) = 343.63, p < .001. Dat verschil zit vrijwel helemaal in de Gentoo’s: die wegen gemiddeld ruim 1300 gram meer dan de Adélies en ruim 1300 gram meer dan de Chinstraps. Adélies en Chinstraps ontlopen elkaar maar 32 gram. Een pinguïn wegen zegt dus veel over welke soort het is, zolang je alleen de Gentoo van de andere twee hoeft te onderscheiden.
Z) Ja, je mag ze geloven. De drie soorten hebben spreidingen van 458.57, 384.34 en 504.12 gram — grootste gedeeld door kleinste is 1.31, ruim onder de 2, dus gelijke spreiding is hier geen probleem. Elke pinguïn is één keer gewogen, dus de waarnemingen zijn onafhankelijk, en de residuen zijn vrijwel symmetrisch (scheefheid 0.18; dat begrip kwam langs in de oefening over de mediaan en de scheefheid).
Vuistregel voor die spreidingen: pas als de grootste meer dan twee keer de kleinste is, gaat gelijke spreiding echt knellen.
En als het antwoord nee was geweest? Dan gaat niet alles hierboven overboord — en dát is waarom deze vraag achteraan mág staan.
- Blijft staan: de gewichten bij f), g) en i). Dat zijn gewoon de verschillen tussen de groepsgemiddelden, en daar verandert geen aanname iets aan.
- Doe je over: de toets bij c), de p’s bij f), en je conclusie bij k). Die rusten alle drie op één gedeelde spreiding.
Let op dat f) in allebei de rijtjes staat: de b van −1375.35 blijft, de p ernaast niet. Eén geval om het te zien — zouden de Gentoo’s twee keer zo wild spreiden als de Chinstraps, dan is die gedeelde spreiding te groot voor de Chinstraps en te klein voor de Gentoo’s, en staat er bij f) een t die niet klopt, terwijl die 1375.35 gram er nog steeds gewoon staat.
Laat R het nu nog eens overdoen, met de woordkolom:
model.soort <- lm(gewicht ~ soort, data = peng)
round(coef(model.soort), 2)
# (Intercept) soortChinstrap soortGentoo
# 3700.66 32.43 1375.35Andere getallen — en niet zomaar andere. Vijf nieuwe vragen, nu op dit model:
a) Welke soort is nu de referentiecategorie, en waaraan zie je dat? (twee aanwijzingen)
b) Bereken de drie gemiddelden opnieuw uit deze vergelijking — kom je op dezelfde drie uit? (rekenmachine erbij, weer twee sommen)
c) Hoeveel procent van het verschil in gewicht tussen deze pinguïns verklaart de soort nu? (rapporteer R²)
d) Eén getal staat in allebei de modellen en één is verdwenen — welk getal is dat, en waar is het gebleven?
e) En het verschil dat je bij j) met de hand uitrekende: staat het hier, en zit er nu wél een toets bij? (rapporteer b, t, df en p)
a) Adélie. Het intercept is 3700.66 en dat is precies het Adélie-gemiddelde; bovendien komt Adélie in geen enkele coëfficiëntnaam voor.
b) Adélie = 3700.66; Chinstrap = 3700.66 + 32.43 = 3733.09; Gentoo = 3700.66 + 1375.35 = 5076.01 — terwijl het echte Gentoo-gemiddelde 5076.02 is. Dezelfde drie gemiddelden als in model 1, dezelfde afrondingskruimel, andere rij.
c) R² = .67 — precies hetzelfde. Tot op vier decimalen gelijk: .6697 in allebei.
d) De 1375.35 staat in allebei, daar met een min ervoor en hier zonder. De −1342.93 is weg en er staat 32.43 voor in de plaats, want dat verschil wordt nu afgemeten aan de Adélies in plaats van aan de Gentoo’s: 3733.09 − 3700.66 = 32.43.
e) Ja, en dit is het paar dat model 1 je niet gaf: b = 32.43 g, t(339) = 0.48, p = .63. Niet significant — met deze data kun je Chinstraps en Adélies op gewicht niet uit elkaar houden. Merk op dat jouw aftreksom bij j) op 32.42 uitkwam en hier 32.43 staat: dat is dezelfde afrondingskruimel als bij b), want je trok twee getallen van elkaar af die al afgerond waren.
En Z)? Die hoeft hier niet. Dit is hetzelfde model als hierboven, alleen vanaf een ander vertrekpunt, dus het draagt precies dezelfde aannames — en die heb je daar al nagelopen.
Twee modellen, dezelfde drie gemiddelden. Het zijn er dus geen twee — het is hetzelfde model, met een ander vertrekpunt.
Je hebt net zelf uitgerekend dat je handgebouwde model en dat van R op dezelfde drie gemiddelden landen. Je kunt het R ook rechtstreeks vragen:
all.equal(fitted(model.dummy), fitted(model.soort)) # → TRUE
deviance(model.dummy) # → 72443483
deviance(model.soort) # → 72443483
summary(model.dummy)$r.squared # → 0.669672
summary(model.soort)$r.squared # → 0.669672Drie vragen, en ze geven alle drie hetzelfde antwoord. fitted() haalt op wat een model voor elke pinguïn voorspelt, en all.equal() legt die twee rijen naast elkaar: TRUE. Elke pinguïn krijgt in allebei de modellen exact hetzelfde voorspelde gewicht. deviance() geeft de restsom — al het verschil in gewicht dat het model níét weet te voorspellen, bij elkaar opgeteld — en die staat allebei de keren op 72443483, tot op de gram gelijk. En het deel van de verschillen in gewicht dat het model wél te pakken krijgt, de R², is in allebei .67.
Twee modellen die elke pinguïn hetzelfde voorspellen, evenveel laten liggen en even veel van de verschillen in gewicht verklaren, zijn niet twee modellen.
anova() staat
Bij de oefening over meervoudige regressie liet je anova() twee modellen tegen elkaar afwegen, en daar deed hij echt werk: het tweede model was het eerste plus de soort. Zulke modellen heten genest — het kleine zit helemaal in het grote, en wat de toets je vertelt is precies wat het grote model erbij kreeg.
Deze twee hier zijn niet genest. Geen van beide zit in de ander, want het is één model dat twee keer is opgeschreven. Doe je het tóch — en de kans is groot dat je het een keer per ongeluk doet, want je hébt die aanroep net geleerd — dan krijg je dit te zien:
> anova(model.dummy, model.soort)
Res.Df RSS Df Sum of Sq F Pr(>F)
1 339 72443483
2 339 72443483 0 2.9802e-08Dit hoef je niet na te doen; het staat er zodat je het herkent. RSS is de restsom die je hierboven met deviance() ophaalde — dezelfde 72443483, twee keer. Kijk nu naar de kolom Df: daar staat een 0, en die telt hoeveel voorspellers het tweede model méér gebruikt dan het eerste. Er kwam er geen bij, dus er valt niets te toetsen — en de kolommen F en Pr(>F), waar anders je toets en je p-waarde staan, blijven dan ook leeg. Wat er bij Sum of Sq staat is machineruis: het kruimeltje dat een computer overhoudt als hij met kommagetallen rekent. Op jouw computer ligt er misschien een andere kruimel. Het is geen uitkomst, en je kunt er niets uit aflezen.
Vandaar de vuistregel voor de rest van dit oefenboek: anova() met twee modellen erin is er voor geneste modellen. Wil je weten of twee schrijfwijzen hetzelfde model beschrijven, kijk dan naar wat ze voorspellen — precies wat je hierboven deed.
En dan de derde gedaante. Draai dezelfde vraag als eenweg ANOVA:
model.anova <- aov(gewicht ~ soort, data = peng)
summary(model.anova)
# soort Df 2 F value 343.6 Pr(>F) <2e-16
summary(model.dummy)$fstatistic # jouw eigen dummy's
# value numdf dendf
# 343.6263 2.0000 339.0000
summary(model.soort)$fstatistic # die van R
# value numdf dendf
# 343.6263 2.0000 339.0000Drie keer hetzelfde getal. De ANOVA-tabel drukt het afgerond af als 343.6, $fstatistic geeft het voluit als 343.6263, en netjes op twee decimalen is dat F(2, 339) = 343.63 — één getal, drie schrijfwijzen. Wat die F precies méét, daar gaat de ANOVA-oefening straks over; kijk er nu alleen naar dat je hem drie keer tegenkomt terwijl je drie verschillende dingen leek te doen.
En let op wat er in de kop van dat kader staat: eenweg ANOVA, en niet alleen “ANOVA”. Soort op gewicht doe je niet met ANOVA, maar met eenweg ANOVA. Ik zeg altijd: alles is regressie, en alles is ook analysis of variance — we doen niet anders.
En dat is geen muggenzifterij. Als álles variantieanalyse is, zegt het woord “ANOVA” niets meer over wát je gedaan hebt; het aantal factoren is dan het enige dat je lezer nog iets vertelt. Hier heb je er één — soort — en dus is het een eenweg ANOVA. Bij twee factoren wordt het een tweeweg ANOVA, en dat is een ander gesprek.
snap het
Wie je weglaat, bepaalt wat alle andere getallen betekenen
referentiecategorie · k − 1 dummy’s · intercept
Eerst het tellen. Hoeveel dummy’s heb je nodig voor een groepsvariabele? Eentje minder dan het aantal groepen. Man en vrouw zijn twee groepen en dus één dummy; drie pinguïnsoorten zijn twee dummy’s. Statistici schrijven k − 1, met k voor het aantal groepen.
Waarom niet gewoon drie? Omdat de derde niets nieuws vertelt. Wie op adelie een 0 heeft staan en op chinstrap ook, is een Gentoo — dat weet je al, zonder ernaar te kijken. De weggelaten groep is de referentiecategorie, en die herken je hieraan: hij scoort nul op álle dummy’s.
En dan verdient hij dat woord ook echt, want alles wat er in de tabel staat wordt vanaf hem gemeten:
- Het intercept is het gemiddelde van de referentiegroep. Vul in de vergelijking voor alle dummy’s een 0 in en er blijft niets over dan b0.
- Elk regressiegewicht is het verschil met de referentiegroep, in de eenheid van de uitkomst. Hier dus in grammen.
Daarmee wordt het verschil tussen je twee tabellen doodgewoon:
| soort | gemiddeld gewicht | Gentoo weggelaten | Adélie weggelaten |
|---|---|---|---|
| Adélie | 3700.66 g | −1375.35 | 3700.66 |
| Chinstrap | 3733.09 g | −1342.93 | 32.43 |
| Gentoo | 5076.02 g | 5076.02 | 1375.35 |
Lees de twee rechterkolommen van boven naar beneden. Het intercept verspringt van de onderste rij naar de bovenste, en alle andere getallen kantelen mee. Er is geen enkel feit over pinguïns veranderd; de drie gemiddelden staan er nog gewoon. Wat veranderde is vanaf wélk gemiddelde je de rest afmeet.
Dus kies dat vertrekpunt met opzet, in plaats van het te laten gebeuren. R kiest namelijk zelf als jij het niet doet, en hij kiest op alfabet — vandaar dat de Adélies in lm(gewicht ~ soort) het intercept kregen. Wil je een andere:
peng$soort <- relevel(factor(peng$soort), ref = "Gentoo")
model.gentoo.ref <- lm(gewicht ~ soort, data = peng)
round(coef(model.gentoo.ref), 2)
# (Intercept) soortAdelie soortChinstrap
# 5076.02 -1375.35 -1342.93Precies je handwerk van daarnet, in één regel. Wat je met de hand deed en wat relevel() doet is dezelfde ingreep: iemand moet het vertrekpunt zijn, en jij zegt wie.
Dat handwerk is trouwens geen achterhaalde omweg. In SPSS moet je de dummy’s zelf bouwen voor je ze in een regressie kunt stoppen; daar is dit gewoon de weg. In R hoeft het niet, maar het is nog steeds de eerlijkste manier om te zien wat er gebeurt. Allebei goed — als je maar weet wie je weggelaten hebt.
Welke groep hoort het te zijn? Drie vuistregels, in die volgorde:
- Is er een controlegroep, dan is dat ’m. Bij een controlegroep en twee behandelingen lees je de twee gewichten meteen als “wat de behandeling oplevert ten opzichte van niets doen”.
- Anders de grootste groep, zodat het intercept op de meeste waarnemingen rust.
- Nooit de meest bijzondere groep, want dan meet je alles af aan de uitzondering.
Bij deze pinguïns is er geen controlegroep. De Gentoo is de uitschieter en de Adélie de gewone maat, en de Adélies zijn er bovendien de meeste: 151 van de 342. Vandaar dat wij in de rest van dit hoofdstuk de Adélies als referentie aanhouden.
En dan de aansluiting, want daar is dit hoofdstuk om begonnen. Kijk terug naar wat je gedaan hebt: je hebt drie groepen vergeleken, en er is geen nieuwe toets aan te pas gekomen. Je gebruikte de regressie die je al kende, met nullen en enen in plaats van karaat.
Dat is geen handigheidje dat toevallig ook werkt. lm(), aov() en straks de t-toets zijn één model met drie deuren, en twee van je getallen komen zo weer terug:
- 1375.35 gram tussen Gentoo’s en Adélies — dat is het verschil dat de t-toets in de volgende oefening onderzoekt.
- F(2, 339) = 343.63 — die kom je in de ANOVA-oefening weer tegen.
Dezelfde getallen, andere jas.
De les van dit blok, in drie zinnen:
- Een groep wordt een getal door hem in nullen en enen te hakken. Voor k groepen heb je k − 1 dummy’s nodig; de groep die je weglaat is de referentiecategorie.
- Het intercept is het gemiddelde van de referentiegroep, elk gewicht is een verschil daarmee. Ruil je de referentie om, dan kantelen alle getallen mee zonder dat er één feit verandert.
- Groepen vergelijken is regressie. De t-toets en de eenweg ANOVA zijn geen aparte machines; het is dezelfde som, anders opgeschreven.
Zou je dit model in één eerlijke zin opschrijven, dan hóórt de referentiecategorie erin — laat je ’m weg, dan weet je lezer niet waar de getallen vandaan geteld zijn. De standaardfouten (SE) staan in dezelfde tabel als de t’s, die je bij vraag f) al opvroeg met round(summary(model.soort)$coefficients, 2):
Om te onderzoeken hoe het gewicht van pinguïns samenhangt met hun soort, voorspelden we hun gewicht uit twee dummyvariabelen (N = 342), met de Adélies als referentiecategorie. Samen verklaarden de twee dummy’s 67% van het verschil in gewicht tussen deze pinguïns, R² = .67, F(2, 339) = 343.63, p < .001. Gentoo’s wogen gemiddeld meer dan Adélies (b = 1375.35 g, 95% CI [1264.91, 1485.80], SE = 56.15, t(339) = 24.50, p < .001). Tussen Chinstraps en Adélies vonden we geen aanwijzing voor een verschil (b = 32.43 g, 95% CI [−100.37, 165.22], SE = 67.51, t(339) = 0.48, p = .63).
Zo streng is de zin met opzet. Alleen de symbolen zijn cursief (N, b, SE, t, p, R², F), de gewone woorden staan rechtop, en de getallen staan op twee decimalen. Alles staat in de verleden tijd, want je vertelt wat je gevónden hebt. En de zin staat op eigen benen: wie hem uit dit blok tilt, weet nog steeds wie er gewogen zijn, hoeveel het er waren, en vanaf wie er geteld is. Elk getal dat een lezer nodig heeft om je na te rekenen, staat erin — dat is de lat.
Waarom er staat wat er staat. Dit stuk hoort bij de les en niet bij de zin: wie de zin overschrijft, laat dit hier liggen. Maar je moet het één keer gezien hebben, want anders schrijf je die zin over zonder te weten waarom hij zo lang is.
- De referentiecategorie staat erin omdat je lezer anders niet weet vanaf waar de 1375.35 geteld is. Met een andere soort als vertrekpunt was precies hetzelfde getal −1375.35 geweest.
- N = 342 staat erin omdat p < .001 zonder steekproefgrootte niets zegt over hoeveel bewijs er ligt.
- De F staat erin naast de R² omdat die twee een paar zijn: R² zegt hoevéél de soort verklaart, de F zegt of dat meer is dan toeval zou opleveren. Een percentage zonder toets is een getal zonder maat.
- R² staat naast de significantie omdat significant “vrij zeker” betekent, en niet “groot” of “belangrijk”. De p zegt hoe zeker je bent dát de soort met het gewicht samenhangt; R² zegt hoeveel van het verschil in gewicht hij verklaart.
- En de Chinstrap blijft erin, juist omdat hij niet significant is. Laat je hem weg, dan rapporteer je alleen wat goed uitkwam — en ziet je lezer niet dat de soort het gewicht goed voorspelt voor de één en nauwelijks voor de ander.
En let op wat er níét staat. Er staat niet dat Chinstraps en Adélies nauwelijks verschillen, want dat weten we niet. Kijk naar het interval: 95% CI [−100.37, 165.22]. Daar past een verschil van 165 gram in, en dat is bijna 4% van wat een pinguïn gemiddeld weegt. Niet-significant betekent hier “we hebben geen verschil kunnen aantonen”, niet “er is er geen”. Dat verschil tussen die twee zinnen is het interval, en daarom staat het erbij.
jouw beurt
Voetnoten
Palmer Penguins — Gorman, Williams & Fraser (2014), Palmer Station Antarctica LTER; via het R-pakket
palmerpenguins(CC0).penguins.csv, 342 pinguïns. Alle getallen op die échte data nagerekend.↩︎