Question 01

Reading the data

library(GAD)
## Loading required package: matrixStats
## Loading required package: R.methodsS3
## R.methodsS3 v1.8.2 (2022-06-13 22:00:14 UTC) successfully loaded. See ?R.methodsS3 for help.
dat <- read.csv("https://raw.githubusercontent.com/tmatis12/datafiles/main/PowderProduction.csv")

Item a

The equation for a three factor analysis is the following:

\[ y_{ijlk}=\mu+\alpha_i+\beta_j+\gamma_l+\alpha\beta_{ij}+\alpha\gamma_{il}+\beta\gamma_{jl}+\alpha\beta\gamma_{ijl}+\epsilon_{ijlk} \]

Item b

Hypothesis

The hypothesis that we are trying to reject are the following (not in order of relevance):

\[ Ho: \alpha\beta\gamma_{ijl}= 0\\ Ha: \alpha\beta\gamma_{ijl}\neq 0 \]

\[ Ho: \alpha\beta_{ij}= 0\\ Ha: \alpha\beta_{ij}\neq 0 \]

\[ Ho: \beta\gamma_{jl}= 0\\ Ha: \beta\gamma_{jl}\neq 0 \]

\[ Ho: \alpha\gamma_{il}= 0\\ Ha: \alpha\gamma_{il}\neq 0 \]

\[ Ho: \alpha_{i}= 0\\ Ha: \alpha_{i}\neq 0 \]

\[ Ho: \beta_{j}= 0\\ Ha: \beta_{j}\neq 0 \]

\[ Ho: \gamma_{l}= 0 \\ Ha: \gamma_{l}\neq 0 \]

Preparing the data

dat$Ammonium <- as.fixed(dat$Ammonium)
dat$StirRate <- as.fixed(dat$StirRate)
dat$Temperature <- as.fixed(dat$Temperature)

Statistical Analysis

To test our hypothesis, let’s perform the following:

Testing the three factor interaction

\[ Ho: \alpha\beta\gamma_{ijl}= 0\\ Ha: \alpha\beta\gamma_{ijl}\neq 0 \]

model.aov <- aov(Density~Ammonium+StirRate+Temperature+Ammonium*StirRate+Ammonium*Temperature+StirRate*Temperature+Ammonium*StirRate*Temperature, data=dat)
summary(model.aov)
##                               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

From the results above, it we cannot reject Ho, since the p-value is greater than \(\alpha\). Therefore, the interaction between the three factors is not significant.

Now, let’s test the least significant effect in order to validate the two factors interaction by removing the three factor interaction.

Testing alpha and gamma interaction

\[ Ho: \alpha\gamma_{il}= 0\\ Ha: \alpha\gamma_{il}\neq 0 \]

model.aov <- aov(Density~Ammonium+StirRate+Temperature+Ammonium*StirRate+Ammonium*Temperature+StirRate*Temperature, data=dat)
summary(model.aov)
##                      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

Since the p-value is greater than \(\alpha\), we cannot reject Ho. Therefore, we need to check the next least significant interaction by removing alpha gamma interation.

Testing beta gamma interaction

\[ Ho: \beta\gamma_{jl}= 0\\ Ha: \beta\gamma_{jl}\neq 0 \]

model.aov <- aov(Density~Ammonium+StirRate+Temperature+Ammonium*StirRate+StirRate*Temperature, data=dat)
summary(model.aov)
##                      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

Since the p-value is greater than \(\alpha\), we cannot reject Ho. Therefore, we need to test the last two factor interaction to see its relevance by removing beta gamma interaction.

Testing alpha beta interaction

\[ Ho: \alpha\beta_{ij}= 0\\ Ha: \alpha\beta_{ij}\neq 0 \]

model.aov <- aov(Density~Ammonium+StirRate+Temperature+Ammonium*StirRate, data=dat)
summary(model.aov)
##                   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

Since the p-value of the last two factor interaction is lesser than \(\alpha\) (\(p-value=0.022<\alpha=0.05\)), we can reject Ho and consider that there is a significant interaction between the Ammonium percentage and the Stir Rate.

Since there is an interaction going on, we can check the plot by

interaction.plot(dat$Ammonium,dat$StirRate,dat$Density)

From the plot, since the curves are not parallel, we can confirm the interaction between the two factors.

Question 02

Part a

\[ Ho: \alpha\beta_{ij}=0 \\ Ha: \alpha\beta_{ij}\neq0 \]

library(GAD)
position <- c(rep(1,9),rep(2,9))
temperature <- c(rep(seq(800,850,25),6))
response <- c(570,  1063,   565,
              565,  1080,   510,
              583,  1043,   590,
              528,  988,    526,
              547,  1026,   538,
              521,  1004,   532)
