Laten we eerst wat data pakken. Het pakket MASS bevat gegevens over 93 auto's die in 1993 in de VS te koop waren. Ze zijn opgeslagen in het object Cars93 en bevatten 27 kenmerken per auto, waarvan sommige categorisch zijn. Laten we dus het pakket MASS laden en kijken naar het type voertuigen in cars93:
library(MASS)
Cars93$Type
## [1] Small Midsize Compact Midsize Midsize Midsize Large Large ## [9] Midsize Large Midsize Compact Compact Sporty Midsize Van ## [17] Van Large Sporty Large Compact Large Small Small ## [25] Compact Van Midsize Sporty Small Large Small Small ## [33] Compact Sporty Sporty Van Midsize Large Small Sporty## [41] Sporty Small Compact Small Small Sporty Midsize Midsize## [49] Midsize Midsize Midsize Large Small Small Compact Van ## [57] Sporty Compact Midsize Sporty Midsize Small Midsize Small ## [65] Compact Van Midsize Compact Midsize Van Large Sporty## [73] Small Compact Sporty Midsize Large Compact Small Small ## [81] Small Compact Small Small Sporty Midsize Van Small ## [89] Van Compact Sporty Compact Midsize## Levels: Compact Large Midsize Small Sporty VanMogelijk gemaakt door
We hebben daar 6 autotypen. De functie table vertelt hoeveel we er van elk type hebben:
table(Cars93$Type)#### Compact Large Midsize Small Sporty Van## 16 11 22 21 14 9Mogelijk gemaakt door
prop.table zet dit om naar fracties:
prop.table(table(Cars93$Type))#### Compact Large Midsize Small Sporty Van## 0.17204301 0.11827957 0.23655914 0.22580645 0.15053763 0.09677419Mogelijk gemaakt door
Hetzelfde met de herkomst van de auto's:
table(Cars93$Origin)#### USA non-USA## 48 45
prop.table(table(Cars93$Origin))#### USA non-USA## 0.516129 0.483871Mogelijk gemaakt door
Hoe maak je een kruistabel
Mooi, we zagen dat onze dataset een vergelijkbaar aantal Amerikaanse en niet-Amerikaanse auto's bevat en dat de meest voorkomende typen Midsize en Small zijn. Maar misschien verschillen de Amerikaanse en niet-Amerikaanse auto's in type?
Laten we naar de autotypen kijken in relatie tot hun herkomst. We kunnen table opnieuw gebruiken, maar nu met twee argumenten. Het eerste wordt de rijvariabele en het tweede de kolomvariabele:
table(Cars93$Type, Cars93$Origin)## ## USA non-USA## Compact 7 9## Large 11 0## Midsize 10 12## Small 7 14## Sporty 8 6## Van 5 4Mogelijk gemaakt door
Nu zagen we wat iedereen weet: Amerikanen houden van grote voertuigen! De bovenstaande tabel laat de gezamenlijke verdeling zien van twee categorische variabelen (Type en Origin). Zulke tabellen noemen we kruistabellen.
Hoe krijg je marginalen uit een kruistabel
De functies rowSums en colSums spreken eigenlijk voor zich
(tab1<-table(Cars93$Type, Cars93$Origin))## ## USA non-USA## Compact 7 9## Large 11 0## Midsize 10 12## Small 7 14## Sporty 8 6## Van 5 4
rowSums(tab1)## Compact Large Midsize Small Sporty Van## 16 11 22 21 14 9
colSums(tab1)## USA non-USA## 48 45Mogelijk gemaakt door
Hoe krijg je percentages uit een kruistabel
prop.table genest met table geeft frequenties
prop.table(table(Cars93$Type, Cars93$Origin))## ## USA non-USA## Compact 0.07526882 0.09677419## Large 0.11827957 0.00000000## Midsize 0.10752688 0.12903226## Small 0.07526882 0.15053763## Sporty 0.08602151 0.06451613## Van 0.05376344 0.04301075Mogelijk gemaakt door
Omzetten naar percentages is simpelweg vermenigvuldigen met 100:
prop.table(table(Cars93$Type, Cars93$Origin))*100## ## USA non-USA## Compact 7.526882 9.677419## Large 11.827957 0.000000## Midsize 10.752688 12.903226## Small 7.526882 15.053763## Sporty 8.602151 6.451613## Van 5.376344 4.301075Mogelijk gemaakt door
Let op: dit is een gezamenlijke kansverdeling, waaruit we bijvoorbeeld zien dat zo'n 7,5% van de auto's klein is en van Amerikaanse herkomst.
Vaker zijn we geïnteresseerd in de verdeling van één variabele binnen groepen die door een andere variabele zijn gevormd. Hier lijkt de verdeling van autotypen onder de Amerikaanse en (afzonderlijk) niet-Amerikaanse auto's interessant. Hiervoor gebruiken we het argument margin van de functie prop.table. Het geeft aan waar de groeperingsvariabele staat: in rijen (margin=1) of in kolommen (margin=2):
prop.table(table(Cars93$Type, Cars93$Origin), margin=2)*100## ## USA non-USA## Compact 14.583333 20.000000## Large 22.916667 0.000000## Midsize 20.833333 26.666667## Small 14.583333 31.111111## Sporty 16.666667 13.333333## Van 10.416667 8.888889Mogelijk gemaakt door
Nu zien we makkelijk dat kleine auto's twee keer zo vaak voorkomen in het niet-Amerikaanse deel van onze dataset als in het Amerikaanse.
Let ook op dat kolommen tot 100 optellen, terwijl in de gezamenlijke verdelingstabel (degene zonder het argument margin) 100 de som van de hele tabel was.
(tab2<-prop.table(table(Cars93$Type, Cars93$Origin), margin=2)*100)## ## USA non-USA## Compact 14.583333 20.000000## Large 22.916667 0.000000## Midsize 20.833333 26.666667## Small 14.583333 31.111111## Sporty 16.666667 13.333333## Van 10.416667 8.888889
colSums(tab2)## USA non-USA## 100 100Mogelijk gemaakt door
Chi-kwadraattoets
De meest voorkomende vraag bij kruistabellen is of de rij- en kolomvariabelen onafhankelijk zijn. De meest basale manier om dit te beantwoorden is een chi-kwadraattoets uitvoeren. Die wordt uitgebreid behandeld in deze tutorial. Laten we controleren of Type en Origin onafhankelijk zijn:
chisq.test(Cars93$Type, Cars93$Origin)## Warning in chisq.test(Cars93$Type, Cars93$Origin): Chi-squared## approximation may be incorrect#### Pearson's Chi-squared test#### data: Cars93$Type and Cars93$Origin## X-squared = 14.08, df = 5, p-value = 0.01511Mogelijk gemaakt door
Blijkbaar niet, maar we kregen ook de waarschuwing Chi-squared approximation may be incorrect. Dit komt doordat de chi-kwadraatstatistiek slechts bij benadering een chi-kwadraatverdeling volgt. Hoe meer observaties we hebben, hoe beter de benadering. De functie chisq.test geeft bovenstaande waarschuwing zodra een van de verwachte aantallen lager is dan 5 (voor wat een "expected count" is, zie de gelinkte tutorial).
Exacte toets van Fisher
De exacte toets van Fisher is een alternatief voor de chi-kwadraattoets, dat vooral wordt gebruikt wanneer de chi-kwadraatbenadering niet toereikend is. Laten we die uitvoeren:
fisher.test(Cars93$Type, Cars93$Origin)#### Fisher's Exact Test for Count Data#### data: Cars93$Type and Cars93$Origin## p-value = 0.007248## alternative hypothesis: two.sidedMogelijk gemaakt door
De resultaten lijken behoorlijk op die van de chi-kwadraat, maar dat is niet altijd zo. Een belangrijk nadeel van de exacte toets van Fisher is dat deze voor grote tabellen (of grote steekproeven) computationeel inefficiënt wordt.
G-toets
Een ander alternatief is de zogeheten G-toets. De statistiek volgt ook bij benadering een chi-kwadraatverdeling, maar voor kleine steekproeven ligt deze benadering dichterbij dan die van de chi-kwadraattoets. Voor de G-toets kunnen we de functie GTest uit het pakket DescTools gebruiken. De resultaten zijn opnieuw vergelijkbaar met de twee vorige toetsen: Type en Origin zijn niet onafhankelijk.
library(DescTools)
GTest(Cars93$Type, Cars93$Origin)#### Log likelihood ratio (G-test) test of independence without## correction#### data: Cars93$Type and Cars93$Origin## G = 18.362, X-squared df = 5, p-value = 0.002526Mogelijk gemaakt door
Correctie van Yates
In 2x2-kruistabellen kan de chi-kwadraattoets worden verbeterd met de continuïteitscorrectie van Yates. Deze trekt simpelweg 0,5 af van elke term | Observed - Expected | in de chi-kwadraatstatistiek. Raadpleeg opnieuw deze tutorial als je nu de draad kwijt bent. Bovendien past R de correctie van Yates automatisch toe wanneer dat nodig is. Laten we kijken naar de beschikbaarheid van uitvoeringen met handgeschakelde transmissie bij Amerikaanse en niet-Amerikaanse auto's:
(tab3<-table(Cars93$Man.trans.avail, Cars93$Origin))## ## USA non-USA## No 26 6## Yes 22 39
chisq.test(tab3)#### Pearson's Chi-squared test with Yates' continuity correction#### data: tab3## X-squared = 15.397, df = 1, p-value = 8.712e-05Mogelijk gemaakt door
De correctie van Yates werd automatisch toegepast, en R gaf dat aan.
3 (of meer) dimensionale tabel
Laten we tegelijk kijken naar Man.trans.avail, Origin en Type. table splitst de dataset op volgens de als derde opgegeven variabele:
table(Cars93$Man.trans.avail, Cars93$Origin, Cars93$Type)## , , = Compact#### ## USA non-USA## No 2 0## Yes 5 9#### , , = Large#### ## USA non-USA## No 11 0## Yes 0 0#### , , = Midsize#### ## USA non-USA## No 9 4## Yes 1 8#### , , = Small#### ## USA non-USA## No 0 0## Yes 7 14#### , , = Sporty#### ## USA non-USA## No 0 0## Yes 8 6#### , , = Van#### ## USA non-USA## No 4 2## Yes 1 2Mogelijk gemaakt door
ftable geeft een compacter overzicht:
ftable(Cars93$Man.trans.avail, Cars93$Origin, Cars93$Type)## Compact Large Midsize Small Sporty Van## ## No USA 2 11 9 0 0 4## non-USA 0 0 4 0 0 2## Yes USA 5 0 1 7 8 1## non-USA 9 0 8 14 6 2Mogelijk gemaakt door
Cochran-Mantel-Haenszel-toets
Uit de bovenstaande table-resultaten zou duidelijk moeten zijn dat de relatie tussen Origin en Man.trans.avail verschilt per Type. Voor kleine en sportieve auto's is er bijvoorbeeld geen verband: elke Amerikaanse, evenals elke niet-Amerikaanse kleine of sportieve auto, is verkrijgbaar in een handgeschakelde versie. Aan de andere kant hebben de meeste middelgrote Amerikaanse auto's geen handgeschakelde versie, terwijl de meeste niet-Amerikaanse auto's van dit type die wel hebben. Om rekening te houden met mogelijk verschillende relaties per laag (stratum), kan de Cochran-Mantel-Haenszel-toets worden gebruikt:
mantelhaen.test(Cars93$Man.trans.avail, Cars93$Origin, Cars93$Type)#### Mantel-Haenszel chi-squared test with continuity correction#### data: Cars93$Man.trans.avail and Cars93$Origin and Cars93$Type## Mantel-Haenszel X-squared = 8.0153, df = 1, p-value = 0.004638## alternative hypothesis: true common odds ratio is not equal to 1## 95 percent confidence interval:## 2.226531 76.891307## sample estimates:## common odds ratio## 13.08438Mogelijk gemaakt door
Het derde argument van mantelhaen.test identificeert de strata. Vergelijk bovenstaande resultaten met die zonder stratificatie (het voorbeeld met de correctie van Yates hierboven). Er is nog steeds een verband, maar het bewijs ervoor is een stuk zwakker.
Associatiematen
Nu we enkele verbanden tussen variabelen hebben ontdekt, is het tijd om hun sterkte te meten. Er zijn tal van mogelijke maten. Veel daarvan worden hier beschreven. Ik focus nu op de twee meest gebruikte.
Cramers V
V is gebaseerd op de chi-kwadraatstatistiek:
$$ V = \sqrt{\frac{\chi^2/N}{\min(C-1, R-1)}}, $$ waarbij:
N de totale som in de kruistabel is (som van alle cellen),
C het aantal kolommen is
R het aantal rijen is.
V ∈ [0; 1]. Hoe groter V, hoe sterker de relatie tussen variabelen. V = 0 kan worden geïnterpreteerd als onafhankelijkheid (omdat V = 0 dan en slechts dan wanneer χ2 = 0). Het belangrijkste nadeel van V is het gebrek aan precieze interpretatie. Is V = 0,6 een sterke, middelmatige of zwakke associatie?
De functie CramerV uit DescTools kan dit voor ons berekenen:
Nogmaals: is dit een sterke, middelmatige of zwakke associatie?
Goodman en Kruskal's lambda
Goodman en Kruskal's lambda is een voorbeeld van een maat gebaseerd op de proportionele reductie in variatie. Dit soort maten probeert R2 — de determinatiecoëfficiënt uit de lineaire regressie — na te bootsen:
ze nemen waarden uit [0, 1]
ze drukken een fractie uit van de variatie die door de onafhankelijke variabele wordt verklaard.
Laten we een kolomvariabele gebruiken als de onafhankelijke. De formule van Goodman en Kruskal's lambda is dan: $$\lambda = \frac{L - \sum_j L_j}{L},$$ waarbij:
Lj de som is van niet-modale frequenties in de j-de kolom,
L de som is van niet-modale frequenties in de "totaal"-kolom
Dit is het best uit te leggen met een voorbeeld. Laten we kijken naar Type en Origin, waarbij de laatste de onafhankelijke is:
table(Cars93$Type, Cars93$Origin)## ## USA non-USA## Compact 7 9## Large 11 0## Midsize 10 12## Small 7 14## Sporty 8 6## Van 5 4Mogelijk gemaakt door
De modus van de VS-kolom is 11, dus L1 = 7 + 10 + 7 + 8 + 5 = 37. De modus van de niet-VS-kolom is 14, dus L2 = 9 + 0 + 12 + 6 + 4 = 31. De "totaal"-kolom is
rowSums(table(Cars93$Type, Cars93$Origin))## Compact Large Midsize Small Sporty Van## 16 11 22 21 14 9Mogelijk gemaakt door
De modus is 22, dus L = 16 + 11 + 21 + 14 + 9 = 71 en $$\lambda = \frac{71 - (37+31)}{71}=0.042$$
We kunnen in plaats daarvan de functie Lambda uit het pakket DescTools gebruiken:
Lambda(Cars93$Type, Cars93$Origin, direction='row')## [1] 0.04225352Mogelijk gemaakt door
De parameter direction bepaalt waar (in row of in column) de afhankelijke variabele staat.
Origin verklaart slechts zo'n 4% van de variatie in Type. Let op dat lambda variatie definieert als een dichotomie tussen wel of niet behoren tot de grootste groep.
Het is het vermelden waard dat lambda nul is wanneer de modale categorie van elke kolom dezelfde is. Kijk naar deze tabel:
(lambdaTab<-cbind(c(0,100), c(49,51)))## [,1] [,2]## [1,] 0 49## [2,] 100 51
chisq.test(lambdaTab)#### Pearson's Chi-squared test with Yates' continuity correction#### data: lambdaTab## X-squared = 62.279, df = 1, p-value = 2.981e-15Mogelijk gemaakt door
De modale categorie is in elke kolom dezelfde (dat is de tweede rij). Dus moet lambda nul zijn ondanks een significant en zichtbaar verband:
Hoe maak je een tabel
Laten we eerst wat data pakken. Het pakket
MASSbevat gegevens over 93 auto's die in 1993 in de VS te koop waren. Ze zijn opgeslagen in het objectCars93en bevatten 27 kenmerken per auto, waarvan sommige categorisch zijn. Laten we dus het pakketMASSladen en kijken naar het type voertuigen incars93:We hebben daar 6 autotypen. De functie
tablevertelt hoeveel we er van elk type hebben:prop.tablezet dit om naar fracties:Hetzelfde met de herkomst van de auto's:
Hoe maak je een kruistabel
Mooi, we zagen dat onze dataset een vergelijkbaar aantal Amerikaanse en niet-Amerikaanse auto's bevat en dat de meest voorkomende typen Midsize en Small zijn. Maar misschien verschillen de Amerikaanse en niet-Amerikaanse auto's in type?
Laten we naar de autotypen kijken in relatie tot hun herkomst. We kunnen
tableopnieuw gebruiken, maar nu met twee argumenten. Het eerste wordt de rijvariabele en het tweede de kolomvariabele:Nu zagen we wat iedereen weet: Amerikanen houden van grote voertuigen! De bovenstaande tabel laat de gezamenlijke verdeling zien van twee categorische variabelen (
TypeenOrigin). Zulke tabellen noemen we kruistabellen.Hoe krijg je marginalen uit een kruistabel
De functies
rowSumsencolSumsspreken eigenlijk voor zichHoe krijg je percentages uit een kruistabel
prop.tablegenest mettablegeeft frequentiesOmzetten naar percentages is simpelweg vermenigvuldigen met 100:
Let op: dit is een gezamenlijke kansverdeling, waaruit we bijvoorbeeld zien dat zo'n 7,5% van de auto's klein is en van Amerikaanse herkomst.
Vaker zijn we geïnteresseerd in de verdeling van één variabele binnen groepen die door een andere variabele zijn gevormd. Hier lijkt de verdeling van autotypen onder de Amerikaanse en (afzonderlijk) niet-Amerikaanse auto's interessant. Hiervoor gebruiken we het argument
marginvan de functieprop.table. Het geeft aan waar de groeperingsvariabele staat: in rijen (margin=1) of in kolommen (margin=2):Nu zien we makkelijk dat kleine auto's twee keer zo vaak voorkomen in het niet-Amerikaanse deel van onze dataset als in het Amerikaanse.
Let ook op dat kolommen tot 100 optellen, terwijl in de gezamenlijke verdelingstabel (degene zonder het argument
margin) 100 de som van de hele tabel was.Chi-kwadraattoets
De meest voorkomende vraag bij kruistabellen is of de rij- en kolomvariabelen onafhankelijk zijn. De meest basale manier om dit te beantwoorden is een chi-kwadraattoets uitvoeren. Die wordt uitgebreid behandeld in deze tutorial. Laten we controleren of
TypeenOriginonafhankelijk zijn:Blijkbaar niet, maar we kregen ook de waarschuwing
Chi-squared approximation may be incorrect. Dit komt doordat de chi-kwadraatstatistiek slechts bij benadering een chi-kwadraatverdeling volgt. Hoe meer observaties we hebben, hoe beter de benadering. De functiechisq.testgeeft bovenstaande waarschuwing zodra een van de verwachte aantallen lager is dan 5 (voor wat een "expected count" is, zie de gelinkte tutorial).Exacte toets van Fisher
De exacte toets van Fisher is een alternatief voor de chi-kwadraattoets, dat vooral wordt gebruikt wanneer de chi-kwadraatbenadering niet toereikend is. Laten we die uitvoeren:
De resultaten lijken behoorlijk op die van de chi-kwadraat, maar dat is niet altijd zo. Een belangrijk nadeel van de exacte toets van Fisher is dat deze voor grote tabellen (of grote steekproeven) computationeel inefficiënt wordt.
G-toets
Een ander alternatief is de zogeheten G-toets. De statistiek volgt ook bij benadering een chi-kwadraatverdeling, maar voor kleine steekproeven ligt deze benadering dichterbij dan die van de chi-kwadraattoets. Voor de G-toets kunnen we de functie
GTestuit het pakketDescToolsgebruiken. De resultaten zijn opnieuw vergelijkbaar met de twee vorige toetsen:TypeenOriginzijn niet onafhankelijk.Correctie van Yates
In 2x2-kruistabellen kan de chi-kwadraattoets worden verbeterd met de continuïteitscorrectie van Yates. Deze trekt simpelweg 0,5 af van elke term | Observed - Expected | in de chi-kwadraatstatistiek. Raadpleeg opnieuw deze tutorial als je nu de draad kwijt bent. Bovendien past
Rde correctie van Yates automatisch toe wanneer dat nodig is. Laten we kijken naar de beschikbaarheid van uitvoeringen met handgeschakelde transmissie bij Amerikaanse en niet-Amerikaanse auto's:De correctie van Yates werd automatisch toegepast, en
Rgaf dat aan.3 (of meer) dimensionale tabel
Laten we tegelijk kijken naar
Man.trans.avail,OriginenType.tablesplitst de dataset op volgens de als derde opgegeven variabele:ftablegeeft een compacter overzicht:Cochran-Mantel-Haenszel-toets
Uit de bovenstaande
table-resultaten zou duidelijk moeten zijn dat de relatie tussenOriginenMan.trans.availverschilt perType. Voor kleine en sportieve auto's is er bijvoorbeeld geen verband: elke Amerikaanse, evenals elke niet-Amerikaanse kleine of sportieve auto, is verkrijgbaar in een handgeschakelde versie. Aan de andere kant hebben de meeste middelgrote Amerikaanse auto's geen handgeschakelde versie, terwijl de meeste niet-Amerikaanse auto's van dit type die wel hebben. Om rekening te houden met mogelijk verschillende relaties per laag (stratum), kan de Cochran-Mantel-Haenszel-toets worden gebruikt:Het derde argument van
mantelhaen.testidentificeert de strata. Vergelijk bovenstaande resultaten met die zonder stratificatie (het voorbeeld met de correctie van Yates hierboven). Er is nog steeds een verband, maar het bewijs ervoor is een stuk zwakker.Associatiematen
Nu we enkele verbanden tussen variabelen hebben ontdekt, is het tijd om hun sterkte te meten. Er zijn tal van mogelijke maten. Veel daarvan worden hier beschreven. Ik focus nu op de twee meest gebruikte.
Cramers V
V is gebaseerd op de chi-kwadraatstatistiek:
$$ V = \sqrt{\frac{\chi^2/N}{\min(C-1, R-1)}}, $$ waarbij:
V ∈ [0; 1]. Hoe groter V, hoe sterker de relatie tussen variabelen. V = 0 kan worden geïnterpreteerd als onafhankelijkheid (omdat V = 0 dan en slechts dan wanneer χ2 = 0). Het belangrijkste nadeel van V is het gebrek aan precieze interpretatie. Is V = 0,6 een sterke, middelmatige of zwakke associatie?
De functie
CramerVuitDescToolskan dit voor ons berekenen:Nogmaals: is dit een sterke, middelmatige of zwakke associatie?
Goodman en Kruskal's lambda
Goodman en Kruskal's lambda is een voorbeeld van een maat gebaseerd op de proportionele reductie in variatie. Dit soort maten probeert R2 — de determinatiecoëfficiënt uit de lineaire regressie — na te bootsen:
Laten we een kolomvariabele gebruiken als de onafhankelijke. De formule van Goodman en Kruskal's lambda is dan: $$\lambda = \frac{L - \sum_j L_j}{L},$$ waarbij:
Dit is het best uit te leggen met een voorbeeld. Laten we kijken naar
TypeenOrigin, waarbij de laatste de onafhankelijke is:De modus van de VS-kolom is 11, dus L1 = 7 + 10 + 7 + 8 + 5 = 37. De modus van de niet-VS-kolom is 14, dus L2 = 9 + 0 + 12 + 6 + 4 = 31. De "totaal"-kolom is
De modus is 22, dus L = 16 + 11 + 21 + 14 + 9 = 71 en $$\lambda = \frac{71 - (37+31)}{71}=0.042$$
We kunnen in plaats daarvan de functie
Lambdauit het pakketDescToolsgebruiken:De parameter
directionbepaalt waar (inrowof incolumn) de afhankelijke variabele staat.Originverklaart slechts zo'n 4% van de variatie inType. Let op dat lambda variatie definieert als een dichotomie tussen wel of niet behoren tot de grootste groep.Het is het vermelden waard dat lambda nul is wanneer de modale categorie van elke kolom dezelfde is. Kijk naar deze tabel:
De modale categorie is in elke kolom dezelfde (dat is de tweede rij). Dus moet lambda nul zijn ondanks een significant en zichtbaar verband:
Wil je meer leren over statistiek in R, volg dan de DataCamp-cursus Statistical Modeling in R (Part 1).