require(nlme)## Loading required package: nlmerequire(lme4)## Loading required package: lme4## Loading required package: Matrix##
## Attaching package: 'lme4'## The following object is masked from 'package:nlme':
##
## lmListQ7
# read data
dta <- read.csv("attitudes2Sci.csv", header = TRUE)
head(dta) State PrivPub school class sex like Class
1 ACT public 1 1 f 8 1.1
2 ACT public 1 1 f 6 1.1
3 ACT public 1 1 f 5 1.1
4 ACT public 1 1 f 2 1.1
5 ACT public 1 1 f 5 1.1
6 ACT public 1 1 f 6 1.1str(dta) 'data.frame': 1385 obs. of 7 variables:
$ State : Factor w/ 2 levels "ACT","NSW": 1 1 1 1 1 1 1 1 1 1 ...
$ PrivPub: Factor w/ 2 levels "private","public": 2 2 2 2 2 2 2 2 2 2 ...
$ school : int 1 1 1 1 1 1 1 1 1 1 ...
$ class : int 1 1 1 1 1 1 1 1 1 1 ...
$ sex : Factor w/ 2 levels "f","m": 1 1 1 1 1 1 1 1 2 2 ...
$ like : int 8 6 5 2 5 6 3 7 6 3 ...
$ Class : num 1.1 1.1 1.1 1.1 1.1 1.1 1.1 1.1 1.1 1.1 ...dta$school<-as.factor(dta$school)
dta$class<-as.factor(dta$class)
dta$Class<-as.factor(dta$Class)
# model
lm <- lmer(data = dta, like ~ State + PrivPub + sex + (1|school) + (1|Class))
summary(lm) Linear mixed model fit by REML ['lmerMod']
Formula: like ~ State + PrivPub + sex + (1 | school) + (1 | Class)
Data: dta
REML criterion at convergence: 5546.3
Scaled residuals:
Min 1Q Median 3Q Max
-2.6941 -0.6893 0.1513 0.6822 2.5651
Random effects:
Groups Name Variance Std.Dev.
Class (Intercept) 0.3281 0.5728
school (Intercept) 0.0000 0.0000
Residual 3.0522 1.7471
Number of obs: 1383, groups: Class, 66; school, 41
Fixed effects:
Estimate Std. Error t value
(Intercept) 4.71753 0.17539 26.898
StateNSW 0.02762 0.45471 0.061
PrivPubpublic 0.41586 0.19714 2.109
sexm 0.18275 0.09827 1.860
Correlation of Fixed Effects:
(Intr) SttNSW PrvPbp
StateNSW -0.360
PrivPubpblc -0.818 0.314
sexm -0.297 0.029 0.021#95% CI
confint(lm, oldNames = FALSE) Computing profile confidence intervals ... 2.5 % 97.5 %
sd_(Intercept)|Class 0.40906315 0.7147376
sd_(Intercept)|school 0.00000000 0.3592051
sigma 1.68180258 1.8152490
(Intercept) 4.37617881 5.0579610
StateNSW -0.85520600 0.9123157
PrivPubpublic 0.03312581 0.7995526
sexm -0.01060332 0.3745478Q5
#read data
dta <- read.table("http://titan.ccunix.ccu.edu.tw/~psycfs/lmm/Data/iq_language.txt", header = TRUE)
head(dta) School Pupil IQ Language Group_size IQ_c School_mean Group_mean
1 1 17001 15.0 46 29 3.1659379 -1.51406 5.9
2 1 17002 14.5 45 29 2.6659379 -1.51406 5.9
3 1 17003 9.5 33 29 -2.3340621 -1.51406 5.9
4 1 17004 11.0 46 29 -0.8340621 -1.51406 5.9
5 1 17005 8.0 20 29 -3.8340621 -1.51406 5.9
6 1 17006 9.5 30 29 -2.3340621 -1.51406 5.9
SES_c
1 -4.811981
2 -17.811981
3 -12.811981
4 -4.811981
5 -17.811981
6 -17.811981str(dta) 'data.frame': 2287 obs. of 9 variables:
$ School : int 1 1 1 1 1 1 1 1 1 1 ...
$ Pupil : int 17001 17002 17003 17004 17005 17006 17007 17008 17009 17010 ...
$ IQ : num 15 14.5 9.5 11 8 9.5 9.5 13 9.5 11 ...
$ Language : int 46 45 33 46 20 30 30 57 36 36 ...
$ Group_size : int 29 29 29 29 29 29 29 29 29 29 ...
$ IQ_c : num 3.166 2.666 -2.334 -0.834 -3.834 ...
$ School_mean: num -1.51 -1.51 -1.51 -1.51 -1.51 ...
$ Group_mean : num 5.9 5.9 5.9 5.9 5.9 5.9 5.9 5.9 5.9 5.9 ...
$ SES_c : num -4.81 -17.81 -12.81 -4.81 -17.81 ...dta$School<-as.factor(dta$School)
dta$Pupil<-as.factor(dta$Pupil)
# model
lm <- lmer(data = dta, Language ~ (1|School) + IQ_c)
summary(lm) Linear mixed model fit by REML ['lmerMod']
Formula: Language ~ (1 | School) + IQ_c
Data: dta
REML criterion at convergence: 15255.8
Scaled residuals:
Min 1Q Median 3Q Max
-4.0939 -0.6375 0.0579 0.7061 3.1448
Random effects:
Groups Name Variance Std.Dev.
School (Intercept) 9.602 3.099
Residual 42.245 6.500
Number of obs: 2287, groups: School, 131
Fixed effects:
Estimate Std. Error t value
(Intercept) 40.60823 0.30819 131.8
IQ_c 2.48759 0.07008 35.5
Correlation of Fixed Effects:
(Intr)
IQ_c 0.018#95% CI
confint(lm, oldNames = FALSE) Computing profile confidence intervals ... 2.5 % 97.5 %
sd_(Intercept)|School 2.619559 3.628170
sigma 6.308665 6.697510
(Intercept) 39.997753 41.212270
IQ_c 2.349851 2.626424