Problem 1

(a)

Let \(A\) denote Ammonium, \(B\) denote Stir Rate, and \(C\) denote Temperature. The fixed-effects \(2^3\) factorial model is

\[Y_{ijkl}=\mu+\alpha_i+\beta_j+\gamma_k+(\alpha\beta)_{ij}+(\alpha\gamma)_{ik}+(\beta\gamma)_{jk}+(\alpha\beta\gamma)_{ijk}+\varepsilon_{ijkl},\]

where \(i,j,k=1,2\), \(l=1,2\), and \(\varepsilon_{ijkl}{\sim}N(0,\sigma^2)\). The standard sum-to-zero constraints apply to the effects.

(b)

powder <- data.frame(Ammonium = factor(c(2,2,30,30,2,2,30,30,2,2,30,30,2,2,30,30)),StirRate = factor(c(100,100,100,100,150,150,150,150,100,100,100,100,150,150,150,150)),Temperature = factor(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))

full_model <- aov(Density ~ Ammonium * StirRate * Temperature, data = powder)
summary(full_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 ' ' 1

At \(\alpha=0.05\), Ammonium (\(p=0.0102\)), Stir Rate (\(p=0.0029\)), and Ammonium \(\times\) Stir Rate (\(p=0.0288\)) are significant in the full model.

Temperature (\(p=0.7812\)), Ammonium \(\times\) Temperature (\(p=0.9428\)), Stir Rate \(\times\) Temperature (\(p=0.1489\)), and the three-way interaction (\(p=0.5534\)) are not significant.

Interaction plots

par(mfrow = c(2,2))
with(powder, interaction.plot(Ammonium, StirRate, Density, xlab="Ammonium (%)", ylab="Mean density", trace.label="Stir rate (RPM)",main="Ammonium x Stir Rate"))

with(powder, interaction.plot(Ammonium, Temperature, Density,xlab="Ammonium (%)", ylab="Mean density",trace.label="Temperature (C)", main="Ammonium x Temperature"))

with(powder, interaction.plot(StirRate, Temperature, Density,xlab="Stir rate (RPM)", ylab="Mean density", trace.label="Temperature (C)",main="Stir Rate x Temperature"))

plot.new()

par(mfrow = c(1,1))

The Ammonium \(\times\) Stir Rate interaction is statistically significant, indicating that the effect of Ammonium depends on the Stir Rate. The other two-way interactions are not statistically significant. The three-way interaction is evaluated in the ANOVA above.

model_no_ABC <- aov(Density ~ (Ammonium + StirRate + Temperature)^2,
                    data = powder)
summary(model_no_ABC)
##                      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
model_no_AC <- aov(Density ~ Ammonium * StirRate +
                     StirRate:Temperature + Temperature, data = powder)
summary(model_no_AC)
##                      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
model_no_BC <- aov(Density ~ Ammonium * StirRate + Temperature,
                    data = powder)
summary(model_no_BC)
##                   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
final_model <- aov(Density ~ Ammonium * StirRate, data = powder)
summary(final_model)
##                   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

The three-way interaction is removed first (\(p=0.5534\) in the full model). Ammonium \(\times\) Temperature is removed next (\(p=0.9405\) in the model without the three-way interaction). Stir Rate \(\times\) Temperature is then removed (\(p=0.1323\)), followed by Temperature (\(p=0.7786\)). Ammonium \(\times\) Stir Rate remains significant, so both constituent main effects must remain in the final model.

The final model is

\[Y_{ijl}=\mu+\alpha_i+\beta_j+(\alpha\beta)_{ij}+\varepsilon_{ijl}.\]

Problem 2

A full \(2^3\) factorial design contains eight treatment combinations. With three replicates arranged in a randomized complete block design, there are \(8\times3=24\) experimental units. Each block contains all eight treatment combinations once.

library(agricolae)

set.seed(12345)

tomato <- design.ab(trt = c(2, 2, 2),r = 3,design = "rcbd",seed = 12345)

tomato$book
##    plots block A B C
## 1    101     1 2 2 1
## 2    102     1 1 2 2
## 3    103     1 2 1 2
## 4    104     1 1 1 2
## 5    105     1 2 1 1
## 6    106     1 1 2 1
## 7    107     1 2 2 2
## 8    108     1 1 1 1
## 9    109     2 1 1 2
## 10   110     2 2 2 2
## 11   111     2 1 2 1
## 12   112     2 1 2 2
## 13   113     2 1 1 1
## 14   114     2 2 1 2
## 15   115     2 2 2 1
## 16   116     2 2 1 1
## 17   117     3 2 1 2
## 18   118     3 2 2 2
## 19   119     3 2 2 1
## 20   120     3 1 2 1
## 21   121     3 2 1 1
## 22   122     3 1 2 2
## 23   123     3 1 1 1
## 24   124     3 1 1 2

Complete R Code

# Problem 1

powder <- data.frame(Ammonium = factor(c(2,2,30,30,2,2,30,30,2,2,30,30,2,2,30,30)),
StirRate = factor(c(100,100,100,100,150,150,150,150,100,100,100,100,150,150,150,150)),Temperature = factor(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))

powder

model1 <- aov(Density ~ Ammonium * StirRate * Temperature, data = powder)
summary(model1)

par(mfrow = c(2,2))
with(powder, interaction.plot(Ammonium, StirRate, Density,main = "Ammonium x Stir Rate"))

with(powder, interaction.plot(Ammonium, Temperature, Density, main = "Ammonium x Temperature"))

with(powder, interaction.plot(StirRate, Temperature, Density,main = "Stir Rate x Temperature"))

plot.new()
par(mfrow = c(1,1))


model2 <- aov(Density ~ (Ammonium + StirRate + Temperature)^2,data = powder)
summary(model2)

model3 <- aov(Density ~ Ammonium * StirRate +StirRate:Temperature + Temperature, data = powder)
summary(model3)

model4 <- aov(Density ~ Ammonium * StirRate + Temperature,data = powder)
summary(model4)


model5 <- aov(Density ~ Ammonium * StirRate, data = powder)
summary(model5)


# Problem 2

library(agricolae)

set.seed(12345)

tomato <- design.ab(trt = c(2,2,2), r = 3, design = "rcbd", seed = 12345)

tomato$book