Replicate Model 2 from Schlueter, Masso, and
Davidov (2019) .
R solution -> dont’ peek to early ;) !
library(tidyverse) # Data manipulation
library(haven) # Read Stata dta files
library(essurvey) # API to ESS
ESS <- import_rounds(rounds = 7, ess_email = YOUR_EMAIL) %>%
recode_missings() %>%
# redefine some variables as factors
mutate(
# Educational attainment
edu_attain = case_when(
edulvlb < 200 ~ 1,
edulvlb > 199 & edulvlb < 300 ~ 2,
edulvlb > 299 & edulvlb < 400 ~ 3,
edulvlb > 399 & edulvlb < 500 ~ 4,
edulvlb > 499 & edulvlb < 600 ~ 5,
edulvlb > 599 & edulvlb < 700 ~ 6,
edulvlb > 699 & edulvlb < 801 ~ 7,
TRUE ~ as.numeric(NA)
),
# Watching TV, but not news on TV
TV_not_news = tvtot - tvpol,
# Immigrant friends
dfegcf = max(dfegcf, na.rm = TRUE) - dfegcf,
brncntr = as_factor(brncntr),
ctzcntr = as_factor(ctzcntr),
rlgblg = as_factor(rlgblg),
rlgdnm = as_factor(rlgdnm),
uempla = as_factor(uempla),
rlgdnm = case_when( # Add atheists to the religion variable
rlgblg == "No" ~ "Atheist",
rlgdnm=="Not applicable"|rlgdnm=="Refusal"|rlgdnm=="No answer" ~ as.character(NA),
TRUE ~ as.character(rlgdnm)
) %>% factor(),
gndr = as_factor(gndr),
gndr = case_when(
gndr == "No answer" ~ as.character(NA),
TRUE ~ as.character(gndr)
) %>% factor()
) %>% zap_labels() %>%
# keep only those born and with citizenship in country of interview,
# non-Muslim and drop Israel
filter(brncntr == "Yes" & ctzcntr == "Yes" & rlgdnm != "Islamic" & cntry != "IL") %>%
inner_join(., read_dta("./../../assets/ESS7/Schlueter.dta") %>%
mutate(
z_muslim_j = (muslim_j - mean(muslim_j, na.rm= TRUE)) / sd(muslim_j, na.rm = TRUE)
), by = "cntry") %>%
# Keep only the following variables
select(cntry, almuslv, agea, gndr, edu_attain, TV_not_news,
uempla, hincfel, rlgdnm, dfegcf, tvtot, tvpol,
muslim_j, z_muslim_j, mipex_j, relstate_j, claims_j) %>%
drop_na() # Casewise deletion
library(lme4) # Multilevel modeling
# Model 2
Model_2 <- lmer(
almuslv ~ gndr + agea + edu_attain + uempla + hincfel +
dfegcf + rlgdnm + TV_not_news + muslim_j + mipex_j +
relstate_j + claims_j + (1 | cntry), data = ESS
)
summary(Model_2)
# Linear mixed model fit by REML ['lmerMod']
# Formula: almuslv ~ gndr + agea + edu_attain + uempla + hincfel + dfegcf +
# rlgdnm + TV_not_news + muslim_j + mipex_j + relstate_j + claims_j + (1 | cntry)
# Data: ESS
#
# REML criterion at convergence: 73799
#
# Scaled residuals:
# Min 1Q Median 3Q Max
# -3.557 -0.694 0.009 0.743 3.020
#
# Random effects:
# Groups Name Variance Std.Dev.
# cntry (Intercept) 0.0344 0.186
# Residual 0.6774 0.823
# Number of obs: 30064, groups: cntry, 20
#
# Fixed effects:
# Estimate Std. Error t value
# (Intercept) 3.529811 0.250494 14.09
# gndrMale -0.001773 0.009618 -0.18
# agea 0.005507 0.000272 20.23
# edu_attain -0.090369 0.002825 -31.99
# uemplaMarked -0.026586 0.025022 -1.06
# hincfel 0.101535 0.006686 15.19
# dfegcf -0.234470 0.007446 -31.49
# rlgdnmEastern Orthodox 0.000852 0.054319 0.02
# rlgdnmEastern religions -0.092259 0.103340 -0.89
# rlgdnmJewish -0.022459 0.171995 -0.13
# rlgdnmOther Christian denomination -0.077413 0.050117 -1.54
# rlgdnmOther non-Christian religions -0.221022 0.100889 -2.19
# rlgdnmProtestant 0.032939 0.015282 2.16
# rlgdnmRoman Catholic 0.057076 0.013169 4.33
# TV_not_news 0.028947 0.002832 10.22
# muslim_j -0.088233 0.018832 -4.69
# mipex_j -0.010364 0.003905 -2.65
# relstate_j -0.017433 0.014196 -1.23
# claims_j 0.122111 0.096072 1.27
Stata solution -> dont’ peek to early ;) !
. use "./../../assets/ESS7/ESS7e02_2.dta", cl. quietly do "./../../assets/ESS7/ESS7e02_2_formats_unicode.do"
.
. * Micro data preparation
. ** Keep only those born in and with citizenship of country of interview,
. ** Non-Muslims, and drop Israel
. keep if brncntr == 1 & ctzcntr == 1 & rlgdnm != 6 & cntry != "IL"
(6,789 observations deleted)
.
. ** Add atheists to the religion variable
. replace rlgdnm = 0 if rlgblg == 2
(15,411 real changes made)
.
. ** Educational attainment
. gen edu_attain = edulvlb
(93 missing values generated)
. recode edu_attain (0 113 129 = 1) (212 213 221 222 223 229 = 2) (311 312 313
> 321 322 323 = 3) (412 413 421 422 423 = 4) (510 520 = 5) (610 620 = 6) (710 7
> 20 800 = 7) (else = .)
(edu_attain: 33396 changes made)
. label var edu_attain "Education"
.
. ** Watching TV, but not news on TV
. gen TV_not_news = tvtot - tvpol
(1,506 missing values generated)
. label var TV_not_news "TV exposure"
.
. ** Immigrant friends
. sum dfegcf
Variable | Obs Mean Std. Dev. Min Max
-------------+---------------------------------------------------------
dfegcf | 33,278 2.426979 .6874448 1 3
. replace dfegcf = `r(max)' - dfegcf
(33,396 real changes made, 118 to missing)
. lab var dfegcf "Friendships with immigrants"
.
. ** Some more labels
. label define uempla 0 "Not unemployed" 1 "Unemployed", replace
. label variable hincfel "Economic deprivation"
.
. ** Use the following variables for the analysis
. keep cntry almuslv agea gndr edu_attain uempla hincfel rlgblg rlgdnm dfegcf T
> V_not_news
. ** Listwise deletion of missing values
. foreach var in cntry almuslv agea gndr edu_attain uempla hincfel rlgblg rlgdn
> m dfegcf TV_not_news {
2. drop if missing(`var')
3. }
(0 observations deleted)
(1,217 observations deleted)
(54 observations deleted)
(20 observations deleted)
(160 observations deleted)
(0 observations deleted)
(242 observations deleted)
(114 observations deleted)
(65 observations deleted)
(75 observations deleted)
(1,385 observations deleted)
.
. ** Save memory
. compress
variable edu_attain was float now byte
variable TV_not_news was float now byte
(180,384 bytes saved)
. ** Make temporary file and save it
. tempfile micro_data
. save "`micro_data'"
file /var/folders/2_/0vc2cx65781cm0qm9tw72fg80000gn/T//St33330.000001 saved
.
. * Macro data
. use "./../../assets/ESS7/Schlueter.dta", clear
. ** z-standardize macro predictor
. egen z_muslim_j = std(muslim_j)
(20 missing values generated)
.
. * Conjoin Schlueter et al. macro data
. merge 1:m cntry using "`micro_data'", keep(match) nogen
Result # of obs.
-----------------------------------------
not matched 0
matched 30,064
-----------------------------------------
.
. * Model 2
. mixed almuslv i.gndr agea edu_attain i.uempla hincfel dfegcf i.rlgdnm TV_not_
> news muslim_j mipex_j relstate_j claims_j || cntry: , reml
Performing EM optimization:
Performing gradient-based optimization:
Iteration 0: log restricted-likelihood = -36899.45
Iteration 1: log restricted-likelihood = -36899.45
Computing standard errors:
Mixed-effects REML regression Number of obs = 30,064
Group variable: cntry Number of groups = 20
Obs per group:
min = 953
avg = 1,503.2
max = 2,478
Wald chi2(18) = 4600.41
Log restricted-likelihood = -36899.45 Prob > chi2 = 0.0000
------------------------------------------------------------------------------
almuslv | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
gndr |
Female | .0017727 .0096184 0.18 0.854 -.017079 .0206244
agea | .0055072 .0002722 20.23 0.000 .0049737 .0060406
edu_attain | -.0903693 .0028245 -31.99 0.000 -.0959053 -.0848333
|
uempla |
Unemployed | -.026586 .0250219 -1.06 0.288 -.0756281 .0224561
hincfel | .1015351 .0066859 15.19 0.000 .088431 .1146391
dfegcf | -.2344702 .0074458 -31.49 0.000 -.2490638 -.2198767
|
rlgdnm |
Roman Cat.. | .0570756 .0131686 4.33 0.000 .0312655 .0828856
Protestant | .032939 .0152821 2.16 0.031 .0029866 .0628914
Eastern O.. | .0008523 .0543188 0.02 0.987 -.1056106 .1073151
Other Chr.. | -.0774134 .050117 -1.54 0.122 -.1756409 .0208141
Jewish | -.0224594 .1719952 -0.13 0.896 -.3595639 .314645
Eastern r.. | -.0922595 .1033398 -0.89 0.372 -.2948018 .1102828
Other non.. | -.2210224 .1008889 -2.19 0.028 -.4187609 -.0232839
|
TV_not_news | .028947 .0028321 10.22 0.000 .0233962 .0344979
muslim_j | -.0882334 .0188324 -4.69 0.000 -.1251441 -.0513226
mipex_j | -.0103644 .0039052 -2.65 0.008 -.0180185 -.0027104
relstate_j | -.0174327 .0141957 -1.23 0.219 -.0452557 .0103903
claims_j | .1221106 .0960718 1.27 0.204 -.0661867 .3104079
_cons | 3.528038 .2504508 14.09 0.000 3.037164 4.018913
------------------------------------------------------------------------------
------------------------------------------------------------------------------
Random-effects Parameters | Estimate Std. Err. [95% Conf. Interval]
-----------------------------+------------------------------------------------
cntry: Identity |
var(_cons) | .0344445 .0127685 .0166564 .0712291
-----------------------------+------------------------------------------------
var(Residual) | .6774208 .0055284 .6666716 .6883433
------------------------------------------------------------------------------
LR test vs. linear model: chibar2(01) = 935.93 Prob >= chibar2 = 0.0000
.
Add a random slope for “TV exposure”; interpret the results by
drawing how you iagine the random slope to look like by hand.
R solution -> dont’ peek to early ;) !
ESS <- ESS %>%
mutate(
z_almuslv = (almuslv - mean(almuslv, na.rm = TRUE)) / sd(almuslv, na.rm = TRUE),
z_TV_not_news = (TV_not_news - mean(TV_not_news, na.rm = TRUE)) / sd(TV_not_news, na.rm = TRUE))
# Model 4
Model_3 <- lmer(
z_almuslv ~ gndr + agea + edu_attain + uempla + hincfel +
dfegcf + rlgdnm + z_TV_not_news + muslim_j + mipex_j +
relstate_j + claims_j + (1 + z_TV_not_news | cntry), data = ESS
)
summary(Model_3)
# Linear mixed model fit by REML ['lmerMod']
# Formula: z_almuslv ~ gndr + agea + edu_attain + uempla + hincfel + dfegcf +
# rlgdnm + z_TV_not_news + muslim_j + mipex_j + relstate_j + claims_j + (1 + z_TV_not_news | cntry)
# Data: ESS
#
# REML criterion at convergence: 75542
#
# Scaled residuals:
# Min 1Q Median 3Q Max
# -3.561 -0.695 0.008 0.742 3.103
#
# Random effects:
# Groups Name Variance Std.Dev. Corr
# cntry (Intercept) 0.03771 0.1942
# z_TV_not_news 0.00148 0.0384 -0.22
# Residual 0.71730 0.8469
# Number of obs: 30064, groups: cntry, 20
#
# Fixed effects:
# Estimate Std. Error t value
# (Intercept) 0.925346 0.258058 3.59
# gndrMale -0.001809 0.009901 -0.18
# agea 0.005674 0.000281 20.23
# edu_attain -0.092874 0.002912 -31.90
# uemplaMarked -0.027001 0.025754 -1.05
# hincfel 0.104613 0.006885 15.19
# dfegcf -0.240715 0.007669 -31.39
# rlgdnmEastern Orthodox 0.001690 0.055905 0.03
# rlgdnmEastern religions -0.097854 0.106378 -0.92
# rlgdnmJewish -0.020189 0.177027 -0.11
# rlgdnmOther Christian denomination -0.076974 0.051589 -1.49
# rlgdnmOther non-Christian religions -0.226086 0.103827 -2.18
# rlgdnmProtestant 0.034904 0.015736 2.22
# rlgdnmRoman Catholic 0.058991 0.013561 4.35
# z_TV_not_news 0.052573 0.010113 5.20
# muslim_j -0.087053 0.019369 -4.49
# mipex_j -0.010417 0.004029 -2.59
# relstate_j -0.019605 0.014609 -1.34
# claims_j 0.121544 0.098840 1.23
Stata solution -> dont’ peek to early ;) !
. quietly do "./../../assets/ESS7/9-exercise..
. * Z-standardization
. egen z_TV_not_news = std(TV_not_news)
. egen z_almuslv = std(almuslv)
.
. ** Model 3
. mixed z_almuslv i.gndr agea edu_attain i.uempla hincfel dfegcf i.rlgdnm z_TV_
> not_news muslim_j mipex_j relstate_j claims_j || cntry: z_TV_not_news, reml
Performing EM optimization:
Performing gradient-based optimization:
Iteration 0: log restricted-likelihood = -37771.034
Iteration 1: log restricted-likelihood = -37771.034
Computing standard errors:
Mixed-effects REML regression Number of obs = 30,064
Group variable: cntry Number of groups = 20
Obs per group:
min = 953
avg = 1,503.2
max = 2,478
Wald chi2(18) = 4222.24
Log restricted-likelihood = -37771.034 Prob > chi2 = 0.0000
------------------------------------------------------------------------------
z_almuslv | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
gndr |
Female | .0018031 .0099012 0.18 0.856 -.017603 .0212091
agea | .0056736 .0002805 20.22 0.000 .0051237 .0062234
edu_attain | -.0928767 .0029118 -31.90 0.000 -.0985837 -.0871696
|
uempla |
Unemployed | -.02694 .0257548 -1.05 0.296 -.0774185 .0235385
hincfel | .104622 .0068853 15.20 0.000 .0911271 .1181169
dfegcf | -.2407303 .0076692 -31.39 0.000 -.2557615 -.225699
|
rlgdnm |
Roman Cat.. | .0588494 .0135636 4.34 0.000 .0322653 .0854335
Protestant | .0347839 .0157364 2.21 0.027 .0039411 .0656268
Eastern O.. | .0019552 .0559046 0.03 0.972 -.1076158 .1115261
Other Chr.. | -.0770574 .0515893 -1.49 0.135 -.1781706 .0240557
Jewish | -.0196228 .1770271 -0.11 0.912 -.3665896 .3273441
Eastern r.. | -.0978406 .1063807 -0.92 0.358 -.306343 .1106617
Other non.. | -.2260544 .1038275 -2.18 0.029 -.4295526 -.0225562
|
z_TV_not_n~s | .0523385 .0101408 5.16 0.000 .0324629 .0722141
muslim_j | -.0918123 .0197136 -4.66 0.000 -.1304502 -.0531745
mipex_j | -.0106394 .004088 -2.60 0.009 -.0186517 -.002627
relstate_j | -.0166677 .0148585 -1.12 0.262 -.0457899 .0124546
claims_j | .1287731 .1005631 1.28 0.200 -.068327 .3258732
_cons | .9210309 .2619374 3.52 0.000 .407643 1.434419
------------------------------------------------------------------------------
------------------------------------------------------------------------------
Random-effects Parameters | Estimate Std. Err. [95% Conf. Interval]
-----------------------------+------------------------------------------------
cntry: Independent |
var(z_TV_n~s) | .0014856 .0006267 .0006498 .0033963
var(_cons) | .037738 .013992 .0182466 .0780506
-----------------------------+------------------------------------------------
var(Residual) | .717295 .0058555 .7059099 .7288638
------------------------------------------------------------------------------
LR test vs. linear model: chi2(2) = 972.91 Prob > chi2 = 0.0000
Note: LR test is conservative and provided only for reference.
.
Can the % Muslims in a country explain part of the random slope?
Explain your results.
Yes, in line with theat theory, watching TV
increases xenophobia particularly where there is many Muslims
Well, in line with theat theory, watching TV increases
xenophobia particularly where there is many Muslims. But the cross-level
interation is not significant. No, the % Mulsim
coefficient ist still significant. No, the TV exposure
coefficient ist still significant.
R solution -> dont’ peek to early ;) !
# Model 4
Model_4 <- lmer(
z_almuslv ~ gndr + agea + edu_attain + uempla + hincfel +
dfegcf + rlgdnm + z_TV_not_news*z_muslim_j + mipex_j +
relstate_j + claims_j + (1 + z_TV_not_news | cntry), data = ESS
)
summary(Model_4)
# Linear mixed model fit by REML ['lmerMod']
# Formula: z_almuslv ~ gndr + agea + edu_attain + uempla + hincfel + dfegcf +
# rlgdnm + z_TV_not_news * z_muslim_j + mipex_j + relstate_j + claims_j + (1 + z_TV_not_news | cntry)
# Data: ESS
#
# REML criterion at convergence: 75545
#
# Scaled residuals:
# Min 1Q Median 3Q Max
# -3.555 -0.695 0.007 0.743 3.106
#
# Random effects:
# Groups Name Variance Std.Dev. Corr
# cntry (Intercept) 0.03761 0.1939
# z_TV_not_news 0.00134 0.0366 -0.22
# Residual 0.71730 0.8469
# Number of obs: 30064, groups: cntry, 20
#
# Fixed effects:
# Estimate Std. Error t value
# (Intercept) 0.663256 0.268713 2.47
# gndrMale -0.001866 0.009901 -0.19
# agea 0.005677 0.000281 20.24
# edu_attain -0.092806 0.002912 -31.87
# uemplaMarked -0.026783 0.025754 -1.04
# hincfel 0.104567 0.006885 15.19
# dfegcf -0.240643 0.007669 -31.38
# rlgdnmEastern Orthodox 0.001761 0.055905 0.03
# rlgdnmEastern religions -0.098644 0.106378 -0.93
# rlgdnmJewish -0.018014 0.177032 -0.10
# rlgdnmOther Christian denomination -0.076695 0.051589 -1.49
# rlgdnmOther non-Christian religions -0.226189 0.103827 -2.18
# rlgdnmProtestant 0.035037 0.015736 2.23
# rlgdnmRoman Catholic 0.059160 0.013561 4.36
# z_TV_not_news 0.052603 0.009764 5.39
# z_muslim_j -0.231661 0.049479 -4.68
# mipex_j -0.010465 0.004028 -2.60
# relstate_j -0.019461 0.014604 -1.33
# claims_j 0.122244 0.098800 1.24
# z_TV_not_news:z_muslim_j 0.015001 0.009830 1.53
# optimizer (nloptwrap) convergence code: 0 (OK)
# Model failed to converge with max|grad| = 0.0024852 (tol = 0.002, component 1)
Stata solution -> dont’ peek to early ;) !
. quietly do "./../../assets/ESS7/9-exercise..
. * Z-standardization
. egen z_TV_not_news = std(TV_not_news)
. egen z_almuslv = std(almuslv)
.
. ** Model 4
. mixed z_almuslv i.gndr agea edu_attain i.uempla hincfel dfegcf i.rlgdnm c.z_T
> V_not_news##c.muslim_j mipex_j relstate_j claims_j || cntry: z_TV_not_news, r
> eml
Performing EM optimization:
Performing gradient-based optimization:
Iteration 0: log restricted-likelihood = -37774.531
Iteration 1: log restricted-likelihood = -37774.531
Computing standard errors:
Mixed-effects REML regression Number of obs = 30,064
Group variable: cntry Number of groups = 20
Obs per group:
min = 953
avg = 1,503.2
max = 2,478
Wald chi2(19) = 4233.82
Log restricted-likelihood = -37774.531 Prob > chi2 = 0.0000
------------------------------------------------------------------------------
z_almuslv | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
gndr |
Female | .0018604 .0099012 0.19 0.851 -.0175456 .0212664
agea | .0056768 .0002805 20.24 0.000 .0051269 .0062266
edu_attain | -.092809 .0029121 -31.87 0.000 -.0985166 -.0871015
|
uempla |
Unemployed | -.0267224 .0257551 -1.04 0.299 -.0772014 .0237566
hincfel | .1045781 .0068852 15.19 0.000 .0910833 .1180728
dfegcf | -.2406598 .0076692 -31.38 0.000 -.2556912 -.2256284
|
rlgdnm |
Roman Cat.. | .0590203 .0135638 4.35 0.000 .0324358 .0856048
Protestant | .0349144 .0157365 2.22 0.027 .0040715 .0657573
Eastern O.. | .0020329 .0559045 0.04 0.971 -.1075379 .1116036
Other Chr.. | -.0767766 .0515894 -1.49 0.137 -.1778899 .0243367
Jewish | -.0174538 .1770322 -0.10 0.921 -.3644304 .3295229
Eastern r.. | -.0986338 .1063808 -0.93 0.354 -.3071363 .1098687
Other non.. | -.2261565 .1038273 -2.18 0.029 -.4296543 -.0226588
|
z_TV_not_n~s | .0344196 .0153439 2.24 0.025 .004346 .0644931
muslim_j | -.0917756 .0197202 -4.65 0.000 -.1304265 -.0531247
|
c. |
z_TV_not_n~s#|
c.muslim_j | .0059433 .0039112 1.52 0.129 -.0017225 .0136091
|
mipex_j | -.0106688 .0040894 -2.61 0.009 -.0186839 -.0026537
relstate_j | -.0165834 .0148638 -1.12 0.265 -.0457158 .0125491
claims_j | .1291909 .1005977 1.28 0.199 -.0679769 .3263588
_cons | .9212703 .2620248 3.52 0.000 .4077111 1.434829
------------------------------------------------------------------------------
------------------------------------------------------------------------------
Random-effects Parameters | Estimate Std. Err. [95% Conf. Interval]
-----------------------------+------------------------------------------------
cntry: Independent |
var(z_TV_n~s) | .0013476 .0006016 .0005618 .0032326
var(_cons) | .0377645 .0140022 .018259 .0781071
-----------------------------+------------------------------------------------
var(Residual) | .7172967 .0058555 .7059115 .7288656
------------------------------------------------------------------------------
LR test vs. linear model: chi2(2) = 972.43 Prob > chi2 = 0.0000
Note: LR test is conservative and provided only for reference.
.
References
Schlueter, Elmar, Anu Masso, and Eldad Davidov. 2019.
“What
Factors Explain Anti-Muslim Prejudice? An
Assessment of the Effects of Muslim Population Size,
Institutional Characteristics and Immigration-Related Media
Claims.” Journal of Ethnic and Migration Studies
0(0):1–16. doi:
10.1080/1369183X.2018.1550160 .