################################################################################################
# 
# SALARY DATA
#
################################################################################################

data    = read.table("salary.txt",header=TRUE)
n       = nrow(data)

y       = data[,1]
x1      = data[,2]
x2      = data[,3]
x3      = data[,4]

pdf(file="salary-regression.pdf",width=20,height=15)
par(mfrow=c(2,4))

plot(y,xlab="observation index",ylab="Salary",pch=16)
abline(h=mean(y),col=2)
title("No regression \n s=12.5")
abline(h=mean(y)-2*sqrt(var(y)),lty=2,col=2)
abline(h=mean(y)+2*sqrt(var(y)),lty=2,col=2)

reg  = lm(y~x1)
coef = reg$coef
s    = sqrt(mean(reg$res^2))
plot(x1,y,xlab="Position",ylab="Salary",pch=16)
abline(coef,col=2)
abline(coef[1]+2*s,coef[2],col=2,lty=2)
abline(coef[1]-2*s,coef[2],col=2,lty=2)
title("Salary = 113.8+5.74*Position \n s=7.2 - R2=0.67")

reg  = lm(y~x2)
coef = reg$coef
s    = sqrt(mean(reg$res^2))
plot(x2,y,xlab="Years of Experience",ylab="Salary",pch=16)
abline(coef,col=2)
abline(coef[1]+2*s,coef[2],col=2,lty=2)
abline(coef[1]-2*s,coef[2],col=2,lty=2)
title("Salary = 135.0+0.75*Years \n s=11.9 - R2=0.09")

reg  = lm(y~x3)
coef = reg$coef
s    = sqrt(mean(reg$res^2))
plot(x3,y,xlab="Gender",ylab="Salary",pch=16)
abline(coef,col=2)
abline(coef[1]+2*s,coef[2],col=2,lty=2)
abline(coef[1]-2*s,coef[2],col=2,lty=2)
title("Salary = 140.5+3.64*Gender \n s=12.4 - R2=0.02")

reg=lm(y~x1+x2)
coef = reg$coef
s    = sqrt(mean(reg$res^2))
plot(reg$fit,y,xlab="Fitted salary",ylab="Salary",ylim=range(reg$fit,y),
     xlim=range(reg$fit,y),pch=16)
abline(0,1,col=2)
abline(-2*s,1,col=2,lty=2)
abline(2*s,1,col=2,lty=2)
title("Salary = 114.9+6.65*Position-0.55*Years \n s=6.7 - R2=0.70")

reg=lm(y~x1+x3)
coef = reg$coef
s    = sqrt(mean(reg$res^2))
plot(reg$fit,y,xlab="Fitted salary",ylab="Salary",ylim=range(reg$fit,y),
     xlim=range(reg$fit,y),pch=16)
abline(0,1,col=2)
abline(-2*s,1,col=2,lty=2)
abline(2*s,1,col=2,lty=2)
title("Salary = 114.6+6.06*Position-3.72*Gender \n s=7.0 - R2=0.69")

reg=lm(y~x2+x3)
coef = reg$coef
s    = sqrt(mean(reg$res^2))
plot(reg$fit,y,xlab="Fitted salary",ylab="Salary",ylim=range(reg$fit,y),
     xlim=range(reg$fit,y),pch=16)
abline(0,1,col=2)
abline(-2*s,1,col=2,lty=2)
abline(2*s,1,col=2,lty=2)
title("Salary = 135.0+0.75*Years-0.03*Gender \n s=11.9 - R2=0.091")

reg=lm(y~x1+x2+x3)
coef = reg$coef
s    = sqrt(mean(reg$res^2))
plot(reg$fit,y,xlab="Fitted salary",ylab="Salary",ylim=range(reg$fit,y),
     xlim=range(reg$fit,y),pch=16)
abline(0,1,col=2)
abline(-2*s,1,col=2,lty=2)
abline(2*s,1,col=2,lty=2)
title("Salary = 115.3+6.71*Position-0.47*Years-2.2*Gender \n s=6.7 - R2=0.71")

dev.off()