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).