mango <- read.table("http://www.stat.ufl.edu/~winner/data/mangowine.dat",header=F, col.names=c("runnum","c_temp","c_pH","c_inoc","ethanol","glycerol","acidity", "one","temperature","pH","inoculum")) attach(mango) mango.1 <- lm(ethanol ~ temperature + pH + inoculum + I(temperature*pH) + I(temperature*inoculum) + I(pH*inoculum) + I(temperature^2) + I(pH^2) + I(inoculum^2)) summary(mango.1) anova(mango.1) Beta1 <- coef(mango.1)[2:4] int.coef <- coef(mango.1)[5:7] quad.coef <- coef(mango.1)[8:10] Beta2 <- matrix(c(quad.coef[1],0.5*int.coef[1:2],0.5*int.coef[1],quad.coef[2], 0.5*int.coef[3],0.5*int.coef[2:3],quad.coef[3]), ncol=3) stat.pnt <- -0.5 * solve(Beta2) %*% Beta1 stat.pnt install.packages("rsm") library(rsm) mango.2 <- rsm(ethanol ~ SO(temperature,pH,inoculum)) summary(mango.2) win.graph(height=5.5, width=7.0) par(mfrow=c(1,3)) contour(mango.2, ~ temperature + pH + inoculum, at=summary(mango.2)$canonical$xs) win.graph(height=5.5, width=7.0) par(mfrow=c(1,3)) persp(mango.2, ~ temperature + pH + inoculum, at=summary(mango.2)$canonical$xs)