Statistical Analysis of Split Plot Experiments using R
NE Milla, Jr.
2026-02-26
Loading libraries
Analysis of Split-Plot Experiments in RCBD
An experiment was carried out to compare the effects of three concentrations of a chemical seed dressing on the yield of oats. Three varieties were used in the trial because it is suspected that the response to the seed treatment would depend on the variety used. A split design was laid out in five randomized blocks. Varieties were assigned at random to the main plots within each block, and the seed treatment and control were assigned at random to the subplots within each main plot.
Importing the data
## # A tibble: 6 × 4
## Block Variety Concentration Yield
## <dbl> <dbl> <dbl> <dbl>
## 1 1 1 1 42.9
## 2 1 1 4 44.4
## 3 1 1 3 49.5
## 4 1 1 2 53.8
## 5 1 2 3 59.8
## 6 1 2 1 53.3ANOVA for split plot experiments
with(spplot, split2.rbd(factor1 = Variety,
factor2 = Concentration,
block = Block,
resp = Yield,
quali = c(TRUE, TRUE),
fac.names = c("Variety", "Concentration"),
mcomp = "tukey",
unfold = "1"))## ------------------------------------------------------------------------
## Legend:
## FACTOR 1 (plot): Variety
## FACTOR 2 (split-plot): Concentration
## ------------------------------------------------------------------------
##
## ------------------------------------------------------------------------
## Analysis of Variance Table
## ------------------------------------------------------------------------
## DF SS MS Fc Pr(>Fc)
## Variety 2 3631.5 1815.77 60.775 1.5e-05 ***
## Block 4 2931.3 732.82 24.528 0.000152 ***
## Error a 8 239.0 29.88
## Concentration 3 350.2 116.73 8.316 0.000248 ***
## Variety*Concentration 6 368.3 61.38 4.372 0.002049 **
## Error b 36 505.4 14.04
## Total 59 8025.7
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## ------------------------------------------------------------------------
## CV 1 = 10.6146 %
## CV 2 = 7.275885 %
##
## No significant interaction: analyzing the simple effects
## ------------------------------------------------------------------------
## Variety
## Tukey's test
## ------------------------------------------------------------------------
## Groups Treatments Means
## a 3 60.255
## b 2 52.88
## c 1 41.35
## ------------------------------------------------------------------------
##
## Concentration
## Tukey's test
## ------------------------------------------------------------------------
## Groups Treatments Means
## a 2 55.45333
## b 3 51.22
## b 4 50.29333
## b 1 49.01333
## ------------------------------------------------------------------------
##
##
##
##
## Significant interaction: analyzing the interaction
## ------------------------------------------------------------------------
##
## Analyzing Variety inside of each level of Concentration
## ------------------------------------------------------------------------
## DF SS MS Fc p.value
## Variety : Concentration 1 2.00000 1573.2013 786.60067 43.70568 0.000000
## Variety : Concentration 2 2.00000 495.6253 247.81267 13.76915 0.000048
## Variety : Concentration 3 2.00000 422.6680 211.33400 11.74229 0.000148
## Variety : Concentration 4 2.00000 1508.3253 754.16267 41.90333 0.000000
## Pooled Error 32.22141 579.9106 17.99768 NA NA
## ------------------------------------------------------------------------
##
##
## Variety inside of Concentration 1
## ------------------------------------------------------------------------
## Tukey's test
## ------------------------------------------------------------------------
## Groups Treatments Means
## a 3 60.84
## b 2 50.34
## c 1 35.86
## ------------------------------------------------------------------------
##
## Variety inside of Concentration 2
## ------------------------------------------------------------------------
## Tukey's test
## ------------------------------------------------------------------------
## Groups Treatments Means
## a 3 62.78
## b 2 54.84
## b 1 48.74
## ------------------------------------------------------------------------
##
## Variety inside of Concentration 3
## ------------------------------------------------------------------------
## Tukey's test
## ------------------------------------------------------------------------
## Groups Treatments Means
## a 3 56.96
## a 2 52.54
## b 1 44.16
## ------------------------------------------------------------------------
##
## Variety inside of Concentration 4
## ------------------------------------------------------------------------
## Tukey's test
## ------------------------------------------------------------------------
## Groups Treatments Means
## a 3 60.44
## b 2 53.8
## c 1 36.64
## ------------------------------------------------------------------------
##
##
## Analyzing Concentration inside of each level of Variety
## ------------------------------------------------------------------------
## DF SS MS Fc p.value
## Concentration : Variety 1 3 574.1620 191.38733 13.633626 0.000004
## Concentration : Variety 2 3 56.2760 18.75867 1.336288 0.277818
## Concentration : Variety 3 3 88.0455 29.34850 2.090663 0.118618
## Error b 36 505.3640 14.03789 NA NA
## ------------------------------------------------------------------------
##
##
## Concentration inside of Variety 1
## ------------------------------------------------------------------------
## Tukey's test
## ------------------------------------------------------------------------
## Groups Treatments Means
## a 2 48.74
## a 3 44.16
## b 4 36.64
## b 1 35.86
## ------------------------------------------------------------------------
## ------------------------------------------------------------------------
##
##
## Concentration inside of Variety 2
## ------------------------------------------------------------------------
## According to F test, the means of this factor are not different.
## ------------------------------------------------------------------------
## Levels Means
## 1 1 50.34
## 2 2 54.84
## 3 3 52.54
## 4 4 53.80
## ------------------------------------------------------------------------
##
## Concentration inside of Variety 3
## ------------------------------------------------------------------------
## According to F test, the means of this factor are not different.
## ------------------------------------------------------------------------
## Levels Means
## 1 1 60.84
## 2 2 62.78
## 3 3 56.96
## 4 4 60.44
## ------------------------------------------------------------------------Another way to generate ANOVA for split plot experiments
#Converts integer codes to nominal codes
spplot <- spplot %>%
dplyr::mutate(Block = as.factor(Block),
Variety = as.factor(Variety),
Concentration = as.factor(Concentration))#Generates ANOVA
spplot.out <- aov(Yield ~ Block + Variety + Block:Variety + Concentration
+ Variety:Concentration,
data = spplot)
anova(spplot.out)## Analysis of Variance Table
##
## Response: Yield
## Df Sum Sq Mean Sq F value Pr(>F)
## Block 4 2931.3 732.82 52.2028 1.691e-14 ***
## Variety 2 3631.5 1815.77 129.3477 < 2.2e-16 ***
## Concentration 3 350.2 116.73 8.3156 0.0002483 ***
## Block:Variety 8 239.0 29.88 2.1283 0.0583461 .
## Variety:Concentration 6 368.3 61.38 4.3725 0.0020493 **
## Residuals 36 505.4 14.04
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1We need to manually correct the F ratio and p-value for Block and Variety
#Extract the Mean Square for *Block:Variety* and use this as *MSE(a)*
MSE_a <- anova(spplot.out)["Mean Sq"][4,1]
MS_Blk <- anova(spplot.out)["Mean Sq"][1,1]
MS_Var <- anova(spplot.out)["Mean Sq"][2,1]
print(c(MSE_a,MS_Blk, MS_Var))## [1] 29.87704 732.81692 1815.76850The corrected F ratio and p-value for Block is:
\[ \begin{aligned} F_{Block} &= \frac{MS_{Block}}{MSE(a)} \notag \\ &= 24.528\notag \\ p-value &= 1.5179\times 10^{-4} \end{aligned} \]
The corrected F ratio and p-value for Variety is:
\[ \begin{aligned} F_{Var} &= \frac{MS_{Var}}{MSE(a)} \notag \\ &= 60.775\notag \\ p-value &= 5.02\times 10^{-6} \end{aligned} \]
Alternative code for analysis of split plot in RCBD.
sp.out <- aov(Yield ~ Block + Variety + Error(Block/Variety)
+ Variety*Concentration,
data = spplot)
summary(sp.out)##
## Error: Block
## Df Sum Sq Mean Sq
## Block 4 2931 732.8
##
## Error: Block:Variety
## Df Sum Sq Mean Sq F value Pr(>F)
## Variety 2 3632 1815.8 60.77 1.45e-05 ***
## Residuals 8 239 29.9
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Error: Within
## Df Sum Sq Mean Sq F value Pr(>F)
## Concentration 3 350.2 116.73 8.316 0.000248 ***
## Variety:Concentration 6 368.3 61.38 4.372 0.002049 **
## Residuals 36 505.4 14.04
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Generate the visualization of the interaction effect
VC.means <- spplot |>
group_by(Variety, Concentration) |>
dplyr::summarize(Mean = mean(Yield),
SD = sd(Yield),
n = length(Yield)) |>
dplyr::mutate(SE = SD/sqrt(n))
VC.means## # A tibble: 12 × 6
## # Groups: Variety [3]
## Variety Concentration Mean SD n SE
## <fct> <fct> <dbl> <dbl> <int> <dbl>
## 1 1 1 35.9 6.27 5 2.80
## 2 1 2 48.7 7.20 5 3.22
## 3 1 3 44.2 7.10 5 3.18
## 4 1 4 36.6 6.46 5 2.89
## 5 2 1 50.3 12.7 5 5.66
## 6 2 2 54.8 9.92 5 4.44
## 7 2 3 52.5 10.1 5 4.52
## 8 2 4 53.8 7.45 5 3.33
## 9 3 1 60.8 9.57 5 4.28
## 10 3 2 62.8 5.67 5 2.54
## 11 3 3 57.0 9.78 5 4.37
## 12 3 4 60.4 10.0 5 4.47#Generate summary statistics for all treatment combinations
VC.means1 <- emmeans(spplot.out, specs=c("Variety", "Concentration"))#Assign the letters to the treatment means
VC.means.cld <- cld(VC.means1,
Letters=letters,
decreasing = TRUE) |>
dplyr::arrange(Variety, Concentration) |>
dplyr::mutate(SE1 = VC.means$SE)#Comparing the levels of Concentration for each level of Variety
ggplot(VC.means.cld, aes(x = Variety,
y = emmean,
fill = Concentration,
label = .group)) +
geom_col(position = "dodge") +
geom_errorbar(aes(ymin = emmean - SE1,
ymax = emmean + SE1),
width = 0.2,
size = 0.7,
position = position_dodge(0.9)) +
geom_text(vjust=-4.0,
position = position_dodge(0.9),
size = 2.5) +
labs(x= "Variety", y= "Mean Yield") +
lims(y = c(0, 75)) +
scale_fill_brewer(palette = "Greens") +
theme_classic()
#Comparing the levels of Variety for each level of Concentration
ggplot(VC.means.cld, aes(x = Concentration,
y = emmean,
fill = Variety,
label = .group)) +
geom_col(position = "dodge") +
geom_errorbar(aes(ymin = emmean - SE1,
ymax = emmean + SE1),
width = 0.2,
size = 0.7,
position = position_dodge(0.9)) +
geom_text(vjust=-4.0,
position = position_dodge(0.9),
size = 2.5) +
labs(x= "Concentration", y= "Mean Yield") +
lims(y = c(0, 75)) +
scale_fill_brewer(palette = "Greens") +
theme_classic()