##############################
#
#  assignment05a.R
#
##############################



####################
# Problem 05

wheat <- read.csv("wheat3.csv", header=TRUE)
names(wheat)


# Boxplots 

boxplot(height~strip, data=wheat, las=1, ylab="Wheat height (cm)", xlab="Strip ID")

par(mar=c(5,6,2,2)+0.1)
boxplot(height~fertilizer, data=wheat, horizontal=TRUE, las=1, xlab="Wheat height (cm)")

par(mar=c(5,7,2,2)+0.1)
boxplot(height~strip+fertilizer, data=wheat, horizontal=TRUE, las=1,
	xlab="Wheat height (cm)", ylab="")


# Assumptions of ANOVA: Normality
with(wheat, shapiro.test(height[strip=="A"]))
with(wheat, shapiro.test(height[strip=="B"]))
with(wheat, shapiro.test(height[strip=="C"]))

with(wheat, shapiro.test(height[fertilizer=="N18P51K20"]))
with(wheat, shapiro.test(height[fertilizer=="N18P18K18"]))
with(wheat, shapiro.test(height[fertilizer=="N13P00K44"]))



# Assumptions of ANOVA: Equal variance
bartlett.test(height~strip*fertilizer, data=wheat)
fligner.test(height~strip*fertilizer, data=wheat)


# Use ANOVA

model.lm.1 <- lm(height~strip*fertilizer, data=wheat)
model.aov.1 <- aov(height~strip*fertilizer, data=wheat)
summary(model.lm.1)
summary(model.aov.1)

# So, interaction is not statistically significant


# So, drop it
model.lm.2 <- lm(height~strip+fertilizer, data=wheat)
model.aov.2 <- aov(height~strip+fertilizer, data=wheat)
summary(model.lm.2)
summary(model.aov.2)
anova(model.lm.2)

model.aov.2$coefficients


# Alternatively,
strips <- c("C","A")
fertilizers <- c("N18P51K20","N18P18K18")
predict(model.lm.2, newdata=data.frame(strip=strips, fertilizer=fertilizers))

# And alternatively,
predict(model.lm.2, newdata=data.frame(strip=strips, fertilizer=fertilizers))[1]-
  predict(model.lm.2, newdata=data.frame(strip=strips, fertilizer=fertilizers))[2]