Homework 8
Eghorieta /
11/2/2021
Question 5.4
depth<-rep(seq(1:4),9)
feedrate<-c(rep(1,12),rep(2,12),rep(3,12))
response<-c(74,79,82,99,64,68,88,104,60,73,92,96,92,98,99,104,86,104,108,110,88,88,95,99,99,104,108,114,98,99,110,111,102,95,99,107)
data.frame(feedrate,depth,response)## feedrate depth response
## 1 1 1 74
## 2 1 2 79
## 3 1 3 82
## 4 1 4 99
## 5 1 1 64
## 6 1 2 68
## 7 1 3 88
## 8 1 4 104
## 9 1 1 60
## 10 1 2 73
## 11 1 3 92
## 12 1 4 96
## 13 2 1 92
## 14 2 2 98
## 15 2 3 99
## 16 2 4 104
## 17 2 1 86
## 18 2 2 104
## 19 2 3 108
## 20 2 4 110
## 21 2 1 88
## 22 2 2 88
## 23 2 3 95
## 24 2 4 99
## 25 3 1 99
## 26 3 2 104
## 27 3 3 108
## 28 3 4 114
## 29 3 1 98
## 30 3 2 99
## 31 3 3 110
## 32 3 4 111
## 33 3 1 102
## 34 3 2 95
## 35 3 3 99
## 36 3 4 107Ho: \(\alpha_{i} = 0\) - Null Hypothesis
Ha: \(\alpha_{i} \ne 0\) - Alternative Hypothesis
Ho: \(\beta_{j} = 0\) - Null Hypothesis
Ha: \(\beta_{j} \ne 0\) - Alternative Hypothesis
Ho: \(\alpha \beta_{ij} = 0\) - Null Hypothesis
Ha: \(\alpha \beta_{ij} \ne 0\) - Alternative Hypothesis
Model Equation
\(y_{ijk} = \mu + \alpha_{i} + \beta_j + \alpha \beta_{ij} + \epsilon_{ijk}\)
library(GAD)## Loading required package: matrixStats## Loading required package: R.methodsS3## R.methodsS3 v1.8.1 (2020-08-26 16:20:06 UTC) successfully loaded. See ?R.methodsS3 for help.depth<-as.factor(depth)
feedrate<-as.factor(feedrate)
model<-aov(response~depth*feedrate)
summary(model)## Df Sum Sq Mean Sq F value Pr(>F)
## depth 3 2125.1 708.4 24.663 1.65e-07 ***
## feedrate 2 3160.5 1580.2 55.018 1.09e-09 ***
## depth:feedrate 6 557.1 92.8 3.232 0.018 *
## Residuals 24 689.3 28.7
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1First we test the interaction hypothesis and if we fail to reject the null hypothesis Ho, then we test the main effects. From the result, the interaction effects /fo/ is 3.232 with a correspondingp-value of 0.018 < \(\alpha\) = 0.05. Hence we reject Ho hypothesis that the means are equal.
Residual plot
plot(model)
There seems to be nothing unusual about the plots. Normal plot seems to be on a straight line with 2 extreme outliers on the tail ends of the data distribution. Other than that, everything is fairly normal.
Point estimates of the mean surface finish at each feed rate
model.tables(model, type="means", cterms="feedrate")## Tables of means
## Grand mean
##
## 94.33333
##
## feedrate
## feedrate
## 1 2 3
## 81.58 97.58 103.83#tapply(response, list(feedrate), mean)P-values
1.65e-07 - depth
1.09e-09 - feedrate
0.018 - depth:feedrate
Question 5.34
depth<-rep(seq(1:4),9)
feedrate<-rep(1:3, each = 12)
block<-c(rep(1:3, each = 4),rep(1:3, each = 4),rep(1:3, each = 4))
response<-c(74,79,82,99,64,68,88,104,60,73,92,96,92,98,99,104,86,104,108,110,88,88,95,99,99,104,108,114,98,99,110,111,102,95,99,107)
data.frame(depth,feedrate,block,response)## depth feedrate block response
## 1 1 1 1 74
## 2 2 1 1 79
## 3 3 1 1 82
## 4 4 1 1 99
## 5 1 1 2 64
## 6 2 1 2 68
## 7 3 1 2 88
## 8 4 1 2 104
## 9 1 1 3 60
## 10 2 1 3 73
## 11 3 1 3 92
## 12 4 1 3 96
## 13 1 2 1 92
## 14 2 2 1 98
## 15 3 2 1 99
## 16 4 2 1 104
## 17 1 2 2 86
## 18 2 2 2 104
## 19 3 2 2 108
## 20 4 2 2 110
## 21 1 2 3 88
## 22 2 2 3 88
## 23 3 2 3 95
## 24 4 2 3 99
## 25 1 3 1 99
## 26 2 3 1 104
## 27 3 3 1 108
## 28 4 3 1 114
## 29 1 3 2 98
## 30 2 3 2 99
## 31 3 3 2 110
## 32 4 3 2 111
## 33 1 3 3 102
## 34 2 3 3 95
## 35 3 3 3 99
## 36 4 3 3 107Ho: \(\alpha_{i} = 0\) - Null Hypothesis
Ha: \(\alpha_{i} \ne 0\) - Alternative Hypothesis
Ho: \(\beta_{j} = 0\) - Null Hypothesis
Ha: \(\beta_{j} \ne 0\) - Alternative Hypothesis
Ho: \(\alpha \beta_{ij} = 0\) - Null Hypothesis
Ha: \(\alpha \beta_{ij} \ne 0\) - Alternative Hypothesis
Model Equation
\(y_{ijkl} = \mu + \alpha_{i} + \beta_{j} + \alpha \beta_{ij} + \gamma_{k} + \epsilon_{ijkl}\)
In the model equation, Gamma_k represents the block in the model and the block does not interact with the factors.
depth<-as.factor(depth)
feedrate<-as.factor(feedrate)
block<-as.factor(block)
bmodel<-aov(response~depth*feedrate+block)
summary(bmodel)## Df Sum Sq Mean Sq F value Pr(>F)
## depth 3 2125.1 708.4 30.637 4.89e-08 ***
## feedrate 2 3160.5 1580.2 68.346 3.64e-10 ***
## block 2 180.7 90.3 3.907 0.03532 *
## depth:feedrate 6 557.1 92.8 4.015 0.00726 **
## Residuals 22 508.7 23.1
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1estimate the variance component for blocks. Does it appear that blocking was useful in this experiment?
mblk<-90.3
mse<-23.1
blk<-(mblk-mse)/(3*4)
Source<-c("block")
data.frame(Source,blk)## Source blk
## 1 block 5.6The block is a very different set of experimental condition at which the experiment is run under and it has no interaction with the factors in the model. Hence it cannot influence the response because its not one of the factors we are studying. What the blocking does is that it reduces the sums of squares and also reduces the degrees of freedom in the error.
To test the hypothesis, First we consider the interaction hypothesis and if we fail to reject the null hypothesis Ho, then we test the main effects. From the result, the interaction effects /fo/ is 4.015 as compared to previous model without blocking was 3.232. The corresponding p-value is still significant at 0.00726 < \(\alpha\) = 0.05. Therefore, we reject Ho hypothesis that the means are equal.
Residual plot
plot(bmodel)
There seems to be no difference in the residual plots as there is nothing unusual about the plots. Normal plot seems to be on a straight line with some extreme outliers on the tail ends of the data distribution. Other than that, everything is fairly normal.
Point estimates of the mean surface finish at each feed rate
model.tables(bmodel, type="means", cterms="feedrate")## Tables of means
## Grand mean
##
## 94.33333
##
## feedrate
## feedrate
## 1 2 3
## 81.58 97.58 103.83P-values
4.89e-08 - depth
3.64e-10 - feedrate
0.00726 - depth:feedrate
The p-values in the model with blocking for the different factors appeared to have increased. The P-value for block is to be ignored because its not a factor.
Question 13.5
temperature<-rep(c(800,825,850),6)
position<-c(rep(1,9),rep(2,9))
readings<-c(570,1063,565,565,1080,510,583,1043,590,528,988,526,547,1026,538,521,1004,532)
data.frame(temperature,position,readings)## temperature position readings
## 1 800 1 570
## 2 825 1 1063
## 3 850 1 565
## 4 800 1 565
## 5 825 1 1080
## 6 850 1 510
## 7 800 1 583
## 8 825 1 1043
## 9 850 1 590
## 10 800 2 528
## 11 825 2 988
## 12 850 2 526
## 13 800 2 547
## 14 825 2 1026
## 15 850 2 538
## 16 800 2 521
## 17 825 2 1004
## 18 850 2 532Ho: \(\alpha_{i} = 0\) - Null Hypothesis
Ha: \(\alpha_{i} \ne 0\) - Alternative Hypothesis
Ho: \(\beta_{j} = 0\) - Null Hypothesis
Ha: \(\beta_{j} \ne 0\) - Alternative Hypothesis
Ho: \(\alpha \beta_{ij} = 0\) - Null Hypothesis
Ha: \(\alpha \beta_{ij} \ne 0\) - Alternative Hypothesis
Model Equation
\(y_{ijk} = \mu + \alpha_{i} + \beta_j + \alpha \beta_{ij} + \epsilon_{ijk}\)
temperature<-as.fixed(temperature)
position<-as.random(position)
baked<-aov(readings~position*temperature)
gad(baked)## Analysis of Variance Table
##
## Response: readings
## Df Sum Sq Mean Sq F value Pr(>F)
## position 1 7160 7160 15.998 0.0017624 **
## 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 ' ' 1Estimate the variance components
tMSA<-7160
pMSB<-472671
bMSAB<-409
bMSE<-448
bVA<-(tMSA-bMSAB)/(3*3)
bVB<-(pMSB-bMSAB)/(2*3)
bVAB<-(bMSAB-bMSE)/(3)
bVE<-448
b<-c(bVA,bVB,bVAB,bVE)
source<-c("A","B","AB","Error")
data.frame(source,b)## source b
## 1 A 750.1111
## 2 B 78710.3333
## 3 AB -13.0000
## 4 Error 448.000013.6
partN<-rep(seq(1:10),6)
operator<-c(rep(1,30),rep(2,30))
measurements<-c(50,52,53,49,48,52,51,52,50,47,49,52,50,51,49,50,51,50,51,46,50,51,50,50,48,50,51,49,50,49,50,51,54,48,48,52,51,53,51,46,48,51,52,50,49,50,50,48,48,47,51,51,51,51,48,50,50,50,49,48)
data.frame(partN,operator,measurements)## partN operator measurements
## 1 1 1 50
## 2 2 1 52
## 3 3 1 53
## 4 4 1 49
## 5 5 1 48
## 6 6 1 52
## 7 7 1 51
## 8 8 1 52
## 9 9 1 50
## 10 10 1 47
## 11 1 1 49
## 12 2 1 52
## 13 3 1 50
## 14 4 1 51
## 15 5 1 49
## 16 6 1 50
## 17 7 1 51
## 18 8 1 50
## 19 9 1 51
## 20 10 1 46
## 21 1 1 50
## 22 2 1 51
## 23 3 1 50
## 24 4 1 50
## 25 5 1 48
## 26 6 1 50
## 27 7 1 51
## 28 8 1 49
## 29 9 1 50
## 30 10 1 49
## 31 1 2 50
## 32 2 2 51
## 33 3 2 54
## 34 4 2 48
## 35 5 2 48
## 36 6 2 52
## 37 7 2 51
## 38 8 2 53
## 39 9 2 51
## 40 10 2 46
## 41 1 2 48
## 42 2 2 51
## 43 3 2 52
## 44 4 2 50
## 45 5 2 49
## 46 6 2 50
## 47 7 2 50
## 48 8 2 48
## 49 9 2 48
## 50 10 2 47
## 51 1 2 51
## 52 2 2 51
## 53 3 2 51
## 54 4 2 51
## 55 5 2 48
## 56 6 2 50
## 57 7 2 50
## 58 8 2 50
## 59 9 2 49
## 60 10 2 48#length(measurements)Ho: \(\alpha_{i} = 0\) - Null Hypothesis
Ha: \(\alpha_{i} \ne 0\) - Alternative Hypothesis
Ho: \(\beta_{j} = 0\) - Null Hypothesis
Ha: \(\beta_{j} \ne 0\) - Alternative Hypothesis
Ho: \(\alpha \beta_{ij} = 0\) - Null Hypothesis
Ha: \(\alpha \beta_{ij} \ne 0\) - Alternative Hypothesis
Model Equation
\(y_{ijk} = \mu + \alpha_{i} + \beta_j + \alpha \beta_{ij} + \epsilon_{ijk}\)
operator<-as.fixed(operator)
partN<-as.random(partN)
MSystem<-aov(measurements~partN*operator)
gad(MSystem)## Analysis of Variance Table
##
## Response: measurements
## Df Sum Sq Mean Sq F value Pr(>F)
## partN 9 99.017 11.0019 7.3346 3.216e-06 ***
## operator 1 0.417 0.4167 0.6923 0.4269
## partN:operator 9 5.417 0.6019 0.4012 0.9270
## Residual 40 60.000 1.5000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Estimate the variance components
#install.packages("lme4")
#library("lme4")
#fm1 <- lmer(measurements ~ partN + (partN | operator), REML = TRUE)
MSA<-11.0019
MSB<-0.4167
MSAB<-0.6019
MSE<-1.5000
VA<-(MSA-MSAB)/(2*3)
VB<-(MSB-MSAB)/(10*3)
VAB<-(MSAB-MSE)/(3)
VE<-1.5000
v<-c(VA,VB,VAB,VE)
source<-c("A","B","AB","Error")
data.frame(source,v)## source v
## 1 A 1.733333333
## 2 B -0.006173333
## 3 AB -0.299366667
## 4 Error 1.500000000#**Question 5.4**
depth<-rep(seq(1:4),9)
feedrate<-c(rep(1,12),rep(2,12),rep(3,12))
response<-c(74,79,82,99,64,68,88,104,60,73,92,96,92,98,99,104,86,104,108,110,88,88,95,99,99,104,108,114,98,99,110,111,102,95,99,107)
data.frame(feedrate,depth,response)
library(GAD)
depth<-as.fixed(depth)
feedrate<-as.fixed(feedrate)
model<-aov(response~depth*feedrate)
summary(model)
#**Residual plot**
plot(model)
#**Point estimates of the mean surface finish at each feed rate**
model.tables(model, type="means", cterms="feedrate")
#**Question 5.34**
depth<-rep(seq(1:4),9)
feedrate<-rep(1:3, each = 12)
block<-c(rep(1:3, each = 4),rep(1:3, each = 4),rep(1:3, each = 4))
response<-c(74,79,82,99,64,68,88,104,60,73,92,96,92,98,99,104,86,104,108,110,88,88,95,99,99,104,108,114,98,99,110,111,102,95,99,107)
data.frame(depth,feedrate,block,response)
depth<-as.factor(depth)
feedrate<-as.factor(feedrate)
block<-as.factor(block)
bmodel<-aov(response~feedrate*depth+block)
summary(bmodel)
mblk<-90.3
mse<-23.1
blk<-(mblk-mse)/(3*4)
Source<-c("block")
data.frame(Source,blk)
#**Residual plot**
plot(bmodel)
#**Point estimates of the mean surface finish at each feed rate**
model.tables(bmodel, type="means", cterms="feedrate")
#**Question 13.5**
temperature<-rep(c(800,825,850),6)
position<-c(rep(1,9),rep(2,9))
readings<-c(570,1063,565,565,1080,510,583,1043,590,528,988,526,547,1026,538,521,1004,532)
data.frame(temperature,position,readings)
temperature<-as.fixed(temperature)
position<-as.random(position)
baked<-aov(readings~position*temperature)
gad(baked)
tMSA<-7160
pMSB<-472671
bMSAB<-409
bMSE<-448
bVA<-(tMSA-bMSAB)/(3*3)
bVB<-(pMSB-bMSAB)/(2*3)
bVAB<-(bMSAB-bMSE)/(3)
bVE<-448
b<-c(bVA,bVB,bVAB,bVE)
source<-c("A","B","AB","Error")
data.frame(source,b)
#Question **13.6**
partN<-rep(seq(1:10),6)
operator<-c(rep(1,30),rep(2,30))
measurements<-c(50,52,53,49,48,52,51,52,50,47,49,52,50,51,49,50,51,50,51,46,50,51,50,50,48,50,51,49,50,49,50,51,54,48,48,52,51,53,51,46,48,51,52,50,49,50,50,48,48,47,51,51,51,51,48,50,50,50,49,48)
h<-data.frame(partN,operator,measurements)
#length(measurements)
operator<-as.fixed(operator)
partN<-as.random(partN)
MSystem<-aov(measurements~partN*operator)
gad(MSystem)
MSA<-11.0019
MSB<-0.4167
MSAB<-0.6019
MSE<-1.5000
VA<-(MSA-MSAB)/(2*3)
VB<-(MSB-MSAB)/(10*3)
VAB<-(MSAB-MSE)/(3)
VE<-1.5000
v<-c(VA,VB,VAB,VE)
source<-c("A","B","AB","Error")
data.frame(source,v)