# Three factor factorial by regression

T <- c(125, 158, 121, 160, 118, 163, 122, 165, 140, 140, 140, 140)
P <- c(41, 40, 82, 80, 39, 40, 80, 83, 60, 60, 60, 60)
C <- c(14, 15, 15, 15, 33, 30, 30, 30, rep(22.5,4))
y <- c(32, 46, 57, 65, 36, 48, 57, 68, 50, 44, 53, 56)

data <- data.frame(y, T, P, C)
data

g <- lm(y~T + P + C + T*P + P*C + I(T^2) + I(C^2)+ I(P^2))
# See help(formula) for an explanation of I().
summary(g) # To get the coefficients

h <- lm(y~T + P + C)
summary(h) # To get the coefficients

anova(g)  # To get SSE
anova(h)  # To get SSE

F0 <- ((102.22-78.75)/5)/26.25
p.value <- 1-pf(F0, 5, 3)

F0 
p.value

h$coeff

resids <- h$resid
fits <- h$fitted

par(mfrow=c(2,3))
plot(fits, resids);abline(h=0)
plot(T, resids);abline(h=0)
plot(P, resids);abline(h=0)
plot(C, resids);abline(h=0)

qqnorm(resids)
qqline(resids)