## Lab 3 - Script 

#Set work directory
setwd("E:/MacEwan/Teaching/Stat378/Lab/R")

# load library -- may have to be installed first
library(car)



# Display data set we are going to use
Prestige
summary(Prestige)
pairs(Prestige)
boxplot(Prestige$prestige~Prestige$type)
boxplot(Prestige$education~Prestige$type)
boxplot(Prestige$income~Prestige$type)
boxplot(Prestige$women~Prestige$type)
boxplot(Prestige$census~Prestige$type)
#identify(census~type)
Prestige[6]
Prestige[c(28,54,96,31,33),]


# remove records with missing values
Pr.complete<-Prestige[!is.na(Prestige$type),] 

# Fit the the multiple linear regression model
model.Prestige.1<-lm(prestige~education+income+women+census+type, data=Pr.complete)
summary(model.Prestige.1)
# test for effect by type
model.Prestige.r1<-lm(prestige~education+income+women+census, data=Pr.complete)
anova(model.Prestige.1,model.Prestige.r1)


# Remove "women"
model.Prestige.2<-lm(prestige~education+income+census+type, data=Pr.complete)
summary(model.Prestige.2)
# test for effect by type
model.Prestige.r2<-lm(prestige~education+income+census, data=Pr.complete)
anova(model.Prestige.2,model.Prestige.r2)



# Remove "census"
model.Prestige.3<-lm(prestige~education+income+type, data=Pr.complete)
summary(model.Prestige.3)
# test for effect by type
model.Prestige.r3<-lm(prestige~education+income, data=Pr.complete)
anova(model.Prestige.3,model.Prestige.r3)


# Add interaction between "type" and income
model.Prestige.4<-lm(prestige~education+income+type+type*income, data=Pr.complete)
summary(model.Prestige)

# test for interaction effect
model.Prestige.r4<-lm(prestige~education+income+type, data=Pr.complete)
anova(model.Prestige.4,model.Prestige.r4)

# Visualize
scatterplot(Prestige$prestige~Prestige$income, groups=Prestige$type,smooth=F, reg.line=T)
scatterplot(Prestige$prestige~Prestige$type)