Split-plot and repeated measures ANOVA

Worked example 1

Our first worked example looks at a long term experiment to assess the effects of nitrogen application and thatch accumulation on chlorophyll content of grass (you can find it analyzed using Minitab by Stephen Arnold, Penn State University))

The experiment was laid out as a split-plot design with two blocks (replications) . Each block contained four main plots each of which contained 3 subplots. The four levels of nitrogen application were randomly allocated to the main plots within each block. The three levels of thatch accumulation (2, 5 and 8 years) were randomly allocated to the subplots within each plot. It is unclear why a split plot design was used for this - Arnold suggests it may be to avoid the frertilizer blowing over onto other plots. Small problem with design - because thatch accumulation depends on time this aspect is esssentially unreplicated - depends on particular conditions during thise periods. Theoretically would be better to use a staggered start - but then of course the experiment would take even longer to complete!!

Effect of nitrogen application and thatch
accumulation on chlorophyll content of grass
Nitrogen Date Blocks
B1 B2
N1 D1 3.8 3.9
D2 5.3 5.4
D3 5.9 4.3
N2 D1 5.2 6.0
D2 5.6 6.1
D3 5.4 6.2
N3 D1 6.0 7.0
D2 5.6 6.4
D3 7.8 7.8
N4 D1 6.8 7.9
D2 8.6 8.6
D3 8.5 8.4

Draw boxplots and assess normality
Plot out data to get a visual assessment of the treatment and block effects, and assess how appropriate (parametric) ANOVA is for the set of data.
Draw interaction plots
If there were no interaction between nitrogen treatment and date, the three lines for the different dates should follow the same trends and be roughly parallel.

Get table of means

Using R
tapply(chlo,nitr,mean) n1 n2 n3 n4 4.766667 5.750000 6.766667 8.133333 > tapply(chlo,date,mean) d1 d2 d3 5.8250 6.4500 6.7875 > tapply(chlo,blck,mean) b1 b2 6.208333 6.500000

Carry out analysis of variance

Sums of squares can be calculated manually as follows:

SS_total = 1017.79 −− 152.5² = 48.7796

24

SS_blocks(S) = 74.5² + 78² − 152.5²

12 12 24
= 0.510416

SS_nitr(A) = 28.6² + 34.5² + 40.6² + 48.8² − 152.5²

6 6 6 6 24
= 37.3246

SS_{Subgrps (A×S)} = 15.0² + 16.2² + ... 24.9² − 152.5²

3 3 3 24
= 39.09292

SS_{Main plot error} = 39.0929 − 37.3246 − 0.5104 = 1.2579

SS_date(B) = 46.6² + 51.6² + 54.3² − 152.5²

8 8 8 24
= 3.8158

SS_{Subgrps (A×B)} = 7.7² + 11.2² + ... 16.9² − 152.5²

2 2 2 24
= 45.2946

SS_{A ×B} = 45.2946 − 37.3246 − 3.8158 = 4.1542

SS_{Subplot error} = 48.7796 − 0.5104 − 37.3246 − 1.2579 − 3.8158 − 4.1542 = 1.7167

ANOVA table
Source of
variation df SS MS F-
ratio P
Blocks (S) 1 0.5104 0.5104
Nitrogen (A) 3 37.3246 12.4415 29.67 0.0099
Main plot error 3 1.2579 0.4193
Dates (B) 2 3.8158 1.9079 8.89 0.0093
A × B 6 4.1542 0.6924 3.23 0.0646
Subplot Error 8 1.7167 0.2146
Total 23 48.7796

Using R
Error: blck Df Sum Sq Mean Sq F value Pr(>F) Residuals 1 0.51042 0.51042 Error: blck:nitr Df Sum Sq Mean Sq F value Pr(>F) nitr 3 37.325 12.442 29.672 0.009896 ** Residuals 3 1.258 0.419 --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Error: blck:nitr:date Df Sum Sq Mean Sq F value Pr(>F) date 2 3.8158 1.9079 8.8913 0.00927 ** nitr:date 6 4.1542 0.6924 3.2265 0.06460 . Residuals 8 1.7167 0.2146 ---

Check diagnostics

