---
title: "STAT 516 Lec 06"
format: html
---
Tensile strength data
```{r}
y <- c(68,63,65,71,66,66,126,128,133,107,110,116,
93,101,98,63,60,59,56,59,57,40,41,44)
agg <- as.factor(rep(c(rep("B",3),rep("S",3)),4))
comp <- as.factor(c(rep("st",6),rep("r",6),
rep("l",6),rep("vl",6)))
```
```{r}
cbind(y,agg,comp)
```
```{r}
lm_out <- lm(y ~ agg + comp + agg:comp)
lm_out
aggregate(y ~ agg + comp, FUN = mean)
```
```{r}
summary(lm_out)
```
```{r}
plot(lm_out, which = 2)
plot(lm_out,which = 1, add.smooth=F)
```
Get the ANOVA table quantities:
```{r}
a <- 2
b <- 4
n <- 3
y... <- mean(y)
yij. <- predict(lm_out) # returns fitted value, in this case the trt group means
SStot <- sum((y - y...)^2)
SSerror <- sum((y - yij.)^2)
SStrt <- SStot - SSerror
MSerror <- SSerror/( a*b*(n-1))
MStrt <- SStrt / (a*b - 1)
Ftest <- MStrt / MSerror
alpha <- 0.05
Fcrit <- qf(1-alpha,a*b - 1, a*b*(n-1))
pval <- 1 - pf(Ftest, a*b - 1, a*b*(n-1))
```
The `summary()` output prints the F test statistic and p-value:
```{r}
summary(lm_out)
```
Use the `anova()` function on the `lm()` output:
```{r}
anova(lm_out)
```
```{r}
interaction.plot(agg,comp,y)
interaction.plot(comp,agg,y)
```
Construct Tukey's CIs for all differences in pairs of means:
```{r}
TukeyHSD(aov(lm_out))
```
Here is how we can do Dunnett's to compare all treatments with a baseline treatment. Say the baseline is basalt aggregate and static compaction.
```{r}
agg_comp <- as.factor(paste(agg,comp,sep="x"))
library(DescTools)
DunnettTest(y ~ agg_comp, control="Bxst", conf.level=0.95)
```
Serum glucose data:
```{r}
y <- c(42.5,43.3,42.9,138.4,144.4,142.7,180.9,180.5,183.0,
39.8,40.3,41.2,132.4,132.4,130.3,176.8,173.6,174.9)
method <- as.factor(c(1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,2,2))
glc <- as.factor(c(1,1,1,2,2,2,3,3,3,1,1,1,2,2,2,3,3,3))
```
```{r}
stripchart(y ~ method + glc, method = "jitter", vertical = T)
```
Look at interaction plots:
```{r}
interaction.plot(method, glc, log(y))
interaction.plot(glc,method, log(y))
```