We will simulate two variations on a balanced design with 5 treatments per block and 20 blocks. In the first blocks have been chosen, perversely, to maximize within plot variance.
Inflate the within block variance so that in most simulations the within block residual variance will exceed the between block variance:
set.seed(31)
## Add variance that is common across all blocks
df0 <- data.frame(block=rep(1:20, rep(5,20)), trt=rep(1:5,20),
y=rnorm(100,4,1))
df0$y <- df0$y + df0$trt*1.01 ## Differences between treatments
## Add extra within block variation
df <- df0
for(i in 1:20)df$y[df$block==i] <- df$y[df$block==i]+rnorm(5,0,sd=4.5)
df$block = factor(df$block); df$trt =factor(df$trt)
y.aov <- aov(y ~ trt +Error(block), data=df)
summary(y.aov)
##
## Error: block
## Df Sum Sq Mean Sq F value Pr(>F)
## Residuals 19 275.5 14.5
##
## Error: Within
## Df Sum Sq Mean Sq F value Pr(>F)
## trt 4 418 104.50 4.606 0.0022
## Residuals 76 1724 22.69
## NB that the between block variance is much smaller than the
## residual within block variance
library(lme4)
## Loading required package: Matrix
y.lmer <- lmer(y~trt + (1|block), data=df)
## boundary (singular) fit: see help('isSingular')
## Add between blocks variance
df <- df0 ## Start with same df0 as above
## Add extra between block variation
blocksvar <- rnorm(20,0,sd=5)
with(df, for(i in 1:20)df$y[block==i]<- y[block==i]+blocksvar[i])
df$block = factor(df$block); df$trt =factor(df$trt)
y.aov <- aov(y ~ trt +Error(block), data=df)
summary(y.aov)
##
## Error: block
## Df Sum Sq Mean Sq F value Pr(>F)
## Residuals 19 23.05 1.213
##
## Error: Within
## Df Sum Sq Mean Sq F value Pr(>F)
## trt 4 184.81 46.2 57.72 <2e-16
## Residuals 76 60.84 0.8
## Between block component of variance is (1.21-0.8)/5 = 0.082
y.lmer <- lmer(y~trt + (1|block), data=df)
y.lmer
## Linear mixed model fit by REML ['lmerMod']
## Formula: y ~ trt + (1 | block)
## Data: df
## REML criterion at convergence: 271.3305
## Random effects:
## Groups Name Std.Dev.
## block (Intercept) 0.2872
## Residual 0.8947
## Number of obs: 100, groups: block, 20
## Fixed Effects:
## (Intercept) trt2 trt3 trt4 trt5
## 4.881 1.382 2.312 3.031 3.946
## 0.287^2 = 0.082
The variance and components of variance estimates vary greatly from one simulation to the next. They are themselves random variables. Repeated simulations of the type shown above have an important role in making clear the extent of the randomness.