library(lme4)
## Loading required package: Matrix
library(nlme)
## 
## Attaching package: 'nlme'
## The following object is masked from 'package:lme4':
## 
##     lmList
Achieve <- read.csv("Achieve.csv")
Model4.1 <- lme( fixed = geread~1, random = ~1|school/class, data = Achieve)
summary(Model4.1)
## Linear mixed-effects model fit by REML
##  Data: Achieve 
##     AIC      BIC logLik
##   46154 46182.97 -23073
## 
## Random effects:
##  Formula: ~1 | school
##         (Intercept)
## StdDev:   0.5583923
## 
##  Formula: ~1 | class %in% school
##         (Intercept) Residual
## StdDev:   0.5221676 2.201589
## 
## Fixed effects: geread ~ 1 
##               Value  Std.Error   DF  t-value p-value
## (Intercept) 4.30806 0.05499164 9752 78.34027       0
## 
## Standardized Within-Group Residuals:
##        Min         Q1        Med         Q3        Max 
## -2.3052003 -0.6289603 -0.2093700  0.3049085  3.8673234 
## 
## Number of Observations: 10320
## Number of Groups: 
##            school class %in% school 
##               160               568
intervals(Model4.1)
## Approximate 95% confidence intervals
## 
##  Fixed effects:
##                lower    est.    upper
## (Intercept) 4.200265 4.30806 4.415855
## attr(,"label")
## [1] "Fixed effects:"
## 
##  Random Effects:
##   Level: school 
##                     lower      est.     upper
## sd((Intercept)) 0.4702517 0.5583923 0.6630533
##   Level: class 
##                     lower      est.     upper
## sd((Intercept)) 0.4545912 0.5221676 0.5997895
## 
##  Within-group standard error:
##    lower     est.    upper 
## 2.170908 2.201589 2.232704
Model4.2 <- lme( fixed = geread~gevocab+clenroll+cenroll, random = ~1|school/class, data = Achieve)
summary(Model4.2)
## Linear mixed-effects model fit by REML
##  Data: Achieve 
##        AIC      BIC    logLik
##   43144.87 43195.56 -21565.43
## 
## Random effects:
##  Formula: ~1 | school
##         (Intercept)
## StdDev:   0.2766194
## 
##  Formula: ~1 | class %in% school
##         (Intercept) Residual
## StdDev:   0.3007871 1.922991
## 
## Fixed effects: geread ~ gevocab + clenroll + cenroll 
##                  Value  Std.Error   DF  t-value p-value
## (Intercept)  1.6751266 0.20809605 9751  8.04978  0.0000
## gevocab      0.5075566 0.00842654 9751 60.23313  0.0000
## clenroll     0.0189860 0.00955860  407  1.98628  0.0477
## cenroll     -0.0000037 0.00000364  158 -1.02193  0.3084
##  Correlation: 
##          (Intr) gevocb clnrll
## gevocab  -0.124              
## clenroll -0.961 -0.062       
## cenroll  -0.134  0.025 -0.007
## 
## Standardized Within-Group Residuals:
##        Min         Q1        Med         Q3        Max 
## -3.2211629 -0.5672782 -0.2079045  0.3183508  4.4736276 
## 
## Number of Observations: 10320
## Number of Groups: 
##            school class %in% school 
##               160               568
Model4.3 <- lme( fixed = geread~gevocab+clenroll+cenroll+gevocab*cenroll, random = ~1|school/class, data = Achieve)
summary(Model4.3)
## Linear mixed-effects model fit by REML
##  Data: Achieve 
##        AIC      BIC    logLik
##   43167.75 43225.69 -21575.88
## 
## Random effects:
##  Formula: ~1 | school
##         (Intercept)
## StdDev:   0.2740961
## 
##  Formula: ~1 | class %in% school
##         (Intercept) Residual
## StdDev:   0.2975919 1.923059
## 
## Fixed effects: geread ~ gevocab + clenroll + cenroll + gevocab * cenroll 
##                      Value  Std.Error   DF  t-value p-value
## (Intercept)      1.7515430 0.20999286 9750  8.34096  0.0000
## gevocab          0.4899998 0.01168332 9750 41.94013  0.0000
## clenroll         0.0188007 0.00951172  407  1.97659  0.0488
## cenroll         -0.0000132 0.00000563  158 -2.33721  0.0207
## gevocab:cenroll  0.0000023 0.00000107 9750  2.18957  0.0286
##  Correlation: 
##                 (Intr) gevocb clnrll cenrll
## gevocab         -0.203                     
## clenroll        -0.949 -0.041              
## cenroll         -0.212  0.542  0.000       
## gevocab:cenroll  0.166 -0.693 -0.007 -0.766
## 
## Standardized Within-Group Residuals:
##        Min         Q1        Med         Q3        Max 
## -3.1901563 -0.5682666 -0.2060729  0.3183307  4.4723839 
## 
## Number of Observations: 10320
## Number of Groups: 
##            school class %in% school 
##               160               568
Model4.4 <- lme( fixed = geread~1, random = ~1|corp/school/class, data = Achieve)
summary(Model4.4)
## Linear mixed-effects model fit by REML
##  Data: Achieve 
##        AIC      BIC    logLik
##   46113.22 46149.43 -23051.61
## 
## Random effects:
##  Formula: ~1 | corp
##         (Intercept)
## StdDev:   0.4209979
## 
##  Formula: ~1 | school %in% corp
##         (Intercept)
## StdDev:    0.295833
## 
##  Formula: ~1 | class %in% school %in% corp
##         (Intercept) Residual
## StdDev:   0.5247746 2.201587
## 
## Fixed effects: geread ~ 1 
##               Value  Std.Error   DF  t-value p-value
## (Intercept) 4.32583 0.07197848 9752 60.09894       0
## 
## Standardized Within-Group Residuals:
##        Min         Q1        Med         Q3        Max 
## -2.2995234 -0.6304742 -0.2130696  0.3028624  3.9448255 
## 
## Number of Observations: 10320
## Number of Groups: 
##                        corp            school %in% corp 
##                          59                         160 
## class %in% school %in% corp 
##                         568
intervals(Model4.4)
## Approximate 95% confidence intervals
## 
##  Fixed effects:
##                lower    est.    upper
## (Intercept) 4.184738 4.32583 4.466923
## attr(,"label")
## [1] "Fixed effects:"
## 
##  Random Effects:
##   Level: corp 
##                    lower      est.     upper
## sd((Intercept)) 0.321723 0.4209979 0.5509065
##   Level: school 
##                     lower     est.     upper
## sd((Intercept)) 0.2003532 0.295833 0.4368144
##   Level: class 
##                     lower      est.     upper
## sd((Intercept)) 0.4578135 0.5247746 0.6015295
## 
##  Within-group standard error:
##    lower     est.    upper 
## 2.170912 2.201587 2.232695
Model4.5 <- lme( fixed = geread~gevocab+gender, random = ~gender|school/class, data = Achieve)
summary(Model4.5)
## Linear mixed-effects model fit by REML
##  Data: Achieve 
##        AIC      BIC    logLik
##   43127.93 43200.35 -21553.97
## 
## Random effects:
##  Formula: ~gender | school
##  Structure: General positive-definite, Log-Cholesky parametrization
##             StdDev    Corr  
## (Intercept) 0.2207298 (Intr)
## gender      0.1101379 -0.019
## 
##  Formula: ~gender | class %in% school
##  Structure: General positive-definite, Log-Cholesky parametrization
##             StdDev      Corr  
## (Intercept) 0.302941367 (Intr)
## gender      0.004914493 -0.022
## Residual    1.922517135       
## 
## Fixed effects: geread ~ gevocab + gender 
##                 Value  Std.Error   DF  t-value p-value
## (Intercept) 2.0150114 0.07570370 9750 26.61708  0.0000
## gevocab     0.5091256 0.00840838 9750 60.54976  0.0000
## gender      0.0175517 0.03929581 9750  0.44666  0.6551
##  Correlation: 
##         (Intr) gevocb
## gevocab -0.527       
## gender  -0.758  0.039
## 
## Standardized Within-Group Residuals:
##        Min         Q1        Med         Q3        Max 
## -3.2117003 -0.5676457 -0.2072293  0.3161412  4.4474380 
## 
## Number of Observations: 10320
## Number of Groups: 
##            school class %in% school 
##               160               568
Model4.6 <- lme( fixed = geread~gevocab+gender, random = list(school = ~1, class = ~gender), data = Achieve)
summary(Model4.6)
## Linear mixed-effects model fit by REML
##  Data: Achieve 
##        AIC      BIC    logLik
##   43125.18 43183.11 -21554.59
## 
## Random effects:
##  Formula: ~1 | school
##         (Intercept)
## StdDev:   0.2737238
## 
##  Formula: ~gender | class %in% school
##  Structure: General positive-definite, Log-Cholesky parametrization
##             StdDev    Corr  
## (Intercept) 0.3623693 (Intr)
## gender      0.1651168 -0.563
## Residual    1.9215122       
## 
## Fixed effects: geread ~ gevocab + gender 
##                 Value  Std.Error   DF  t-value p-value
## (Intercept) 2.0128841 0.07726503 9750 26.05168  0.0000
## gevocab     0.5090473 0.00841459 9750 60.49582  0.0000
## gender      0.0190565 0.03880626 9750  0.49107  0.6234
##  Correlation: 
##         (Intr) gevocb
## gevocab -0.516       
## gender  -0.768  0.039
## 
## Standardized Within-Group Residuals:
##        Min         Q1        Med         Q3        Max 
## -3.2117225 -0.5676842 -0.2072084  0.3182780  4.4324392 
## 
## Number of Observations: 10320
## Number of Groups: 
##            school class %in% school 
##               160               568
Model4.7 <- lmer(geread~1+(1|school/class), data = Achieve)
summary(Model4.7)
## Linear mixed model fit by REML ['lmerMod']
## Formula: geread ~ 1 + (1 | school/class)
##    Data: Achieve
## 
## REML criterion at convergence: 46146
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.3052 -0.6290 -0.2094  0.3049  3.8673 
## 
## Random effects:
##  Groups       Name        Variance Std.Dev.
##  class:school (Intercept) 0.2727   0.5222  
##  school       (Intercept) 0.3118   0.5584  
##  Residual                 4.8470   2.2016  
## Number of obs: 10320, groups:  class:school, 568; school, 160
## 
## Fixed effects:
##             Estimate Std. Error t value
## (Intercept)  4.30806    0.05499   78.34
Model4.8 <- lmer( geread~gevocab+clenroll+cenroll+(1|school/class), data = Achieve)
## Warning: Some predictor variables are on very different scales: consider
## rescaling
summary(Model4.8)
## Linear mixed model fit by REML ['lmerMod']
## Formula: geread ~ gevocab + clenroll + cenroll + (1 | school/class)
##    Data: Achieve
## 
## REML criterion at convergence: 43130.9
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.2212 -0.5673 -0.2079  0.3184  4.4736 
## 
## Random effects:
##  Groups       Name        Variance Std.Dev.
##  class:school (Intercept) 0.09047  0.3008  
##  school       (Intercept) 0.07652  0.2766  
##  Residual                 3.69790  1.9230  
## Number of obs: 10320, groups:  class:school, 568; school, 160
## 
## Fixed effects:
##               Estimate Std. Error t value
## (Intercept)  1.675e+00  2.081e-01    8.05
## gevocab      5.076e-01  8.427e-03   60.23
## clenroll     1.899e-02  9.559e-03    1.99
## cenroll     -3.721e-06  3.642e-06   -1.02
## 
## Correlation of Fixed Effects:
##          (Intr) gevocb clnrll
## gevocab  -0.124              
## clenroll -0.961 -0.062       
## cenroll  -0.134  0.025 -0.007
## fit warnings:
## Some predictor variables are on very different scales: consider rescaling
Model4.9 <- lmer( geread~gevocab+clenroll+cenroll+gevocab*cenroll+(1|school/class), data = Achieve)
## Warning: Some predictor variables are on very different scales: consider
## rescaling
summary(Model4.9)
## Linear mixed model fit by REML ['lmerMod']
## Formula: geread ~ gevocab + clenroll + cenroll + gevocab * cenroll + (1 |  
##     school/class)
##    Data: Achieve
## 
## REML criterion at convergence: 43151.8
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.1902 -0.5683 -0.2061  0.3183  4.4724 
## 
## Random effects:
##  Groups       Name        Variance Std.Dev.
##  class:school (Intercept) 0.08856  0.2976  
##  school       (Intercept) 0.07513  0.2741  
##  Residual                 3.69816  1.9231  
## Number of obs: 10320, groups:  class:school, 568; school, 160
## 
## Fixed effects:
##                   Estimate Std. Error t value
## (Intercept)      1.752e+00  2.100e-01    8.34
## gevocab          4.900e-01  1.168e-02   41.94
## clenroll         1.880e-02  9.512e-03    1.98
## cenroll         -1.316e-05  5.628e-06   -2.34
## gevocab:cenroll  2.340e-06  1.069e-06    2.19
## 
## Correlation of Fixed Effects:
##             (Intr) gevocb clnrll cenrll
## gevocab     -0.203                     
## clenroll    -0.949 -0.041              
## cenroll     -0.212  0.542  0.000       
## gevcb:cnrll  0.166 -0.693 -0.007 -0.766
## fit warnings:
## Some predictor variables are on very different scales: consider rescaling
Model4.10 <- lmer( geread~1+(1|corp/school/class), data = Achieve)
summary(Model4.10)
## Linear mixed model fit by REML ['lmerMod']
## Formula: geread ~ 1 + (1 | corp/school/class)
##    Data: Achieve
## 
## REML criterion at convergence: 46103.2
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.2995 -0.6305 -0.2131  0.3029  3.9448 
## 
## Random effects:
##  Groups              Name        Variance Std.Dev.
##  class:(school:corp) (Intercept) 0.27539  0.5248  
##  school:corp         (Intercept) 0.08748  0.2958  
##  corp                (Intercept) 0.17726  0.4210  
##  Residual                        4.84699  2.2016  
## Number of obs: 10320, groups:  
## class:(school:corp), 568; school:corp, 160; corp, 59
## 
## Fixed effects:
##             Estimate Std. Error t value
## (Intercept)  4.32583    0.07198    60.1
Model4.11 <- lmer( geread~gevocab+gender+(gender|school/class), data = Achieve)
summary(Model4.11)
## Linear mixed model fit by REML ['lmerMod']
## Formula: geread ~ gevocab + gender + (gender | school/class)
##    Data: Achieve
## 
## REML criterion at convergence: 43107.8
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.2043 -0.5680 -0.2069  0.3171  4.4455 
## 
## Random effects:
##  Groups       Name        Variance Std.Dev. Corr 
##  class:school (Intercept) 0.148959 0.38595       
##               gender      0.019818 0.14078  -0.62
##  school       (Intercept) 0.033157 0.18209       
##               gender      0.006802 0.08247  0.58 
##  Residual                 3.692434 1.92157       
## Number of obs: 10320, groups:  class:school, 568; school, 160
## 
## Fixed effects:
##             Estimate Std. Error t value
## (Intercept) 2.015564   0.075629   26.65
## gevocab     0.509097   0.008408   60.55
## gender      0.017237   0.039243    0.44
## 
## Correlation of Fixed Effects:
##         (Intr) gevocb
## gevocab -0.527       
## gender  -0.757  0.039