library(reshape2)
library(gplots)
library(multcomp)

DF = read.csv('../../Data/ASCII_Comma/Chapter_12/tablets2.txt', quote='\'')
DF_melt = melt(DF, id.vars=c(), variable.name='Laboratory', value.name='Measurement')

fit = aov(Measurement ~ Laboratory, data=DF_melt)
print(summary(fit))

#postscript("../../WriteUp/Graphics/Chapter12/prob_22_plotmeans.eps", onefile=FALSE, horizontal=FALSE)
plotmeans(Measurement ~ Laboratory, data=DF_melt)
#dev.off()

# A Tukey comparison of means test:
#
TukeyHSD(fit)

par(las=2)
par(mar=c(5, 8, 4, 2))
plot(TukeyHSD(fit))

# Use the multcomp package:
#
par(las=2)
par(mar=c(5, 4, 6, 2))
tuk = glht(fit, linfct=mcp(Laboratory='Tukey'))
plot(cld(tuk, level=0.05), col='lightgrey')

# A nonparametric test for one-way ANOVA:
#
print(kruskal.test(Measurement ~ Laboratory, data=DF_melt))


# Lets see if there are differences between the manufacture of these tablets and the manufacturer from before:
#
DF_melt$Manufacture = 2 # let this be the second manufacturer

DF_orig = read.csv('../../Data/ASCII_Comma/Chapter_12/tablets1.txt', quote='\'')
DF_orig_melt = melt(DF, id.vars=c(), variable.name='Laboratory', value.name='Measurement')
DF_orig_melt$Manufacture = 1 # these will be samples from the first manufacturer

DF_all_melt = rbind(DF_orig_melt, DF_melt)
fit = aov(Measurement ~ Manufacture, data=DF_all_melt)
print(summary(fit))