Tweede beweging · Verschillen koppelen
Oefening 3.1 · Correlatie en scatterplot
Samenhang tussen twee variabelen · correlation & scatterplot
In deze oefening
Data bij deze oefening — diamantjes.csv de twaalf steentjes uit het boek · penguins.csv 342 pinguïns · bigfive.csv 2634 ingevulde vragenlijsten · wat elke kolom betekent
prikkel
Waar r vandaan komt — de twaalf steentjes eerst
covariantie · sum of products
Twee kolommen, twaalf steentjes: karaat en glans, netjes naast elkaar. En jij wilt één getal dat zegt hoe die twee samen bewegen. Het halve werk heb je in het boek al gedaan — van elk steentje ken je de afwijking van het gemiddelde, de blauwe streepjes. Alleen heeft elk steentje er hier twéé: eentje op karaat (het gemiddelde is 1.60) en eentje op glans (gemiddelde 50.00).
Dan de zet die je uit het boek kent: vermenigvuldig die twee afwijkingen per steentje met elkaar. Wijzen ze dezelfde kant op — allebei erboven, of allebei eronder — dan komt er plus uit; spreken ze elkaar tegen, dan min.
Doe er twee met de hand, dan zie je het gebeuren. De twaalf steentjes hebben in het boek namen, en twee ervan komen hier goed uit. Anna weegt 1.0 karaat en glanst 20, dus haar afwijkingen zijn −0.6 karaat en −30 glanspunten: allebei onder het gemiddelde, en min keer min is plus. Haar bijdrage is +18. Kees weegt 2.1 karaat en glanst 30 — zwaar, maar dof. Zijn afwijkingen zijn +0.5 en −20, die spreken elkaar tegen, en zijn bijdrage is −10.
| afwijking op karaat | afwijking op glans | die twee maal elkaar | |
|---|---|---|---|
| Anna | −0.6 | −30 | +18 |
| Kees | +0.5 | −20 | −10 |
Bij onze twaalf duwen acht steentjes de som omhoog, drie liggen precies op een gemiddelde en tellen dus nul, en alleen Kees telt tegen. Tel alles op en er staat 54.00 — die som heet de sum of products. Deel door n − 1 = 11, en daar is hij: de covariantie (Sxy), 4.91. Het gemiddelde rechthoekje.
En nu de vraag die je bij élk getal hoort te stellen: 4.91 wát? Reken de eenheid maar gewoon mee — een afwijking in karaat keer een afwijking in glans is… karaat-glanzen. Is 4.91 karaat-glanzen veel samenhang? Weet je niet. Niemand denkt in karaat-glanzen. Ruw is ruk. De covariantie zit nog op de schaal van zijn twee variabelen; je kunt er de ríchting aan aflezen — plus, dus zwaarder gaat samen met glanzender — maar niet de sterkte.
De reparatie is één deling. Deel de covariantie door de twee standaardafwijkingen — 0.3908 karaat en 21.7423 glanspunten — en kijk wat er met de eenheden gebeurt: karaat streept weg tegen karaat, glans tegen glans. Wat overblijft is een kaal getal zonder eenheid, de gestandaardiseerde versie van dezelfde samenhang, en dat getal kan nooit buiten −1 en +1 komen: r = .58. Dát kun je in één keer lezen, bij diamantjes net zo goed als bij pinguïns.
d <- read.csv("diamantjes.csv")
# de ruwe samenhangsmaat: nog op de schaal van karaat en glans
cov(d$karaat, d$glans) # → 4.91
# deel door de twee standaardafwijkingen — de eenheid valt weg
cov(d$karaat, d$glans) / (sd(d$karaat) * sd(d$glans)) # → 0.5777
# en dat is precies wat cor() in één keer doet
cor(d$karaat, d$glans, method = "pearson") # → 0.5777 ("pearson" is de standaard)
# en welke lijn hoort daarbij
model.karaat <- lm(glans ~ karaat, data = d) # → -1.43 + 32.14 · karaat
plot(d$karaat, d$glans)
abline(model.karaat)library(dplyr); library(ggplot2)
d <- read.csv("diamantjes.csv")
d |> summarise(S_xy = cov(karaat, glans),
r = cor(karaat, glans, method = "pearson"))
lm(glans ~ karaat, data = d)
# in ggplot: eerst de assen, dan de lagen erop
ggplot(d, aes(karaat, glans)) +
geom_point() +
geom_smooth(method = "lm", se = FALSE)cor() is dus geen zwarte doos: het is jouw covariantie, gedeeld door jouw twee standaardafwijkingen. R doet alleen het strepen sneller.
Let op het verschil tussen wat je scherm zegt en wat je zelf opschrijft. R drukt 0.5777 af, met nul en al; in de tekst hierboven staat r = .58. Dat is geen slordigheid van een van de twee. Een uitvoerregel citeert je scherm; jouw eigen zin volgt APA, en daar laat je de nul weg omdat het ontbréken ervan iets zegt: dit getal kan niet voorbij de 1 komen. Alle # →-regels in dit boekje staan daarom in schermvorm, zodat je kunt nakijken of jij hetzelfde voor je hebt. Kwadrateer die .58 en je hebt r² = .33 — een derde van de verschillen in glans loopt langs het gewicht, de rest niet.
Twee vragen voor je verder leest, en dan mag je los op iets nieuws.
a) Een karaat is gewoon een gewicht: 200 milligram. Weeg de steentjes dus eens in milligram (d$karaat * 200) en vraag de covariantie opnieuw op. (één getal, en één zin over wat er wél en niet veranderd is aan de twaalf steentjes)
b) Deel die nieuwe covariantie door de twee standaardafwijkingen die er nu bij horen. (één getal, en zeg wat je herkent)
a) 981.82 milligram-glanzen — tweehonderd keer de covariantie van daarnet. (Reken je 200 × 4.91 met de hand, dan krijg je 982.00; R vermenigvuldigt de onafgeronde 4.9091. Zelfde soort scheefte als overal waar je met een gedrukt getal verder rekent.) En er is niets aan de steentjes veranderd: geen enkel gewicht, geen enkele glans, geen enkel punt in de wolk. Alleen de maatstok waarmee je het gewicht opschrijft is anders. Een getal dat tweehonderd keer zo groot wordt zodra je een andere maat kiest, kan niet zeggen hoe sterk een verband is.
cov(d$karaat * 200, d$glans) # → 981.82b) .58 — hetzelfde getal als daarnet. Delen door de standaardafwijkingen haalt de maatstok eruit, en dat is precies waarom de correlatie bestaat.
cov(d$karaat * 200, d$glans) / (sd(d$karaat * 200) * sd(d$glans)) # → 0.5777Nu jij — bij de pinguïns, die je inmiddels kent. Nieuwe vraag aan bekende beesten: hangt de lengte van een snavel samen met zijn diepte?
Haal de pinguïns binnen en gooi ze in een plaatje. Elke pinguïn één stip: snavellengte op de ene as, snaveldiepte op de andere.
peng <- read.csv("penguins.csv")
plot(peng$snavellengte, peng$snaveldiepte) # 342 stippen, een wolk
cor(peng$snavellengte, peng$snaveldiepte, method = "pearson") # → -0.235Daar is-ie: een wolk die van linksboven naar rechtsonder hangt, en het getal −0.235 is één cijfer voor “hoe strak kantelt die wolk omláág?”. De schaal loopt van −1 (een perfecte dalende lijn) via 0 (geen rechtlijnig verband) naar +1 (een perfecte stijgende lijn). −.24 is een flauwe daling met veel rommel eromheen — zwak, maar het staat er. Langere snavel, ondiepere snavel. Zegt het getal.
snap het
Drie soorten in één wolk
confounding · derde speler · omgeklapt teken · r²
De verleiding is nu te schrijven: “bij pinguïns gaat een langere snavel samen met een ondiepere.” Het getal is netjes uitgerekend, de zin lijkt te volgen. En toch is-ie mis. Kleur de wolk eens per soort:
plot(peng$snavellengte, peng$snaveldiepte, col = factor(peng$soort), pch = 16)
# en reken r nu apart uit binnen één soort
gentoo <- peng[peng$soort == "Gentoo", ]
cor(gentoo$snavellengte, gentoo$snaveldiepte, method = "pearson") # → +0.643Twee handelingen in vier regels: één soort uit de tabel snijden, en de wolk per soort kleuren. Het uitsnijden staat bij rijen — filteren en sorteren, het kleuren bij kleuren en vormen — met daar ook de vraag die hier onder tafel blijft: wanneer geef je de soorten kleuren, en wanneer elk hun eigen vakje?
Ineens zie je het. Er zijn drie soorten in die wolk, elk als een eigen kluitje. En binnen elk kluitje loopt de samenhang de andere kant op: langere snavel, juist diepere snavel. Reken maar na — Adélie +.39, Chinstrap +.65, Gentoo +.64. Alle drie positief. Het teken is omgeklapt.
Wat gebeurde er? De Gentoo’s (rechtsonder) hebben lange, ondiepe snavels; de andere twee soorten (linksboven) korte, diepe. Gooi je ze op één hoop, dan trekt dat soortverschil de hele wolk schuin naar beneden — en dáár komt die −.24 vandaan. Niet uit een verband tussen snavellengte en -diepte, maar uit het verschil tússen de soorten. De soort is een derde speler die zich als het verband vermomt. Het boek noemt dat straks confounding; jij zag het hier gebeuren, in de kleuren van een scatterplot.
En dat is de dubbele les van dit blok, mooier dan één van de twee alleen:
- De rekenmachine zei: over deze 342 dieren is r precies −.24. Daar valt niets op af te dingen; reken het maar na.
- De scatterplot zei: klopt — maar niet wát je denkt. De richting is een luchtspiegeling van de soort.
Uitgerekend en waar zijn dus twee verschillende dingen. Een correlatie kan tot achter de komma kloppen en je tóch de verkeerde kant op wijzen. Daarom lees je r nooit zonder de wolk ernaast — en zeker niet zonder je af te vragen wie er nog meer in die wolk zit. Zou je het in één eerlijke zin moeten opschrijven, dan niet “langere snavels zijn ondieper”, maar: over drie soorten heen hangt de wolk omlaag, maar binnen elke soort omhoog — het min-teken is een soorteffect.
(En dat kleine getal r² dat je weleens ziet — hier .06 — is het “aandeel gedeelde variantie”: van alle verschil in snaveldiepte deelt maar een paar procent iets met snavellengte. Ook dat is geen “verklaring”, alleen boekhouding.)
jouw beurt
Van vin naar vragenlijst
schaalscore · dezelfde route, ander onderwerp
Pinguïns hebben één groot voordeel: je kunt ze wegen. Bij mensen meet je zelden iets wat op een weegschaal past, en dan bouw je je variabele zélf uit een vragenlijst.
In bigfive.csv staan 2634 mensen die vijfentwintig uitspraken over zichzelf hebben beoordeeld op een schaal van 1 tot 6. Je gebruikt er hier tien: vijf gaan over neuroticisme — hoe snel raak je van slag — en vijf over openheid, de smaak voor nieuwe dingen en nieuwe ideeën. Eén rij is één mens.
Tien losse antwoorden zijn nog geen variabele. Je maakt er twee getallen van door per persoon het gemiddelde van zijn vijf antwoorden te nemen; dan staat er voor iedereen één neuroticisme-score en één openheid-score, allebei nog steeds op die schaal van 1 tot 6.
En dan loop je precies hetzelfde rondje als hierboven. Dat is de bedoeling: de vinlengte is nu de openheid, en het gewicht is het neuroticisme, maar de route is regel voor regel dezelfde.
b <- read.csv("bigfive.csv")
# per persoon het gemiddelde van zijn vijf antwoorden
b$neuroticisme <- rowMeans(b[, c("N1", "N2", "N3", "N4", "N5")])
b$openheid <- rowMeans(b[, c("O1", "O2", "O3", "O4", "O5")])Kijk of je twee kolommen erbij hebt en of ze binnen de schaal blijven: range(b$openheid) hoort iets tussen 1 en 6 te geven. Staat er een 0 of een 7, dan heb je een verkeerde kolom te pakken.
Vier stappen, en dan de zin waar dit hoofdstuk op uitloopt.
a) Vraag van allebei de schalen het gemiddelde en de standaardafwijking op. (vier getallen, met de schaal erbij waarop ze staan)
b) Reken de covariantie uit tussen openheid en neuroticisme. Doe het daarna nog eens met de rijsóm in plaats van het rijgemiddelde (rowSums). (twee getallen, en zeg hoeveel keer ze schelen en waarom)
c) Deel allebei door de bijbehorende standaardafwijkingen. (twee getallen, en één zin over wat je hier herkent van de pinguïns)
d) Teken de wolk, en kwadrateer de correlatie. (rapporteer r², en één zin over hoe zo’n wolk eruitziet vergeleken met die van vin en gewicht)
e) Schrijf op wat je gevonden hebt, in twee zinnen en in deze vaste vorm:
Binnen deze n mensen geldt: … (r = …). Of dat ook buiten deze n opgaat, weet je nog niet.
(vul de n, de bevinding in gewone taal en de r in — en laat die tweede zin staan zoals hij er staat)
a) Neuroticisme: M = 3.16, SD = 1.19 schaalpunten. Openheid: M = 3.87, SD = 0.55 schaalpunten. Allebei op de schaal van 1 tot 6 waarop de vragen beantwoord zijn — dus de gemiddelde deelnemer zit bij neuroticisme net onder het midden en bij openheid er duidelijk boven.
mean(b$neuroticisme) # → 3.16
sd(b$neuroticisme) # → 1.1946
mean(b$openheid) # → 3.87
sd(b$openheid) # → 0.5455b) Met de rijgemiddelden is de covariantie 0.1050, met de rijsommen 2.6253. Dat scheelt precies een factor 25, en niet toevallig: een som is vijf keer een gemiddelde, en in een covariantie zitten twee variabelen, dus vijf maal vijf. Weer hetzelfde als bij de kilo’s van de pinguïns — je hebt de maat veranderd en niet de mensen.
cov(b$neuroticisme, b$openheid) # → 0.105
N_som <- rowSums(b[, c("N1", "N2", "N3", "N4", "N5")])
O_som <- rowSums(b[, c("O1", "O2", "O3", "O4", "O5")])
cov(N_som, O_som) # → 2.6253c) Allebei .16. Dat is wat je bij de pinguïns ook zag: zodra je door de twee standaardafwijkingen deelt, vallen de eenheden weg en maakt het niet meer uit of je in gemiddelden of in sommen rekende.
cov(b$neuroticisme, b$openheid) / (sd(b$neuroticisme) * sd(b$openheid)) # → 0.1612
N_som <- rowSums(b[, c("N1", "N2", "N3", "N4", "N5")])
O_som <- rowSums(b[, c("O1", "O2", "O3", "O4", "O5")])
cov(N_som, O_som) / (sd(N_som) * sd(O_som)) # → 0.1612d) r² = .03. En de wolk zíet er ook zo uit: een dikke, ronde vlek zonder duidelijke richting, waar je met moeite een flauwe helling in ontdekt. Zet dat naast de vin-en-gewicht-wolk van hiervoor — daar lagen de punten in een strakke sliert. Hier delen de twee maten nog geen drie procent van hun verschillen.
plot(b$openheid, b$neuroticisme)
cor(b$openheid, b$neuroticisme)^2 # → 0.026Dat de punten in blokjes lijken te vallen is geen fout: een gemiddelde van vijf hele getallen kan maar een beperkt aantal waarden aannemen, dus veel mensen landen op precies hetzelfde plekje.
e) Bijvoorbeeld:
Binnen deze 2634 mensen geldt: wie hoger scoort op openheid, scoort ook iets hoger op neuroticisme (r = .16). Of dat ook buiten deze 2634 opgaat, weet je nog niet.
Die tweede zin is geen slag om de arm en geen bescheidenheid. Alles wat je in dit hoofdstuk hebt uitgerekend — de covariantie, de correlatie, de r² — beschrijft de mensen die je vóór je hebt, en verder niemand. Wat er nodig is om die stap naar buiten wél te zetten, is precies waar de volgende oefening mee begint.
◠ Zelfde vorm in het boek
Hier: één getal voor de samenhang, en de wolk die verklapt wat het getal verzwijgt.
In het boek: W4 — Hangen ze samen?
Daar begint het bij de wolk, niet bij het getal — en de rechte lijn komt pas als je hebt gezien wie er allemaal in de wolk zitten.
Voetnoten
Palmer Penguins — Gorman, Williams & Fraser (2014), Palmer Station Antarctica LTER; via het R-pakket
palmerpenguins(CC0). Alle getallen in dit blok zijn op die échte data nagerekend.↩︎