## Health data
health = read.table("/home/gilles/biogillesg/teaching/02418/2013/week5/exercises/data/mlr07.txt",sep=";",header=TRUE)

## EDA
pairs(health)

cor(health)


## Fitting full model
full = lm(data=health,formula=death.rate ~ 1 + doctor + hospital + income  + pop.density)
full

summary(full)

## Fitting null model
null = lm(data=health,formula=death.rate ~ 1)

null

summary(null)

## Variable selection
step.back = step(full,direction='backward',
                 scope=list(lower=(death.rate ~ 1),
                            upper=(death.rate ~ .)))
summary(step.back)
step.back$anova

step.for = step(null,direction='forward',
                scope=list(lower=(death.rate ~ 1),
                           upper=(death.rate ~ .)))
step.for$anova