Let A = Ammonium (\(i = 1,2\)), B = Stir Rate (\(j = 1,2\)), C = Temperature (\(k = 1,2\)), with \(n = 2\) replicates (\(l = 1,2\)). All effects are fixed.
\[y_{ijkl} = \mu + \alpha_i + \beta_j + \gamma_k + (\alpha\beta)_{ij} + (\alpha\gamma)_{ik} + (\beta\gamma)_{jk} + (\alpha\beta\gamma)_{ijk} + \epsilon_{ijkl}, \quad \epsilon_{ijkl} \sim N(0,\sigma^2)\]
Ammonium <- c(2,2,30,30,2,2,30,30,2,2,30,30,2,2,30,30)
StirRate <- c(rep(100,4), rep(150,4), rep(100,4), rep(150,4))
Temperature <- c(rep(8,8), rep(40,8))
Density <- c(14.68,15.18,15.12,17.48,7.54,6.66,12.46,12.62,
10.95,17.68,12.65,15.96,8.03,8.84,14.96,14.96)
data <- data.frame(Ammonium, StirRate, Temperature, Density)
Ammonium <- as.factor(Ammonium)
StirRate <- as.factor(StirRate)
Temperature <- as.factor(Temperature)
Step 1: full model (all main effects and all interactions)
model1 <- aov(Density ~ Ammonium*StirRate*Temperature)
summary(model1)
## Df Sum Sq Mean Sq F value Pr(>F)
## Ammonium 1 44.39 44.39 11.180 0.01018 *
## StirRate 1 70.69 70.69 17.804 0.00292 **
## Temperature 1 0.33 0.33 0.083 0.78117
## Ammonium:StirRate 1 28.12 28.12 7.082 0.02875 *
## Ammonium:Temperature 1 0.02 0.02 0.005 0.94281
## StirRate:Temperature 1 10.13 10.13 2.551 0.14889
## Ammonium:StirRate:Temperature 1 1.52 1.52 0.383 0.55341
## Residuals 8 31.76 3.97
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
The three-factor interaction Ammonium:StirRate:Temperature has the largest p-value (p = 0.553), so it is dropped: there is no evidence that the two-factor interactions change with the level of the third factor.
Step 2: main effects and all two-factor interactions
model2 <- aov(Density ~ Ammonium + StirRate + Temperature + Ammonium:StirRate + Ammonium:Temperature + StirRate:Temperature)
summary(model2)
## Df Sum Sq Mean Sq F value Pr(>F)
## Ammonium 1 44.39 44.39 12.004 0.00711 **
## StirRate 1 70.69 70.69 19.115 0.00179 **
## Temperature 1 0.33 0.33 0.089 0.77268
## Ammonium:StirRate 1 28.12 28.12 7.603 0.02221 *
## Ammonium:Temperature 1 0.02 0.02 0.006 0.94054
## StirRate:Temperature 1 10.13 10.13 2.739 0.13232
## Residuals 9 33.28 3.70
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Ammonium:Temperature has the largest p-value (p = 0.941), so it is dropped.
Step 3: remove Ammonium:Temperature
model3 <- aov(Density ~ Ammonium + StirRate + Temperature + Ammonium:StirRate + StirRate:Temperature)
summary(model3)
## Df Sum Sq Mean Sq F value Pr(>F)
## Ammonium 1 44.39 44.39 13.329 0.00446 **
## StirRate 1 70.69 70.69 21.225 0.00097 ***
## Temperature 1 0.33 0.33 0.098 0.76019
## Ammonium:StirRate 1 28.12 28.12 8.443 0.01568 *
## StirRate:Temperature 1 10.13 10.13 3.041 0.11178
## Residuals 10 33.30 3.33
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
StirRate:Temperature (p = 0.112) is not significant at alpha = 0.05, so it is dropped.
Step 4: remove StirRate:Temperature
model4 <- aov(Density ~ Ammonium + StirRate + Temperature + Ammonium:StirRate)
summary(model4)
## Df Sum Sq Mean Sq F value Pr(>F)
## Ammonium 1 44.39 44.39 11.242 0.00644 **
## StirRate 1 70.69 70.69 17.903 0.00141 **
## Temperature 1 0.33 0.33 0.083 0.77861
## Ammonium:StirRate 1 28.12 28.12 7.121 0.02185 *
## Residuals 11 43.43 3.95
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Temperature (p = 0.779) is not significant and is no longer part of any interaction in the model, so it is dropped as well.
Final model
model5 <- aov(Density ~ Ammonium*StirRate)
summary(model5)
## Df Sum Sq Mean Sq F value Pr(>F)
## Ammonium 1 44.39 44.39 12.17 0.004472 **
## StirRate 1 70.69 70.69 19.38 0.000861 ***
## Ammonium:StirRate 1 28.12 28.12 7.71 0.016751 *
## Residuals 12 43.76 3.65
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
interaction.plot(Ammonium, StirRate, Density)
Summary of dropped terms
| Step | Term dropped | p-value | Reason |
|---|---|---|---|
| 1 | Ammonium:StirRate:Temperature | 0.553 | not significant (largest p-value in the full model) |
| 2 | Ammonium:Temperature | 0.941 | not significant |
| 3 | StirRate:Temperature | 0.112 | not significant |
| 4 | Temperature | 0.779 | not significant and not in any remaining interaction |
Conclusion: In the final model, Ammonium (p = 0.00447), Stir Rate (p = 0.00086) and the interaction Ammonium:StirRate (p = 0.01675) are significant at alpha = 0.05. Because the interaction is significant, the two factors must be interpreted together. The interaction plot shows that at a stir rate of 100 RPM the density is high and almost does not change with Ammonium (about 14.6 to 15.3), while at 150 RPM the density is much lower at 2% Ammonium and increases strongly at 30% (about 7.8 to 13.8). Temperature does not affect the density.
#install.packages("agricolae")
library(agricolae)
# 2^3 factorial: three factors with two levels each, 3 replications, RCBD
trt <- c(2, 2, 2)
design <- design.ab(trt, r = 3, design = "rcbd", seed = 5342)
design$book
## plots block A B C
## 1 101 1 1 1 2
## 2 102 1 2 2 2
## 3 103 1 1 2 2
## 4 104 1 2 2 1
## 5 105 1 2 1 1
## 6 106 1 1 1 1
## 7 107 1 2 1 2
## 8 108 1 1 2 1
## 9 109 2 1 1 2
## 10 110 2 2 1 2
## 11 111 2 1 1 1
## 12 112 2 2 1 1
## 13 113 2 2 2 1
## 14 114 2 1 2 2
## 15 115 2 1 2 1
## 16 116 2 2 2 2
## 17 117 3 1 1 1
## 18 118 3 1 2 2
## 19 119 3 2 2 2
## 20 120 3 2 1 1
## 21 121 3 2 2 1
## 22 122 3 1 2 1
## 23 123 3 1 1 2
## 24 124 3 2 1 2
Factor 1 is A (Fertilizer Type: Organic compost, Chemical
fertilizer), Factor 2 is B (Irrigation Frequency: once or twice per day)
and Factor 3 is C (Temperature Setting: 22 C, 28 C). The
block column of the layout identifies the 3 blocks
(greenhouse benches), and each block contains all 8 treatment
combinations in random order, so the book has 3 x 8 = 24 rows.