Analysis of a Split Plot Experiment
Savanna Torres
5/12/2020
The Experiment
A meat science researcher designed a study to investigate the impact of varying the proportion of grain in the daily ration for cattle on the tenderness of beef steaks obtained from the cattle. Twelve steers of the same breed, age, and weight were selected for the study. Four of the steers were randomly assigned to each of the following three rations (factor A, Ration)
- A1, 75% grain, 25% hay
- A2, 50% grain, 50% hay
- A3, 25% grain, 75% hay
After being fed the Ration for 90 days, the steers were butchered and four sirloin steaks were obtained from each carcass. The steaks were then randomly assigned to one of four aging times (factor B, Age): 1, 7, 14, or 21 days. After being stored at 1â—¦C. for 90 days, the steaks were then thawed and cooked to an internal temperature of 70â—¦C. Next, 68 cores (1.27cm diameter) were removed parallel to fiber orientation from each steak, and the peak shear force was measured on each core using a Warner-Bratzler shearing device. The mean shear force (kg) for each steak is the response variable.
Data
data <- read.delim("http://faculty.business.utsa.edu/manderso/STA4723/EXAMS/2020-Spring/Tenderness12.txt")
View(data)
head(data)## Ration Steer Age ShearForce carcass
## 1 A1 1 1 3.1 1
## 2 A1 1 7 2.9 1
## 3 A1 1 14 2.4 1
## 4 A1 1 21 2.1 1
## 5 A1 2 1 3.2 2
## 6 A1 2 7 2.1 2str(data)## 'data.frame': 48 obs. of 5 variables:
## $ Ration : Factor w/ 3 levels "A1","A2","A3": 1 1 1 1 1 1 1 1 1 1 ...
## $ Steer : int 1 1 1 1 2 2 2 2 3 3 ...
## $ Age : int 1 7 14 21 1 7 14 21 1 7 ...
## $ ShearForce: num 3.1 2.9 2.4 2.1 3.2 2.1 2.5 2.1 4.9 4.1 ...
## $ carcass : int 1 1 1 1 2 2 2 2 3 3 ...data$Age <- as.factor(data$Age)
data$Steer <- as.numeric(data$Steer)
str(data)## 'data.frame': 48 obs. of 5 variables:
## $ Ration : Factor w/ 3 levels "A1","A2","A3": 1 1 1 1 1 1 1 1 1 1 ...
## $ Steer : num 1 1 1 1 2 2 2 2 3 3 ...
## $ Age : Factor w/ 4 levels "1","7","14","21": 1 2 3 4 1 2 3 4 1 2 ...
## $ ShearForce: num 3.1 2.9 2.4 2.1 3.2 2.1 2.5 2.1 4.9 4.1 ...
## $ carcass : int 1 1 1 1 2 2 2 2 3 3 ...Analysis
(a) Which is the whole plot factor, and which is the split plot factor? Ration is the whole plot factor since we are investigating the impact of varying the proportion of grain in the daily ration for cattle on the tenderness of beef steaks obtained from the cattle. Four of the steers were randomly assigned to each of of the following three rations.
Age is the split plot factor. After being fed the Ration for 90 days, the steers were butchered and four sirloin steaks were obtained from each carcass. The steaks were then randomly assigned to one of four aging times.
(b) What are the mean shear forces for each ration-aging combination?
x <- aov(ShearForce ~ Ration:Age, data = data)
model.tables(x, "mean")## Tables of means
## Grand mean
##
## 4.19375
##
## Ration:Age
## Age
## Ration 1 7 14 21
## A1 4.300 3.575 3.350 3.225
## A2 4.625 4.550 4.325 3.850
## A3 4.925 4.775 4.550 4.275input = ("Ration 1 7 14 21
A1 4.300 3.575 3.350 3.225
A2 4.625 4.550 4.325 3.850
A3 4.925 4.775 4.550 4.275")
Means = as.matrix(read.table(textConnection(input), header = TRUE, row.names = 1))
Means## X1 X7 X14 X21
## A1 4.300 3.575 3.350 3.225
## A2 4.625 4.550 4.325 3.850
## A3 4.925 4.775 4.550 4.275barplot(Means, beside = TRUE, legend = TRUE, ylim = c(0,10), xlab = "Ration:Age", ylab = "Shear Force")
(c) Show a Haase diagram for this design.
(d) Use ANOVA to analyze this design, and determine which factors are statistically significant.
Hints: Which variable identifies an experimental unit, Steer or carcass? Treat Age as a CLASS or factor variable to get the correct degrees of freedom.
library(agricolae)## Warning: package 'agricolae' was built under R version 3.6.3model <- with(data, sp.plot(Steer, Ration, Age, ShearForce))##
## ANALYSIS SPLIT PLOT: ShearForce
## Class level information
##
## Ration : A1 A2 A3
## Age : 1 7 14 21
## Steer : 1 2 3 4
##
## Number of observations: 48
##
## Analysis of Variance Table
##
## Response: ShearForce
## Df Sum Sq Mean Sq F value Pr(>F)
## Steer 3 7.4323 2.4774 0.4822 0.70670
## Ration 2 8.7988 4.3994 0.8562 0.47085
## Ea 6 30.8296 5.1383
## Age 3 4.4723 1.4908 27.3697 2.425e-08 ***
## Ration:Age 6 0.7246 0.1208 2.2172 0.07223 .
## Eb 27 1.4706 0.0545
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## cv(a) = 54.1 %, cv(b) = 5.6 %, Mean = 4.19375gla <- model$gl.a
glb <- model$gl.b
Ea <- model$Ea
Eb <- model$Eb
Out1 <- with(data,LSD.test(ShearForce, Ration, gla, Ea, console = TRUE))##
## Study: ShearForce ~ Ration
##
## LSD t Test for ShearForce
##
## Mean Square Error: 5.138264
##
## Ration, means and individual ( 95 %) CI
##
## ShearForce std r LCL UCL Min Max
## A1 3.61250 1.2862737 16 2.225851 4.999149 2.1 6.0
## A2 4.33750 1.0812801 16 2.950851 5.724149 2.1 5.9
## A3 4.63125 0.4142765 16 3.244601 6.017899 3.8 5.3
##
## Alpha: 0.05 ; DF Error: 6
## Critical Value of t: 2.446912
##
## least Significant Difference: 1.961018
##
## Treatments with the same letter are not significantly different.
##
## ShearForce groups
## A3 4.63125 a
## A2 4.33750 a
## A1 3.61250 abar.err(Out1$means, variation = "SE", ylim = c(0,10))
title(cex.main = 0.8, main = "Shear Force based on Ration", xlab = "Ration", ylab = "Shear Force")
Out2 <- with(data,LSD.test(ShearForce, Age, glb, Eb, console = TRUE))##
## Study: ShearForce ~ Age
##
## LSD t Test for ShearForce
##
## Mean Square Error: 0.05446759
##
## Age, means and individual ( 95 %) CI
##
## ShearForce std r LCL UCL Min Max
## 1 4.616667 0.9962049 12 4.478431 4.754902 3.1 6.0
## 14 4.075000 0.9891456 12 3.936764 4.213236 2.4 5.1
## 21 3.783333 1.1043989 12 3.645098 3.921569 2.1 5.0
## 7 4.300000 1.1345804 12 4.161764 4.438236 2.1 5.8
##
## Alpha: 0.05 ; DF Error: 27
## Critical Value of t: 2.051831
##
## least Significant Difference: 0.1954947
##
## Treatments with the same letter are not significantly different.
##
## ShearForce groups
## 1 4.616667 a
## 7 4.300000 b
## 14 4.075000 c
## 21 3.783333 dbar.err(Out2$means, variation = "SE", ylim = c(0,10))
title(cex.main = 0.8, main = "Shear Force by Age", xlab = "Age", ylab = "Shear Force")
Out3 <- with(data,LSD.test(ShearForce, Ration:Age, glb, Eb, console = TRUE))##
## Study: ShearForce ~ Ration:Age
##
## LSD t Test for ShearForce
##
## Mean Square Error: 0.05446759
##
## Ration:Age, means and individual ( 95 %) CI
##
## ShearForce std r LCL UCL Min Max
## A1:1 4.300 1.4023789 4 4.060569 4.539431 3.1 6.0
## A1:14 3.350 1.2503333 4 3.110569 3.589431 2.4 5.1
## A1:21 3.225 1.4032700 4 2.985569 3.464431 2.1 5.0
## A1:7 3.575 1.3598407 4 3.335569 3.814431 2.1 5.2
## A2:1 4.625 1.1586630 4 4.385569 4.864431 3.1 5.9
## A2:14 4.325 0.8655441 4 4.085569 4.564431 3.1 5.1
## A2:21 3.850 1.2556539 4 3.610569 4.089431 2.1 4.9
## A2:7 4.550 1.3000000 4 4.310569 4.789431 2.8 5.8
## A3:1 4.925 0.2629956 4 4.685569 5.164431 4.7 5.3
## A3:14 4.550 0.4358899 4 4.310569 4.789431 4.2 5.1
## A3:21 4.275 0.4272002 4 4.035569 4.514431 3.8 4.8
## A3:7 4.775 0.3095696 4 4.535569 5.014431 4.5 5.2
##
## Alpha: 0.05 ; DF Error: 27
## Critical Value of t: 2.051831
##
## least Significant Difference: 0.3386067
##
## Treatments with the same letter are not significantly different.
##
## ShearForce groups
## A3:1 4.925 a
## A3:7 4.775 ab
## A2:1 4.625 abc
## A2:7 4.550 bcd
## A3:14 4.550 bcd
## A2:14 4.325 cd
## A1:1 4.300 cd
## A3:21 4.275 d
## A2:21 3.850 e
## A1:7 3.575 ef
## A1:14 3.350 fg
## A1:21 3.225 gbar.err(Out3$means, variation = "SE", ylim = c(0,10))
title(cex.main = 0.8, main = "Tenderness", xlab = "Ration:Age", ylab = "Shear Force")
(e) Lower mean shear force indicates a more tender steak. Which combination of Ration and Age yields the most tender steak?
We find that the Ration A1 that consists of 75% grain and 25% hay with the Aging Time of 21 days produces the most tender steak.