Broiler Chicken Drumstick Weight
# Read in data and place in memory
dswt <- read.table("http://www.stat.ufl.edu/~winner/data/drumstick.dat",
header=F,col.names=c("diet.ds","base.ds","meth.ds","weight.ds"))
attach(dswt)
head(dswt)
## diet.ds base.ds meth.ds weight.ds
## 1 1 1 1 116.33
## 2 1 1 1 99.43
## 3 1 1 1 106.58
## 4 1 1 1 109.64
## 5 1 1 1 78.58
## 6 1 1 1 93.18
# Create factor variables for base and methionine factors
base.ds.f <- factor(base.ds, levels=1:2, labels=c("sorg","corn"))
meth.ds.f <- factor(meth.ds, levels=1:2, labels=c("m.absent","m.present"))
# Mean weights by base, by meth, by base/meth
tapply(weight.ds, base.ds.f, mean)
## sorg corn
## 99.8750 104.9998
tapply(weight.ds, meth.ds.f, mean)
## m.absent m.present
## 103.6248 101.2499
tapply(weight.ds, list(base.ds.f, meth.ds.f), mean)
## m.absent m.present
## sorg 106.0798 93.67017
## corn 101.1698 108.82967
# Interaction plot (X-axis, groups, y)
interaction.plot(base.ds.f, meth.ds.f, weight.ds)
# Set options so that factor effects sum to 0
options(contrasts=c("contr.sum", "contr.poly"))
# Additive model with main effects for base, meth (no interaction)
aov.add <- aov(weight.ds ~ base.ds.f + meth.ds.f)
summary.lm(aov.add)
##
## Call:
## aov(formula = weight.ds ~ base.ds.f + meth.ds.f)
##
## Residuals:
## Min 1Q Median 3Q Max
## -54.142 -10.980 0.065 10.731 62.228
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 102.437 1.115 91.896 <2e-16 ***
## base.ds.f1 -2.562 1.115 -2.299 0.0224 *
## meth.ds.f1 1.187 1.115 1.065 0.2878
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 17.27 on 237 degrees of freedom
## Multiple R-squared: 0.02637, Adjusted R-squared: 0.01815
## F-statistic: 3.209 on 2 and 237 DF, p-value: 0.04214
anova(aov.add)
## Analysis of Variance Table
##
## Response: weight.ds
## Df Sum Sq Mean Sq F value Pr(>F)
## base.ds.f 1 1576 1575.78 5.2840 0.02239 *
## meth.ds.f 1 338 338.41 1.1348 0.28784
## Residuals 237 70678 298.22
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Interaction model with main effects and interaction and Tukey Pairwise comparisons
aov.int <- aov(weight.ds ~ base.ds.f * meth.ds.f)
summary.lm(aov.int)
##
## Call:
## aov(formula = weight.ds ~ base.ds.f * meth.ds.f)
##
## Residuals:
## Min 1Q Median 3Q Max
## -59.16 -10.94 0.78 9.58 57.21
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 102.437 1.068 95.892 < 2e-16 ***
## base.ds.f1 -2.562 1.068 -2.399 0.0172 *
## meth.ds.f1 1.187 1.068 1.112 0.2674
## base.ds.f1:meth.ds.f1 5.017 1.068 4.697 4.49e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 16.55 on 236 degrees of freedom
## Multiple R-squared: 0.1096, Adjusted R-squared: 0.09828
## F-statistic: 9.683 on 3 and 236 DF, p-value: 4.731e-06
anova(aov.int)
## Analysis of Variance Table
##
## Response: weight.ds
## Df Sum Sq Mean Sq F value Pr(>F)
## base.ds.f 1 1576 1575.8 5.7535 0.01723 *
## meth.ds.f 1 338 338.4 1.2356 0.26745
## base.ds.f:meth.ds.f 1 6042 6041.8 22.0596 4.488e-06 ***
## Residuals 236 64636 273.9
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
TukeyHSD(aov.int)
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = weight.ds ~ base.ds.f * meth.ds.f)
##
## $base.ds.f
## diff lwr upr p adj
## corn-sorg 5.12475 0.9156597 9.33384 0.0172332
##
## $meth.ds.f
## diff lwr upr p adj
## m.present-m.absent -2.374917 -6.584007 1.834174 0.2674495
##
## $`base.ds.f:meth.ds.f`
## diff lwr upr p adj
## corn:m.absent-sorg:m.absent -4.910000 -12.727911 2.9079114 0.3665304
## sorg:m.present-sorg:m.absent -12.409667 -20.227578 -4.5917553 0.0003206
## corn:m.present-sorg:m.absent 2.749833 -5.068078 10.5677447 0.7994898
## sorg:m.present-corn:m.absent -7.499667 -15.317578 0.3182447 0.0653388
## corn:m.present-corn:m.absent 7.659833 -0.158078 15.4777447 0.0571899
## corn:m.present-sorg:m.present 15.159500 7.341589 22.9774114 0.0000061
E-reader Reading Times
# Read in data (Table form)
eread <- read.table("http://www.stat.ufl.edu/~winner/data/ereader1.dat",
header=F,col.names=c("device","illum","readtime"))
attach(eread)
device <- factor(device)
illum <- factor(illum)
# Transform readtime so sums of squares aren't huge
readtime <- readtime/100
# Interaction Plot (X-axis, groups, Y)
interaction.plot(illum, device, readtime)
# Additive Model
eread.mod1 <- aov(readtime ~ device + illum)
anova(eread.mod1)
## Analysis of Variance Table
##
## Response: readtime
## Df Sum Sq Mean Sq F value Pr(>F)
## device 2 70.70 35.348 5.1987 0.0086140 **
## illum 3 148.11 49.369 7.2606 0.0003531 ***
## Residuals 54 367.17 6.800
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
TukeyHSD(eread.mod1,"device")
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = readtime ~ device + illum)
##
## $device
## diff lwr upr p adj
## 2-1 -2.206260 -4.193515 -0.219005 0.0262395
## 3-1 -2.388265 -4.375520 -0.401010 0.0148045
## 3-2 -0.182005 -2.169260 1.805250 0.9735138
TukeyHSD(eread.mod1,"illum")
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = readtime ~ device + illum)
##
## $illum
## diff lwr upr p adj
## 2-1 -1.119987 -3.644038 1.4040644 0.6442676
## 3-1 -3.361987 -5.886038 -0.8379356 0.0046176
## 4-1 -3.806987 -6.331038 -1.2829356 0.0010910
## 3-2 -2.242000 -4.766051 0.2820510 0.0984819
## 4-2 -2.687000 -5.211051 -0.1629490 0.0327612
## 4-3 -0.445000 -2.969051 2.0790510 0.9658741
# Interaction Model
eread.mod2 <- aov(readtime ~ device * illum)
anova(eread.mod2)
## Analysis of Variance Table
##
## Response: readtime
## Df Sum Sq Mean Sq F value Pr(>F)
## device 2 70.70 35.348 4.6483 0.0142790 *
## illum 3 148.11 49.369 6.4920 0.0008906 ***
## device:illum 6 2.15 0.359 0.0472 0.9995253
## Residuals 48 365.02 7.605
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
TukeyHSD(eread.mod2,"device")
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = readtime ~ device * illum)
##
## $device
## diff lwr upr p adj
## 2-1 -2.206260 -4.315285 -0.09723488 0.0384849
## 3-1 -2.388265 -4.497290 -0.27923988 0.0230557
## 3-2 -0.182005 -2.291030 1.92702012 0.9762840
TukeyHSD(eread.mod2,"illum")
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = readtime ~ device * illum)
##
## $illum
## diff lwr upr p adj
## 2-1 -1.119987 -3.799852 1.559878183 0.6838249
## 3-1 -3.361987 -6.041852 -0.682121817 0.0085849
## 4-1 -3.806987 -6.486852 -1.127121817 0.0023697
## 3-2 -2.242000 -4.921865 0.437864849 0.1306514
## 4-2 -2.687000 -5.366865 -0.007135151 0.0491639
## 4-3 -0.445000 -3.124865 2.234864849 0.9708523
Air Permeability of Jeggings
jeg <- read.csv("http://www.stat.ufl.edu/~winner/data/jegging_comfort.csv")
attach(jeg); names(jeg)
## [1] "fleeceMat" "yarnCnt" "airPerm"
fleeceMat <- factor(fleeceMat)
yarnCnt <- factor(yarnCnt)
interaction.plot(yarnCnt, fleeceMat, airPerm)
jeg.mod1 <- aov(airPerm ~ yarnCnt + fleeceMat)
anova(jeg.mod1)
## Analysis of Variance Table
##
## Response: airPerm
## Df Sum Sq Mean Sq F value Pr(>F)
## yarnCnt 1 171125 171125 97.459 < 2.2e-16 ***
## fleeceMat 8 322082 40260 22.929 < 2.2e-16 ***
## Residuals 170 298497 1756
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
TukeyHSD(jeg.mod1, "fleeceMat")
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = airPerm ~ yarnCnt + fleeceMat)
##
## $fleeceMat
## diff lwr upr p adj
## 2-1 -78.4990 -120.135919 -36.862081 0.0000006
## 3-1 49.5000 7.863081 91.136919 0.0076755
## 4-1 46.5010 4.864081 88.137919 0.0163788
## 5-1 -50.9995 -92.636419 -9.362581 0.0051578
## 6-1 47.0000 5.363081 88.636919 0.0144893
## 7-1 13.0010 -28.635919 54.637919 0.9871443
## 8-1 22.0005 -19.636419 63.637419 0.7695528
## 9-1 19.5000 -22.136919 61.136919 0.8671064
## 3-2 127.9990 86.362081 169.635919 0.0000000
## 4-2 125.0000 83.363081 166.636919 0.0000000
## 5-2 27.4995 -14.137419 69.136419 0.4936647
## 6-2 125.4990 83.862081 167.135919 0.0000000
## 7-2 91.5000 49.863081 133.136919 0.0000000
## 8-2 100.4995 58.862581 142.136419 0.0000000
## 9-2 97.9990 56.362081 139.635919 0.0000000
## 4-3 -2.9990 -44.635919 38.637919 0.9999998
## 5-3 -100.4995 -142.136419 -58.862581 0.0000000
## 6-3 -2.5000 -44.136919 39.136919 0.9999999
## 7-3 -36.4990 -78.135919 5.137919 0.1372018
## 8-3 -27.4995 -69.136419 14.137419 0.4936647
## 9-3 -30.0000 -71.636919 11.636919 0.3699858
## 5-4 -97.5005 -139.137419 -55.863581 0.0000000
## 6-4 0.4990 -41.137919 42.135919 1.0000000
## 7-4 -33.5000 -75.136919 8.136919 0.2261513
## 8-4 -24.5005 -66.137419 17.136419 0.6494738
## 9-4 -27.0010 -68.637919 14.635919 0.5195065
## 6-5 97.9995 56.362581 139.636419 0.0000000
## 7-5 64.0005 22.363581 105.637419 0.0001041
## 8-5 73.0000 31.363081 114.636919 0.0000046
## 9-5 70.4995 28.862581 112.136419 0.0000113
## 7-6 -33.9990 -75.635919 7.637919 0.2091189
## 8-6 -24.9995 -66.636419 16.637419 0.6238691
## 9-6 -27.5000 -69.136919 14.136919 0.4936389
## 8-7 8.9995 -32.637419 50.636419 0.9989845
## 9-7 6.4990 -35.137919 48.135919 0.9999098
## 9-8 -2.5005 -44.137419 39.136419 0.9999999
jeg.mod2 <- aov(airPerm ~ yarnCnt * fleeceMat)
anova(jeg.mod2)
## Analysis of Variance Table
##
## Response: airPerm
## Df Sum Sq Mean Sq F value Pr(>F)
## yarnCnt 1 171125 171125 226.676 < 2.2e-16 ***
## fleeceMat 8 322082 40260 53.330 < 2.2e-16 ***
## yarnCnt:fleeceMat 8 176198 22025 29.174 < 2.2e-16 ***
## Residuals 162 122299 755
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
TukeyHSD(jeg.mod2, "fleeceMat")
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = airPerm ~ yarnCnt * fleeceMat)
##
## $fleeceMat
## diff lwr upr p adj
## 2-1 -78.4990 -105.8180152 -51.1799848 0.0000000
## 3-1 49.5000 22.1809848 76.8190152 0.0000020
## 4-1 46.5010 19.1819848 73.8200152 0.0000103
## 5-1 -50.9995 -78.3185152 -23.6804848 0.0000009
## 6-1 47.0000 19.6809848 74.3190152 0.0000079
## 7-1 13.0010 -14.3180152 40.3200152 0.8558525
## 8-1 22.0005 -5.3185152 49.3195152 0.2246999
## 9-1 19.5000 -7.8190152 46.8190152 0.3825169
## 3-2 127.9990 100.6799848 155.3180152 0.0000000
## 4-2 125.0000 97.6809848 152.3190152 0.0000000
## 5-2 27.4995 0.1804848 54.8185152 0.0471336
## 6-2 125.4990 98.1799848 152.8180152 0.0000000
## 7-2 91.5000 64.1809848 118.8190152 0.0000000
## 8-2 100.4995 73.1804848 127.8185152 0.0000000
## 9-2 97.9990 70.6799848 125.3180152 0.0000000
## 4-3 -2.9990 -30.3180152 24.3200152 0.9999940
## 5-3 -100.4995 -127.8185152 -73.1804848 0.0000000
## 6-3 -2.5000 -29.8190152 24.8190152 0.9999985
## 7-3 -36.4990 -63.8180152 -9.1799848 0.0014200
## 8-3 -27.4995 -54.8185152 -0.1804848 0.0471336
## 9-3 -30.0000 -57.3190152 -2.6809848 0.0198075
## 5-4 -97.5005 -124.8195152 -70.1814848 0.0000000
## 6-4 0.4990 -26.8200152 27.8180152 1.0000000
## 7-4 -33.5000 -60.8190152 -6.1809848 0.0051063
## 8-4 -24.5005 -51.8195152 2.8185152 0.1176444
## 9-4 -27.0010 -54.3200152 0.3180152 0.0554146
## 6-5 97.9995 70.6804848 125.3185152 0.0000000
## 7-5 64.0005 36.6814848 91.3195152 0.0000000
## 8-5 73.0000 45.6809848 100.3190152 0.0000000
## 9-5 70.4995 43.1804848 97.8185152 0.0000000
## 7-6 -33.9990 -61.3180152 -6.6799848 0.0041568
## 8-6 -24.9995 -52.3185152 2.3195152 0.1020593
## 9-6 -27.5000 -54.8190152 -0.1809848 0.0471258
## 8-7 8.9995 -18.3195152 36.3185152 0.9817851
## 9-7 6.4990 -20.8200152 33.8180152 0.9979669
## 9-8 -2.5005 -29.8195152 24.8185152 0.9999985