#Packages
library(rio) #laden van data
library(sjPlot) #overzichten van data objecten
library(tidyverse) #datamanipulatie en grafieken
library(marginaleffects) #voorspelde waarden berekenen
#Data
ESS9NL <- import("ESS9e03, Netherlands.sav")11 Voorspelde kansen
In het vorige hoofdstuk lag de focus op marginale effecten: hoeveel verandert de probabiliteit dat Y=1 gemiddeld genomen als X 1 eenheid stijgt. In dit hoofdstuk gebruiken we de predictions() functie uit het marginaleffects package om de onderliggende kansen te bestuderen. We gebruiken deze functie voor drie soorten voorspelde kansen:
- De voorspelde kans voor elke observatie in het model
- De verwachte kansen als een bepaalde predictor verschillende waarden aanneemt
- De verwachte kans als alle predictoren specifieke waarden aannemen
Zie Hoofdstuk 5 om na te gaan hoe je deze functie gebruikt voor lineaire modellen.
We gebruiken deze packages en data:
We maken gebruik van hetzelfde model als in het vorige hoofdstuk. We voorspelden daarin de kans dat iemand ging stemmen op basis van gender, leeftijd, vertrouwen in politici en links-rechtsideologie. We doen hier eerst de nodige data management, meer informatie over deze stappen vind je in vorige hoofdstukken:
#Datamanagement
ESS9NL <- ESS9NL |>
#Factor maken van categorische variabelen
mutate(gndr_F = factorize(gndr),
vote_F = factorize(vote)) |>
#Not Eligible op missing zetten
mutate(vote_F = na_if(vote_F,"Not eligible to vote")) |>
#Relevel van variabelen
mutate(vote_F = relevel(vote_F, "No"),
gndr_F = relevel(gndr_F, "Female"))
#Het model
Vote_model_mp <- glm(vote_F ~ gndr_F + agea + trstplt + lrscale,
data = ESS9NL, family = "binomial")
#Resultaten printen
summary(Vote_model_mp)
Call:
glm(formula = vote_F ~ gndr_F + agea + trstplt + lrscale, family = "binomial",
data = ESS9NL)
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -0.284194 0.380455 -0.747 0.455
gndr_FMale 0.043281 0.154201 0.281 0.779
agea 0.018349 0.004503 4.075 4.61e-05 ***
trstplt 0.195020 0.038706 5.039 4.69e-07 ***
lrscale 0.029257 0.039306 0.744 0.457
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 1173.9 on 1424 degrees of freedom
Residual deviance: 1135.3 on 1420 degrees of freedom
(248 observations deleted due to missingness)
AIC: 1145.3
Number of Fisher Scoring iterations: 411.1 Voorspelde kans voor individuele observaties
Om op basis van het logistische regressiemodel de kans dat de afhankelijke variabele Y gelijk is aan 1 (hier: dat een respondent heeft gestemd) te voorspellen voor elke observatie in het model maken we gebruik van de predictions functie. De resultaten worden altijd opgeslaan in een nieuwe dataset, die je een naam geeft naar keuze (hier: Vote_pred).
#Resultaten opslaan in nieuw object
Vote_pred <- predictions(Vote_model_mp,
conf_level = 0.95,
newdata = ESS9NL)
#tibble() gebruiken voor overzicht
tibble(Vote_pred)# A tibble: 1,673 × 580
rowid estimate p.value s.value conf.low conf.high name essround edition
<int> <dbl> <dbl> <dbl> <dbl> <dbl> <chr> <dbl> <chr>
1 1 0.835 6.13e-35 114. 0.796 0.867 ESS9e03 9 3
2 2 0.910 1.71e-56 185. 0.884 0.931 ESS9e03 9 3
3 3 0.904 4.26e-45 147. 0.874 0.928 ESS9e03 9 3
4 4 NA NA NA NA NA ESS9e03 9 3
5 5 0.864 9.71e-39 126. 0.828 0.894 ESS9e03 9 3
6 6 0.912 7.81e-53 173. 0.884 0.933 ESS9e03 9 3
7 7 0.800 1.32e-12 39.5 0.731 0.854 ESS9e03 9 3
8 8 0.914 9.43e-31 99.7 0.877 0.941 ESS9e03 9 3
9 9 0.877 1.88e-47 155. 0.845 0.903 ESS9e03 9 3
10 10 0.944 2.64e-39 128. 0.917 0.963 ESS9e03 9 3
# ℹ 1,663 more rows
# ℹ 571 more variables: proddate <chr>, idno <dbl>, cntry <chr>, nwspol <dbl>,
# netusoft <dbl>, netustm <dbl>, ppltrst <dbl>, pplfair <dbl>, pplhlp <dbl>,
# polintr <dbl>, psppsgva <dbl>, actrolga <dbl>, psppipla <dbl>,
# cptppola <dbl>, trstprl <dbl>, trstlgl <dbl>, trstplc <dbl>, trstplt <dbl>,
# trstprt <dbl>, trstep <dbl>, trstun <dbl>, vote <dbl>, prtvtcat <dbl>,
# prtvtdbe <dbl>, prtvtdbg <dbl>, prtvtgch <dbl>, prtvtbcy <dbl>, …Dit is de syntax-uitleg:
Vote_pred <--
Nieuw data object met voorspelde kansen.
predictions(Vote_model_mp,-
We voeren de predictions functie uit op het model tussen haakjes.
conf_level = 0.95,-
Standaard betrouwbaarheidsniveau. De waarde kan veranderd worden (bv.
conf_level = 0.99). newdata = ESS9NL)-
We kopiëren de variabelen uit de originele dataset. Dit gedeelte kan weggelaten worden als je dit niet nodig acht.
11.2 Gemiddelde voorspelde kansen
We kunnen de predictions() functie ook gebruiken om de gemiddelde voorspelde kans dat Y=1 te berekenen voor specifieke waarden van een onafhankelijke variabele. De andere onafhankelijke variabelen worden constant gehouden op hun gemiddelde (continue variabelen) of modus (factor variabelen). Deze voorspellingen kunnen we ook weergeven in een figuur zoals besproken in Paragraaf 14.5 .
11.2.1 Continue onafhankelijke variabele
De volgende code gebruiken we als de predictor die ons interesseert continu is. Hier berekenen we de gemiddelde voorspelde kans om te stemmen als vertrouwen in politici (trstplt) verandert.
ESS9NL |>
select(trstplt) |>
view_df()| ID | Name | Label | Values | Value Labels |
| 1 | trstplt | Trust in politicians | 0 1 2 3 4 5 6 7 8 9 10 77 88 99 |
No trust at all 1 2 3 4 5 6 7 8 9 Complete trust Refusal Don't know No answer |
table(ESS9NL$trstplt)
0 1 2 3 4 5 6 7 8 9 10
43 41 68 90 173 292 457 375 87 16 8 De variabele loopt van 0 tot 10 dus deze waarden gebruiken we als minimum en maximum. We berekenen verder de kans per interval van 2 eenheden (missing waarden zijn reeds op NA gezet). We zouden ook voor elke eenheid de kans kunnen berekenen maar dit geeft vrij veel output, wellicht meer dan we nodig hebben. We zouden in plaats van deze intervallen ook minimum, 1ste kwartiel, mediaan, 3rde kwartiel en maximum van de variabele kunnen gebruiken (zie Paragraaf 5.3.1 voor het proces om deze waarden te verkrijgen).
#Voorspellingen in nieuw object
Pred_conts <- predictions(Vote_model_mp,
newdata = datagrid(trstplt = seq(from = 0, to = 10, by = 2))) newdata = datagrid(trstplt-
Alle predictoren in het model worden op hun gemiddelde/modus gehouden behalve de predictor die tussen haakjes staat.
= seq(from=0,to=10,by=2)))-
We vragen hier voorspellingen voor een sequentie (
seq) van waarden: van (from) het minimum tot (to) het maximum met tussenstappen (by) van 2. We zouden als alternatief deze code kunnen gebruiken:trstplt = c(0,2,4,6,8,10)).
Laten we de voorspellingen bekijken:
tibble(Pred_conts)- 1
-
tibble()wordt gebruikt om de onderliggende data beter te kunnen zien.
# A tibble: 6 × 11
rowid estimate p.value s.value conf.low conf.high gndr_F agea lrscale
<int> <dbl> <dbl> <dbl> <dbl> <dbl> <fct> <dbl> <dbl>
1 1 0.699 9.37e- 5 13.4 0.603 0.779 Male 50.7 5.15
2 2 0.774 1.09e-15 49.7 0.717 0.822 Male 50.7 5.15
3 3 0.835 1.61e-46 152. 0.802 0.863 Male 50.7 5.15
4 4 0.882 1.29e-64 212. 0.856 0.904 Male 50.7 5.15
5 5 0.917 6.30e-48 157. 0.889 0.938 Male 50.7 5.15
6 6 0.942 3.65e-34 111. 0.912 0.962 Male 50.7 5.15
# ℹ 2 more variables: trstplt <dbl>, vote_F <fct>
Output uitleg
De output is gelijkaardig aan die voor voorspellingen voor lineaire regressiemodellen (zie Hoofdstuk 5):
- De
estimatekolom bevat de voorspelde kans. We doen dit getal *100 om de kans in percent uit te drukken. bv. Voor iemand met een score van ‘0’ op de schaaltrstpltverwachten we een kans om te stemmen van 69,86%. - De
p.valuet.e.m.conf.highkolommen geven de onzekerheid van de schatting weer. - We kunnen ook de kolommen zien voor de andere onafhankelijke variabelen in het model (
gndr_F,agea,lrscale). Deze kolommen tonen de waarde waarop deze variabelen constant worden gehouden. Depredictions()functie houdt continue variabelen op hun gemiddelde en factor variabelen op hun modus als jenewdata = datagrid()gebruikt zoals we hierboven gedaan hebben. - De laatste 2 kolommen tonen
trstplt, met de waarden gebruikt om de voorspelling te berekenen en een kolom (niet zichtbaar hier) met de Y (vote_F) die toont welke categorie voorspeld wordt.
11.2.2 Factor onafhankelijke variabele
De code voor categorische variabelen is licht anders. We gebruiken hier de by= optie. In dit voorbeeld berekenen we de gemiddelde voorspelde kans voor mannen en vrouwen (met andere predictoren constant gehouden op hun gemiddelde).
#voorspellingen in nieuw object
Pred_cat <- predictions(Vote_model_mp,
by="gndr_F",
newdata = "mean")
#tibble voor overzicht
tibble(Pred_cat)# A tibble: 2 × 13
rowid gndr_F estimate std.error statistic p.value s.value conf.low conf.high
<int> <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 1 Female 0.863 0.0133 65.1 0 Inf 0.837 0.889
2 2 Male 0.868 0.0127 68.2 0 Inf 0.843 0.893
# ℹ 4 more variables: agea <dbl>, trstplt <dbl>, lrscale <dbl>,
# rowid_dedup <int>by="gndr_F"-
Met deze optie duiden we aan dat we voor elk niveau van de factor variabele een voorspelde kans willen berekenen.
newdata = "mean"-
Deze optie moeten we toevoegen om duidelijk te maken dat andere predictoren op hun gemiddelde/modus gehouden moeten worden. Deze optie moet samen met
by=gebruikt worden.
11.3 Voorspelde kansen voor specifieke waarden van predictoren
Ten slotte kunnen we de voorspelde kans op Y berekenen als een observatie bepaalde, hypothetische waarden zou aannemen.
Bijvoorbeeld, hier berekenen we de voorspelde kans om te stemmen voor een man (gndr), die 33 jaar oud is (agea), een score van 2 heeft voor vertrouwen in politici (trstplt) en een score van 8 heeft op de links-rechts schaal (lrscale). We moeten hiervoor de waarden voor alle predictoren verduidelijken tussen haakjes bij newdata=datagrid.
#Berekenen en opslaan in object
Pred_specific <- predictions(Vote_model_mp,
newdata = datagrid(gndr_F=c("Male"),
agea=c(33),
trstplt=c(2),
lrscale=c(8)))
#bekijken
Pred_specific- 1
- We gebruiken haakjes omdat dit een factor variabele is met labels voor categorieën.
gndr_F agea trstplt lrscale Estimate Pr(>|z|) S 2.5 % 97.5 %
Male 33 2 8 0.729 <0.001 20.2 0.645 0.799
Type: invlink(link)