Module 11: Cross-Classified Linear Mixed Models
Jake Reynolds - November 11, 2020
Part 1: Running Separate and Combined Null Models
1.Using reading score (p7read) as the outcome, run and interpret three separate null models, one for school (schid), one for neighborhood (neighid), and one additive null model combining school and neighborhood.
The AIC and BIC for the three null models are as follows: school = 18558.8 and 18576.0, neighborhood= 18606.3 and 18623, and combining school and neighborhood = 18524.9 and 18547.8. The AIC and BIC both improve when you combine school and neightborhood. The third null model is the best between the three.
Part 2: Add Predictors!
- sing the model-building approach of your choice, add a series of predictors to the model: at the student level, add momed and daded; at the school level, add dadocc; and at the neighborhood level, add deprive. Interpret the results and use AIC and BIC to select the best-fitting model. Based on AIC and BIC, the best model below is Model 5 ( p7read ~ momed.fac + daded.fac + deprive + dadocc + momed.fac:daded.fac + (1|schid) + (1|neighid)) with an AIC and BIC of 18247.3 and 18299.0, respectively.
Part 3: Test an Interaction Effect
- Test a model that includes an interaction between momed and daded. Interpret the results (you don’t need to visualize). What research question is this interaction testing?
The model is testing to see if the value of one variable impacts the estimate of the other variable. In this model it appears that there is an interaction between momed.fac and daded.fac. The test shows that this interaction is significant.
Load in Our MVP Packages
suppressPackageStartupMessages(library(tidyverse))
library(lme4)
library(Hmisc)Part 1: Loading in the Data, and then Recoding Variables
Load in the Data
scotland <- haven::read_dta("scotland.dta")
glimpse(scotland)Rows: 2,310
Columns: 11
$ neighid <dbl> 26, 26, 26, 26, 26, 27, 29, 29, 29, 29, 29, 29, 29, 29, 29, 30, 30, 31, 31, 3...
$ schid <dbl> 20, 20, 20, 20, 20, 20, 18, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, 18, 20, 1...
$ attain <dbl> 1.5177, -1.3276, 0.5610, 1.5177, -1.3276, -0.1325, 0.0293, 0.5610, 0.1582, -0...
$ p7vrq <dbl> 17.972, -10.028, 2.972, 1.972, -1.028, 3.972, 8.972, -0.028, 4.972, -8.028, -...
$ p7read <dbl> 17.134, -27.866, 6.134, 11.134, -0.866, -0.866, 6.134, -5.866, 11.134, -13.86...
$ dadocc <dbl> 16.196, -3.454, 2.316, -9.094, -3.454, -3.454, 16.196, -3.454, -11.494, -3.45...
$ dadunemp <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,...
$ daded <dbl> 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0,...
$ momed <dbl> 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0,...
$ male <dbl> 1, 0, 0, 0, 1, 0, 1, 0, 1, 0, 0, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0, 1, 0, 1,...
$ deprive <dbl> -0.551, -0.551, -0.551, -0.551, -0.551, 0.147, -0.083, -0.083, -0.083, -0.083...The describe function from Hmisc is a great codebook…
describe(scotland)scotland
11 Variables 2310 Observations
---------------------------------------------------------------------------------------------------
neighid Format:%12.0g
n missing distinct Info Mean Gmd .05 .10 .25 .50 .75
2310 0 524 1 495.3 305.5 61.0 143.9 240.2 530.0 707.0
.90 .95
808.0 861.0
lowest : 26 27 29 30 31, highest: 1092 1095 1096 1097 1098
---------------------------------------------------------------------------------------------------
schid Format:%12.0g
n missing distinct Info Mean Gmd .05 .10 .25 .50 .75
2310 0 17 0.995 10.01 7.174 0 2 5 9 16
.90 .95
19 20
lowest : 0 1 2 3 5, highest: 16 17 18 19 20
Value 0 1 2 3 5 6 7 8 9 10 13 15 16 17
Frequency 146 22 146 159 155 101 286 112 136 133 92 190 111 154
Proportion 0.063 0.010 0.063 0.069 0.067 0.044 0.124 0.048 0.059 0.058 0.040 0.082 0.048 0.067
Value 18 19 20
Frequency 91 102 174
Proportion 0.039 0.044 0.075
---------------------------------------------------------------------------------------------------
attain Format:%12.0g
n missing distinct Info Mean Gmd .05 .10 .25 .50 .75
2310 0 14 0.987 0.0934 1.124 -1.3276 -1.3276 -0.5812 0.1582 0.7350
.90 .95
1.5177 1.5177
lowest : -1.3276 -0.5812 -0.3600 -0.1325 0.0293, highest: 0.7350 0.9127 1.1405 1.5177 2.4151
---------------------------------------------------------------------------------------------------
p7vrq Format:%12.0g
n missing distinct Info Mean Gmd .05 .10 .25 .50 .75
2310 0 68 0.999 0.5058 11.95 -17.028 -13.028 -7.028 -0.028 7.972
.90 .95
13.972 16.972
lowest : -27.028 -26.028 -25.028 -24.028 -23.028, highest: 35.972 36.972 40.972 41.972 42.972
---------------------------------------------------------------------------------------------------
p7read Format:%12.0g
n missing distinct Info Mean Gmd .05 .10 .25 .50 .75
2310 0 61 1 -0.04435 15.8 -23.866 -17.866 -9.866 -0.866 9.134
.90 .95
19.134 24.134
lowest : -31.866 -30.866 -29.866 -28.866 -27.866, highest: 24.134 25.134 26.134 27.134 28.134
---------------------------------------------------------------------------------------------------
dadocc Format:%12.0g
n missing distinct Info Mean Gmd
2310 0 7 0.933 -0.4642 12.33
lowest : -23.454 -11.494 -9.094 -3.454 2.316, highest: -9.094 -3.454 2.316 16.196 29.226
Value -23.454 -11.494 -9.094 -3.454 2.316 16.196 29.226
Frequency 91 285 303 884 242 397 108
Proportion 0.039 0.123 0.131 0.383 0.105 0.172 0.047
---------------------------------------------------------------------------------------------------
dadunemp Format:%12.0g
n missing distinct Info Sum Mean Gmd
2310 0 2 0.292 252 0.1091 0.1945
---------------------------------------------------------------------------------------------------
daded Format:%12.0g
n missing distinct Info Sum Mean Gmd
2310 0 2 0.507 497 0.2152 0.3379
---------------------------------------------------------------------------------------------------
momed Format:%12.0g
n missing distinct Info Sum Mean Gmd
2310 0 2 0.56 574 0.2485 0.3736
---------------------------------------------------------------------------------------------------
male Format:%12.0g
n missing distinct Info Sum Mean Gmd
2310 0 2 0.749 1109 0.4801 0.4994
---------------------------------------------------------------------------------------------------
deprive Format:%12.0g
n missing distinct Info Mean Gmd .05 .10 .25 .50 .75
2310 0 458 1 0.02167 0.6664 -0.8250 -0.6930 -0.3970 -0.0620 0.2957
.90 .95
0.8410 1.1400
lowest : -1.082 -1.048 -1.030 -0.983 -0.975, highest: 2.330 2.419 2.438 2.498 2.959
---------------------------------------------------------------------------------------------------Part One: Data Exploration and Coding/ Labelling
Start with variable labelling using the label function
label(scotland$neighid)="Neighborhood ID"
label(scotland$schid)="School ID"
label(scotland$attain)="Student Measure of Educational Attainment"
label(scotland$p7vrq)="Primary 7 Verbal Reasoning Quotient"
label(scotland$p7read)="Primary 7 Reading Test Scores"
label(scotland$dadocc)="School Mean of Dad's Occupation Score on Hope-Goldthorpe Scale"
label(scotland$dadunemp)="Dad Currently Unemployed"
label(scotland$daded)="Dad School After Age 15"
label(scotland$momed)="Mom School After Age 15"
label(scotland$male)="Male"
label(scotland$deprive)="Neighborhood Deprivation Score (Poverty, Health, and Housing)"
scot.clean <- scotland %>%
mutate(.,
dadunemp.fac = as_factor(dadunemp),
daded.fac = as_factor(daded),
momed.fac = as_factor(momed),
male.fac = as_factor(male))
levels(scot.clean$dadunemp.fac) = c("No","Yes")
levels(scot.clean$daded.fac) = c("No","Yes")
levels(scot.clean$momed.fac) = c("No","Yes")
levels(scot.clean$male.fac) = c("No","Yes")
glimpse(scot.clean)Rows: 2,310
Columns: 15
$ neighid <labelled> 26, 26, 26, 26, 26, 27, 29, 29, 29, 29, 29, 29, 29, 29, 29, 30, 30, ...
$ schid <labelled> 20, 20, 20, 20, 20, 20, 18, 20, 20, 20, 20, 20, 20, 20, 20, 20, 20, ...
$ attain <labelled> 1.5177, -1.3276, 0.5610, 1.5177, -1.3276, -0.1325, 0.0293, 0.5610, 0...
$ p7vrq <labelled> 17.972, -10.028, 2.972, 1.972, -1.028, 3.972, 8.972, -0.028, 4.972, ...
$ p7read <labelled> 17.134, -27.866, 6.134, 11.134, -0.866, -0.866, 6.134, -5.866, 11.13...
$ dadocc <labelled> 16.196, -3.454, 2.316, -9.094, -3.454, -3.454, 16.196, -3.454, -11.4...
$ dadunemp <labelled> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,...
$ daded <labelled> 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0,...
$ momed <labelled> 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0,...
$ male <labelled> 1, 0, 0, 0, 1, 0, 1, 0, 1, 0, 0, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0,...
$ deprive <labelled> -0.551, -0.551, -0.551, -0.551, -0.551, 0.147, -0.083, -0.083, -0.08...
$ dadunemp.fac <fct> No, No, No, No, No, No, No, No, No, No, No, No, No, No, No, No, No, No, N...
$ daded.fac <fct> No, No, No, No, No, No, Yes, No, No, No, No, No, No, No, No, Yes, No, Yes...
$ momed.fac <fct> No, No, No, No, No, Yes, Yes, No, No, No, No, No, No, No, No, No, No, Yes...
$ male.fac <fct> Yes, No, No, No, Yes, No, Yes, No, Yes, No, No, No, Yes, No, Yes, Yes, Ye...describe(scot.clean)scot.clean
15 Variables 2310 Observations
---------------------------------------------------------------------------------------------------
neighid : Neighborhood ID Format:%12.0g
n missing distinct Info Mean Gmd .05 .10 .25 .50 .75
2310 0 524 1 495.3 305.5 61.0 143.9 240.2 530.0 707.0
.90 .95
808.0 861.0
lowest : 26 27 29 30 31, highest: 1092 1095 1096 1097 1098
---------------------------------------------------------------------------------------------------
schid : School ID Format:%12.0g
n missing distinct Info Mean Gmd .05 .10 .25 .50 .75
2310 0 17 0.995 10.01 7.174 0 2 5 9 16
.90 .95
19 20
lowest : 0 1 2 3 5, highest: 16 17 18 19 20
Value 0 1 2 3 5 6 7 8 9 10 13 15 16 17
Frequency 146 22 146 159 155 101 286 112 136 133 92 190 111 154
Proportion 0.063 0.010 0.063 0.069 0.067 0.044 0.124 0.048 0.059 0.058 0.040 0.082 0.048 0.067
Value 18 19 20
Frequency 91 102 174
Proportion 0.039 0.044 0.075
---------------------------------------------------------------------------------------------------
attain : Student Measure of Educational Attainment Format:%12.0g
n missing distinct Info Mean Gmd .05 .10 .25 .50 .75
2310 0 14 0.987 0.0934 1.124 -1.3276 -1.3276 -0.5812 0.1582 0.7350
.90 .95
1.5177 1.5177
lowest : -1.3276 -0.5812 -0.3600 -0.1325 0.0293, highest: 0.7350 0.9127 1.1405 1.5177 2.4151
---------------------------------------------------------------------------------------------------
p7vrq : Primary 7 Verbal Reasoning Quotient Format:%12.0g
n missing distinct Info Mean Gmd .05 .10 .25 .50 .75
2310 0 68 0.999 0.5058 11.95 -17.028 -13.028 -7.028 -0.028 7.972
.90 .95
13.972 16.972
lowest : -27.028 -26.028 -25.028 -24.028 -23.028, highest: 35.972 36.972 40.972 41.972 42.972
---------------------------------------------------------------------------------------------------
p7read : Primary 7 Reading Test Scores Format:%12.0g
n missing distinct Info Mean Gmd .05 .10 .25 .50 .75
2310 0 61 1 -0.04435 15.8 -23.866 -17.866 -9.866 -0.866 9.134
.90 .95
19.134 24.134
lowest : -31.866 -30.866 -29.866 -28.866 -27.866, highest: 24.134 25.134 26.134 27.134 28.134
---------------------------------------------------------------------------------------------------
dadocc : School Mean of Dad's Occupation Score on Hope-Goldthorpe Scale Format:%12.0g
n missing distinct Info Mean Gmd
2310 0 7 0.933 -0.4642 12.33
lowest : -23.454 -11.494 -9.094 -3.454 2.316, highest: -9.094 -3.454 2.316 16.196 29.226
Value -23.454 -11.494 -9.094 -3.454 2.316 16.196 29.226
Frequency 91 285 303 884 242 397 108
Proportion 0.039 0.123 0.131 0.383 0.105 0.172 0.047
---------------------------------------------------------------------------------------------------
dadunemp : Dad Currently Unemployed Format:%12.0g
n missing distinct Info Sum Mean Gmd
2310 0 2 0.292 252 0.1091 0.1945
---------------------------------------------------------------------------------------------------
daded : Dad School After Age 15 Format:%12.0g
n missing distinct Info Sum Mean Gmd
2310 0 2 0.507 497 0.2152 0.3379
---------------------------------------------------------------------------------------------------
momed : Mom School After Age 15 Format:%12.0g
n missing distinct Info Sum Mean Gmd
2310 0 2 0.56 574 0.2485 0.3736
---------------------------------------------------------------------------------------------------
male : Male Format:%12.0g
n missing distinct Info Sum Mean Gmd
2310 0 2 0.749 1109 0.4801 0.4994
---------------------------------------------------------------------------------------------------
deprive : Neighborhood Deprivation Score (Poverty, Health, and Housing) Format:%12.0g
n missing distinct Info Mean Gmd .05 .10 .25 .50 .75
2310 0 458 1 0.02167 0.6664 -0.8250 -0.6930 -0.3970 -0.0620 0.2957
.90 .95
0.8410 1.1400
lowest : -1.082 -1.048 -1.030 -0.983 -0.975, highest: 2.330 2.419 2.438 2.498 2.959
---------------------------------------------------------------------------------------------------
dadunemp.fac : Dad Currently Unemployed
n missing distinct
2310 0 2
Value No Yes
Frequency 2058 252
Proportion 0.891 0.109
---------------------------------------------------------------------------------------------------
daded.fac : Dad School After Age 15
n missing distinct
2310 0 2
Value No Yes
Frequency 1813 497
Proportion 0.785 0.215
---------------------------------------------------------------------------------------------------
momed.fac : Mom School After Age 15
n missing distinct
2310 0 2
Value No Yes
Frequency 1736 574
Proportion 0.752 0.248
---------------------------------------------------------------------------------------------------
male.fac : Male
n missing distinct
2310 0 2
Value No Yes
Frequency 1201 1109
Proportion 0.52 0.48
---------------------------------------------------------------------------------------------------Explore our Outcome, attain
ggplot(data = scot.clean, mapping = aes(x = attain)) +
geom_histogram(bins = 20) +
labs(title = "Distribution of Educational Attainment Scores",
x = "Standardized Educational Attainment Score") +
theme_minimal()
Explore a few other predictors
ggplot(data = scot.clean, mapping = aes(x = dadocc)) +
geom_histogram(bins = 30) +
labs(title = "Distribution of (School Mean) of Dad's Occupation Scores",
x = "Standardized Educational Attainment Score") +
theme_minimal()
ggplot(data = scot.clean, mapping = aes(x = deprive)) +
geom_histogram(bins = 50) +
labs(title = "Distribution of Neighborhood Context Scores",
x = "Standardized Neighborhood Deprivation Score") +
theme_minimal()
Part Two: Modeling
First, run two separate models with school and neighborhood and the L2 clusters…
model.null.school <- lmer(p7read ~ (1|schid), REML = FALSE, data = scot.clean)
summary(model.null.school)Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's method [lmerModLmerTest]
Formula: p7read ~ (1 | schid)
Data: scot.clean
AIC BIC logLik deviance df.resid
18558.8 18576.0 -9276.4 18552.8 2307
Scaled residuals:
Min 1Q Median 3Q Max
-2.7611 -0.6895 -0.0027 0.6854 2.5664
Random effects:
Groups Name Variance Std.Dev.
schid (Intercept) 18.75 4.33
Residual 176.64 13.29
Number of obs: 2310, groups: schid, 17
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) -0.2912 1.0951 16.1869 -0.266 0.794Calculate School ICC
icc.school <- 18.75/(18.75+176.64)
icc.school[1] 0.09596192model.null.neigh <- lmer(p7read ~ (1|neighid), REML = FALSE, data = scot.clean)
summary(model.null.neigh)Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's method [lmerModLmerTest]
Formula: p7read ~ (1 | neighid)
Data: scot.clean
AIC BIC logLik deviance df.resid
18606.3 18623.5 -9300.1 18600.3 2307
Scaled residuals:
Min 1Q Median 3Q Max
-2.65885 -0.67581 -0.01035 0.67042 2.51256
Random effects:
Groups Name Variance Std.Dev.
neighid (Intercept) 32.89 5.735
Residual 160.37 12.664
Number of obs: 2310, groups: neighid, 524
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) -0.2407 0.3816 454.9417 -0.631 0.529Calculate Neighborhood ICC
icc.neigh <- 32.89/(32.89+160.37)
icc.neigh[1] 0.1701852Now, a null additive cross-classified model…
model.null.crossed <- lmer(p7read ~ (1|schid) + (1|neighid), REML = FALSE, data = scot.clean)
summary(model.null.crossed)Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's method [lmerModLmerTest]
Formula: p7read ~ (1 | schid) + (1 | neighid)
Data: scot.clean
AIC BIC logLik deviance df.resid
18524.9 18547.8 -9258.4 18516.9 2306
Scaled residuals:
Min 1Q Median 3Q Max
-2.61239 -0.66099 0.01158 0.66532 2.55698
Random effects:
Groups Name Variance Std.Dev.
neighid (Intercept) 17.21 4.148
schid (Intercept) 16.59 4.074
Residual 160.36 12.663
Number of obs: 2310, groups: neighid, 524; schid, 17
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) -0.3701 1.0522 16.2037 -0.352 0.73Calculate Neighborhood ICC in Crossed Model
icc.neigh.crossed <- 17.21/(17.21+16.59+160.36)
icc.neigh.crossed[1] 0.08863824Calculate School ICC in Crossed Model
icc.school.crossed <- 16.59/(16.59+17.21+160.36)
icc.school.crossed[1] 0.08544499Add an L1 Predictor, Male…
model.1 <- lmer(p7read ~ momed.fac + (1|schid) + (1|neighid), REML = FALSE, data = scot.clean)
summary(model.1)Linear mixed model fit by maximum likelihood . t-tests use
Satterthwaite's method [lmerModLmerTest]
Formula: p7read ~ momed.fac + (1 | schid) + (1 | neighid)
Data: scot.clean
AIC BIC logLik deviance df.resid
18472.8 18501.5 -9231.4 18462.8 2305
Scaled residuals:
Min 1Q Median 3Q Max
-2.81338 -0.67650 0.01877 0.65461 2.60452
Random effects:
Groups Name Variance Std.Dev.
neighid (Intercept) 15.17 3.895
schid (Intercept) 14.54 3.813
Residual 157.87 12.564
Number of obs: 2310, groups: neighid, 524; schid, 17
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) -1.5528 1.0027 16.9863 -1.549 0.14
momed.facYes 4.7026 0.6339 2286.4349 7.418 1.67e-13
(Intercept)
momed.facYes ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
momed.facYs -0.161L1 Plus School-Level Predictor Only (dadocc)…
model.2 <- lmer(p7read ~ momed.fac + daded.fac + (1|schid) + (1|neighid), REML = FALSE, data = scot.clean)
summary(model.2)Linear mixed model fit by maximum likelihood . t-tests use
Satterthwaite's method [lmerModLmerTest]
Formula:
p7read ~ momed.fac + daded.fac + (1 | schid) + (1 | neighid)
Data: scot.clean
AIC BIC logLik deviance df.resid
18424.4 18458.9 -9206.2 18412.4 2304
Scaled residuals:
Min 1Q Median 3Q Max
-2.99961 -0.66409 0.01792 0.66336 2.63998
Random effects:
Groups Name Variance Std.Dev.
neighid (Intercept) 13.24 3.639
schid (Intercept) 13.06 3.614
Residual 155.66 12.476
Number of obs: 2310, groups: neighid, 524; schid, 17
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) -2.1132 0.9580 17.2968 -2.206 0.041199
momed.facYes 2.4043 0.7069 2272.7606 3.401 0.000683
daded.facYes 5.3704 0.7504 2286.9432 7.157 1.11e-12
(Intercept) *
momed.facYes ***
daded.facYes ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) mmd.fY
momed.facYs -0.110
daded.facYs -0.083 -0.460L1 Plus Neighborhood-Level Predictors Only (deprive)…
model.3 <- lmer(p7read ~ momed.fac + daded.fac + deprive + (1|schid) + (1|neighid), REML = FALSE, data = scot.clean)
summary(model.3)Linear mixed model fit by maximum likelihood .
t-tests use Satterthwaite's method [lmerModLmerTest]
Formula:
p7read ~ momed.fac + daded.fac + deprive + (1 | schid) + (1 |
neighid)
Data: scot.clean
AIC BIC logLik deviance df.resid
18344.5 18384.7 -9165.3 18330.5 2303
Scaled residuals:
Min 1Q Median 3Q Max
-3.1962 -0.6647 0.0125 0.6529 3.2538
Random effects:
Groups Name Variance Std.Dev.
neighid (Intercept) 5.699 2.387
schid (Intercept) 7.762 2.786
Residual 156.133 12.495
Number of obs: 2310, groups: neighid, 524; schid, 17
Fixed effects:
Estimate Std. Error df t value
(Intercept) -1.7545 0.7651 16.9213 -2.293
momed.facYes 2.1601 0.6980 2292.3330 3.095
daded.facYes 4.9723 0.7419 2300.8818 6.702
deprive -4.7584 0.4970 394.3489 -9.574
Pr(>|t|)
(Intercept) 0.03492 *
momed.facYes 0.00199 **
daded.facYes 2.58e-11 ***
deprive < 2e-16 ***
---
Signif. codes:
0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) mmd.fY ddd.fY
momed.facYs -0.137
daded.facYs -0.105 -0.457
deprive -0.048 0.047 0.086L1 Plus Both School and Neighborhood Predictors…
model.4 <- lmer(p7read ~ momed.fac + momed.fac + deprive + dadocc + (1|schid) + (1|neighid), REML = FALSE, data = scot.clean)
summary(model.4)Linear mixed model fit by maximum likelihood .
t-tests use Satterthwaite's method [lmerModLmerTest]
Formula:
p7read ~ momed.fac + momed.fac + deprive + dadocc + (1 | schid) +
(1 | neighid)
Data: scot.clean
AIC BIC logLik deviance df.resid
18271.5 18311.7 -9128.7 18257.5 2303
Scaled residuals:
Min 1Q Median 3Q Max
-3.2922 -0.6614 0.0065 0.6538 3.1382
Random effects:
Groups Name Variance Std.Dev.
neighid (Intercept) 4.027 2.007
schid (Intercept) 6.620 2.573
Residual 152.662 12.356
Number of obs: 2310, groups: neighid, 524; schid, 17
Fixed effects:
Estimate Std. Error df t value
(Intercept) -0.76145 0.71112 17.07699 -1.071
momed.facYes 2.96492 0.62353 2296.61991 4.755
deprive -3.92655 0.48878 414.83943 -8.033
dadocc 0.26231 0.02379 2282.96420 11.025
Pr(>|t|)
(Intercept) 0.299
momed.facYes 2.11e-06 ***
deprive 9.93e-15 ***
dadocc < 2e-16 ***
---
Signif. codes:
0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) mmd.fY depriv
momed.facYs -0.227
deprive -0.028 0.053
dadocc 0.056 -0.197 0.211Finally, test a model with school-neighborhood interaction…
model.5 <- lmer(p7read ~ momed.fac + daded.fac + deprive + dadocc + momed.fac:daded.fac + (1|schid) + (1|neighid), REML = FALSE, data = scot.clean)
summary(model.5)Linear mixed model fit by maximum likelihood .
t-tests use Satterthwaite's method [lmerModLmerTest]
Formula:
p7read ~ momed.fac + daded.fac + deprive + dadocc + momed.fac:daded.fac +
(1 | schid) + (1 | neighid)
Data: scot.clean
AIC BIC logLik deviance df.resid
18247.3 18299.0 -9114.6 18229.3 2301
Scaled residuals:
Min 1Q Median 3Q Max
-3.3143 -0.6715 0.0026 0.6539 3.1707
Random effects:
Groups Name Variance Std.Dev.
neighid (Intercept) 3.513 1.874
schid (Intercept) 6.506 2.551
Residual 151.222 12.297
Number of obs: 2310, groups: neighid, 524; schid, 17
Fixed effects:
Estimate Std. Error
(Intercept) -1.37684 0.71429
momed.facYes 3.12965 0.84626
daded.facYes 5.34180 1.01514
deprive -3.75830 0.48407
dadocc 0.23530 0.02432
momed.facYes:daded.facYes -4.24203 1.43861
df t value Pr(>|t|)
(Intercept) 18.13613 -1.928 0.069719
momed.facYes 2292.70713 3.698 0.000222
daded.facYes 2296.98834 5.262 1.56e-07
deprive 410.53776 -7.764 6.62e-14
dadocc 2284.28187 9.676 < 2e-16
momed.facYes:daded.facYes 2272.83105 -2.949 0.003224
(Intercept) .
momed.facYes ***
daded.facYes ***
deprive ***
dadocc ***
momed.facYes:daded.facYes **
---
Signif. codes:
0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) mmd.fY ddd.fY depriv dadocc
momed.facYs -0.181
daded.facYs -0.163 0.143
deprive -0.039 0.059 0.067
dadocc 0.087 -0.079 -0.192 0.195
mmd.fcYs:.Y 0.103 -0.587 -0.678 -0.057 0.032Step 6: Compare Original “Final” Model, Robust Model, and Trimmed Model for Differences (Sensitivity Analysis)
library(modelsummary)
library(broom.mixed)
models <- list(model.null.neigh, model.null.school, model.null.crossed, model.1, model.2, model.3, model.4, model.5)
modelsummary(models, output = "markdown")| Model 1 | Model 2 | Model 3 | Model 4 | Model 5 | Model 6 | Model 7 | Model 8 | |
|---|---|---|---|---|---|---|---|---|
| (Intercept) | -0.241 | -0.291 | -0.370 | -1.553 | -2.113 | -1.754 | -0.761 | -1.377 |
| (0.382) | (1.095) | (1.052) | (1.003) | (0.958) | (0.765) | (0.711) | (0.714) | |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.639 | 2.387 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.639 | 2.387 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.639 | 2.387 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.639 | 2.387 | 2.573 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.639 | 2.786 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.639 | 2.786 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.639 | 2.786 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.639 | 2.786 | 2.573 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.614 | 2.387 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.614 | 2.387 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.614 | 2.387 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.614 | 2.387 | 2.573 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.614 | 2.786 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.614 | 2.786 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.614 | 2.786 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.895 | 3.614 | 2.786 | 2.573 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.639 | 2.387 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.639 | 2.387 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.639 | 2.387 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.639 | 2.387 | 2.573 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.639 | 2.786 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.639 | 2.786 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.639 | 2.786 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.639 | 2.786 | 2.573 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.614 | 2.387 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.614 | 2.387 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.614 | 2.387 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.614 | 2.387 | 2.573 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.614 | 2.786 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.614 | 2.786 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.614 | 2.786 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.148 | 3.813 | 3.614 | 2.786 | 2.573 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.639 | 2.387 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.639 | 2.387 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.639 | 2.387 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.639 | 2.387 | 2.573 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.639 | 2.786 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.639 | 2.786 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.639 | 2.786 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.639 | 2.786 | 2.573 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.614 | 2.387 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.614 | 2.387 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.614 | 2.387 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.614 | 2.387 | 2.573 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.614 | 2.786 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.614 | 2.786 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.614 | 2.786 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.895 | 3.614 | 2.786 | 2.573 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.639 | 2.387 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.639 | 2.387 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.639 | 2.387 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.639 | 2.387 | 2.573 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.639 | 2.786 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.639 | 2.786 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.639 | 2.786 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.639 | 2.786 | 2.573 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.614 | 2.387 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.614 | 2.387 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.614 | 2.387 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.614 | 2.387 | 2.573 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.614 | 2.786 | 2.007 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.614 | 2.786 | 2.007 | 2.551 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.614 | 2.786 | 2.573 | 1.874 |
| sd__(Intercept) | 5.735 | 4.330 | 4.074 | 3.813 | 3.614 | 2.786 | 2.573 | 2.551 |
| sd__Observation | 12.664 | 13.291 | 12.663 | 12.564 | 12.476 | 12.495 | 12.356 | 12.297 |
| momed.facYes | 4.703 | 2.404 | 2.160 | 2.965 | 3.130 | |||
| (0.634) | (0.707) | (0.698) | (0.624) | (0.846) | ||||
| daded.facYes | 5.370 | 4.972 | 5.342 | |||||
| (0.750) | (0.742) | (1.015) | ||||||
| deprive | -4.758 | -3.927 | -3.758 | |||||
| (0.497) | (0.489) | (0.484) | ||||||
| dadocc | 0.262 | 0.235 | ||||||
| (0.024) | (0.024) | |||||||
| momed.facYes × daded.facYes | -4.242 | |||||||
| (1.439) | ||||||||
| AIC | 18606.3 | 18558.8 | 18524.9 | 18472.8 | 18424.4 | 18344.5 | 18271.5 | 18247.3 |
| BIC | 18623.5 | 18576.0 | 18547.8 | 18501.5 | 18458.9 | 18384.7 | 18311.7 | 18299.0 |
| Log.Lik. | -9300.147 | -9276.398 | -9258.426 | -9231.392 | -9206.212 | -9165.253 | -9128.727 | -9114.635 |
LS0tDQp0aXRsZTogJ01vZHVsZSAxMTogQ3Jvc3MtQ2xhc3NpZmllZCBMaW5lYXIgTWl4ZWQgTW9kZWxzJw0KYXV0aG9yOiAnSmFrZSBSZXlub2xkcyAtIE5vdmVtYmVyIDExLCAyMDIwJw0Kb3V0cHV0OiBodG1sX25vdGVib29rDQotLS0NClBhcnQgMTogUnVubmluZyBTZXBhcmF0ZSBhbmQgQ29tYmluZWQgTnVsbCBNb2RlbHMNCg0KMS5Vc2luZyByZWFkaW5nIHNjb3JlIChwN3JlYWQpIGFzIHRoZSBvdXRjb21lLCBydW4gYW5kIGludGVycHJldCB0aHJlZSBzZXBhcmF0ZSBudWxsIG1vZGVscywgb25lIGZvciBzY2hvb2wgKHNjaGlkKSwgb25lIGZvciBuZWlnaGJvcmhvb2QgKG5laWdoaWQpLCBhbmQgb25lIGFkZGl0aXZlIG51bGwgbW9kZWwgY29tYmluaW5nIHNjaG9vbCBhbmQgbmVpZ2hib3Job29kLg0KDQpUaGUgQUlDIGFuZCBCSUMgZm9yIHRoZSB0aHJlZSBudWxsIG1vZGVscyBhcmUgYXMgZm9sbG93czogc2Nob29sID0gMTg1NTguOCBhbmQgMTg1NzYuMCwgbmVpZ2hib3Job29kPSAxODYwNi4zIGFuZCAxODYyMywgYW5kIGNvbWJpbmluZyBzY2hvb2wgYW5kIG5laWdoYm9yaG9vZCA9IDE4NTI0LjkgYW5kIDE4NTQ3LjguIFRoZSBBSUMgYW5kIEJJQyBib3RoIGltcHJvdmUgd2hlbiB5b3UgY29tYmluZSBzY2hvb2wgYW5kIG5laWdodGJvcmhvb2QuIFRoZSB0aGlyZCBudWxsIG1vZGVsIGlzIHRoZSBiZXN0IGJldHdlZW4gdGhlIHRocmVlLiANCg0KUGFydCAyOiBBZGQgUHJlZGljdG9ycyENCg0KMi4gc2luZyB0aGUgbW9kZWwtYnVpbGRpbmcgYXBwcm9hY2ggb2YgeW91ciBjaG9pY2UsIGFkZCBhIHNlcmllcyBvZiBwcmVkaWN0b3JzIHRvIHRoZSBtb2RlbDogYXQgdGhlIHN0dWRlbnQgbGV2ZWwsIGFkZCBtb21lZCBhbmQgZGFkZWQ7IGF0IHRoZSBzY2hvb2wgbGV2ZWwsIGFkZCBkYWRvY2M7IGFuZCBhdCB0aGUgbmVpZ2hib3Job29kIGxldmVsLCBhZGQgZGVwcml2ZS4gSW50ZXJwcmV0IHRoZSByZXN1bHRzIGFuZCB1c2UgQUlDIGFuZCBCSUMgdG8gc2VsZWN0IHRoZSBiZXN0LWZpdHRpbmcgbW9kZWwuDQpCYXNlZCBvbiBBSUMgYW5kIEJJQywgdGhlIGJlc3QgbW9kZWwgYmVsb3cgaXMgTW9kZWwgNSAoIHA3cmVhZCB+IG1vbWVkLmZhYyArIGRhZGVkLmZhYyArIGRlcHJpdmUgKyBkYWRvY2MgKyBtb21lZC5mYWM6ZGFkZWQuZmFjICsgKDF8c2NoaWQpICsgKDF8bmVpZ2hpZCkpIHdpdGggYW4gQUlDIGFuZCBCSUMgb2YgMTgyNDcuMyBhbmQgMTgyOTkuMCwgcmVzcGVjdGl2ZWx5LiANCg0KDQpQYXJ0IDM6IFRlc3QgYW4gSW50ZXJhY3Rpb24gRWZmZWN0DQoNCjMuICBUZXN0IGEgbW9kZWwgdGhhdCBpbmNsdWRlcyBhbiBpbnRlcmFjdGlvbiBiZXR3ZWVuIG1vbWVkIGFuZCBkYWRlZC4gSW50ZXJwcmV0IHRoZSByZXN1bHRzICh5b3UgZG9u4oCZdCBuZWVkIHRvIHZpc3VhbGl6ZSkuIFdoYXQgcmVzZWFyY2ggcXVlc3Rpb24gaXMgdGhpcyBpbnRlcmFjdGlvbiB0ZXN0aW5nPw0KDQpUaGUgbW9kZWwgaXMgdGVzdGluZyB0byBzZWUgaWYgdGhlIHZhbHVlIG9mIG9uZSB2YXJpYWJsZSBpbXBhY3RzIHRoZSBlc3RpbWF0ZSBvZiB0aGUgb3RoZXIgdmFyaWFibGUuIEluIHRoaXMgbW9kZWwgaXQgYXBwZWFycyB0aGF0IHRoZXJlIGlzIGFuIGludGVyYWN0aW9uIGJldHdlZW4gbW9tZWQuZmFjIGFuZCBkYWRlZC5mYWMuIFRoZSB0ZXN0IHNob3dzIHRoYXQgdGhpcyBpbnRlcmFjdGlvbiBpcyBzaWduaWZpY2FudC4gDQoNCg0KDQoNCiMgTG9hZCBpbiBPdXIgTVZQIFBhY2thZ2VzDQpgYGB7cn0NCnN1cHByZXNzUGFja2FnZVN0YXJ0dXBNZXNzYWdlcyhsaWJyYXJ5KHRpZHl2ZXJzZSkpDQpsaWJyYXJ5KGxtZTQpDQpsaWJyYXJ5KEhtaXNjKQ0KYGBgDQoNCiMgUGFydCAxOiBMb2FkaW5nIGluIHRoZSBEYXRhLCBhbmQgdGhlbiBSZWNvZGluZyBWYXJpYWJsZXMNCiMjIExvYWQgaW4gdGhlIERhdGENCmBgYHtyfQ0KDQpzY290bGFuZCA8LSBoYXZlbjo6cmVhZF9kdGEoInNjb3RsYW5kLmR0YSIpDQoNCmdsaW1wc2Uoc2NvdGxhbmQpDQpgYGANCiMjIFRoZSBgZGVzY3JpYmVgIGZ1bmN0aW9uIGZyb20gYEhtaXNjYCBpcyBhIGdyZWF0IGNvZGVib29rLi4uDQpgYGB7cn0NCmRlc2NyaWJlKHNjb3RsYW5kKQ0KYGBgDQojIFBhcnQgT25lOiBEYXRhIEV4cGxvcmF0aW9uIGFuZCBDb2RpbmcvIExhYmVsbGluZw0KIyMgU3RhcnQgd2l0aCB2YXJpYWJsZSBsYWJlbGxpbmcgdXNpbmcgdGhlIGBsYWJlbGAgZnVuY3Rpb24NCmBgYHtyfQ0KbGFiZWwoc2NvdGxhbmQkbmVpZ2hpZCk9Ik5laWdoYm9yaG9vZCBJRCINCmxhYmVsKHNjb3RsYW5kJHNjaGlkKT0iU2Nob29sIElEIg0KbGFiZWwoc2NvdGxhbmQkYXR0YWluKT0iU3R1ZGVudCBNZWFzdXJlIG9mIEVkdWNhdGlvbmFsIEF0dGFpbm1lbnQiDQpsYWJlbChzY290bGFuZCRwN3ZycSk9IlByaW1hcnkgNyBWZXJiYWwgUmVhc29uaW5nIFF1b3RpZW50Ig0KbGFiZWwoc2NvdGxhbmQkcDdyZWFkKT0iUHJpbWFyeSA3IFJlYWRpbmcgVGVzdCBTY29yZXMiDQpsYWJlbChzY290bGFuZCRkYWRvY2MpPSJTY2hvb2wgTWVhbiBvZiBEYWQncyBPY2N1cGF0aW9uIFNjb3JlIG9uIEhvcGUtR29sZHRob3JwZSBTY2FsZSINCmxhYmVsKHNjb3RsYW5kJGRhZHVuZW1wKT0iRGFkIEN1cnJlbnRseSBVbmVtcGxveWVkIg0KbGFiZWwoc2NvdGxhbmQkZGFkZWQpPSJEYWQgU2Nob29sIEFmdGVyIEFnZSAxNSINCmxhYmVsKHNjb3RsYW5kJG1vbWVkKT0iTW9tIFNjaG9vbCBBZnRlciBBZ2UgMTUiDQpsYWJlbChzY290bGFuZCRtYWxlKT0iTWFsZSINCmxhYmVsKHNjb3RsYW5kJGRlcHJpdmUpPSJOZWlnaGJvcmhvb2QgRGVwcml2YXRpb24gU2NvcmUgKFBvdmVydHksIEhlYWx0aCwgYW5kIEhvdXNpbmcpIg0KDQpgYGANCg0KDQpgYGB7cn0NCnNjb3QuY2xlYW4gPC0gc2NvdGxhbmQgJT4lDQogIG11dGF0ZSguLCANCiAgICAgICAgIGRhZHVuZW1wLmZhYyA9IGFzX2ZhY3RvcihkYWR1bmVtcCksDQogICAgICAgICBkYWRlZC5mYWMgPSBhc19mYWN0b3IoZGFkZWQpLA0KICAgICAgICAgbW9tZWQuZmFjID0gYXNfZmFjdG9yKG1vbWVkKSwNCiAgICAgICAgIG1hbGUuZmFjID0gYXNfZmFjdG9yKG1hbGUpKQ0KDQpsZXZlbHMoc2NvdC5jbGVhbiRkYWR1bmVtcC5mYWMpID0gYygiTm8iLCJZZXMiKQ0KbGV2ZWxzKHNjb3QuY2xlYW4kZGFkZWQuZmFjKSA9IGMoIk5vIiwiWWVzIikNCmxldmVscyhzY290LmNsZWFuJG1vbWVkLmZhYykgPSBjKCJObyIsIlllcyIpDQpsZXZlbHMoc2NvdC5jbGVhbiRtYWxlLmZhYykgPSBjKCJObyIsIlllcyIpDQoNCmdsaW1wc2Uoc2NvdC5jbGVhbikNCmBgYA0KDQpgYGB7cn0NCmRlc2NyaWJlKHNjb3QuY2xlYW4pDQpgYGANCg0KIyMgRXhwbG9yZSBvdXIgT3V0Y29tZSwgYGF0dGFpbmANCmBgYHtyfQ0KZ2dwbG90KGRhdGEgPSBzY290LmNsZWFuLCBtYXBwaW5nID0gYWVzKHggPSBhdHRhaW4pKSArDQogIGdlb21faGlzdG9ncmFtKGJpbnMgPSAyMCkgKyANCiAgbGFicyh0aXRsZSA9ICJEaXN0cmlidXRpb24gb2YgRWR1Y2F0aW9uYWwgQXR0YWlubWVudCBTY29yZXMiLA0KICAgICAgIHggPSAiU3RhbmRhcmRpemVkIEVkdWNhdGlvbmFsIEF0dGFpbm1lbnQgU2NvcmUiKSArDQogIHRoZW1lX21pbmltYWwoKQ0KYGBgDQoNCiMjIEV4cGxvcmUgYSBmZXcgb3RoZXIgcHJlZGljdG9ycw0KDQpgYGB7cn0NCmdncGxvdChkYXRhID0gc2NvdC5jbGVhbiwgbWFwcGluZyA9IGFlcyh4ID0gZGFkb2NjKSkgKw0KICBnZW9tX2hpc3RvZ3JhbShiaW5zID0gMzApICsgDQogIGxhYnModGl0bGUgPSAiRGlzdHJpYnV0aW9uIG9mIChTY2hvb2wgTWVhbikgb2YgRGFkJ3MgT2NjdXBhdGlvbiBTY29yZXMiLA0KICAgICAgIHggPSAiU3RhbmRhcmRpemVkIEVkdWNhdGlvbmFsIEF0dGFpbm1lbnQgU2NvcmUiKSArDQogIHRoZW1lX21pbmltYWwoKQ0KYGBgDQoNCmBgYHtyfQ0KZ2dwbG90KGRhdGEgPSBzY290LmNsZWFuLCBtYXBwaW5nID0gYWVzKHggPSBkZXByaXZlKSkgKw0KICBnZW9tX2hpc3RvZ3JhbShiaW5zID0gNTApICsgDQogIGxhYnModGl0bGUgPSAiRGlzdHJpYnV0aW9uIG9mIE5laWdoYm9yaG9vZCBDb250ZXh0IFNjb3JlcyIsDQogICAgICAgeCA9ICJTdGFuZGFyZGl6ZWQgTmVpZ2hib3Job29kIERlcHJpdmF0aW9uIFNjb3JlIikgKw0KICB0aGVtZV9taW5pbWFsKCkNCmBgYA0KDQojIFBhcnQgVHdvOiBNb2RlbGluZw0KIyMgRmlyc3QsIHJ1biB0d28gc2VwYXJhdGUgbW9kZWxzIHdpdGggc2Nob29sIGFuZCBuZWlnaGJvcmhvb2QgYW5kIHRoZSBMMiBjbHVzdGVycy4uLg0KDQpgYGB7cn0NCm1vZGVsLm51bGwuc2Nob29sIDwtIGxtZXIocDdyZWFkIH4gKDF8c2NoaWQpLCBSRU1MID0gRkFMU0UsIGRhdGEgPSBzY290LmNsZWFuKQ0Kc3VtbWFyeShtb2RlbC5udWxsLnNjaG9vbCkNCmBgYA0KDQojIyBDYWxjdWxhdGUgU2Nob29sIElDQw0KYGBge3J9DQppY2Muc2Nob29sIDwtIDE4Ljc1LygxOC43NSsxNzYuNjQpDQppY2Muc2Nob29sDQpgYGANCg0KYGBge3J9DQptb2RlbC5udWxsLm5laWdoIDwtIGxtZXIocDdyZWFkIH4gKDF8bmVpZ2hpZCksIFJFTUwgPSBGQUxTRSwgZGF0YSA9IHNjb3QuY2xlYW4pDQpzdW1tYXJ5KG1vZGVsLm51bGwubmVpZ2gpDQpgYGANCg0KIyMgQ2FsY3VsYXRlIE5laWdoYm9yaG9vZCBJQ0MNCmBgYHtyfQ0KaWNjLm5laWdoIDwtIDMyLjg5LygzMi44OSsxNjAuMzcpDQppY2MubmVpZ2gNCmBgYA0KDQojIyBOb3csIGEgbnVsbCBhZGRpdGl2ZSBjcm9zcy1jbGFzc2lmaWVkIG1vZGVsLi4uDQoNCmBgYHtyfQ0KbW9kZWwubnVsbC5jcm9zc2VkIDwtIGxtZXIocDdyZWFkIH4gKDF8c2NoaWQpICsgKDF8bmVpZ2hpZCksIFJFTUwgPSBGQUxTRSwgZGF0YSA9IHNjb3QuY2xlYW4pDQpzdW1tYXJ5KG1vZGVsLm51bGwuY3Jvc3NlZCkNCmBgYA0KDQojIyBDYWxjdWxhdGUgTmVpZ2hib3Job29kIElDQyBpbiBDcm9zc2VkIE1vZGVsDQpgYGB7cn0NCmljYy5uZWlnaC5jcm9zc2VkIDwtIDE3LjIxLygxNy4yMSsxNi41OSsxNjAuMzYpDQppY2MubmVpZ2guY3Jvc3NlZA0KYGBgDQoNCiMjIENhbGN1bGF0ZSBTY2hvb2wgSUNDIGluIENyb3NzZWQgTW9kZWwNCmBgYHtyfQ0KaWNjLnNjaG9vbC5jcm9zc2VkIDwtIDE2LjU5LygxNi41OSsxNy4yMSsxNjAuMzYpDQppY2Muc2Nob29sLmNyb3NzZWQNCmBgYA0KDQojIyBBZGQgYW4gTDEgUHJlZGljdG9yLCBNYWxlLi4uDQoNCmBgYHtyfQ0KbW9kZWwuMSA8LSBsbWVyKHA3cmVhZCB+IG1vbWVkLmZhYyArICgxfHNjaGlkKSArICgxfG5laWdoaWQpLCBSRU1MID0gRkFMU0UsIGRhdGEgPSBzY290LmNsZWFuKQ0Kc3VtbWFyeShtb2RlbC4xKQ0KYGBgDQoNCiMjIEwxIFBsdXMgU2Nob29sLUxldmVsIFByZWRpY3RvciBPbmx5IChgZGFkb2NjYCkuLi4NCg0KYGBge3J9DQptb2RlbC4yIDwtIGxtZXIocDdyZWFkIH4gbW9tZWQuZmFjICsgZGFkZWQuZmFjICsgKDF8c2NoaWQpICsgKDF8bmVpZ2hpZCksIFJFTUwgPSBGQUxTRSwgZGF0YSA9IHNjb3QuY2xlYW4pDQpzdW1tYXJ5KG1vZGVsLjIpDQpgYGANCg0KIyMgTDEgUGx1cyBOZWlnaGJvcmhvb2QtTGV2ZWwgUHJlZGljdG9ycyBPbmx5IChgZGVwcml2ZWApLi4uDQoNCmBgYHtyfQ0KbW9kZWwuMyA8LSBsbWVyKHA3cmVhZCB+IG1vbWVkLmZhYyArIGRhZGVkLmZhYyArIGRlcHJpdmUgKyAoMXxzY2hpZCkgKyAoMXxuZWlnaGlkKSwgUkVNTCA9IEZBTFNFLCBkYXRhID0gc2NvdC5jbGVhbikNCnN1bW1hcnkobW9kZWwuMykNCmBgYA0KDQojIyBMMSBQbHVzIEJvdGggU2Nob29sIGFuZCBOZWlnaGJvcmhvb2QgUHJlZGljdG9ycy4uLg0KDQpgYGB7cn0NCm1vZGVsLjQgPC0gbG1lcihwN3JlYWQgfiBtb21lZC5mYWMgKyBtb21lZC5mYWMgKyBkZXByaXZlICsgZGFkb2NjICsgKDF8c2NoaWQpICsgKDF8bmVpZ2hpZCksIFJFTUwgPSBGQUxTRSwgZGF0YSA9IHNjb3QuY2xlYW4pDQpzdW1tYXJ5KG1vZGVsLjQpDQpgYGANCg0KIyMgRmluYWxseSwgdGVzdCBhIG1vZGVsIHdpdGggc2Nob29sLW5laWdoYm9yaG9vZCBpbnRlcmFjdGlvbi4uLg0KYGBge3J9DQptb2RlbC41IDwtIGxtZXIocDdyZWFkIH4gbW9tZWQuZmFjICsgZGFkZWQuZmFjICsgZGVwcml2ZSArIGRhZG9jYyArIG1vbWVkLmZhYzpkYWRlZC5mYWMgKyAoMXxzY2hpZCkgKyAoMXxuZWlnaGlkKSwgUkVNTCA9IEZBTFNFLCBkYXRhID0gc2NvdC5jbGVhbikNCnN1bW1hcnkobW9kZWwuNSkNCmBgYA0KDQojIFN0ZXAgNjogQ29tcGFyZSBPcmlnaW5hbCAiRmluYWwiIE1vZGVsLCBSb2J1c3QgTW9kZWwsIGFuZCBUcmltbWVkIE1vZGVsIGZvciBEaWZmZXJlbmNlcyAoU2Vuc2l0aXZpdHkgQW5hbHlzaXMpDQpgYGB7cn0NCmxpYnJhcnkobW9kZWxzdW1tYXJ5KQ0KbGlicmFyeShicm9vbS5taXhlZCkNCm1vZGVscyA8LSBsaXN0KG1vZGVsLm51bGwubmVpZ2gsIG1vZGVsLm51bGwuc2Nob29sLCBtb2RlbC5udWxsLmNyb3NzZWQsIG1vZGVsLjEsIG1vZGVsLjIsIG1vZGVsLjMsIG1vZGVsLjQsIG1vZGVsLjUpDQptb2RlbHN1bW1hcnkobW9kZWxzLCBvdXRwdXQgPSAibWFya2Rvd24iKQ0KYGBgDQoNCg==