Modelos Econometricos aplicados no R

Prévia do material em texto

regressaolinear_sigaa
# R e Modelo de Regressão Linear - Econometria II
#Limpa todo o workspace
rm(list=ls())
dados <- read.table(“/caminho para a pasta/cps1988.txt", header=T)
#dimensão dos dados
dim(dados)
#Nome das variáveis
names(dados)
#criando uma nova coluna experiencia ao quadrado e coluna de 1s
dados$exp2 <- dados$experience^2
dados$x1 <- rep(1,nrow(dados))
#para se referir às variáveis sem precisar chamar a matriz
attach(dados)
#mostrar primeiras duas linhas da base
head(dados, n=2)
dados[1:2,]
#numero de linhas na base
nrow(dados)
#mostrar ultimas 5 linhas
tail(dados, n=5)
#mostrar linha 50
dados[50,]
#missing values na coluna smsa
sum(is.na(dados$smsa))
#média da coluna education
mean(dados$education, na.rm=TRUE)
#subconjunto dos dados para education maior que 9 e experience menor que 30, média 
de wage
sub_dados <- subset(dados, dados$education>9 & dados$experience<30)
mean(sub_dados$wage)
#media de "wage" quando "education" é igual a 6
mean(dados$wage[dados$education==6])
aggregate(wage ~ education, data=dados, mean)
aggregate(wage ~ region, data=dados, mean)
#Transformando em vetores e matrizes
y <- log(dados$wage)
x <- dados[c("x1", "education", "experience", "exp2")]
#Estatísticas
mean(y)
Página 1
regressaolinear_sigaa
median(y)
var(y)
sd(y)
summary(y)
summary(x)
#Histograma
hist(education)
hist(y)
hist(y, freq=F, ylim=c(0, 0.8)) #densidade
lines(density(y), col=4)
#Variáveis de categoria
summary(region)
prop.table(table(region))
barplot(table(region))
pie(table(region))
# 2 variáveis
xtabs(~parttime + region)
plot(parttime~ region)
plot(parttime~ ethnicity)
#Regressão
reg <- lm(log(wage) ~ education + experience + exp2)
sreg<- summary(reg)
class(sreg) # retorna a classe do objeto
names(sreg)
# Construindo matriz xlx
xlx <- t(x) %*% x
# Transformando x em matriz
x <- as.matrix(x)
xlx <- t(x) %*% x
xlxi <- solve(xlx)
beta <- xlxi %*% t(x) %*% y
# residuos: reg$res
# estimador da variância do erro: e'e/(n-k)
estsig2 <- (reg$res %*% reg$res)/(nrow(x)-ncol(x))
estsig2 <- diag(estsig2)
# erro padrão
estsig <- sqrt(estsig2)
#Variancia dos betas estimados
varbeta <- estsig2*xlxi
sebeta <- sqrt(diag(varbeta))
summary(reg)$coef[,2]
Página 2
regressaolinear_sigaa
#Calculando R2
r2 <- 1-sum(reg$res^2)/sum((y-mean(y))^2)
r2alt <- 1-(reg$res %*% reg$res)/((y-mean(y)) %*% (y-mean(y)))
#Teste F
ftest <- (r2/(ncol(x)-1))/((1-r2)/(nrow(x)-ncol(x)))
#Intervalos de Confiança
leftbeta <- beta - qt(0.975, df=nrow(x)-ncol(x))*sebeta
rightbeta <- beta + qt(0.975, df=nrow(x)-ncol(x))*sebeta
confint(reg)
icbeta <- cbind(leftbeta, rightbeta)
icbeta
Página 3