setwd("/Users/hlopes/Desktop/TEACHING/2016-2to6-Econometria-graduacaoeconomia/aulas/aula0316")


data = read.table("temco.txt",header=TRUE)


attach(data)

n = nrow(data)

# Criando uma variavel dummy (1=caracteristica de interesse; 0=qualquer outra caracteristicas)
engenharia = rep(0,n)
engenharia[DEPT==2]=1




par(mfrow=c(2,4))

plot(EDUC,SALARIO)

plot(ANOSEMP,SALARIO)

plot(EXPPREV,SALARIO)

plot(SUPER,SALARIO)

plot(DEPT,SALARIO)

plot(VENDAS,SALARIO)

plot(engenharia,SALARIO)

plot(SEXO,SALARIO)


# Algumas regressoes
reg1 = lm(SALARIO~EDUC)
reg2 = lm(SALARIO~ANOSEMP)
reg3 = lm(SALARIO~EDUC+ANOSEMP)
reg4 = lm(SALARIO~EDUC+SEXO+engenharia)
reg5 = lm(SALARIO~ANOSEMP+SEXO+engenharia)
reg6 = lm(SALARIO~EDUC+ANOSEMP+SEXO+engenharia)



par(mfrow=c(2,3))
plot(reg1$fit,reg1$res/summary(reg1)$sigma,xlab="Valores ajustados",ylab="residuos padron",main="Reg 1")
abline(h=2,lty=2);abline(h=-2,lty=2)

plot(reg2$fit,reg2$res/summary(reg2)$sigma,xlab="Valores ajustados",ylab="residuos padron",main="Reg 2")
abline(h=2,lty=2);abline(h=-2,lty=2)

plot(reg3$fit,reg3$res/summary(reg3)$sigma,xlab="Valores ajustados",ylab="residuos padron",main="Reg 2")
abline(h=2,lty=2);abline(h=-2,lty=2)

plot(reg4$fit,reg4$res/summary(reg4)$sigma,xlab="Valores ajustados",ylab="residuos padron",main="Reg 4")
abline(h=2,lty=2);abline(h=-2,lty=2)

plot(reg5$fit,reg5$res/summary(reg5)$sigma,xlab="Valores ajustados",ylab="residuos padron",main="Reg 5")
abline(h=2,lty=2);abline(h=-2,lty=2)

plot(reg6$fit,reg6$res/summary(reg6)$sigma,xlab="Valores ajustados",ylab="residuos padron",main="Reg 6")
abline(h=2,lty=2);abline(h=-2,lty=2)


# Mais algumas regressoes
reg7 = lm(SALARIO~ANOSEMP+SEXO+engenharia+EXPPREV+SUPER)
reg8 = lm(SALARIO~ANOSEMP+engenharia+EXPPREV+SUPER)
reg9 = lm(SALARIO~ANOSEMP+engenharia+SUPER)

par(mfrow=c(1,3))
plot(reg7$fit,reg7$res/summary(reg7)$sigma,xlab="Valores ajustados",ylab="residuos padron",main="Reg 7")
abline(h=2,lty=2);abline(h=-2,lty=2)

plot(reg8$fit,reg8$res/summary(reg9)$sigma,xlab="Valores ajustados",ylab="residuos padron",main="Reg 8")
abline(h=2,lty=2);abline(h=-2,lty=2)

plot(reg9$fit,reg9$res/summary(reg9)$sigma,xlab="Valores ajustados",ylab="residuos padron",main="Reg 9")
abline(h=2,lty=2);abline(h=-2,lty=2)



# Mais 2 regressoes

super = rep(0,n)
super[SUPER>=2]=1

reg10 = lm(SALARIO~ANOSEMP+engenharia+super)
reg11 = lm(SALARIO~ANOSEMP+engenharia)

plot(reg10$fit,reg10$res/summary(reg10)$sigma,xlab="Valores ajustados",ylab="residuos padron",main="Reg 10")
abline(h=2,lty=2);abline(h=-2,lty=2)

plot(reg11$fit,reg11$res/summary(reg11)$sigma,xlab="Valores ajustados",ylab="residuos padron",main="Reg 11")
abline(h=2,lty=2);abline(h=-2,lty=2)


r2.ajustado = c(summary(reg1)$adj.r,
summary(reg2)$adj.r,
summary(reg3)$adj.r,
summary(reg4)$adj.r,
summary(reg5)$adj.r,
summary(reg6)$adj.r,
summary(reg7)$adj.r,
summary(reg8)$adj.r,
summary(reg9)$adj.r,
summary(reg10)$adj.r,
summary(reg11)$adj.r)


r2 = c(summary(reg1)$r.sq,
summary(reg2)$r.sq,
summary(reg3)$r.sq,
summary(reg4)$r.sq,
summary(reg5)$r.sq,
summary(reg6)$r.sq,
summary(reg7)$r.sq,
summary(reg8)$r.sq,
summary(reg9)$r.sq,
summary(reg10)$r.sq,
summary(reg11)$r.sq)

cbind(r2,r2.ajustado)



par(mfrow=c(1,1))
plot(r2,r2.ajustado,xlim=range(r2,r2.ajustado),ylim=range(r2,r2.ajustado),col=0)
text(r2,r2.ajustado,1:11)
abline(0,1)
legend("topleft",legend=c("1: EDUC",
"2: ANOSEMP",
"3: EDUC+ANOSEMP",
"4: EDUC+SEXO+engenharia",
"5: ANOSEMP+SEXO+engenharia",
"6: EDUC+ANOSEMP+SEXO+engenharia",
"7: ANOSEMP+SEXO+engenharia+EXPPREV+SUPER",
"8: ANOSEMP+engenharia+EXPPREV+SUPER",
"9: ANOSEMP+engenharia+SUPER",
"10: ANOSEMP+engenharia+super",
"11: ANOSEMP+engenharia"))


# http://www.statmethods.net/stats/rdiagnostics.html
# ncvTest
# Computes a score test of the hypothesis of constant error variance against the alternative that 
# the error variance changes with the level of the response (fitted values), or with a linear 
# combination of predictors.
install.packages("car")
library(car)

res = round(c(ncvTest(reg1)$p,
ncvTest(reg2)$p,
ncvTest(reg3)$p,
ncvTest(reg4)$p,
ncvTest(reg5)$p,
ncvTest(reg6)$p,
ncvTest(reg7)$p,
ncvTest(reg8)$p,
ncvTest(reg9)$p,
ncvTest(reg10)$p,
ncvTest(reg11)$p),4)


cbind(r2,r2.ajustado,res)