library(mosaic) library(mosaicData) library(rgl) ## for 3d plots library(RColorBrewer) ## for colors data(SAT) ## Fit a model to predict SAT using ## (1) % taking test, and ## (2) expenditure / pupil main.effects.model <- lm(sat ~ frac + expend, data = SAT) main.effects.model ## Create a prediction function for this model f.hat <- makeFun(main.effects.model) ## A 3D plot of the three variables ## (This is just for demo purposes: ## Not expecting you will do this kind of thing yourself!) with(SAT, plot3d(x = frac, y = expend, z = sat, col = brewer.pal(7, "Blues"))) ## plotting the "regression plane" persp3d(f.hat, xlim = range(~frac, data = SAT), ylim = range(~expend, data = SAT), zlim = range(~sat, data = SAT), col = rainbow, add = TRUE) summary(main.effects.model) plot(main.effects.model, which = 1) ## Fit a model to predict SAT using % taking and expenditure / pupil ## with an interaction term interaction.model <- lm(sat ~ frac + expend + frac:expend, data = SAT) interaction.model ## Create a function object for this model f.hat <- makeFun(interaction.model) ## A 3D plot of the three variables with(SAT, plot3d(x = frac, y = expend, z = sat, col = brewer.pal(7, "Blues"))) ## the regression plane persp3d(f.hat, xlim = range(~frac, data = SAT), ylim = range(~expend, data = SAT), zlim = range(~sat, data = SAT), col = rainbow, add = TRUE) summary(interaction.model) plot(interaction.model, which = 1) ## Fitting a model which is quadratic in fraction taking SAT quad.linear.model <- lm(sat ~ frac + I(frac^2) + expend, data = SAT) quad.linear.model f.hat <- makeFun(quad.linear.model) with(SAT, plot3d(x = frac, y = expend, z = sat, col = brewer.pal(7, "Blues"))) persp3d(f.hat, xlim = range(~frac, data = SAT), ylim = range(~expend, data = SAT), zlim = range(~sat, data = SAT), col = rainbow, add = TRUE) summary(quad.linear.model) plot(quad.linear.model, which = 1) ## Fitting a model with quadratic terms for both variables ## and an interaction second.order.model <- lm(sat ~ frac + expend + I(frac^2) + expend:frac + I(expend^2), data = SAT) second.order.model f.hat <- makeFun(second.order.model) with(SAT, plot3d(x = frac, y = expend, z = sat, col = brewer.pal(7, "Blues"))) persp3d(f.hat, xlim = range(~frac, data = SAT), ylim = range(~expend, data = SAT), zlim = range(~sat, data = SAT), col = rainbow, add = TRUE) summary(second.order.model) plot(second.order.model, which = 1) quad.frac <- lm(sat ~ frac + I(frac^2), data = SAT) anova(quad.frac, second.order.model) anova(quad.linear.model, second.order.model) anova(main.effects.model, interaction.model) anova(interaction.model, second.order.model)