if (!suppressWarnings(require("GFD"))) {
    install.packages("GFD")
    library("GFD")
}

# two-sample problem
data("weightgain", package = "HSAUR")
weightgain2 <- subset(weightgain, type == "High")
set.seed(123)
two_sample <- GFD(weightgain ~ source, data = weightgain2, nperm = 10000,
    alpha = 0.05)
pdf("plot_two-sample.pdf", width = 5)
plot(two_sample, main = "Two-sample test", cex.axis = 1.5, cex.lab = 1.5, cex.main = 1.5, lwd = 2)
dev.off()
two_sample

t.test(weightgain ~ source, data = weightgain2)

# one-way layout
data("startup", package = "GFD")
set.seed(456)
model1 <- GFD(Costs ~ company, data = startup, nperm = 10000, alpha = 0.05)
summary(model1)
pdf("Startup.pdf", width = 10)
plot(model1, main = "Startup Costs", cex.axis = 1.5, cex.lab = 1.5, cex.main = 1.5, lwd = 2)
dev.off()

# two-way layout
set.seed(789)
model2 <- GFD(weightgain ~ source * type, data = weightgain)
summary(model2)
pdf("weightgain_interaction.pdf", width = 5)
plot(model2, factor = "source:type", main = "Interaction", xlab = "Type", cex.axis = 1.5, cex.lab = 1.5, cex.main = 1.5)
dev.off()
pdf("weightgain_source.pdf", width = 5)
plot(model2, factor = "source", main = "Mean weight gain", xlab = "source",  cex.axis = 1.5, cex.lab = 1.5, cex.main = 1.5)
dev.off()

# three-way layout
data("pizza", package = "GFD")
set.seed(1234)
model3 <- GFD(Delivery ~ Crust * Coke * Bread, data = pizza)
summary(model3)
pdf("pizza_threeway.pdf", width = 10)
plot(model3, factor = "Crust:Coke:Bread", legendpos = "center", main = "Delivery time of pizza", xlab = "Bread",  cex.axis = 1.5, cex.lab = 1.5, cex.main = 1.5, lwd = 2)
dev.off()
pdf("pizza_interaction.pdf", width = 10)
plot(model3, factor = "Crust:Coke", legendpos = "topleft", main = "Two-way interaction", xlab = "Coke",  cex.axis = 1.5, cex.lab = 1.5, cex.main = 1.5, lwd = 2)
dev.off()

# nested design
data("curdies", package = "GFD")
set.seed(987)
nested <- GFD(dugesia ~ season + season:site, data = curdies)
summary(nested)
pdf("nested_interaction.pdf", width = 10)
plot(nested, factor="season:site", xlab = "site",  cex.axis = 1.5, cex.lab = 1.5, cex.main = 1.5, lwd = 2)
dev.off()