This website is a replication package for the paper “Exposure to Asylum Seekers and Changing Support for the Radical Right: A Natural Experiment in the Netherlands” by Tolsma, Laméris, and Savelkoul (n.d.). It contains R code to replicate all Tables in the manuscript.
Use the top menu to navigate to the section of interest. The section replicates the Appendix of the manuscript.
To copy the code click the button in the upper right corner of the code-chunks.
We start with a datafile in long format we obtained after performing our dataprep, see here.
rm(list = ls())
library(haven)
dataset <- read_stata("data/evazc v22122017.dta")
# dependent variable. delete not allowed to vote
table(dataset$vote, useNA = "always")
sum(table(dataset$vote, useNA = "always"))
dataset <- dataset[dataset$vote != 16, ]
dataset <- dataset[dataset$vote != 17, ]
# add 'VNL' to other party
dataset$vote[dataset$vote == 12] <- 18
#'1'='vvd', '2'='pvda', '3'='pvv', '4'='cda', '5'='sp'
#'6'='d66', '7'='gl', '8'='cu', '9'='sgp', '10'='dier'
#'11'='50plus', '12'='vnl','13='blanco','14'='dontknow','15='no vote', '18'='other'
# define dichotomous dependent variable
dataset$y_pvv <- ifelse(dataset$vote == 3, 1, 0)
# define multinomial dependent variable
dataset$vote_rwsp <- NA
dataset$vote_rwsp[dataset$vote == 2 | dataset$vote == 6 | dataset$vote == 7 | dataset$vote == 8 | dataset$vote ==
10 | dataset$vote == 11 | dataset$vote == 12 | dataset$vote == 18] <- 0 #other
dataset$vote_rwsp[dataset$vote == 1 | dataset$vote == 4 | dataset$vote == 9] <- 1 # right-wing
dataset$vote_rwsp[dataset$vote == 5] <- 2 # anti-establishment
dataset$vote_rwsp[dataset$vote == 13 | dataset$vote == 14 | dataset$vote == 15] <- 3 #demobilized
dataset$vote_rwsp[dataset$vote == 3] <- 4
# 1: right-wing; 2: anti-establishment; 3: demobilized; 0: other; 4: PVV
# delete respondents with missing data on zipcode
sum(is.na(dataset$pc4))
dataset <- dataset[!is.na(dataset$pc4), ]
# delete respondents with only one observation.
dataset$one <- 1
test <- aggregate(dataset$one, by = list(dataset$PanelistIdQuestion), FUN = sum)
names(test) <- c("PanelistIdQuestion", "ncases")
dataset <- merge(dataset, test)
table(dataset$ncases)
dataset <- dataset[dataset$ncases != 1, ]
# recode tijd/time variable
table(dataset$tijd, useNA = "always")
dataset$tijd <- dataset$tijd - 1
# capacity of ASC s1 is regular, s2=temporary, s3=crisis, s123=total
dataset$s1pc4 <- ifelse(dataset$tijd == 0, dataset$s1cap_pc4t1rel, dataset$s1cap_pc4t2rel)
dataset$s2pc4 <- ifelse(dataset$tijd == 0, dataset$s2cap_pc4t1rel, dataset$s2cap_pc4t2rel)
dataset$s3pc4 <- ifelse(dataset$tijd == 0, dataset$s3cap_pc4t1rel, dataset$s3cap_pc4t2rel)
dataset$s123pc4 <- dataset$s1pc4 + dataset$s2pc4 + dataset$s3pc4
# compute variables for hybrid model: within individual means (between vars) and deviations from mean
# (within vars).
test <- aggregate(dataset[, c("contactnw", "threat", "s123pc4", "s1pc4", "s2pc4", "s3pc4", "y_pvv", "vote_rwsp")],
by = list(dataset$PanelistIdQuestion), FUN = mean)
names(test) <- c("PanelistIdQuestion", "contactnw_mean", "threat_mean", "s123pc4_mean", "s1pc4_mean",
"s2pc4_mean", "s3pc4_mean", "y_pvv_mean", "vote_rwsp_mean")
dataset <- merge(dataset, test)
dataset$contactnw_afw = dataset$contactnw - dataset$contactnw_mean
dataset$threat_afw = dataset$threat - dataset$threat_mean
dataset$s123pc4_afw = dataset$s123pc4 - dataset$s123pc4_mean
dataset$s1pc4_afw = dataset$s1pc4 - dataset$s1pc4_mean
dataset$s2pc4_afw = dataset$s2pc4 - dataset$s2pc4_mean
dataset$s3pc4_afw = dataset$s3pc4 - dataset$s3pc4_mean
# create filter vars for analyses: fe is fixed effects, fem is fixed effects multinomial
dataset$sel_fe = (dataset$y_pvv != dataset$y_pvv_mean)
dataset$sel_fem = (dataset$vote_rwsp != dataset$vote_rwsp_mean)
dataset$sel_s123pc4 = (dataset$s123pc4 != dataset$s123pc4_mean)
# match additional data on crisis ASC these data have been added in a later phase
duur <- read.csv2("data/Adressen alle AZCs_17022017.csv")
names(duur)
duur <- duur[, c("PC4", "Soort", "Cap_Gem", "Duur.tot.4.11")]
names(duur) <- c("pc4", "soort", "cap", "duur")
duur <- duur[duur$soort == 3, ]
sort(unique(duur$pc4[duur$cap > 0]))
# we are multiplying length of stay with number of AS
duur$s3_altop <- duur$duur * duur$cap
duur2 <- aggregate(duur[, c("cap", "duur", "s3_altop")], by = list(duur$pc4), FUN = sum)
# duur2
names(duur2)[1] <- c("pc4")
dataset <- merge(dataset, duur2, all.x = T)
dataset$duur[is.na(dataset$duur)] <- 0
dataset$s3_altop[is.na(dataset$s3_altop)] <- 0
dataset$cap[is.na(dataset$cap)] <- 0
# asylum seekers in crisis centers per 1000 inhabitants weighted by length of stay
dataset$s3_altop2 <- 1000 * dataset$s3_altop/(dataset$inw2014_pc4 * dataset$pauto2014_pc4)
dataset$cap2 <- 1000 * dataset$cap/(dataset$inw2014_pc4 * dataset$pauto2014_pc4)
dataset$s3_altop <- 1000 * dataset$s3_altop/dataset$inw2014_pc4
dataset$cap <- 1000 * dataset$cap/dataset$inw2014_pc4In our manuscript we exploit the panel design of our sample and analyze native Dutch respondents who participated in wave 1 and wave 2. (In our fixed effects models we further select on respondents who changed voting intention). Following the suggestion of reviewer#1 We investigated the possibility that PVV supporters in wave 1 are more motivated to accept the invitation to the 2nd survey than non-PVV supporters. We would like to point out that if this type of selection occurred our results will be biased against our findings, because we are then less likely to observe an increase in PVV support between wave 1 and wave 2.
# this is our data in long format, before we select on panel members
data_attrition <- haven::read_dta("data\\evazc attrition v12042017.dta")
table(data_attrition$tijd)#>
#> 1 2
#> 26064 28434# we check whether respondents of wave1 also participated in wave2
data_attrition$w2_participation <- NA
data_attrition$w2_participation[data_attrition$tijd == 1] <- data_attrition$PanelistIdQuestion[data_attrition$tijd ==
1] %in% data_attrition$PanelistIdQuestion[data_attrition$tijd == 2]
# we check whether zipcodes will experience increase in ASC
data_attrition$s1cap_pc4difrel <- data_attrition$s1cap_pc4t2rel - data_attrition$s1cap_pc4t1rel
data_attrition$s2cap_pc4difrel <- data_attrition$s2cap_pc4t2rel - data_attrition$s2cap_pc4t1rel
data_attrition$s3cap_pc4difrel <- data_attrition$s3cap_pc4t2rel - data_attrition$s3cap_pc4t1rel
data_attrition$s123cap_pc4difrel <- (data_attrition$s1cap_pc4t2rel + data_attrition$s2cap_pc4t2rel +
data_attrition$s3cap_pc4t2rel) - (data_attrition$s1cap_pc4t1rel + data_attrition$s2cap_pc4t1rel +
data_attrition$s3cap_pc4t1rel)
data_at_t1 <- data_attrition[data_attrition$tijd == 1, ]
data_at_t2 <- data_attrition[data_attrition$tijd == 2, ]Is support for PVV in wave1 related to participation in wave 2?
m1 <- glm(w2_participation ~ PVV, family = binomial, data = data_at_t1)
summary(m1)#>
#> Call:
#> glm(formula = w2_participation ~ PVV, family = binomial, data = data_at_t1)
#>
#> Deviance Residuals:
#> Min 1Q Median 3Q Max
#> -1.7826 0.6758 0.7385 0.7385 0.7385
#>
#> Coefficients:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 1.15997 0.01586 73.126 < 2e-16 ***
#> PVV 0.20058 0.04147 4.837 1.32e-06 ***
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> (Dispersion parameter for binomial family taken to be 1)
#>
#> Null deviance: 28306 on 26063 degrees of freedom
#> Residual deviance: 28282 on 26062 degrees of freedom
#> AIC: 28286
#>
#> Number of Fisher Scoring iterations: 4logLik(m1)#> 'log Lik.' -14141.17 (df=2)The results of this test indeed shows that PVV supporters in wave 1 are more likely to participate in wave 2 than non-PVV supporters in wave 1. We hence conclude that we are likely to underestimate the positive relationship between an increase in local exposure to asylum seekers and support for the PVV.
Perhaps more importantly, we would like to investigate whether respondents of Wave 1 who are going to experience an increase in local experience are more likely to participate in Wave 2.
m2 <- glm(w2_participation ~ s1cap_pc4difrel + s2cap_pc4difrel + s3cap_pc4difrel, family = binomial,
data = data_at_t1)
summary(m2)#>
#> Call:
#> glm(formula = w2_participation ~ s1cap_pc4difrel + s2cap_pc4difrel +
#> s3cap_pc4difrel, family = binomial, data = data_at_t1)
#>
#> Deviance Residuals:
#> Min 1Q Median 3Q Max
#> -1.9392 0.7175 0.7277 0.7277 0.8132
#>
#> Coefficients:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 1.193718 0.015086 79.125 <2e-16 ***
#> s1cap_pc4difrel -0.000445 0.001347 -0.330 0.741
#> s2cap_pc4difrel -0.001026 0.001192 -0.861 0.389
#> s3cap_pc4difrel 0.001967 0.002707 0.727 0.467
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> (Dispersion parameter for binomial family taken to be 1)
#>
#> Null deviance: 27552 on 25404 degrees of freedom
#> Residual deviance: 27550 on 25401 degrees of freedom
#> (659 observations deleted due to missingness)
#> AIC: 27558
#>
#> Number of Fisher Scoring iterations: 4logLik(m2)#> 'log Lik.' -13775.13 (df=4)The answer is NO. There is no selectivity with respect to our ‘treatment’, increased exposure to asylum seekers. Not in general, and not among the PVV supporters in wave 1.
And in particular if this is the case for PVV supporters in Wave 1.
m3 <- glm(w2_participation ~ s1cap_pc4difrel * PVV + s2cap_pc4difrel * PVV + s3cap_pc4difrel * PVV, family = binomial,
data = data_at_t1)
summary(m3)#>
#> Call:
#> glm(formula = w2_participation ~ s1cap_pc4difrel * PVV + s2cap_pc4difrel *
#> PVV + s3cap_pc4difrel * PVV, family = binomial, data = data_at_t1)
#>
#> Deviance Residuals:
#> Min 1Q Median 3Q Max
#> -1.9803 0.6781 0.7369 0.7369 0.8256
#>
#> Coefficients:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 1.1647791 0.0163214 71.365 < 2e-16 ***
#> s1cap_pc4difrel -0.0004059 0.0013610 -0.298 0.766
#> PVV 0.1881509 0.0429038 4.385 1.16e-05 ***
#> s2cap_pc4difrel -0.0010531 0.0012648 -0.833 0.405
#> s3cap_pc4difrel 0.0010852 0.0028119 0.386 0.700
#> s1cap_pc4difrel:PVV -0.0006326 0.0088767 -0.071 0.943
#> PVV:s2cap_pc4difrel 0.0005045 0.0038195 0.132 0.895
#> PVV:s3cap_pc4difrel 0.0084162 0.0092573 0.909 0.363
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> (Dispersion parameter for binomial family taken to be 1)
#>
#> Null deviance: 27552 on 25404 degrees of freedom
#> Residual deviance: 27527 on 25397 degrees of freedom
#> (659 observations deleted due to missingness)
#> AIC: 27543
#>
#> Number of Fisher Scoring iterations: 4logLik(m3)#> 'log Lik.' -13763.69 (df=8)Unfortunately, in R it is not possible to estimate a fixed effects multinomial logit model. we have to switch to STATA. I will use R to call Stata.
To replicate this part, please make sure to have Stata installed and within stata the package femlogit.
Party labels:
Set your stata path with chooseStataBin()
require(RStata)
dataset_sel <- dataset[dataset$sel_fe == 1, ]
dataset_selfem <- dataset[dataset$sel_fem == 1, ]
stata_src <- "
femlogit vote_rwsp s1pc4 s2pc4 s3pc4 tijd threat contactnw, group(PanelistIdQuestion) base(4)
"
stata(stata_src, data.in = dataset_selfem, stata.version = 15, stata.path = "\"C:\\Program Files (x86)\\Stata15\\StataSE-64\"")#> .
#> . femlogit vote_rwsp s1pc4 s2pc4 s3pc4 tijd threat contactnw, group(PanelistIdQ
#> > uestion) base(4)
#>
#> Iteration 0: log likelihood = -2624.4138
#> Iteration 1: log likelihood = -2562.8586
#> Iteration 2: log likelihood = -2545.3
#> Iteration 3: log likelihood = -2545.2376
#> Iteration 4: log likelihood = -2545.2362
#> Iteration 5: log likelihood = -2545.2362
#>
#> Fixed-effects multinomial logistic regression Number of obs = 8,466
#> LR chi2(24) = 777.71
#> Prob > chi2 = 0.0000
#> Log likelihood = -2545.2362 Pseudo R2 = 0.1325
#>
#> ------------------------------------------------------------------------------
#> vote_rwsp | Coef. Std. Err. z P>|z| [95% Conf. Interval]
#> -------------+----------------------------------------------------------------
#> _0 |
#> s1pc4 | -.0183711 .0168844 -1.09 0.277 -.0514639 .0147217
#> s2pc4 | -.0654326 .0706959 -0.93 0.355 -.203994 .0731287
#> s3pc4 | -.0338523 .0187038 -1.81 0.070 -.0705111 .0028065
#> tijd | -.8722234 .0734109 -11.88 0.000 -1.016106 -.7283406
#> threat | -.3237684 .0663498 -4.88 0.000 -.4538117 -.1937252
#> contactnw | -.0562277 .0335763 -1.67 0.094 -.122036 .0095806
#> -------------+----------------------------------------------------------------
#> _1 |
#> s1pc4 | -.003719 .01596 -0.23 0.816 -.035 .027562
#> s2pc4 | -.0801347 .0769849 -1.04 0.298 -.2310223 .0707529
#> s3pc4 | -.0372031 .0196402 -1.89 0.058 -.0756972 .001291
#> tijd | -.9072487 .0814028 -11.15 0.000 -1.066795 -.7477021
#> threat | -.2071291 .0732672 -2.83 0.005 -.3507303 -.063528
#> contactnw | .002693 .0380905 0.07 0.944 -.0719631 .0773491
#> -------------+----------------------------------------------------------------
#> _2 |
#> s1pc4 | -.0049259 .0157157 -0.31 0.754 -.0357281 .0258762
#> s2pc4 | -.07914 .0727013 -1.09 0.276 -.221632 .063352
#> s3pc4 | -.0274022 .0191666 -1.43 0.153 -.064968 .0101637
#> tijd | -2.069922 .0922435 -22.44 0.000 -2.250716 -1.889128
#> threat | -.3650243 .0876528 -4.16 0.000 -.5368208 -.1932279
#> contactnw | -.0619441 .0445424 -1.39 0.164 -.1492456 .0253574
#> -------------+----------------------------------------------------------------
#> _3 |
#> s1pc4 | -.0073686 .0157294 -0.47 0.639 -.0381977 .0234605
#> s2pc4 | -.0644874 .0706524 -0.91 0.361 -.2029636 .0739888
#> s3pc4 | -.0223896 .0187887 -1.19 0.233 -.0592148 .0144356
#> tijd | -.8945652 .0730327 -12.25 0.000 -1.037707 -.7514237
#> threat | -.1908516 .0655231 -2.91 0.004 -.3192744 -.0624287
#> contactnw | -.0690404 .0331943 -2.08 0.038 -.1341 -.0039808
#> -------------+----------------------------------------------------------------
#> _4 | (base outcome)
#> ------------------------------------------------------------------------------We tested whether effects were significantly different across choices with setting constraints. An example of how we did this is given below.
stata_src <- "
cons 7 [0=1]: s2pc4
cons 8 [1=2]: s2pc4
cons 9 [2=3]: s2pc4
femlogit vote_rwsp s1pc4 s2pc4 s3pc4 tijd threat contactnw, group(PanelistIdQuestion) base(4)
estimates store m1
femlogit vote_rwsp s1pc4 s2pc4 s3pc4 tijd threat contactnw, group(PanelistIdQuestion) base(4) const(7/11)
estimates store m2
lrtest m2 m1
"
stata(stata_src, data.in = dataset_selfem, stata.version = 15, stata.path = "\"C:\\Program Files (x86)\\Stata15\\StataSE-64\"")require(survival)
ev_sel <- dataset[dataset$sel_fe == 1, ]
ev_sel2 <- ev_sel[ev_sel$tijd == 1, ]
# need to multiply _afw with 2 (thus becomes t2-t1 (First difference design), instead of t2 -
# mean(t2,t1) (fixed effects design)
ev_sel2$s123pc4_afw <- 2 * ev_sel2$s123pc4_afw
ev_sel2$s1pc4_afw <- 2 * ev_sel2$s1pc4_afw
ev_sel2$s2pc4_afw <- 2 * ev_sel2$s2pc4_afw
ev_sel2$s3pc4_afw <- 2 * ev_sel2$s3pc4_afw
ev_sel2$threat_afw <- 2 * ev_sel2$threat_afw
ev_sel2$contactnw_afw <- 2 * ev_sel2$contactnw_afw
# robustness check with newly matched crisis data.
m1_alt1e <- glm(y_pvv ~ s123pc4_afw, data = ev_sel2, family = binomial)
m1_alt2e <- glm(y_pvv ~ s1pc4_afw + s2pc4_afw + cap, data = ev_sel2, family = binomial)
m1_alt3e <- glm(y_pvv ~ s1pc4_afw + s2pc4_afw + cap + contactnw_afw + threat_afw, data = ev_sel2, family = binomial)
# summary(m1_alt1e) summary(m1_alt2e)
summary(m1_alt3e)
# very small changes, due to 3 additional matched zipcodes#>
#> Call:
#> glm(formula = y_pvv ~ s1pc4_afw + s2pc4_afw + cap + contactnw_afw +
#> threat_afw, family = binomial, data = ev_sel2)
#>
#> Deviance Residuals:
#> Min 1Q Median 3Q Max
#> -2.1826 0.4360 0.6895 0.7633 1.1679
#>
#> Coefficients:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 1.084221 0.064792 16.734 < 2e-16 ***
#> s1pc4_afw 0.007373 0.013654 0.540 0.5892
#> s2pc4_afw 0.061345 0.066965 0.916 0.3596
#> cap 0.027678 0.016781 1.649 0.0991 .
#> contactnw_afw 0.034183 0.029355 1.164 0.2442
#> threat_afw 0.265453 0.058505 4.537 5.7e-06 ***
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> (Dispersion parameter for binomial family taken to be 1)
#>
#> Null deviance: 1534.7 on 1388 degrees of freedom
#> Residual deviance: 1504.6 on 1383 degrees of freedom
#> AIC: 1516.6
#>
#> Number of Fisher Scoring iterations: 7# robustness check with time heterogeneity
m1_alt1c <- glm(y_pvv ~ s123pc4_afw, data = ev_sel2, family = binomial)
m1_alt2c <- glm(y_pvv ~ s1pc4_afw + s2pc4_afw + s3_altop, data = ev_sel2, family = binomial)
m1_alt3c <- glm(y_pvv ~ s1pc4_afw + s2pc4_afw + s3_altop + contactnw_afw + threat_afw, data = ev_sel2,
family = binomial)
# summary(m1_alt1c) summary(m1_alt2c)
summary(m1_alt3c)#>
#> Call:
#> glm(formula = y_pvv ~ s1pc4_afw + s2pc4_afw + s3_altop + contactnw_afw +
#> threat_afw, family = binomial, data = ev_sel2)
#>
#> Deviance Residuals:
#> Min 1Q Median 3Q Max
#> -2.0320 0.4773 0.7002 0.7580 1.1608
#>
#> Coefficients:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 1.100081 0.064642 17.018 < 2e-16 ***
#> s1pc4_afw 0.007528 0.013889 0.542 0.588
#> s2pc4_afw 0.060322 0.066525 0.907 0.365
#> s3_altop 0.001513 0.002318 0.653 0.514
#> contactnw_afw 0.032986 0.029300 1.126 0.260
#> threat_afw 0.265221 0.058424 4.540 5.64e-06 ***
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> (Dispersion parameter for binomial family taken to be 1)
#>
#> Null deviance: 1534.7 on 1388 degrees of freedom
#> Residual deviance: 1507.6 on 1383 degrees of freedom
#> AIC: 1519.6
#>
#> Number of Fisher Scoring iterations: 7# robustness check with time only, not shown in manuscript
m1_alt3c2 <- glm(y_pvv ~ s1pc4_afw + s2pc4_afw + duur + contactnw_afw + threat_afw, data = ev_sel2, family = binomial)
# controlling for days crisis ASC were used.
m1_alt4c2 <- glm(y_pvv ~ s1pc4_afw + s2pc4_afw + s3pc4_afw + duur + contactnw_afw + threat_afw, data = ev_sel2,
family = binomial)
# summary(m1_alt3c2)
summary(m1_alt4c2)#>
#> Call:
#> glm(formula = y_pvv ~ s1pc4_afw + s2pc4_afw + s3pc4_afw + duur +
#> contactnw_afw + threat_afw, family = binomial, data = ev_sel2)
#>
#> Deviance Residuals:
#> Min 1Q Median 3Q Max
#> -2.3398 0.3127 0.6930 0.7584 1.1591
#>
#> Coefficients:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 1.09907 0.06525 16.843 < 2e-16 ***
#> s1pc4_afw 0.01014 0.01625 0.624 0.53264
#> s2pc4_afw 0.06014 0.06639 0.906 0.36501
#> s3pc4_afw 0.06390 0.02468 2.589 0.00962 **
#> duur -0.11347 0.04277 -2.653 0.00797 **
#> contactnw_afw 0.02939 0.02954 0.995 0.31982
#> threat_afw 0.26393 0.05862 4.502 6.73e-06 ***
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> (Dispersion parameter for binomial family taken to be 1)
#>
#> Null deviance: 1534.7 on 1388 degrees of freedom
#> Residual deviance: 1495.7 on 1382 degrees of freedom
#> AIC: 1509.7
#>
#> Number of Fisher Scoring iterations: 7ev_sel2$treatment <- ev_sel2$s123pc4_afw > 0
ev_sel2$treatments1 <- ev_sel2$s1pc4_afw > 0
ev_sel2$treatments2 <- ev_sel2$s2pc4_afw > 0
ev_sel2$treatments3 <- ev_sel2$s3pc4_afw > 0
m1_alt1b <- glm(y_pvv ~ treatment, data = ev_sel2, family = binomial)
m1_alt2b <- glm(y_pvv ~ treatments1 + treatments2 + treatments3, data = ev_sel2, family = binomial)
m1_alt3b <- glm(y_pvv ~ treatments1 + treatments2 + treatments3 + contactnw_afw + threat_afw, data = ev_sel2,
family = binomial)
# summary(m1_alt1b) summary(m1_alt2b)
summary(m1_alt3b)#>
#> Call:
#> glm(formula = y_pvv ~ treatments1 + treatments2 + treatments3 +
#> contactnw_afw + threat_afw, family = binomial, data = ev_sel2)
#>
#> Deviance Residuals:
#> Min 1Q Median 3Q Max
#> -2.0392 0.5090 0.7018 0.7542 1.1588
#>
#> Coefficients:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 1.11183 0.06616 16.805 < 2e-16 ***
#> treatments1TRUE -0.41549 0.51173 -0.812 0.417
#> treatments2TRUE 0.50267 0.64165 0.783 0.433
#> treatments3TRUE 0.12041 0.28623 0.421 0.674
#> contactnw_afw 0.03279 0.02936 1.117 0.264
#> threat_afw 0.26699 0.05843 4.569 4.9e-06 ***
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> (Dispersion parameter for binomial family taken to be 1)
#>
#> Null deviance: 1534.7 on 1388 degrees of freedom
#> Residual deviance: 1509.9 on 1383 degrees of freedom
#> AIC: 1521.9
#>
#> Number of Fisher Scoring iterations: 4This one is not shown in manuscript. There is a discussion in the literature at which geographical scale we could expect ‘context effects’. The crisis centers were a very local phenomenon. We thus do not expect effects at municipality level. Most residents will not be aware of crisis centers. At the municipality level the discussion was more focused on new temporary centers.
# robustness check with exposure on municipality level
ev_sel2$s123gc_afw <- 1000 * (ev_sel2$s1cap_gct2 + ev_sel2$s2cap_gct2 + ev_sel2$s3cap_gct2 - ev_sel2$s1cap_gct1 -
ev_sel2$s2cap_gct1 - ev_sel2$s3cap_gct1)/ev_sel2$inw2014_gc
ev_sel2$s1gc_afw <- 1000 * (ev_sel2$s1cap_gct2 - ev_sel2$s1cap_gct1)/ev_sel2$inw2014_gc
ev_sel2$s2gc_afw <- 1000 * (ev_sel2$s2cap_gct2 - ev_sel2$s2cap_gct1)/ev_sel2$inw2014_gc
ev_sel2$s3gc_afw <- 1000 * (ev_sel2$s3cap_gct2 - ev_sel2$s3cap_gct1)/ev_sel2$inw2014_gc
m1_alt1d <- glm(y_pvv ~ s123gc_afw, data = ev_sel2, family = binomial)
m1_alt2d <- glm(y_pvv ~ s1gc_afw + s2gc_afw + s3gc_afw, data = ev_sel2, family = binomial)
m1_alt3d <- glm(y_pvv ~ s1gc_afw + s2gc_afw + s3gc_afw + contactnw_afw + threat_afw, data = ev_sel2,
family = binomial)
# summary(m1_alt1d) summary(m1_alt2d)
summary(m1_alt3d)#>
#> Call:
#> glm(formula = y_pvv ~ s1gc_afw + s2gc_afw + s3gc_afw + contactnw_afw +
#> threat_afw, family = binomial, data = ev_sel2)
#>
#> Deviance Residuals:
#> Min 1Q Median 3Q Max
#> -2.1231 0.4447 0.6912 0.7671 1.1763
#>
#> Coefficients:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 1.07278 0.07449 14.403 < 2e-16 ***
#> s1gc_afw -0.01145 0.03202 -0.358 0.7207
#> s2gc_afw 0.08222 0.04624 1.778 0.0754 .
#> s3gc_afw 0.01125 0.03037 0.370 0.7111
#> contactnw_afw 0.03025 0.02923 1.035 0.3008
#> threat_afw 0.26756 0.05860 4.566 4.98e-06 ***
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> (Dispersion parameter for binomial family taken to be 1)
#>
#> Null deviance: 1534.7 on 1388 degrees of freedom
#> Residual deviance: 1506.8 on 1383 degrees of freedom
#> AIC: 1518.8
#>
#> Number of Fisher Scoring iterations: 5Let us construct sample weights for the fixed effects sample. We show descriptive statistics pre- and post-weighting in Appendix 4.
require(anesrake)
ev_sel <- dataset[dataset$sel_fe == 1, ]
ev_sel2 <- ev_sel[ev_sel$tijd == 1, ]
# need to multiply _afw with 2 (thus becomes t2-t1 (First difference design), instead of t2 -
# mean(t2,t1) (fixed effects design)
ev_sel2$s123pc4_afw <- 2 * ev_sel2$s123pc4_afw
ev_sel2$s1pc4_afw <- 2 * ev_sel2$s1pc4_afw
ev_sel2$s2pc4_afw <- 2 * ev_sel2$s2pc4_afw
ev_sel2$s3pc4_afw <- 2 * ev_sel2$s3pc4_afw
ev_sel2$threat_afw <- 2 * ev_sel2$threat_afw
ev_sel2$contactnw_afw <- 2 * ev_sel2$contactnw_afw
# ADD WEIGHTS TO Data make caseidn
ev_sel2$caseidn <- 1:length(ev_sel2[, 1])
N <- length(ev_sel2[, 1])
# MAKE TARGET VALUES Geslacht
gender_w <- c(0.503, 0.497) # vrouw - man
# Age: 3cats
age_w <- c(0.238, 0.344, 0.418) #18-35,35-55,55+
names(age_w) <- c("age1", "age2", "age3")
# Opl: 3cats
educ_w <- c(0.323, 0.393, 0.284) #basis/vmbo/mbo1,mbo2-4/HV/VWO,HBO/WO
names(educ_w) <- c("opl1", "opl2", "opl3")
targets <- list(gender_w, age_w, educ_w)
names(targets) <- c("gender_w", "age_w", "educ_w")
# check variables in data table(ev_sel2$gender, useNA='always') #man=1
ev_sel2$gender_w <- as.logical(ev_sel2$gender)
# table(ev_sel2$gender_w, useNA='always')
# age: 3cats table(ev_sel2$age, useNA='always') #18-35,35-55,55+
ev_sel2$age_w <- NA
ev_sel2$age_w[ev_sel2$age <= 35] <- 1
ev_sel2$age_w[ev_sel2$age > 35 & ev_sel2$age <= 55] <- 2
ev_sel2$age_w[ev_sel2$age > 55] <- 3
ev_sel2$age_w <- as.factor(ev_sel2$age_w)
levels(ev_sel2$age_w) <- c("age1", "age2", "age3")
# table(ev_sel2$age_w, useNA='always')
# Opl: 3cats table(ev_sel2$educ, useNA='always')
ev_sel2$educ_w <- NA
ev_sel2$educ_w[ev_sel2$educ <= 8] <- 1
ev_sel2$educ_w[ev_sel2$educ > 8 & ev_sel2$educ <= 12] <- 2
ev_sel2$educ_w[ev_sel2$educ > 12] <- 3
ev_sel2$educ_w <- as.factor(ev_sel2$educ_w)
levels(ev_sel2$educ_w) <- c("opl1", "opl2", "opl3")
# table(ev_sel2$educ_w, useNA='always')
# Education (years of education that constitute the shortest route to obtain a university degree.
# (After this coding, all values are subtracted from the maximum years of schooling necessary to
# obtain a university grade)) (opl1) 5 lagere school not finished (4yrs) /lagere school (6yrs) (opl1)
# 7 lbo, vmbo-kb/bbl (6,5yrs) /mavo, vmbo-tl (8yrs) (opl2) 11 havo (10yrs) / vwo/gymnasium (12yrs)
# (opl2) 9,5 mbo-kort (kmbo) (8,5yrs) /mbo-tussen/lang (mbo) (10,5yrs) (opl3) 14 hbo (14yrs) (opl3)
# 16,5 universiteit (bachelormaster, doctoraal)(16,5)
# make weights
weighted <- anesrake(inputter = targets, dataframe = ev_sel2, choosemethod = "total", caseid = ev_sel2$caseidn,
cap = 3, pctlim = 5, verbose = FALSE)
test <- weightassess(targets, ev_sel2, weighted$weightvec)
# test
# table(ev_sel2$educ) oplp <- test$educ_w[c(1:3),1]*N oplw <- test$educ_w[c(1:3),4] opls <-
# table(ev_sel2$educ_w) test2 <- matrix(c(oplp, oplw, opls), byrow=T, nrow=3) dimnames(test2) <-
# list(weightype=c('population','weighted sample', 'original sample'),opl=c('low','medium','high'))
# test2 chisq.test(test2[1:2,]) agep <- test$age_w[c(1:3),1]*N agew <- test$age_w[c(1:3),4] test2 <-
# matrix(c(agep, agew), byrow=T, nrow=2) dimnames(test2) <-
# list(weightype=c('p','w'),age=c('1','2','3')) test2 chisq.test(test2) geslp <-
# test$gender_w[c(1:2),1]*N geslw <- test$gender_w[c(1:2),4] test2 <- matrix(c(geslp, geslw),
# byrow=T, nrow=2) dimnames(test2) <- list(weightype=c('p','w'),age=c('v','m')) test2
# chisq.test(test2)
ev_sel2 <- cbind(ev_sel2, weighted$weightvec)
names(ev_sel2)[names(ev_sel2) == "weighted$weightvec"] <- "weightvec3"#> [1] "raking achieved only partial convergence, please check the results to ensure that sufficient convergence was achieved."
#> [1] "no improvement was apparent after 85 iterations"
#> [1] "current total change in the iteration is: 0.0396851187449571 average change per weight is: 2.85709998163838e-05"And let us perform a robustness check on the weighted sample.
require(survey)
des2 <- svydesign(id = ~1, weights = ~weightvec3, data = ev_sel2)
# check if weighing gives correct mean svymean(~gender, des2) yes it does.
glm.sampling.weights1 <- svyglm(y_pvv ~ s123pc4_afw, family = quasibinomial(), design = des2)
glm.sampling.weights2 <- svyglm(y_pvv ~ s1pc4_afw + s2pc4_afw + s3pc4_afw, family = quasibinomial(),
design = des2)
glm.sampling.weights3 <- svyglm(y_pvv ~ s1pc4_afw + s2pc4_afw + s3pc4_afw + contactnw_afw + threat_afw,
family = quasibinomial(), design = des2)
summary(glm.sampling.weights3)#>
#> Call:
#> svyglm(formula = y_pvv ~ s1pc4_afw + s2pc4_afw + s3pc4_afw +
#> contactnw_afw + threat_afw, design = des2, family = quasibinomial())
#>
#> Survey design:
#> svydesign(id = ~1, weights = ~weightvec3, data = ev_sel2)
#>
#> Coefficients:
#> Estimate Std. Error t value Pr(>|t|)
#> (Intercept) 1.15674 0.10111 11.440 < 2e-16 ***
#> s1pc4_afw 0.03683 0.05262 0.700 0.4841
#> s2pc4_afw 0.15215 0.08602 1.769 0.0771 .
#> s3pc4_afw 0.03675 0.01986 1.850 0.0645 .
#> contactnw_afw 0.09457 0.04758 1.988 0.0470 *
#> threat_afw 0.42119 0.08623 4.885 1.16e-06 ***
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> (Dispersion parameter for quasibinomial family taken to be 0.9809453)
#>
#> Number of Fisher Scoring iterations: 9We analyze whether persons are more likely to have voted in Time1 for PVV than at Time0. Although we use a fixed effects model (thus focus on within person changes) and thereby control for time-stable heterogeneity, we cannot control for time-varying heterogeneity. Because of possible pre-treatment differences among individuals (both measured, and unmeasured), the causal effect may also differ across individuals.
The idea with pre-processing our data with a matching procedure is that our treated group is as similar as possible to the control group after matching and that the treatment is closer to being independent to (unmeasured time-varying) covariates. After matching we can be more confidant that we have ‘no ommitted variable bias’.
“Recall that under the usual econometric conditions for omitted variable bias, a variable Xi must be controlled for if it is causally prior to Ti, empirically related to Ti, and affects Yi conditional on Ti. If instead one or more of the three conditions do not hold, then Xi may be omitted without any resulting bias (although the variance may increase).”
Ho et al. (2007)
We define the treatment as zipcodes that experienced in increase in asylum seekers in ASC. Although we define a binary treatment variable, in the model we use the continuous ‘treatment’ (i.e. the increase in asylum seekers per 1000 neigbhorhood residents) as predictor.
require("MatchIt")
require(optmatch)
ev_sel2 <- ev_sel[ev_sel$tijd == 1, ]
# start preprocessing data (matching) because we also want to match on ethnic density of NB and SES
# of NB need to remove missing values. table(ev_sel2$pnwal2014_pc4, useNA='always')
ev_sel2 <- ev_sel2[!is.na(ev_sel2$pnwal2014_pc4), ]
# lost 3 respondents, no treated table(ev_sel2$woz2012_pc4, useNA='always')
ev_sel2 <- ev_sel2[!(ev_sel2$woz2012_pc4 == 0), ]
# lost 1 respondents, no treated take log of NBses
ev_sel2$LNwoz <- log(ev_sel2$woz2012_pc4)
summary(ev_sel2$LNwoz)
# define treatment variable: whether there has been an increase in total asylum seekers.
ev_sel2$treatment <- ev_sel2$s123pc4_afw > 0
table(ev_sel2$treatment, useNA = "always")
# 111 treated
# define treatment variable2: whether there has been an increase in asylum seekers in crisis centers.
ev_sel2$treatment2 <- ev_sel2$s3pc4_afw > 0
table(ev_sel2$treatment2, useNA = "always")
# 75 treated
ev_sel2$threat_t1 <- ev_sel2$threat_mean - ev_sel2$threat_afw
ev_sel2$contactnw_t1 <- ev_sel2$contactnw_mean - ev_sel2$contactnw_afw
# let use smaller dataset (to remove missing values in non relevant vars)
ev_sel2 <- ev_sel2[, c("y_pvv", "treatment", "treatment2", "s123pc4_afw", "s1pc4_afw", "s2pc4_afw", "s3pc4_afw",
"gender", "age", "educ", "contactnw_afw", "threat_afw", "contactnw_t1", "threat_t1", "pnwal2014_pc4",
"woz2012_pc4", "LNwoz")]
m.out <- matchit(treatment ~ gender + scale(age) + scale(educ) + scale(contactnw_t1) + scale(threat_t1) +
scale(pnwal2014_pc4) + scale(LNwoz), data = ev_sel2)
data_matched <- match.data(m.out)
# quickly see effect of matching
prop.table(table(ev_sel2$treatment, ev_sel2$y_pvv), margin = 1)
prop.table(table(data_matched$treatment, data_matched$y_pvv), margin = 1)
# conclusion: treated areas less likely to see increase in support for PVV
summary(glm(y_pvv ~ treatment, data = data_matched, family = binomial))
# but not significant
# m1_alt_matched1 <- glm(y_pvv ~ s123pc4_afw , data=data_matched, family=binomial) m1_alt_matched2 <-
# glm(y_pvv ~ s1pc4_afw + s2pc4_afw + s3pc4_afw , data=data_matched, family=binomial)
m1_alt_matched3 <- glm(y_pvv ~ s1pc4_afw + s2pc4_afw + s3pc4_afw + contactnw_afw + threat_afw, data = data_matched,
family = binomial)
# summary(m1_alt_matched1) summary(m1_alt_matched2)
summary(m1_alt_matched3)#> Min. 1st Qu. Median Mean 3rd Qu. Max.
#> 3.746 5.092 5.281 5.270 5.444 6.849
#>
#> FALSE TRUE <NA>
#> 1274 111 0
#>
#> FALSE TRUE <NA>
#> 1310 75 0
#>
#> 0 1
#> FALSE 0.2417582 0.7582418
#> TRUE 0.2342342 0.7657658
#>
#> 0 1
#> FALSE 0.2252252 0.7747748
#> TRUE 0.2342342 0.7657658
#>
#> Call:
#> glm(formula = y_pvv ~ treatment, family = binomial, data = data_matched)
#>
#> Deviance Residuals:
#> Min 1Q Median 3Q Max
#> -1.7267 0.7144 0.7144 0.7306 0.7306
#>
#> Coefficients:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 1.23547 0.22722 5.437 5.41e-08 ***
#> treatmentTRUE -0.05092 0.31915 -0.160 0.873
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> (Dispersion parameter for binomial family taken to be 1)
#>
#> Null deviance: 239.29 on 221 degrees of freedom
#> Residual deviance: 239.27 on 220 degrees of freedom
#> AIC: 243.27
#>
#> Number of Fisher Scoring iterations: 4
#>
#>
#> Call:
#> glm(formula = y_pvv ~ s1pc4_afw + s2pc4_afw + s3pc4_afw + contactnw_afw +
#> threat_afw, family = binomial, data = data_matched)
#>
#> Deviance Residuals:
#> Min 1Q Median 3Q Max
#> -2.1647 0.1668 0.6346 0.7872 1.4076
#>
#> Coefficients:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 1.012934 0.192592 5.259 1.44e-07 ***
#> s1pc4_afw 0.007049 0.015234 0.463 0.6436
#> s2pc4_afw 0.129256 0.139715 0.925 0.3549
#> s3pc4_afw 0.065414 0.038936 1.680 0.0929 .
#> contactnw_afw -0.006188 0.147258 -0.042 0.9665
#> threat_afw 0.769756 0.304061 2.532 0.0114 *
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> (Dispersion parameter for binomial family taken to be 1)
#>
#> Null deviance: 239.29 on 221 degrees of freedom
#> Residual deviance: 225.25 on 216 degrees of freedom
#> AIC: 237.25
#>
#> Number of Fisher Scoring iterations: 8This one is not shown in manuscript. We now define a binary treatment based on an increase in crisis centers.
# ## lets match on crisis centers
m.out2 <- matchit(treatment2 ~ scale(gender) + scale(age) + scale(educ) + scale(contactnw_t1) + scale(threat_t1) +
scale(pnwal2014_pc4) + scale(LNwoz), data = ev_sel2)
summary(m.out2)
data_matched2 <- match.data(m.out2)
# #quickly see effect of matching
prop.table(table(ev_sel2$treatment2, ev_sel2$y_pvv), margin = 1)
prop.table(table(data_matched2$treatment2, data_matched2$y_pvv), margin = 1)
# #conclusion: in treated zipcodes stronger increase in support for PVV
summary(glm(y_pvv ~ treatment2, data = data_matched2, family = binomial))
# but not significant.
m1_alt_matched1b <- glm(y_pvv ~ s123pc4_afw, data = data_matched2, family = binomial)
m1_alt_matched2b <- glm(y_pvv ~ s3pc4_afw, data = data_matched2, family = binomial)
m1_alt_matched3b <- glm(y_pvv ~ s3pc4_afw + contactnw_afw + threat_afw, data = data_matched2, family = binomial)
# summary(m1_alt_matched2b)
summary(m1_alt_matched3b)
# #conclusion: matching makes effect stronger#>
#> Call:
#> matchit(formula = treatment2 ~ scale(gender) + scale(age) + scale(educ) +
#> scale(contactnw_t1) + scale(threat_t1) + scale(pnwal2014_pc4) +
#> scale(LNwoz), data = ev_sel2)
#>
#> Summary of Balance for All Data:
#> Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean eCDF Max
#> distance 0.0609 0.0538 0.3617 0.8503 0.1193 0.2312
#> scale(gender) -0.0629 0.0036 -0.0633 1.1095 0.0133 0.0267
#> scale(age) 0.0385 -0.0022 0.0404 1.0201 0.0157 0.0674
#> scale(educ) 0.1023 -0.0059 0.1107 0.9528 0.0290 0.0754
#> scale(contactnw_t1) -0.0031 0.0002 -0.0031 1.1137 0.0180 0.0453
#> scale(threat_t1) -0.0582 0.0033 -0.0647 0.8985 0.0279 0.0583
#> scale(pnwal2014_pc4) -0.0846 0.0048 -0.1009 0.7762 0.0656 0.1367
#> scale(LNwoz) 0.3466 -0.0198 0.4196 0.7574 0.1192 0.2375
#>
#>
#> Summary of Balance for Matched Data:
#> Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean eCDF Max
#> distance 0.0609 0.0608 0.0023 1.0078 0.0008 0.0267
#> scale(gender) -0.0629 -0.1626 0.0949 0.8964 0.0200 0.0400
#> scale(age) 0.0385 -0.0332 0.0711 1.2220 0.0236 0.0800
#> scale(educ) 0.1023 0.0347 0.0693 0.9196 0.0400 0.1333
#> scale(contactnw_t1) -0.0031 -0.0217 0.0176 1.0471 0.0362 0.0667
#> scale(threat_t1) -0.0582 0.0934 -0.1595 0.8572 0.0560 0.1200
#> scale(pnwal2014_pc4) -0.0846 0.0037 -0.0997 0.6341 0.0482 0.1333
#> scale(LNwoz) 0.3466 0.3600 -0.0154 0.7971 0.0312 0.0800
#> Std. Pair Dist.
#> distance 0.0035
#> scale(gender) 0.8541
#> scale(age) 0.9673
#> scale(educ) 1.1253
#> scale(contactnw_t1) 1.0993
#> scale(threat_t1) 1.1411
#> scale(pnwal2014_pc4) 0.9312
#> scale(LNwoz) 0.3998
#>
#> Percent Balance Improvement:
#> Std. Mean Diff. Var. Ratio eCDF Mean eCDF Max
#> distance 99.4 95.2 99.4 88.5
#> scale(gender) -50.0 -5.3 -50.0 -50.0
#> scale(age) -76.1 -906.7 -49.8 -18.7
#> scale(educ) 37.4 -73.3 -38.1 -76.8
#> scale(contactnw_t1) -461.4 57.3 -101.5 -47.2
#> scale(threat_t1) -146.5 -44.0 -100.7 -105.8
#> scale(pnwal2014_pc4) 1.2 -79.8 26.5 2.5
#> scale(LNwoz) 96.3 18.4 73.8 66.3
#>
#> Sample Sizes:
#> Control Treated
#> All 1310 75
#> Matched 75 75
#> Unmatched 1235 0
#> Discarded 0 0
#>
#>
#> 0 1
#> FALSE 0.2419847 0.7580153
#> TRUE 0.2266667 0.7733333
#>
#> 0 1
#> FALSE 0.2666667 0.7333333
#> TRUE 0.2266667 0.7733333
#>
#> Call:
#> glm(formula = y_pvv ~ treatment2, family = binomial, data = data_matched2)
#>
#> Deviance Residuals:
#> Min 1Q Median 3Q Max
#> -1.7229 0.7170 0.7170 0.7876 0.7876
#>
#> Coefficients:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 1.0116 0.2611 3.874 0.000107 ***
#> treatment2TRUE 0.2156 0.3798 0.568 0.570208
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> (Dispersion parameter for binomial family taken to be 1)
#>
#> Null deviance: 167.59 on 149 degrees of freedom
#> Residual deviance: 167.27 on 148 degrees of freedom
#> AIC: 171.27
#>
#> Number of Fisher Scoring iterations: 4
#>
#>
#> Call:
#> glm(formula = y_pvv ~ s3pc4_afw + contactnw_afw + threat_afw,
#> family = binomial, data = data_matched2)
#>
#> Deviance Residuals:
#> Min 1Q Median 3Q Max
#> -2.5076 0.1407 0.6372 0.7994 1.2257
#>
#> Coefficients:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 0.77855 0.23509 3.312 0.000927 ***
#> s3pc4_afw 0.09363 0.04389 2.133 0.032915 *
#> contactnw_afw 0.29595 0.18900 1.566 0.117382
#> threat_afw 0.39704 0.35489 1.119 0.263236
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> (Dispersion parameter for binomial family taken to be 1)
#>
#> Null deviance: 167.59 on 149 degrees of freedom
#> Residual deviance: 157.77 on 146 degrees of freedom
#> AIC: 165.77
#>
#> Number of Fisher Scoring iterations: 5The weighing statistics.
test#> $gender_w
#> Target Unweighted N Unweighted % Wtd N Wtd % Change in % Resid. Disc. Orig. Disc.
#> FALSE 0.497 279 0.2008639 690.2131 0.4969137 0.2960497 8.632599e-05 0.2961361
#> TRUE 0.503 1110 0.7991361 698.7869 0.5030863 -0.2960497 -8.632599e-05 -0.2961361
#> Total 1.000 1389 1.0000000 1389.0000 1.0000000 0.5920995 1.726520e-04 0.5922721
#>
#> $age_w
#> Target Unweighted N Unweighted % Wtd N Wtd % Change in % Resid. Disc. Orig. Disc.
#> age1 0.238 39 0.02807775 117.0038 0.08423599 0.05615824 0.15376401 0.2099222
#> age2 0.344 266 0.19150468 574.1759 0.41337356 0.22186888 -0.06937356 0.1524953
#> age3 0.418 1084 0.78041757 697.8203 0.50239045 -0.27802712 -0.08439045 -0.3624176
#> Total 1.000 1389 1.00000000 1389.0000 1.00000000 0.55605424 0.30752801 0.7248351
#>
#> $educ_w
#> Target Unweighted N Unweighted % Wtd N Wtd % Change in % Resid. Disc. Orig. Disc.
#> opl1 0.323 300 0.2159827 524.7036 0.3777564 0.16177365 -0.05475637 0.10701728
#> opl2 0.393 520 0.3743701 577.0288 0.4154275 0.04105744 -0.02242749 0.01862995
#> opl3 0.284 569 0.4096472 287.2676 0.2068161 -0.20283109 0.07718386 -0.12564723
#> Total 1.000 1389 1.0000000 1389.0000 1.0000000 0.40566218 0.15436772 0.25129446The balance statistics.
summary(m.out)#>
#> Call:
#> matchit(formula = treatment ~ gender + scale(age) + scale(educ) +
#> scale(contactnw_t1) + scale(threat_t1) + scale(pnwal2014_pc4) +
#> scale(LNwoz), data = ev_sel2)
#>
#> Summary of Balance for All Data:
#> Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean eCDF Max
#> distance 0.0847 0.0797 0.2521 1.0158 0.0809 0.1476
#> gender 0.7838 0.7998 -0.0390 . 0.0161 0.0161
#> scale(age) 0.1081 -0.0094 0.1254 0.8700 0.0210 0.0697
#> scale(educ) 0.0630 -0.0055 0.0709 0.9291 0.0178 0.0388
#> scale(contactnw_t1) 0.0638 -0.0056 0.0681 1.0399 0.0230 0.0701
#> scale(threat_t1) -0.0402 0.0035 -0.0467 0.8687 0.0274 0.0596
#> scale(pnwal2014_pc4) -0.0243 0.0021 -0.0308 0.7150 0.0524 0.1506
#> scale(LNwoz) 0.1810 -0.0158 0.2063 0.9054 0.0539 0.1211
#>
#>
#> Summary of Balance for Matched Data:
#> Means Treated Means Control Std. Mean Diff. Var. Ratio eCDF Mean eCDF Max
#> distance 0.0847 0.0847 0.0007 1.0037 0.0008 0.0180
#> gender 0.7838 0.8649 -0.1970 . 0.0811 0.0811
#> scale(age) 0.1081 0.0514 0.0605 0.9344 0.0175 0.0811
#> scale(educ) 0.0630 0.0730 -0.0103 1.0861 0.0165 0.0450
#> scale(contactnw_t1) 0.0638 0.0638 0.0000 1.0433 0.0180 0.0631
#> scale(threat_t1) -0.0402 0.0938 -0.1430 0.9087 0.0523 0.1081
#> scale(pnwal2014_pc4) -0.0243 -0.0924 0.0796 0.9119 0.0460 0.1351
#> scale(LNwoz) 0.1810 0.3225 -0.1482 0.9139 0.0501 0.1712
#> Std. Pair Dist.
#> distance 0.0029
#> gender 0.5909
#> scale(age) 0.9495
#> scale(educ) 0.9978
#> scale(contactnw_t1) 0.9616
#> scale(threat_t1) 1.1185
#> scale(pnwal2014_pc4) 0.9692
#> scale(LNwoz) 0.7542
#>
#> Percent Balance Improvement:
#> Std. Mean Diff. Var. Ratio eCDF Mean eCDF Max
#> distance 99.7 76.2 99.1 87.8
#> gender -404.9 . -404.9 -404.9
#> scale(age) 51.8 51.2 16.9 -16.4
#> scale(educ) 85.4 -12.3 7.0 -16.1
#> scale(contactnw_t1) 100.0 -8.5 21.5 10.0
#> scale(threat_t1) -206.4 32.0 -90.9 -81.5
#> scale(pnwal2014_pc4) -158.2 72.5 12.2 10.3
#> scale(LNwoz) 28.2 9.4 7.0 -41.4
#>
#> Sample Sizes:
#> Control Treated
#> All 1274 111
#> Matched 111 111
#> Unmatched 1163 0
#> Discarded 0 0To obtain DiD estimators, we estimate our dichotomous outcome variable with a Linear Probability Model. We use “identity” as link function set the variance to “mu(1-mu)” (as in the binomial distribution) and we allow for (possible) overdispersion. We report cluster corrected standard errors, with individuals as clusters. The DiD model is estimated on the ‘hybrid sample’.
require(sandwich)
require(lmtest)
ev <- dataset
# Because we include ethnic density of NB and SES of NB need as time-constant covariates we have to
# remove missing values. Similar to our analyses reported in the main text.
sum(is.na(ev$pnwal2014_pc4))#> [1] 90ev_sel4 <- ev[!is.na(ev$pnwal2014_pc4), ]
# lost 90 respondents
ev_sel4 <- ev_sel4[!(ev_sel4$woz2012_pc4 == 0), ]
# lost 22 respondents
ev_sel4 <- ev_sel4[!is.na(ev_sel4$woz2012_pc4), ]
# lost 3 respondents
# center variables!
ev_sel4$age_c <- ev_sel4$age - mean(ev_sel4$age)
ev_sel4$educ_c <- ev_sel4$educ - mean(ev_sel4$educ)
ev_sel4$pnwal2014_pc4_c <- ev_sel4$pnwal2014_pc4 - mean(ev_sel4$pnwal2014_pc4)
ev_sel4$woz2012_pc4_c <- ev_sel4$woz2012_pc4 - mean(ev_sel4$woz2012_pc4)
# divide to help estimation
ev_sel4$woz2012_pc4_c <- ev_sel4$woz2012_pc4_c/100
ev <- ev_sel4
# define treatments! for continuous treatment vars we multiple the deviations from mean exposure with
# 2
# treatment (total - binary)
ev$treatment_did1b <- ifelse((ev$s123pc4_afw > 0 & ev$tijd == 1) | (ev$s123pc4_afw < 0 & ev$tijd == 0),
1, 0)
# treatment (crisis - continuous)
ev$treatment_did1c <- 2 * ifelse(ev$tijd == 1, ev$s123pc4_afw, -ev$s123pc4_afw)
# treatment (crisis - binary)
ev$treatment_did2b <- ifelse((ev$s3pc4_afw > 0 & ev$tijd == 1) | (ev$s3pc4_afw < 0 & ev$tijd == 0), 1,
0)
# treatment (crisis - continuous)
ev$treatment_did2c <- 2 * ifelse(ev$tijd == 1, ev$s3pc4_afw, -ev$s3pc4_afw)# summary(ev$woz2012_pc4_c)
did_lm <- lm(y_pvv ~ treatment_did1b + tijd + tijd:treatment_did1b + gender + age_c + educ_c + pnwal2014_pc4_c +
woz2012_pc4_c, data = ev)
# summary(did_lm) round(coeftest(did_lm, vcov = vcovCL(did_lm, cluster = ~ PanelistIdQuestion)), 3)
# coef(did_lm)
# did_lpm <- glm(y_pvv ~ treatment_did1b + tijd + tijd:treatment_did1b + gender + age_c + educ_c +
# pnwal2014_pc4_c + woz2012_pc4_c, family=binomial(link='identity'), data=ev,
# start=round(coef(did_lm),2)) summary(did_lpm) round(coeftest(did_lpm_het, vcov =
# vcovCL(did_lpm_het, cluster = ~ PanelistIdQuestion)), 3)
did_lpm_het <- glm(y_pvv ~ treatment_did1b + tijd + tijd:treatment_did1b + gender + age_c + educ_c +
pnwal2014_pc4_c + woz2012_pc4_c, family = quasi(link = "identity", variance = "mu(1-mu)"), data = ev,
start = round(coef(did_lm), 2))
# summary(did_lpm_het)
round(coeftest(did_lpm_het, vcov = vcovCL(did_lpm_het, cluster = ~PanelistIdQuestion)), 3)#>
#> z test of coefficients:
#>
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 0.117 0.004 31.028 <2e-16 ***
#> treatment_did1b 0.010 0.011 0.933 0.351
#> tijd 0.035 0.002 15.849 <2e-16 ***
#> gender 0.073 0.005 16.124 <2e-16 ***
#> age_c -0.001 0.000 -4.268 <2e-16 ***
#> educ_c -0.020 0.001 -28.850 <2e-16 ***
#> pnwal2014_pc4_c 0.085 0.033 2.568 0.010 **
#> woz2012_pc4_c -0.001 0.000 -6.190 <2e-16 ***
#> treatment_did1b:tijd 0.000 0.009 -0.042 0.966
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1did_lm <- lm(y_pvv ~ treatment_did1b + tijd + tijd:treatment_did1b, data = ev)
# summary(did_lm) round(coeftest(did_lm, vcov = vcovCL(did_lm, cluster = ~ PanelistIdQuestion)), 3)
# coef(did_lm)
# did_lpm <- glm(y_pvv ~ treatment_did1b + tijd + tijd:treatment_did1b + gender + age_c + educ_c +
# pnwal2014_pc4_c + woz2012_pc4_c, family=binomial(link='identity'), data=ev,
# start=round(coef(did_lm),2)) summary(did_lpm) round(coeftest(did_lpm_het, vcov =
# vcovCL(did_lpm_het, cluster = ~ PanelistIdQuestion)), 3)
did_lpm_het <- glm(y_pvv ~ treatment_did1b + tijd + tijd:treatment_did1b, family = quasi(link = "identity",
variance = "mu(1-mu)"), data = ev, start = round(coef(did_lm), 2))
# summary(did_lpm_het)
round(coeftest(did_lpm_het, vcov = vcovCL(did_lpm_het, cluster = ~PanelistIdQuestion)), 3)#>
#> z test of coefficients:
#>
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 0.167 0.003 59.562 <2e-16 ***
#> treatment_did1b 0.001 0.010 0.142 0.887
#> tijd 0.037 0.002 18.614 <2e-16 ***
#> treatment_did1b:tijd 0.005 0.008 0.620 0.535
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1did_lm <- lm(y_pvv ~ treatment_did1c + tijd + tijd:treatment_did1c + gender + age_c + educ_c + pnwal2014_pc4_c +
woz2012_pc4_c, data = ev)
# summary(did_lm) round(coeftest(did_lm, vcov = vcovCL(did_lm, cluster = ~ PanelistIdQuestion)), 3)
# coef(did_lm)
# did_lpm <- glm(y_pvv ~ treatment_did1b + tijd + tijd:treatment_did1b + gender + age_c + educ_c +
# pnwal2014_pc4_c + woz2012_pc4_c, family=binomial(link='identity'), data=ev,
# start=round(coef(did_lm),2)) summary(did_lpm) round(coeftest(did_lpm_het, vcov =
# vcovCL(did_lpm_het, cluster = ~ PanelistIdQuestion)), 3)
did_lpm_het <- glm(y_pvv ~ treatment_did1c + tijd + tijd:treatment_did1c + gender + age_c + educ_c +
pnwal2014_pc4_c + woz2012_pc4_c, family = quasi(link = "identity", variance = "mu(1-mu)"), data = ev,
start = round(coef(did_lm), 2))
summary(did_lpm_het)#>
#> Call:
#> glm(formula = y_pvv ~ treatment_did1c + tijd + tijd:treatment_did1c +
#> gender + age_c + educ_c + pnwal2014_pc4_c + woz2012_pc4_c,
#> family = quasi(link = "identity", variance = "mu(1-mu)"),
#> data = ev, start = round(coef(did_lm), 2))
#>
#> Deviance Residuals:
#> Min 1Q Median 3Q Max
#> -1.0041 -0.7094 -0.5682 -0.4093 3.3616
#>
#> Coefficients:
#> Estimate Std. Error t value Pr(>|t|)
#> (Intercept) 1.201e-01 3.151e-03 38.103 < 2e-16 ***
#> treatment_did1c 1.709e-05 1.467e-04 0.117 0.9073
#> tijd 3.286e-02 3.620e-03 9.077 < 2e-16 ***
#> gender 7.005e-02 3.563e-03 19.661 < 2e-16 ***
#> age_c -8.311e-04 1.377e-04 -6.035 1.61e-09 ***
#> educ_c -2.054e-02 5.740e-04 -35.776 < 2e-16 ***
#> pnwal2014_pc4_c 7.678e-02 2.596e-02 2.957 0.0031 **
#> woz2012_pc4_c -4.122e-04 1.812e-03 -0.228 0.8200
#> treatment_did1c:tijd 6.818e-05 2.264e-04 0.301 0.7633
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> (Dispersion parameter for quasi family taken to be 1.003036)
#>
#> Null deviance: 36681 on 38181 degrees of freedom
#> Residual deviance: 34985 on 38173 degrees of freedom
#> AIC: NA
#>
#> Number of Fisher Scoring iterations: 25round(coeftest(did_lpm_het, vcov = vcovCL(did_lpm_het, cluster = ~PanelistIdQuestion)), 3)#>
#> z test of coefficients:
#>
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 0.120 0.004 29.624 <2e-16 ***
#> treatment_did1c 0.000 0.000 0.105 0.916
#> tijd 0.033 0.002 15.833 <2e-16 ***
#> gender 0.070 0.005 14.428 <2e-16 ***
#> age_c -0.001 0.000 -4.308 <2e-16 ***
#> educ_c -0.021 0.001 -26.257 <2e-16 ***
#> pnwal2014_pc4_c 0.077 0.035 2.203 0.028 *
#> woz2012_pc4_c 0.000 0.003 -0.142 0.887
#> treatment_did1c:tijd 0.000 0.000 0.530 0.596
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1did_lm <- lm(y_pvv ~ treatment_did1c + tijd + tijd:treatment_did1c, data = ev)
# summary(did_lm) round(coeftest(did_lm, vcov = vcovCL(did_lm, cluster = ~ PanelistIdQuestion)), 3)
# coef(did_lm)
# did_lpm <- glm(y_pvv ~ treatment_did1b + tijd + tijd:treatment_did1b + gender + age_c + educ_c +
# pnwal2014_pc4_c + woz2012_pc4_c, family=binomial(link='identity'), data=ev,
# start=round(coef(did_lm),2)) summary(did_lpm) round(coeftest(did_lpm_het, vcov =
# vcovCL(did_lpm_het, cluster = ~ PanelistIdQuestion)), 3)
did_lpm_het <- glm(y_pvv ~ treatment_did1c + tijd + tijd:treatment_did1c, family = quasi(link = "identity",
variance = "mu(1-mu)"), data = ev, start = round(coef(did_lm), 2))
# summary(did_lpm_het)
round(coeftest(did_lpm_het, vcov = vcovCL(did_lpm_het, cluster = ~PanelistIdQuestion)), 4)#>
#> z test of coefficients:
#>
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 0.1675 0.0027 61.6763 <2e-16 ***
#> treatment_did1c -0.0002 0.0001 -1.3256 0.1850
#> tijd 0.0371 0.0019 19.0460 <2e-16 ***
#> treatment_did1c:tijd 0.0002 0.0001 1.9531 0.0508 .
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1did_lm <- lm(y_pvv ~ treatment_did2b + tijd + tijd:treatment_did2b + gender + age_c + educ_c + pnwal2014_pc4_c +
woz2012_pc4_c, data = ev)
# summary(did_lm) round(coeftest(did_lm, vcov = vcovCL(did_lm, cluster = ~ PanelistIdQuestion)), 3)
# coef(did_lm)
# did_lpm <- glm(y_pvv ~ treatment_did1b + tijd + tijd:treatment_did1b + gender + age_c + educ_c +
# pnwal2014_pc4_c + woz2012_pc4_c, family=binomial(link='identity'), data=ev,
# start=round(coef(did_lm),2)) summary(did_lpm) round(coeftest(did_lpm_het, vcov =
# vcovCL(did_lpm_het, cluster = ~ PanelistIdQuestion)), 3)
did_lpm_het <- glm(y_pvv ~ treatment_did2b + tijd + tijd:treatment_did2b + gender + age_c + educ_c +
pnwal2014_pc4_c + woz2012_pc4_c, family = quasi(link = "identity", variance = "mu(1-mu)"), data = ev,
start = round(coef(did_lm), 2))
# summary(did_lpm_het)
round(coeftest(did_lpm_het, vcov = vcovCL(did_lpm_het, cluster = ~PanelistIdQuestion)), 3)#>
#> z test of coefficients:
#>
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 0.120 0.004 28.740 <2e-16 ***
#> treatment_did2b 0.002 0.013 0.120 0.904
#> tijd 0.033 0.002 15.406 <2e-16 ***
#> gender 0.070 0.005 14.161 <2e-16 ***
#> age_c -0.001 0.000 -3.897 <2e-16 ***
#> educ_c -0.021 0.001 -25.575 <2e-16 ***
#> pnwal2014_pc4_c 0.077 0.035 2.187 0.029 *
#> woz2012_pc4_c -0.001 0.003 -0.155 0.877
#> treatment_did2b:tijd 0.008 0.011 0.757 0.449
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1did_lm <- lm(y_pvv ~ treatment_did2b + tijd + tijd:treatment_did2b, data = ev)
# summary(did_lm) round(coeftest(did_lm, vcov = vcovCL(did_lm, cluster = ~ PanelistIdQuestion)), 3)
# coef(did_lm)
# did_lpm <- glm(y_pvv ~ treatment_did1b + tijd + tijd:treatment_did1b + gender + age_c + educ_c +
# pnwal2014_pc4_c + woz2012_pc4_c, family=binomial(link='identity'), data=ev,
# start=round(coef(did_lm),2)) summary(did_lpm) round(coeftest(did_lpm_het, vcov =
# vcovCL(did_lpm_het, cluster = ~ PanelistIdQuestion)), 3)
did_lpm_het <- glm(y_pvv ~ treatment_did2b + tijd + tijd:treatment_did2b, family = quasi(link = "identity",
variance = "mu(1-mu)"), data = ev, start = round(coef(did_lm), 2))
# summary(did_lpm_het)
round(coeftest(did_lpm_het, vcov = vcovCL(did_lpm_het, cluster = ~PanelistIdQuestion)), 3)#>
#> z test of coefficients:
#>
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 0.167 0.003 60.452 <2e-16 ***
#> treatment_did2b -0.004 0.013 -0.289 0.772
#> tijd 0.037 0.002 18.858 <2e-16 ***
#> treatment_did2b:tijd 0.008 0.010 0.820 0.412
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1did_lm <- lm(y_pvv ~ treatment_did2c + tijd + tijd:treatment_did2c + gender + age_c + educ_c + pnwal2014_pc4_c +
woz2012_pc4_c, data = ev)
# summary(did_lm) round(coeftest(did_lm, vcov = vcovCL(did_lm, cluster = ~ PanelistIdQuestion)), 3)
# coef(did_lm)
# did_lpm <- glm(y_pvv ~ treatment_did1b + tijd + tijd:treatment_did1b + gender + age_c + educ_c +
# pnwal2014_pc4_c + woz2012_pc4_c, family=binomial(link='identity'), data=ev,
# start=round(coef(did_lm),2)) summary(did_lpm) round(coeftest(did_lpm_het, vcov =
# vcovCL(did_lpm_het, cluster = ~ PanelistIdQuestion)), 3)
did_lpm_het <- glm(y_pvv ~ treatment_did2c + tijd + tijd:treatment_did2c + gender + age_c + educ_c +
pnwal2014_pc4_c + woz2012_pc4_c, family = quasi(link = "identity", variance = "mu(1-mu)"), data = ev,
start = round(coef(did_lm), 2))
# summary(did_lpm_het)
round(coeftest(did_lpm_het, vcov = vcovCL(did_lpm_het, cluster = ~PanelistIdQuestion)), 3)#>
#> z test of coefficients:
#>
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 0.118 0.004 29.423 <2e-16 ***
#> treatment_did2c 0.000 0.001 0.358 0.720
#> tijd 0.034 0.002 16.506 <2e-16 ***
#> gender 0.072 0.005 15.112 <2e-16 ***
#> age_c -0.001 0.000 -3.795 <2e-16 ***
#> educ_c -0.020 0.001 -26.682 <2e-16 ***
#> pnwal2014_pc4_c 0.083 0.034 2.423 0.015 *
#> woz2012_pc4_c -0.001 0.003 -0.327 0.744
#> treatment_did2c:tijd 0.000 0.000 0.503 0.615
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1did_lm <- lm(y_pvv ~ treatment_did2c + tijd + tijd:treatment_did2c, data = ev)
# summary(did_lm) round(coeftest(did_lm, vcov = vcovCL(did_lm, cluster = ~ PanelistIdQuestion)), 3)
# coef(did_lm)
# did_lpm <- glm(y_pvv ~ treatment_did1b + tijd + tijd:treatment_did1b + gender + age_c + educ_c +
# pnwal2014_pc4_c + woz2012_pc4_c, family=binomial(link='identity'), data=ev,
# start=round(coef(did_lm),2)) summary(did_lpm) round(coeftest(did_lpm_het, vcov =
# vcovCL(did_lpm_het, cluster = ~ PanelistIdQuestion)), 3)
did_lpm_het <- glm(y_pvv ~ treatment_did2c + tijd + tijd:treatment_did2c, family = quasi(link = "identity",
variance = "mu(1-mu)"), data = ev, start = round(coef(did_lm), 2))
# summary(did_lpm_het)
round(coeftest(did_lpm_het, vcov = vcovCL(did_lpm_het, cluster = ~PanelistIdQuestion)), 3)#>
#> z test of coefficients:
#>
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 0.167 0.003 61.175 <2e-16 ***
#> treatment_did2c 0.000 0.000 -0.017 0.986
#> tijd 0.037 0.002 19.011 <2e-16 ***
#> treatment_did2c:tijd 0.001 0.000 1.744 0.081 .
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Ho, Daniel E., Kosuke Imai, Gary King, and Elizabeth A. Stuart. 2007. “Matching as Nonparametric Preprocessing for Reducing Model Dependence in Parametric Causal Inference.” Political Analysis 15 (3): 199–236. https://doi.org/10.1093/pan/mpl013.
Tolsma, Jochem, Joran Laméris, and Michael Savelkoul. n.d. “Exposure to Asylum Seekers and Changing Support for the Radical Right: A Natural Experiment in the Netherlands.” PLOS ONE - (-): –. https://journals.plos.org/plosone/.
Copyright © 2020 Jochem Tolsma