UN <- read.table("http://socserv.mcmaster.ca/jfox/Books/Applied-Regression-2E/datasets/UnitedNations.txt", header=TRUE)
UN <- na.omit(UN[,c("tfr", "GDPperCapita", "illiteracyFemale", "contraception", "region")])

mod1 <- lm(tfr ~ GDPperCapita + illiteracyFemale + contraception, data=UN)
summary(mod1)

library(car)
qq.plot(mod1)

rstud <- rstudent(mod1)
plot(density(rstud))
boxplot(rstud)
identify(rep(1, length(rstud)), rstud, names(rstud))

influencePlot(mod1)
av.plots(mod1, ask=FALSE)

cr.plots(mod1, ask=FALSE)

mod2 <- lm(tfr ~ log(GDPperCapita) + logit(illiteracyFemale) + contraception, data=UN)
summary(mod2)

cr.plots(mod2, ask=FALSE)

qq.plot(mod2)

rstud <- rstudent(mod2)
plot(density(rstud))

boxplot(rstud)
identify(rep(1, length(rstud)), rstud, names(rstud))

influencePlot(mod2)

av.plots(mod2, ask=FALSE)

mod3 <- lm(tfr ~ log(GDPperCapita) + logit(illiteracyFemale) + contraception + region, data=UN)
summary(mod3)
Anova(mod3)

cr.plots(mod3, ask=FALSE)

influencePlot(mod3)

qq.plot(mod3)

rstud <- rstudent(mod3)
plot(density(rstud))

boxplot(rstudent(mod3))
identify(rep(1, length(rstud)), rstud, names(rstud))

mod4 <- lm(tfr ~ (log(GDPperCapita) + logit(illiteracyFemale) + contraception)*region, data=UN)
Anova(mod4)

mod5 <- lm(tfr ~ log(GDPperCapita) + logit(illiteracyFemale) + contraception*region, data=UN)
Anova(mod5)

library(effects)
plot(allEffects(mod5, default.levels=100), ask=FALSE)