1 Input Data

# input data
dta1 <- read.csv('C:/Users/HANK/Desktop/HOMEWORK/Chicano.csv', stringsAsFactors = TRUE)
# load package
pacman::p_load(tidyverse, VCA, lme4, nlme)

2 variance component plot

VCA::varPlot(Score ~ Trt/Class/Pupil, 
             Data=dta1,
             YLabel=list(text="Score", 
                         side=2, 
                         cex=1),
             MeanLine=list(var=c("Trt", "Class"),
                           col=c("darkred", "salmon"), 
                           lwd=c(1, 2)))

summary(m1 <- aov(Score ~ Trt + Error(Class), data=dta1))
## 
## Error: Class
##           Df Sum Sq Mean Sq F value Pr(>F)  
## Trt        1    216     216   9.818 0.0351 *
## Residuals  4     88      22                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Error: Within
##           Df Sum Sq Mean Sq F value Pr(>F)
## Residuals 18    198      11
faraway::sumary(m2 <- lmer(Score ~ Trt + (1 | Class), data=dta1))
## Fixed Effects:
##             coef.est coef.se
## (Intercept) 4.00     1.35   
## TrtT        6.00     1.91   
## 
## Random Effects:
##  Groups   Name        Std.Dev.
##  Class    (Intercept) 1.66    
##  Residual             3.32    
## ---
## number of obs: 24, groups: Class, 6
## AIC = 130.9, DIC = 131.8
## deviance = 127.4
confint(m2, method="boot")
## Computing bootstrap confidence intervals ...
## 
## 138 message(s): boundary (singular) fit: see ?isSingular
##                2.5 %   97.5 %
## .sig01      0.000000 3.650900
## .sigma      2.238982 4.307391
## (Intercept) 1.493823 6.565783
## TrtT        2.069215 9.889960

3 The end