getwd()
## [1] "C:/Users/ASUS/Desktop"
Data management
# path to working directory
setwd("C:/Users/ASUS/Desktop")
# package management
library(pacman)
# load them
pacman::p_load(mlmRev, tidyverse, lme4, merTools)
# input data
dta <- read.csv("C:/Users/ASUS/Desktop/data/pts.csv")
#
head(dta)
## School Teacher Pupil Read_0 Read_1
## 1 1 11 18 -0.16407000 4.2426
## 2 1 11 19 -0.98126000 1.0000
## 3 1 11 20 -1.25370000 2.2361
## 4 1 11 15 -0.87230000 3.1623
## 5 1 11 17 -0.00063104 3.4641
## 6 1 11 16 -0.92678000 4.0000
# coerce variables to factor type and
# compute mean Read_0 by school
dta <- dta %>%
mutate(School = factor(School), Teacher = factor(Teacher),
Pupil = factor(Pupil)) %>%
group_by( School ) %>%
mutate(msRead_0 = mean(Read_0))
head(dta)
## # A tibble: 6 x 6
## # Groups: School [1]
## School Teacher Pupil Read_0 Read_1 msRead_0
## <fct> <fct> <fct> <dbl> <dbl> <dbl>
## 1 1 11 18 -0.164 4.24 -0.314
## 2 1 11 19 -0.981 1 -0.314
## 3 1 11 20 -1.25 2.24 -0.314
## 4 1 11 15 -0.872 3.16 -0.314
## 5 1 11 17 -0.000631 3.46 -0.314
## 6 1 11 16 -0.927 4 -0.314
# compute mean Read_0 by teacher and
# centered teacher mean of Read_0 (from respective school means)
#
dta <- dta %>%
group_by( Teacher ) %>%
mutate(mtRead_0 = mean(Read_0), ctRead_0 = mean(Read_0) - msRead_0 )
Analysis
m0
# no predictor
summary(m0 <- lmer(Read_1 ~ (1 | School) + (1 | Teacher), data = dta))
## Linear mixed model fit by REML ['lmerMod']
## Formula: Read_1 ~ (1 | School) + (1 | Teacher)
## Data: dta
##
## REML criterion at convergence: 2486.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.05438 -0.65686 0.06497 0.62523 2.47233
##
## Random effects:
## Groups Name Variance Std.Dev.
## Teacher (Intercept) 0.1370 0.3701
## School (Intercept) 0.1277 0.3574
## Residual 1.3250 1.1511
## Number of obs: 777, groups: Teacher, 46; School, 20
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 3.3755 0.1065 31.69
# m0 examined (the variance of) Read_1 being explained by School and Teacher intercepts (?)
m0_s
# add a school-level predictor
summary(m0_s <- update(m0, . ~ . + msRead_0))
## boundary (singular) fit: see ?isSingular
## Linear mixed model fit by REML ['lmerMod']
## Formula: Read_1 ~ (1 | School) + (1 | Teacher) + msRead_0
## Data: dta
##
## REML criterion at convergence: 2467.4
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.06110 -0.65597 0.06822 0.62772 2.45391
##
## Random effects:
## Groups Name Variance Std.Dev.
## Teacher (Intercept) 0.1234 0.3513
## School (Intercept) 0.0000 0.0000
## Residual 1.3224 1.1500
## Number of obs: 777, groups: Teacher, 46; School, 20
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 3.39822 0.06808 49.918
## msRead_0 1.00159 0.17807 5.625
##
## Correlation of Fixed Effects:
## (Intr)
## msRead_0 0.086
## convergence code: 0
## boundary (singular) fit: see ?isSingular
# m0_s extended m0 to include an examination of the predicting effect of Read_0 per school
m0_t
# add a teacher-level predictor
summary(m0_t <- update(m0, . ~ . + mtRead_0))
## boundary (singular) fit: see ?isSingular
## Linear mixed model fit by REML ['lmerMod']
## Formula: Read_1 ~ (1 | School) + (1 | Teacher) + mtRead_0
## Data: dta
##
## REML criterion at convergence: 2434.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.07539 -0.65979 0.06707 0.63259 2.49973
##
## Random effects:
## Groups Name Variance Std.Dev.
## Teacher (Intercept) 3.296e-15 5.741e-08
## School (Intercept) 4.471e-02 2.114e-01
## Residual 1.309e+00 1.144e+00
## Number of obs: 777, groups: Teacher, 46; School, 20
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 3.40649 0.06307 54.010
## mtRead_0 1.07319 0.11495 9.336
##
## Correlation of Fixed Effects:
## (Intr)
## mtRead_0 0.024
## convergence code: 0
## boundary (singular) fit: see ?isSingular
# m0_t extended m0 to include an examination of the predicting effect of Read_0 per teacher
m0_ct
# add a teacher-level predictor away from respective school means
summary(m0_ct <- update(m0, . ~ . + ctRead_0))
## boundary (singular) fit: see ?isSingular
## Linear mixed model fit by REML ['lmerMod']
## Formula: Read_1 ~ (1 | School) + (1 | Teacher) + ctRead_0
## Data: dta
##
## REML criterion at convergence: 2454.1
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.1804 -0.6461 0.0791 0.6420 2.5272
##
## Random effects:
## Groups Name Variance Std.Dev.
## Teacher (Intercept) 2.462e-10 1.569e-05
## School (Intercept) 1.764e-01 4.200e-01
## Residual 1.311e+00 1.145e+00
## Number of obs: 777, groups: Teacher, 46; School, 20
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 3.3896 0.1029 32.932
## ctRead_0 1.1742 0.1581 7.427
##
## Correlation of Fixed Effects:
## (Intr)
## ctRead_0 0.000
## convergence code: 0
## boundary (singular) fit: see ?isSingular
# m0_ct extended m0 to include an examination of the predicting effect of ctRead_0
# which indicated the Read_0 score per teacher after subtracting Read_0 mean score
# per school (mean centering)
# The End