boston <- read.table("boston.dat", col.names=c("CRIM","NOX","AGE","DIS",
                                               "PTRATIO","LSTAT","MEDV"))

set.seed(seed)
chosen.rows <- sample(dim(boston)[1],30)
boston <- boston[chosen.rows,]



library(leaps)  # functions for the leaps-and-bounds algorithm;
                # may need to install this using:  install.packages("leaps")

boston.model <- regsubsets(MEDV ~ CRIM + NOX + AGE + DIS + PTRATIO + LSTAT,
                           data=boston, nbest=20)
summary(boston.model)
plot(boston.model, scale="r2")


full.model <- lm(MEDV ~ CRIM + NOX + AGE + DIS + PTRATIO + LSTAT, data=boston)
start.model <- lm(MEDV ~ 1, data=boston)

step(start.model, MEDV ~ CRIM + NOX + AGE + DIS + PTRATIO + LSTAT,
     direction="both")

step(start.model, MEDV ~ CRIM + NOX + AGE + DIS + PTRATIO + LSTAT,
     direction="forward")

step(full.model, direction="backward")