### R Textbook Examples Introduction to Multilevel Modeling by Kreft and de Leeuw Chapter 4: Analyses

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)

attach(imm23)

library(Hmisc)

describe(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   <2e-16 ***
homework      2.3354     0.2684   8.701   <2e-16 ***
parented      3.0040     0.2765  10.864   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 8.724 on 516 degrees of freedom
Multiple R-squared: 0.3389,	Adjusted R-squared: 0.3363
F-statistic: 132.2 on 2 and 516 DF,  p-value: < 2.2e-16
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).

The content of this web site should not be construed as an endorsement of any particular web site, book, or software product by the University of California.