An immediate difficulty arises with checking diagnostics of this split-plot design - if you try plotting out the diagnostics using R, you will simply get the word NULL. This is not because you have done anything wrong! It is simply the result of having only one observation for each of the block × nitrogen × date combinations. The model you have fitted is known as the saturated model.

You can get an overview of the situation by fitting the general linear model and assuming that interactions with blocks are non-existent. The sums of squares for all main effects and the A × B interaction will still be correct, and you can assess the various diagnostics. Note, however, that all the F-ratios (and associated P-values) are incorrect because they are no longer using the correct error term.

Using R
Response: chlo Df Sum Sq Mean Sq F value Pr(>F) nitr 3 37.325 12.442 46.0087 1.627e-06 *** date 2 3.816 1.908 7.0555 0.01068 * blck 1 0.510 0.510 1.8875 0.19683 nitr:date 6 4.154 0.692 2.5603 0.08396 . Residuals 11 2.975 0.270 ----

A better approach to diagnostics with a split-plot design is to consider diagnostics separately for assessment of factor A (between mainplots) and for factor B and A × B (within mainplots).

Between mainplots

Using R
Response: chlo Df Sum Sq Mean Sq F value Pr(>F) nitr 3 12.4539 4.1513 29.5239 0.009967 ** blck 1 0.1691 0.1691 1.2024 0.352972 Residuals 3 0.4218 0.1406 ---

Worked example 2

We take our second example from Ogata & Takeuchi et al (2001) on a trial of a feline pheromone analogue to reduce the frequency of urine marking by cats. We previously compared the number of markings pre-treatment and one week post-treatment using the non-parametric Wilcoxon's matched pairs signed ranks test. But the authors also used repeated measures analysis of variance to examine the different courses of urine marking over time relative to aggression status.

No. urine markings pre and post treatment
No. Sex Agg Pre W1 W2 W3 W4 W5 W6 W7 W8
01 m n 7 0 0 0 0 0 1 3 5
02 m n 12 10 9 9 9 9 9 9 9
07 m n 2 1 2 2 0 0 0 1 0
08 m y 63 8 8 3 5 6 8 10 6
09 m y 10 7 4 4 3 3 4 3 3
10 m y 16 14 8 8 6 5 9 9 11
11 m n 6 6 3 2 2 2 1 2 3
13 m y 18 9 5 2 8 9 7 7 16
14 m y 18 9 3 1 0 0 0 0 0
16 m y 9 7 5 7 4 5 6 5 6
20 m y 13 14 14 10 9 9 10 13 14
21 m y 7 2 4 0 2 1 1 0 0
22 m y 28 2 0 1 0 11 5 3 2
26 f n 3 0 0 0 0 0 0 0 0
27 f n 7 0 0 0 0 0 3 3 5
28 f n 8 9 4 2 1 1 1 2 0
29 f n 6 6 3 2 1 1 1 3 4
32 f y 13 10 9 7 11 8 9 9 11
33 f y 1 0 2 1 3 3 1 0 0
34 f n 3 3 2 2 1 2 2 1 1
35 f n 17 14 10 7 7 3 4 4 3

The data for multicat households (where aggression can be assessed) are given below (complete data sets only).

One's first reaction may be (or perhaps should be) that non-parametric analysis would be a much wiser approach given the patently non-normal distribution of the response variable. However, we will attempt an analysis after a transformation.

Draw boxplots and assess normality
Plot out data to get a visual assessment of the treatment and block effects, and assess how appropriate (parametric) ANOVA is for the set of data.

As expected for count data, the distribution of the raw data within groups does not approximate to normal - instead the distributions are right skewed. We try a square root transformation (or to be more precise a √(Y + 0.5) transform given the large number of zeros) as a possible normalizing function for small whole numbered counts.

This looks more hopeful - most of the groups have more or less symetrical distributions, albeit still with a few high outliers. A log transformation brought the high outliers down a little more, but at the cost of making several distributions left-skewed. Hence we proceed with the analysis on the square root transformed data bearing in mind we need to examine diagnostics carefully after model fitting. The interaction plot for aggression and week suggests similar trends over time for both aggressive and non-aggressive cats.