UCLA Academic Technology Services HomeServicesClassesContactJobs
Help the Stat Consulting Group by giving a gift             
Loading

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)

imm23 <- read.dta(file="~/rdata/imm23.dta")

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

How to cite this page

Report an error on this page or leave a comment

UCLA Researchers are invited to our Statistical Consulting Services
We recommend others to our list of Other Resources for Statistical Computing Help
These pages are Copyrighted (c) by UCLA Academic Technology Services


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