# Repeated measures ANOVA in R (split-split-plot design)

I am trying to fit a repeated measures ANOVA from an experiment with a split-split-plot design and several measures over time.

The experimental design is as follows:

I have 9 blocks in the field. Inside each block I have two subplots representing a split factor (named trt, "with" or "without" a specific treatment). Inside each subplot (with or without) I have another split factor with two quadrats (named sp, "species 1", "Species 2"). In each quadrat I have two seedlings of each species in the study (each seedling identified by a unique id). Finally I have monitored a given response variable for each seedling in the experiment along 4 repeated measures in time (weeks).

Therefore, I have 9 blocks, 2 treatments inside each block, 2 species inside each treatment and 2 seedlings from each species. This was monitored for 4 weeks.

I want to understand if time * trt * sp affect my response variable.

Considering my experimental design, is the following code a correct specification of the Error term for fitting an aov repeated measures split-split-plot model?

fit <- aov(response ~ time * sp * trt + Error(block/trt/sp/id), data = d3)
summary(fit)

Error: block
Df Sum Sq Mean Sq F value Pr(>F)
Residuals  1  11.29   11.29

Error: block:trt
Df Sum Sq Mean Sq
trt  1  0.114   0.114

Error: block:trt:sp
Df Sum Sq Mean Sq
sp      1  61.14   61.14
sp:trt  1  10.27   10.27

Error: block:trt:sp:id
Df Sum Sq Mean Sq F value Pr(>F)
sp         1   7.16   7.159   2.299  0.141
trt        1   1.07   1.072   0.344  0.562
sp:trt     1   2.18   2.181   0.701  0.410
Residuals 28  87.17   3.113

Error: Within
Df Sum Sq Mean Sq F value Pr(>F)
time          1   0.38  0.3781   0.237  0.627
time:sp       1   0.93  0.9317   0.585  0.446
time:trt      1   1.91  1.9117   1.201  0.276
time:sp:trt   1   2.73  2.7257   1.712  0.194
Residuals   104 165.59  1.5922

This code results in the following Warning message:

Warning message:
In aov(response ~ time * sp * trt + Error(block/trt/sp/id), data = d3) :
Error() model is singular

I really appreciate any help on this issue. Thank you very much! I am happy to provide any further detail, if need it.

## Data and plots

To illustrate the results

Data for reproducibility is presented below (dput print):

Edit 1 (updated datset):

