Factorial Designs using R
Yi Da Lin
Problem 1
a) Model equation for the full factorial model
Let A = Ammonium (i = 1,2), B = Stir Rate (j = 1,2), C = Temperature (k = 1,2), with n = 2 replicates (l = 1,2).
\[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)\]
All effects are fixed.
b) Significant factors (alpha = 0.05)
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)
model <- aov(Density ~ Ammonium*StirRate*Temperature)
summary(model)## 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 ' ' 1Starting with the three-factor interaction: it is not significant (p = 0.553), so it is removed. Next, the two-factor interactions Ammonium:Temperature (p = 0.943) and StirRate:Temperature (p = 0.149) are not significant. Ammonium:StirRate is significant (p = 0.0288), so it stays in the model, along with both of its main effects (Ammonium, p = 0.0102 and StirRate, p = 0.0029). Temperature (p = 0.781) is not significant.
Refitting without the non-significant terms:
model2 <- aov(Density ~ Ammonium*StirRate + Temperature)
summary(model2)## 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 ' ' 1interaction.plot(Ammonium, StirRate, Density)
Conclusion: In the final model, Ammonium (p = 0.00644), Stir Rate (p = 0.00141) and their interaction Ammonium:StirRate (p = 0.02185) are significant at alpha = 0.05. Because the interaction is significant, the two factors should be interpreted together using the interaction plot: the effect of Ammonium is much larger at a stir rate of 150 RPM than at 100 RPM. Temperature (p = 0.779) is not significant.
Problem 2
Temperature <- as.factor(c(rep(800,3), rep(825,3), rep(850,3), rep(800,3), rep(825,3), rep(850,3)))
Position <- as.factor(c(rep(1,9), rep(2,9)))
Density <- c(570,565,583, 1063,1080,1043, 565,510,590,
528,547,521, 988,1026,1004, 526,538,532)
data2 <- data.frame(Position, Temperature, Density)a) Temperature and Position both fixed
Position <- as.fixed(Position)
Temperature <- as.fixed(Temperature)
model <- aov(Density ~ Position + Temperature + Position*Temperature)
GAD::gad(model)## $anova
## Analysis of Variance Table
##
## Response: Density
## Df Sum Sq Mean Sq F value Pr(>F)
## Position 1 7160 7160 15.998 0.001762 **
## Temperature 2 945342 472671 1056.117 3.25e-14 ***
## Position:Temperature 2 818 409 0.914 0.427110
## Residuals 12 5371 448
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1b) Temperature and Position both random
Position <- as.random(Position)
Temperature <- as.random(Temperature)
model <- aov(Density ~ Position + Temperature + Position*Temperature)
GAD::gad(model)## $anova
## Analysis of Variance Table
##
## Response: Density
## Df Sum Sq Mean Sq F value Pr(>F)
## Position 1 7160 7160 17.504 0.0526583 .
## Temperature 2 945342 472671 1155.518 0.0008647 ***
## Position:Temperature 2 818 409 0.914 0.4271101
## Residuals 12 5371 448
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1c) Position fixed, Temperature random
Position <- as.fixed(Position)
Temperature <- as.random(Temperature)
model <- aov(Density ~ Position + Temperature + Position*Temperature)
GAD::gad(model)## $anova
## Analysis of Variance Table
##
## Response: Density
## Df Sum Sq Mean Sq F value Pr(>F)
## Position 1 7160 7160 17.504 0.05266 .
## Temperature 2 945342 472671 1056.117 3.25e-14 ***
## Position:Temperature 2 818 409 0.914 0.42711
## Residuals 12 5371 448
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1d) Comments
| Term | a) Both fixed | b) Both random | c) Position fixed, Temperature random |
|---|---|---|---|
| Position | 0.00176 | 0.0527 | 0.0527 |
| Temperature | 3.25e-14 | 0.000865 | 3.25e-14 |
| Position:Temperature | 0.427 | 0.427 | 0.427 |
- The sums of squares and mean squares are the same in all three models. Only the denominators of the F tests change.
- Position:Temperature is always tested against MSE, so its p-value (0.427, not significant) is identical in a, b and c.
- Position is significant in (a) (tested against MSE), but not significant at 0.05 in (b) and (c), because there it is tested against the interaction MS (MSAB), which has only 2 denominator degrees of freedom.
- Temperature is highly significant in every case. In (a) and (c) it is tested against MSE (identical p-value), and in (b) against MSAB, which gives a larger but still very small p-value.
- The conclusion about Position depends on whether it is treated as fixed or random, so how the factors are modeled matters.