Note: This page is designed to show the how multilevel model can be done using R and to be able to compare the results with those in the book.
On this page we will use the lmer function which is found in the lme4 package. There are several other possible choices but we will go with lmer.
The data were downloaded in Stata format from here and imported into R using the foreign library from a directory called rdata on the local computer.
library(foreign)
imm23
imm23
18 Variables 519 Observations
---------------------------------------------------------------------------
schid
n missing unique Mean .05 .10 .25 .50 .75 .90 .95
519 0 23 35489 6053 6439 7801 25642 62821 72080 72292
lowest : 6053 6327 6467 7194 7472, highest: 68448 68493 72080 72292 72991
---------------------------------------------------------------------------
stuid
n missing unique Mean .05 .10 .25 .50 .75 .90 .95
519 0 100 49.86 6.0 10.0 24.5 51.0 74.0 89.0 94.0
lowest : 0 1 2 3 4, highest: 95 96 97 98 99
---------------------------------------------------------------------------
ses
n missing unique Mean .05 .10 .25 .50 .75 .90 .95
519 0 248 -0.001272 -1.331 -1.092 -0.620 -0.120 0.730 1.230 1.391
lowest : -2.41 -2.23 -2.08 -2.04 -1.96, highest: 1.68 1.71 1.80 1.82 1.85
---------------------------------------------------------------------------
meanses
n missing unique Mean .05 .10 .25 .50 .75 .90 .95
519 0 23 -0.001272 -0.8480 -0.7009 -0.3990 -0.1968 0.6575 1.0446 1.0446
lowest : -1.0685 -0.8480 -0.7009 -0.5045 -0.4826, highest: 0.1784 0.6575 0.6998 1.0446 1.1762
---------------------------------------------------------------------------
homework
n missing unique Mean
519 0 8 1.971
0 1 2 3 4 5 6 7
Frequency 42 225 111 47 47 38 6 3
% 8 43 21 9 9 7 1 1
---------------------------------------------------------------------------
white
n missing unique Sum Mean
519 0 2 390 0.7514
---------------------------------------------------------------------------
parented
n missing unique Mean
519 0 6 3.289
1 2 3 4 5 6
Frequency 53 103 176 65 72 50
% 10 20 34 13 14 10
---------------------------------------------------------------------------
public
n missing unique Sum Mean
519 0 2 317 0.6108
---------------------------------------------------------------------------
ratio
n missing unique Mean .05 .10 .25 .50 .75 .90 .95
519 0 14 16.76 10 10 13 18 20 22 25
10 11 12 13 14 17 18 19 20 21 22 23 25 28
Frequency 83 20 22 31 64 20 113 28 20 24 43 8 20 23
% 16 4 4 6 12 4 22 5 4 5 8 2 4 4
---------------------------------------------------------------------------
percmin
n missing unique Mean
519 0 7 2.501
0 1 2 3 5 6 7
Frequency 75 149 36 158 42 20 39
% 14 29 7 30 8 4 8
---------------------------------------------------------------------------
math
n missing unique Mean .05 .10 .25 .50 .75 .90 .95
519 0 42 51.72 35.0 38.0 43.0 51.0 61.5 67.0 69.0
lowest : 30 31 32 33 34, highest: 67 68 69 70 71
---------------------------------------------------------------------------
sex
n missing unique Mean
519 0 2 1.520
1 (249, 48%), 2 (270, 52%)
---------------------------------------------------------------------------
race
n missing unique Mean
519 0 5 3.615
1 2 3 4 5
Frequency 20 39 66 390 4
% 4 8 13 75 1
---------------------------------------------------------------------------
sctype
n missing unique Mean
519 0 4 1.977
1 (317, 61%), 2 (43, 8%), 3 (13, 3%), 4 (146, 28%)
---------------------------------------------------------------------------
cstr
n missing unique Mean
519 0 4 3.819
2 (23, 4%), 3 (148, 29%), 4 (248, 48%), 5 (100, 19%)
---------------------------------------------------------------------------
scsize
n missing unique Mean
519 0 7 3.229
1 2 3 4 5 6 7
Frequency 13 181 154 89 23 45 14
% 3 35 30 17 4 9 3
---------------------------------------------------------------------------
urban
n missing unique Mean
519 0 3 1.832
1 (242, 47%), 2 (122, 24%), 3 (155, 30%)
---------------------------------------------------------------------------
region
n missing unique Mean
519 0 3 2.125
1 (100, 19%), 2 (254, 49%), 3 (165, 32%)
---------------------------------------------------------------------------
Page 64, 4.2.2, The Null Model, Model 0.
library(lme4)
lmer(math~(1|schid))
Linear mixed-effects model fit by REML
Formula: math ~ (1 | schid)
AIC BIC logLik MLdeviance REMLdeviance
3803 3811 -1899 3801 3799
Random effects:
Groups Name Variance Std.Dev.
schid (Intercept) 26.124 5.1112
Residual 81.244 9.0135
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 50.759 1.151 44.09
Page 65, 4.2.3 ‘Homework and ‘MathAchievement’, Model 1.
lmer(math~homework + (1|schid))
Linear mixed-effects model fit by REML
Formula: math ~ homework + (1 | schid)
AIC BIC logLik MLdeviance REMLdeviance
3735 3748 -1865 3731 3729
Random effects:
Groups Name Variance Std.Dev.
schid (Intercept) 21.342 4.6197
Residual 71.284 8.4430
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 46.3558 1.1628 39.87
homework 2.3999 0.2772 8.66
Correlation of Fixed Effects:
(Intr)
homework -0.437
Pages 66 and 67, 4.2.4 Random slope for ‘Homework’, model 2.
lmer(math~homework + (homework|schid))
Linear mixed-effects model fit by REML
Formula: math ~ homework + (homework | schid)
AIC BIC logLik MLdeviance REMLdeviance
3646 3667 -1818 3639 3636
Random effects:
Groups Name Variance Std.Dev. Corr
schid (Intercept) 62.393 7.8989
homework 17.703 4.2075 -0.829
Residual 53.297 7.3005
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 46.3255 1.7585 26.343
homework 1.9804 0.9279 2.134
Correlation of Fixed Effects:
(Intr)
homework -0.824
Page 69, 4.2.5 Adding ‘ParentEducation’, model 3.
lmer(math~homework + parented + (homework|schid))
Linear mixed-effects model fit by REML
Formula: math ~ homework + parented + (homework | schid)
AIC BIC logLik MLdeviance REMLdeviance
3612 3638 -1800 3602 3600
Random effects:
Groups Name Variance Std.Dev. Corr
schid (Intercept) 47.864 6.9184
homework 13.881 3.7257 -0.852
Residual 50.779 7.1259
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 40.8546 1.7901 22.822
homework 1.8817 0.8303 2.266
parented 1.8414 0.2959 6.223
Correlation of Fixed Effects:
(Intr) homwrk
homework -0.722
parented -0.489 -0.026
Page 70 and 71, 4.2.6 Traditional regression model.
regmod<-lm(math ~ homework + parented)
summary(regmod)
Call:
lm(formula = math ~ homework + parented)
Residuals:
Min 1Q Median 3Q Max
-22.25720 -6.91785 0.07817 6.40755 21.41751
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 37.2392 0.9963 37.378
Page 73/74, 4.3.2 A model with ‘SchoolSize’ (Model 2).
lmer(math~homework + scsize + (homework|schid))
Linear mixed-effects model fit by REML
Formula: math ~ homework + scsize + (homework | schid)
AIC BIC logLik MLdeviance REMLdeviance
3646 3672 -1817 3639 3634
Random effects:
Groups Name Variance Std.Dev. Corr
schid (Intercept) 66.357 8.1460
homework 17.768 4.2153 -0.836
Residual 53.310 7.3013
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 44.9724 2.7097 16.597
homework 1.9813 0.9295 2.132
scsize 0.4258 0.6407 0.665
Correlation of Fixed Effects:
(Intr) homwrk
homework -0.543
scsize -0.745 -0.014
Page 74/75, 4.3.3 Changing ‘SchoolSize’ to ‘Public’ (Model 3).
lmer(math~homework + public + (homework|schid))
Linear mixed-effects model fit by REML
Formula: math ~ homework + public + (homework | schid)
AIC BIC logLik MLdeviance REMLdeviance
3640 3666 -1814 3635 3628
Random effects:
Groups Name Variance Std.Dev. Corr
schid (Intercept) 60.351 7.7686
homework 17.283 4.1573 -0.851
Residual 53.345 7.3038
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 49.0509 2.1870 22.428
homework 1.9762 0.9178 2.153
public -4.0633 1.9787 -2.053
Correlation of Fixed Effects:
(Intr) homwrk
homework -0.673
public -0.610 0.009
Page 77, 4.3.4 Adding a cross level interaction with ‘Public’, (Model 4).
lmer(math~homework + public + homework:public + (homework|schid))
Linear mixed-effects model fit by REML
Formula: math ~ homework + public + homework:public + (homework | schid)
AIC BIC logLik MLdeviance REMLdeviance
3639 3669 -1813 3635 3625
Random effects:
Groups Name Variance Std.Dev. Corr
schid (Intercept) 62.550 7.9089
homework 18.249 4.2718 -0.857
Residual 53.333 7.3029
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 48.5289 3.0159 16.091
homework 2.2928 1.5912 1.441
public -3.2619 3.7145 -0.878
homework:public -0.4957 1.9726 -0.251
Correlation of Fixed Effects:
(Intr) homwrk public
homework -0.844
public -0.812 0.685
homwrk:pblc 0.681 -0.807 -0.846
Page 80, 4.3.5 Model 4 will full NELS-88 data (we don’t have these data, so this is omitted). Page 80/82, 4.3.6 Deleting ‘HomePublic’ and adding ‘White’ (Model 5).
lmer(math~homework + public + white + (homework|schid))
Linear mixed-effects model fit by REML
Formula: math ~ homework + public + white + (homework | schid)
AIC BIC logLik MLdeviance REMLdeviance
3629 3659 -1808 3623 3615
Random effects:
Groups Name Variance Std.Dev. Corr
schid (Intercept) 56.045 7.4863
homework 16.783 4.0967 -0.875
Residual 52.727 7.2614
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 46.6283 2.1923 21.269
homework 1.8991 0.9052 2.098
public -3.8823 1.7987 -2.158
white 3.3094 0.9698 3.412
Correlation of Fixed Effects:
(Intr) homwrk public
homework -0.657
public -0.567 0.009
white -0.327 -0.027 0.035
Page 82/83, 4.3.7 Adding a random part for ‘White’ (Model 6).
lmer(math~homework + public + white + (homework + white|schid))
Linear mixed-effects model fit by REML
Formula: math ~ homework + public + white + (homework + white | schid)
AIC BIC logLik MLdeviance REMLdeviance
3631 3673 -1805 3619 3611
Random effects:
Groups Name Variance Std.Dev. Corr
schid (Intercept) 69.117 8.3137
homework 16.656 4.0811 -0.842
white 26.395 5.1376 -0.515 0.139
Residual 51.153 7.1521
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 48.222 2.340 20.608
homework 1.939 0.901 2.152
public -4.925 1.652 -2.980
white 2.599 1.552 1.674
Correlation of Fixed Effects:
(Intr) homwrk public
homework -0.657
public -0.483 0.010
white -0.536 0.083 0.006
Page 85, 4.3.8 Making the coefficient of ‘White’ fixed and adding ‘MeanSES’ (Model 7).
lmer(math~homework + public + white + meanses + (homework|schid))
Linear mixed-effects model fit by REML
Formula: math ~ homework + public + white + meanses + (homework | schid)
AIC BIC logLik MLdeviance REMLdeviance
3623 3657 -1803 3617 3607
Random effects:
Groups Name Variance Std.Dev. Corr
schid (Intercept) 54.124 7.3569
homework 16.443 4.0550 -0.906
Residual 52.794 7.2660
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 44.6159 2.2350 19.962
homework 1.9241 0.8964 2.147
public 0.1604 2.2711 0.071
white 3.0969 0.9631 3.216
meanses 4.9771 1.9510 2.551
Correlation of Fixed Effects:
(Intr) homwrk public white
homework -0.657
public -0.596 0.011
white -0.262 -0.027 -0.075
meanses -0.345 0.004 0.711 -0.141
Page 86, 4.3.9 Deleting the school characteristic ‘Public’ (Model 8).
lmer(math~homework + white + meanses + (homework|schid))
Linear mixed-effects model fit by REML
Formula: math ~ homework + white + meanses + (homework | schid)
AIC BIC logLik MLdeviance REMLdeviance
3624 3654 -1805 3617 3610
Random effects:
Groups Name Variance Std.Dev. Corr
schid (Intercept) 53.613 7.3221
homework 16.409 4.0508 -0.911
Residual 52.788 7.2656
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 44.7022 1.7877 25.005
homework 1.9251 0.8954 2.150
white 3.1148 0.9570 3.255
meanses 4.8925 1.3408 3.649
Correlation of Fixed Effects:
(Intr) homwrk white
homework -0.813
white -0.384 -0.026
meanses 0.139 -0.006 -0.126
Page 87/88, 4.3.10 Adding an interaction between ‘HomeWork’ and ‘MeanSES’ (Model 9).
lmer(math~homework + white + meanses + homework:meanses + (homework|schid))
Linear mixed-effects model fit by REML
Formula: math ~ homework + white + meanses + homework:meanses + (homework | schid)
AIC BIC logLik MLdeviance REMLdeviance
3623 3657 -1804 3617 3607
Random effects:
Groups Name Variance Std.Dev. Corr
schid (Intercept) 55.809 7.4705
homework 17.240 4.1520 -0.915
Residual 52.782 7.2651
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 44.6117 1.8351 24.310
homework 1.9744 0.9292 2.125
white 3.1150 0.9571 3.255
meanses 3.9781 3.0209 1.317
homework:meanses 0.5526 1.6458 0.336
Correlation of Fixed Effects:
(Intr) homwrk white meanss
homework -0.824
white -0.375 -0.023
meanses 0.195 -0.156 -0.066
homwrk:mnss -0.150 0.171 0.011 -0.896
Page 88/89, 4.3.11 Adding another student-level variable (Model 10).
lmer(math~homework + white + meanses + ses + (homework|schid))
Linear mixed-effects model fit by REML
Formula: math ~ homework + white + meanses + ses + (homework | schid)
AIC BIC logLik MLdeviance REMLdeviance
3609 3643 -1796 3600 3593
Random effects:
Groups Name Variance Std.Dev. Corr
schid (Intercept) 49.803 7.0571
homework 14.618 3.8233 -0.901
Residual 51.300 7.1624
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 45.6750 1.7509 26.087
homework 1.8257 0.8497 2.149
white 2.1706 0.9733 2.230
meanses 2.9448 1.4284 2.062
ses 2.2060 0.5356 4.118
Correlation of Fixed Effects:
(Intr) homwrk white meanss
homework -0.798
white -0.408 -0.019
meanses 0.088 0.004 -0.036
ses 0.135 -0.031 -0.235 -0.333
Page 88/89, 4.3.12: Analyses with NELS-88 (we don’t have these data, so these analyses are omitted).
Page 91, 4.4.1 ‘SES’ as a student-level explanatory variable (Model 1).
lmer(math~ses + (1|schid))
Linear mixed-effects model fit by REML
Formula: math ~ ses + (1 | schid)
AIC BIC logLik MLdeviance REMLdeviance
3752 3765 -1873 3748 3746
Random effects:
Groups Name Variance Std.Dev.
schid (Intercept) 12.632 3.5541
Residual 75.328 8.6792
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 51.2009 0.8507 60.18
ses 4.3323 0.5663 7.65
Correlation of Fixed Effects:
(Intr)
ses 0.069
Page 92, 4.4.2 Adding a random slope (Model 2). Model does not converge properly.
lmer(math~ses + (ses|schid))
Linear mixed-effects model fit by REML
Formula: math ~ ses + (ses | schid)
AIC BIC logLik MLdeviance REMLdeviance
3756 3777 -1873 3748 3746
Random effects:
Groups Name Variance Std.Dev. Corr
schid (Intercept) 1.2633e+01 3.55426291
ses 3.7664e-08 0.00019407 0.000
Residual 7.5328e+01 8.67914725
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 51.2009 0.8508 60.18
ses 4.3323 0.5663 7.65
Correlation of Fixed Effects:
(Intr)
ses 0.069
Warning message:
In .local(x, ..., value) :
Estimated variance-covariance for factor 'schid' is singular
Page 93, 4.4.3 Adding ‘PercentMinorities’ (Model 3), SAS Program.
lmer(math~ses + percmin + (1|schid))
Linear mixed-effects model fit by REML
Formula: math ~ ses + percmin + (1 | schid)
AIC BIC logLik MLdeviance REMLdeviance
3750 3767 -1871 3743 3742
Random effects:
Groups Name Variance Std.Dev.
schid (Intercept) 10.698 3.2708
Residual 75.173 8.6702
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 53.1267 1.1781 45.10
ses 4.2988 0.5618 7.65
percmin -0.8094 0.3647 -2.22
Correlation of Fixed Effects:
(Intr) ses
ses -0.003
percmin -0.735 0.071
Page 95, 4.4.4 Adding ‘MeanSES’ (Model 4).
lmer(math~ses + percmin + meanses + (1|schid))
Linear mixed-effects model fit by REML
Formula: math ~ ses + percmin + meanses + (1 | schid)
AIC BIC logLik MLdeviance REMLdeviance
3746 3767 -1868 3740 3736
Random effects:
Groups Name Variance Std.Dev.
schid (Intercept) 8.8468 2.9743
Residual 75.2405 8.6741
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 53.0969 1.1027 48.15
ses 3.8848 0.6098 6.37
percmin -0.6922 0.3457 -2.00
meanses 2.8051 1.4794 1.90
Correlation of Fixed Effects:
(Intr) ses percmn
ses 0.000
percmin -0.728 0.000
meanses -0.007 -0.412 0.162
Page 95, 4.4.5 Analyses with NELS-88, models 2 and 3 (we do not have these data, so these analyses are omitted).
Page 99, 4.5.1 Analysis with class size and a cross level interaction (Model 1).
lmer(math~homework + ratio + (homework|schid))
Linear mixed-effects model fit by REML
Formula: math ~ homework + ratio + (homework | schid)
AIC BIC logLik MLdeviance REMLdeviance
3649 3674 -1818 3639 3637
Random effects:
Groups Name Variance Std.Dev. Corr
schid (Intercept) 63.607 7.9754
homework 17.761 4.2144 -0.826
Residual 53.302 7.3008
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 47.96224 4.08369 11.745
homework 1.97941 0.92932 2.130
ratio -0.09448 0.21237 -0.445
Correlation of Fixed Effects:
(Intr) homwrk
homework -0.359
ratio -0.901 0.002
Page 100, 4.5.2 Interaction between ‘Ratio’ and ‘HomeWork’ (Model 2), SAS Program.
lmer(math~homework + homework:ratio + (homework|schid))
Linear mixed-effects model fit by REML
Formula: math ~ homework + homework:ratio + (homework | schid)
AIC BIC logLik MLdeviance REMLdeviance
3650 3675 -1819 3639 3638
Random effects:
Groups Name Variance Std.Dev. Corr
schid (Intercept) 62.557 7.9093
homework 18.099 4.2543 -0.826
Residual 53.298 7.3005
number of obs: 519, groups: schid, 23
Fixed effects:
Estimate Std. Error t value
(Intercept) 46.33064 1.76061 26.315
homework 2.88416 2.15681 1.337
homework:ratio -0.05259 0.11228 -0.468
Correlation of Fixed Effects:
(Intr) homwrk
homework -0.357
homework:rt 0.000 -0.901
Page 100, 4.5.3 Reporting the modeling session with NELS-88 (we do not have these data, so these analyses are omitted).