# --------------------------------------------------------------- # Problem 10.1 # --------------------------------------------------------------- data <- scan(what = list(verapamil=0,crf=0,calcium=0,response=0), multi.line=F) 0 0 0 1.73 0 0 0 1.57 0 0 0 1.53 0 0 0 2.10 0 0 0 1.31 0 0 0 1.45 0 0 0 1.55 0 0 0 1.75 100 0 0 2.14 100 0 0 2.24 100 0 0 2.15 100 0 0 1.87 0 100 0 4.72 0 100 0 2.82 0 100 0 2.76 0 100 0 4.44 100 100 0 4.36 100 100 0 4.05 100 100 0 6.08 100 100 0 4.58 0 0 2 3.53 0 0 2 3.13 0 0 2 3.47 0 0 2 2.99 100 0 2 3.22 100 0 2 2.89 100 0 2 3.32 100 0 2 3.56 0 100 2 13.18 0 100 2 14.26 0 100 2 15.24 0 100 2 11.18 100 100 2 19.53 100 100 2 16.46 100 100 2 17.89 100 100 2 14.69 data Prb10.1 <- data.frame(data) Prb10.1$verapamil <- as.factor(Prb10.1$verapamil) Prb10.1$calcium <- as.factor(Prb10.1$calcium) Prb10.1$crf <- as.factor(Prb10.1$crf) Prb10.1.lm1 <- lm(response ~ verapamil*crf*calcium, data = Prb10.1) # check model assumptions par(mfrow=c(2,2)) plot(Prb10.1.lm1) par(mfrow=c(1,1)) Prb10.1$lresp <- log(Prb10.1$response) Prb10.1.lm2 <- lm(lresp ~ verapamil*crf*calcium, data = Prb10.1) library(car) # this produces the results using Type II SS Anova(Prb10.1.lm2,type="II") # check model assumptions par(mfrow=c(2,2)) plot(Prb10.1.lm2) par(mfrow=c(1,1)) # examine the interaction interaction.plot(Prb10.1$crf,Prb10.1$calcium, Prb10.1$response,type="b") # --------------------------------------------------------------- # Exercise 11.4 # --------------------------------------------------------------- data <- scan(what = list(head=0,weight=0), multi.line=F) 1 15.70 2 15.69 3 15.75 4 15.68 5 15.65 1 15.68 2 15.71 3 15.82 4 15.66 5 15.60 1 15.64 3 15.75 4 15.59 1 15.60 3 15.71 3 15.84 data Ex11.4 <- data.frame(data) Ex11.4$head <- as.factor(Ex11.4$head) # to get variance component estimates library(nlme) Ex11.4.lme <- lme(weight ~ 1, data = Ex11.4, random = ~ 1 | head ) # to get the ANOVA table anova(lm(weight ~ head, data = Ex11.4)) par(mfrow=c(2,2)) plot(lm(weight ~ head, data = Ex11.4)) par(mfrow=c(1,1)) # --------------------------------------------------------------- # Problem 12.2 # --------------------------------------------------------------- # first read the data directly from Dr. Oehlerts website Prb12.2 <- read.table("http://www.stat.umn.edu/~gary/book/fcdae.data/pr12.2",header=TRUE) Prb12.2$alloy <- as.factor(Prb12.2$alloy) Prb12.2$method <- as.factor(Prb12.2$method) Prb12.2$dentist <- as.factor(Prb12.2$dentist) # get the ANOVA table aov(hardness ~ alloy*method*dentist, data = Prb12.2) # check model assumptions (assuming fixed effects, but still useful) par(mfrow=c(2,2)) plot(lm(hardness ~ alloy*method*dentist, data = Prb12.2)) par(mfrow=c(1,1)) # The afex library calls lmer and other functions to give results fairly # similar to those in Proc MIXED library(afex) Prb12.2.afex1 <- mixed(hardness ~ alloy +method +alloy:method +(1|dentist) +(1|dentist:alloy) +(1|dentist:method), data = Prb12.2) summary(Prb12.2.afex1) nice(Prb12.2.afex1) lsmeans(Prb12.2.afex1, ~method) # examine the interaction interaction.plot(Prb12.2$crf,Prb12.2$dentist, Prb12.2$method,type="b") # --------------------------------------------------------------- # Exercise 13.2 # --------------------------------------------------------------- # first read the data directly from Dr. Oehlerts website Ex13.2 <- read.table("http://www.stat.umn.edu/~gary/book/fcdae.data/ex13.2",header=TRUE) Ex13.2 Ex13.2$batch <- as.factor(Ex13.2$batch) Ex13.2$temp <- as.factor(Ex13.2$temp) Ex13.2$pH <- as.factor(Ex13.2$pH) Ex13.2$sugar <- as.factor(Ex13.2$sugar) # perform RCB analysis with a cRF treatment structure Ex13.2.lm1 <- lm(y ~ batch +temp*pH*sugar, data=Ex13.2) anova(Ex13.2.lm1) # check model assumptions par(mfrow=c(2,2)) plot(Ex13.2.lm1) par(mfrow=c(1,1)) # --------------------------------------------------------------- # Exercise 13.3 # --------------------------------------------------------------- # first read the data directly from Dr. Oehlerts website Ex13.3 <- read.table("http://www.stat.umn.edu/~gary/book/fcdae.data/ex13.3",header=TRUE) Ex13.3 Ex13.3$row <- as.factor(Ex13.3$row) Ex13.3$col <- as.factor(Ex13.3$col) Ex13.3$trt <- as.factor(Ex13.3$trt) # perform Latin Square analysis Ex13.3.lm1 <- lm(y ~ row +col +trt, data=Ex13.3) anova(Ex13.3.lm1) # check model assumptions par(mfrow=c(2,2)) plot(Ex13.3.lm1) par(mfrow=c(1,1)) # Tukey tests library(multcomp) # Tukey's multiple comparison method trt.Tukey <- glht(Ex13.3.lm1,linfct=mcp(trt="Tukey")) summary(trt.Tukey) plot(trt.Tukey,sub="Tukeys method")