File:RbdBlock.png

Summary

set.seed(1234)
block <- as.factor(rep(1:5, each=6))
treatment <- rep(c("A","B","C"),5)
block_shift <- rnorm(5, mean = 0, sd = 2)
treatment_shift <- c(A=0, B=4, C=2)
random_effect <- rnorm(30, mean = 0, sd = 1)

yield <- rnorm(30, mean = 10, sd = 2) + 
  treatment_shift[as.integer(factor(treatment))] + 
  block_shift[as.numeric(block)] +
  random_effect

data <- data.frame(block, treatment, yield)
summary(fm1 <- aov(yield ~ treatment + block, data = data))
#             Df Sum Sq Mean Sq F value   Pr(>F)    
# treatment    2  92.09   46.05   9.555 0.000954 ***
# block        4 140.45   35.11   7.286 0.000607 ***
# Residuals   23 110.84    4.82     

interaction.plot(x.factor = data$block,
                 trace.factor = data$treatment, 
                 response = data$yield, fun = median,
                 col=c("green","red", "blue"),
                 xlab="Block",
                 ylab="Yield", trace.label="Treatment",
                 lty=1:3, lwd=2.5 )