position <- as.fixed(position)
temperature <- as.fixed(temperature)
dat <- data.frame(position,temperature,response)
dat
##    position temperature response
## 1         1         800      570
## 2         1         825     1063
## 3         1         850      565
## 4         1         800      565
## 5         1         825     1080
## 6         1         850      510
## 7         1         800      583
## 8         1         825     1043
## 9         1         850      590
## 10        2         800      528
## 11        2         825      988
## 12        2         850      526
## 13        2         800      547
## 14        2         825     1026
## 15        2         850      538
## 16        2         800      521
## 17        2         825     1004
## 18        2         850      532
model1 <- aov(response~position+temperature+position*temperature,data=dat)
gad(model1)
## Analysis of Variance Table
## 
## Response: response
##                      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    
## Residual             12   5371     448                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

P value of the interaction is 0.427110 which is larger than \(\alpha = 0.05\) level of significance, so we don’t reject H0 hence there is no two factors interaction present. Also p value for the position and temperature is 0.001762 and 3.25e-14 respectively.

Part b

\[ Ho: \alpha\beta_{ij}=0 \\ Ha: \alpha\beta_{ij}\neq0 \]

library(GAD)
position <- c(rep(1,9),rep(2,9))
temperature <- c(rep(seq(800,850,25),6))
response <- c(570,  1063,   565,
              565,  1080,   510,
              583,  1043,   590,
              528,  988,    526,
              547,  1026,   538,
              521,  1004,   532)
position <- as.random(position)
temperature <- as.random(temperature)
dat2 <- data.frame(position,temperature,response)
dat2
##    position temperature response
## 1         1         800      570
## 2         1         825     1063
## 3         1         850      565
## 4         1         800      565
## 5         1         825     1080
## 6         1         850      510
## 7         1         800      583
## 8         1         825     1043
## 9         1         850      590
## 10        2         800      528
## 11        2         825      988
## 12        2         850      526
## 13        2         800      547
## 14        2         825     1026
## 15        2         850      538
## 16        2         800      521
## 17        2         825     1004
## 18        2         850      532
model2 <- aov(response~position+temperature+position*temperature,data=dat2)
gad(model2)
## Analysis of Variance Table
## 
## Response: response
##                      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    
## Residual             12   5371     448                       
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

P value of the interaction is 0.4271101 which is larger than \(\alpha = 0.05\) level of significance, so we don’t reject H0 hence there is no two factors interaction present. Also p value for the position and temperature is 0.0526583 and 0.0008647 respectively.

Part c

\[ Ho: \alpha\beta_{ij}=0 \\ Ha: \alpha\beta_{ij}\neq0 \]

library(GAD)
position <- c(rep(1,9),rep(2,9))
temperature <- c(rep(seq(800,850,25),6))
response <- c(570,  1063,   565,
              565,  1080,   510,
              583,  1043,   590,
              528,  988,    526,
              547,  1026,   538,
              521,  1004,   532)
position <- as.fixed(position)
temperature <- as.random(temperature)
dat3 <- data.frame(position,temperature,response)
dat3
##    position temperature response
## 1         1         800      570
## 2         1         825     1063
## 3         1         850      565
## 4         1         800      565
## 5         1         825     1080
## 6         1         850      510
## 7         1         800      583
## 8         1         825     1043
## 9         1         850      590
## 10        2         800      528
## 11        2         825      988
## 12        2         850      526
## 13        2         800      547
## 14        2         825     1026
## 15        2         850      538
## 16        2         800      521
## 17        2         825     1004
## 18        2         850      532
model3 <- aov(response~position+temperature+position*temperature,data=dat3)
gad(model3)
## Analysis of Variance Table
## 
## Response: response
##                      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    
## Residual             12   5371     448                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

P value of the interaction is 0.427110 which is larger than \(\alpha = 0.05\) level of significance, so we don’t reject H0 hence there is no two factors interaction present. Also p value for the position and temperature is 0.05266 and 3.25e-14 respectively.

Part d

Comments on p values:

The interaction term is not significant for all the part A, B and C.

When both factors are fixed, position and temperature both are significant at alpha= 0.01 and 0.001 level of significance.

When both factors are random, position and temperature both are significant at alpha= 0.1 and 0.001 level of significance.

When the factors have mixed effect, position and temperature both are significant at alpha= 0.1 and 0.001 level of significance.

But, it is possible to see a difference in the p-values of the main effects:

For fixed models, we are getting different p-values when compared with the random models. However, for mixed models, we are getting smaller p-values for the main effects.

Temperature p-value is the same for Fixed and Mixed Models, but different in the random

Position p-value is the same for mixed and Random but different for fixed