structure(list(block = c(1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L,
3L, 3L, 3L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 7L,
7L, 7L, 7L, 8L, 8L, 8L, 8L, 9L, 9L, 9L, 9L, 1L, 1L, 1L, 1L, 2L,
2L, 2L, 2L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 6L,
6L, 6L, 6L, 7L, 7L, 7L, 7L, 8L, 8L, 8L, 8L, 9L, 9L, 9L, 9L, 1L,
1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 5L,
5L, 5L, 5L, 6L, 6L, 6L, 6L, 7L, 7L, 7L, 7L, 8L, 8L, 8L, 8L, 9L,
9L, 9L, 9L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 4L,
4L, 4L, 4L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 7L, 7L, 7L, 7L, 8L,
8L, 8L, 8L, 9L, 9L, 9L, 9L), trt = c("without", "without", "with",
"with", "with", "with", "without", "without", "with", "with",
"without", "without", "with", "with", "without", "without", "with",
"with", "without", "without", "with", "with", "without", "without",
"without", "without", "with", "with", "with", "with", "without",
"without", "without", "without", "with", "with", "without", "without",
"with", "with", "with", "with", "without", "without", "with",
"with", "without", "without", "with", "with", "without", "without",
"with", "with", "without", "without", "with", "with", "without",
"without", "without", "without", "with", "with", "with", "with",
"without", "without", "without", "without", "with", "with", "without",
"without", "with", "with", "with", "with", "without", "without",
"with", "with", "without", "without", "with", "with", "without",
"without", "with", "with", "without", "without", "with", "with",
"without", "without", "without", "without", "with", "with", "with",
"with", "without", "without", "without", "without", "with", "with",
"without", "without", "with", "with", "with", "with", "without",
"without", "with", "with", "without", "without", "with", "with",
"without", "without", "with", "with", "without", "without", "with",
"with", "without", "without", "without", "without", "with", "with",
"with", "with", "without", "without", "without", "without", "with",
"with"), sp = c("species 1", "species 2", "species 2", "species 1",
"species 2", "species 1", "species 2", "species 1", "species 1",
"species 2", "species 2", "species 1", "species 2", "species 1",
"species 1", "species 2", "species 1", "species 2", "species 2",
"species 1", "species 2", "species 1", "species 1", "species 2",
"species 1", "species 2", "species 1", "species 2", "species 1",
"species 2", "species 2", "species 1", "species 2", "species 1",
"species 2", "species 1", "species 1", "species 2", "species 2",
"species 1", "species 2", "species 1", "species 2", "species 1",
"species 1", "species 2", "species 2", "species 1", "species 2",
"species 1", "species 1", "species 2", "species 1", "species 2",
"species 2", "species 1", "species 2", "species 1", "species 1",
"species 2", "species 1", "species 2", "species 1", "species 2",
"species 1", "species 2", "species 2", "species 1", "species 2",
"species 1", "species 2", "species 1", "species 1", "species 2",
"species 2", "species 1", "species 2", "species 1", "species 2",
"species 1", "species 1", "species 2", "species 2", "species 1",
"species 2", "species 1", "species 1", "species 2", "species 1",
"species 2", "species 2", "species 1", "species 2", "species 1",
"species 1", "species 2", "species 1", "species 2", "species 1",
"species 2", "species 1", "species 2", "species 2", "species 1",
"species 2", "species 1", "species 2", "species 1", "species 1",
"species 2", "species 2", "species 1", "species 2", "species 1",
"species 2", "species 1", "species 1", "species 2", "species 2",
"species 1", "species 2", "species 1", "species 1", "species 2",
"species 1", "species 2", "species 2", "species 1", "species 2",
"species 1", "species 1", "species 2", "species 1", "species 2",
"species 1", "species 2", "species 1", "species 2", "species 2",
"species 1", "species 2", "species 1", "species 2", "species 1"
), id = structure(c(1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L,
11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 22L, 23L,
24L, 25L, 26L, 27L, 28L, 29L, 30L, 31L, 32L, 33L, 34L, 35L, 36L,
1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 13L, 14L,
15L, 16L, 17L, 18L, 19L, 20L, 21L, 22L, 23L, 24L, 25L, 26L, 27L,
28L, 29L, 30L, 31L, 32L, 33L, 34L, 35L, 36L, 1L, 2L, 3L, 4L,
5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L,
19L, 20L, 21L, 22L, 23L, 24L, 25L, 26L, 27L, 28L, 29L, 30L, 31L,
32L, 33L, 34L, 35L, 36L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L,
10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 22L,
23L, 24L, 25L, 26L, 27L, 28L, 29L, 30L, 31L, 32L, 33L, 34L, 35L,
36L), .Label = c("1", "2", "3", "4", "5", "6", "7", "8", "9",
"10", "11", "12", "13", "14", "15", "16", "17", "18", "19", "20",
"21", "22", "23", "24", "25", "26", "27", "28", "29", "30", "31",
"32", "33", "34", "35", "36"), class = "factor"), time = c(1L,
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L,
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L,
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L,
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L,
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L),
response = c(3.1, 5.7, 4.8, 2.9, 5, 3.9, 5.7, 4.2, 3.6, 4.4,
3.9, 2.9, 3.2, 7.5, 4.3, 4.8, 3, 4.9, 5.6, 3.9, 4.1, 4.2,
2.8, 3.9, 3.9, 7.5, 4, 4.3, 3.1, 7.1, 5.8, 2.5, 6.4, 4.5,
5, 3.6, 3.1, 5.2, 3.6, 2.9, 5.2, 4.6, 4.7, 4.3, 3.9, 4.4,
4.2, 3.6, 3.2, 2.7, 3.4, 5.6, 2.8, 6, 5.1, 3.7, 4.1, 3.4,
3, 4.1, 3.2, 6.7, 3.1, 3.8, 2.9, 6.9, 5.6, 2.1, 5.6, 4.8,
4.8, 2.7, 3, 5.5, 3.4, 3.1, 5.1, 5, 5, 4.8, 4, 4, 4, 2.6,
3, 3, 3.9, 6, 3, 7, 5, 3.5, 4, 4, 3, 4, 3, 6.5, 4, 5, 4,
8, 6, 2.2, 5.9, 4, 6, 3, 3, 5, 3.5, 3, 5, 4, 4, 2, 6.5, 4,
5, 2, 3, 3, 3, 5.5, 2, 5, 6, 2.5, 5, 2.5, 3, 5, 3, 5.5, 3,
2, 3, 6, 5, 5, 5, 3, 4, 15)), row.names = c(NA, -144L), class = "data.frame")