# R commands for random effect model (CRD with a single factor)

y <- c(98, 91, 96, 95, 97, 90, 95, 96, 99, 93, 97, 99, 96, 92, 95, 98)
looms <- as.factor(rep(1:4,times=4))

qwe <- rbind(y[looms==1],y[looms==2],y[looms==3],y[looms==4])
dimnames(qwe) <- list(c("loom1", "loom2", "loom3","loom4"),NULL)
qwe

g <- lm(y~looms)
anova(g)

# 95% CI on the intraclass correlation coefficient
u <- qf(.975,3,12); l <- qf(.025,3,12)
F0 <- 15.681
L <- (F0/u-1)/4; U <- (F0/l-1)/4
L/(1+L); U/(1+U)

# 95% CI on sigma.sqd.epsilon

u <- qchisq(.975,12); l <- qchisq(.025,12)
df.MSE <- 12*1.896
L <- df.MSE/u; U <- df.MSE/l
L; U