In [1]:
# rm(list=ls()) 
options(OutDec = ",")
#==============================================================================
# Exemplo
#     Mercado de gasolina nos EUA de 1953 a 2004.
#     O consumo (G) de gasolina e' calculado como gastos com gasolina dividido
#     pelo preco da gasolina.
#     O consumo de gasolina per capita e' dado pelo consumo de gasolina 
#     dividido pela tamanho da populacao.
#==============================================================================
dataxy   <- read.csv("../dados/d04_mercado_gasolina.txt",header=T)
G        <- dataxy[,2]/dataxy[,4] # consumo de gasolina
lnG.Pop  <- log(G / dataxy[,3])   # consumo per capita
lnR.Pop  <- log(dataxy[,5])       # renda per capita
lnpreco  <- log(dataxy[,4])       # preco
lnPCN    <- log(dataxy[,6])       # preco de carros novos
lnPCU    <- log(dataxy[,7])       # preco de carros usados
In [2]:
#==============================================================================
# Note que teremos uma observacao a menos por conta da defasagem
#==============================================================================
n          <- length(lnG.Pop)  
lnG.Popt   <- lnG.Pop[2:n]     
lnG.Poptm1 <- lnG.Pop[1:(n-1)]
lnR.Popt   <- lnR.Pop[2:n]
lnprecot   <- lnpreco[2:n]
lnPCNt     <- lnPCN[2:n]
lnPCUt     <- lnPCU[2:n]
yfitAR1    <- lm(lnG.Popt~lnG.Poptm1+lnprecot+lnR.Popt+lnPCNt+lnPCUt)
In [3]:
#==============================================================================
# Resumo da estimacao
#==============================================================================
print(summary(yfitAR1))
Call:
lm(formula = lnG.Popt ~ lnG.Poptm1 + lnprecot + lnR.Popt + lnPCNt + 
    lnPCUt)

Residuals:
      Min        1Q    Median        3Q       Max 
-0,054371 -0,009792  0,000249  0,012239  0,034238 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -3,12320    0,99583  -3,136 0,003013 ** 
lnG.Poptm1   0,83097    0,04576  18,158  < 2e-16 ***
lnprecot    -0,06953    0,01473  -4,720 2,33e-05 ***
lnR.Popt     0,16405    0,05503   2,981 0,004620 ** 
lnPCNt      -0,17839    0,05517  -3,233 0,002293 ** 
lnPCUt       0,12701    0,03577   3,551 0,000914 ***
---
Signif. codes:  0 ‘***’ 0,001 ‘**’ 0,01 ‘*’ 0,05 ‘.’ 0,1 ‘ ’ 1

Residual standard error: 0,01682 on 45 degrees of freedom
Multiple R-squared:  0,9951,	Adjusted R-squared:  0,9946 
F-statistic:  1831 on 5 and 45 DF,  p-value: < 2,2e-16

In [4]:
#==============================================================================
# Analise de variancia
#==============================================================================
print(anova(yfitAR1))
Analysis of Variance Table

Response: lnG.Popt
           Df  Sum Sq Mean Sq   F value    Pr(>F)    
lnG.Poptm1  1 2,57602 2,57602 9102,3908 < 2,2e-16 ***
lnprecot    1 0,00237 0,00237    8,3619  0,005882 ** 
lnR.Popt    1 0,00856 0,00856   30,2636 1,706e-06 ***
lnPCNt      1 0,00005 0,00005    0,1655  0,686107    
lnPCUt      1 0,00357 0,00357   12,6072  0,000914 ***
Residuals  45 0,01274 0,00028                        
---
Signif. codes:  0 ‘***’ 0,001 ‘**’ 0,01 ‘*’ 0,05 ‘.’ 0,1 ‘ ’ 1
In [5]:
#==============================================================================
# Estimativa da matriz de covariancia do estimador dos coeficientes
#==============================================================================
options(digits=4)
print(vcov(yfitAR1))
            (Intercept) lnG.Poptm1   lnprecot   lnR.Popt     lnPCNt     lnPCUt
(Intercept)    0,991678  0,0439149 -1,209e-03 -5,260e-02  0,0051016  9,167e-03
lnG.Poptm1     0,043915  0,0020943 -1,109e-04 -2,188e-03  0,0006812  8,570e-05
lnprecot      -0,001209 -0,0001109  2,171e-04  1,622e-05 -0,0002171 -4,055e-05
lnR.Popt      -0,052602 -0,0021881  1,622e-05  3,028e-03 -0,0002471 -6,062e-04
lnPCNt         0,005102  0,0006812 -2,171e-04 -2,471e-04  0,0030440 -1,678e-03
lnPCUt         0,009167  0,0000857 -4,055e-05 -6,062e-04 -0,0016782  1,280e-03
In [6]:
#==============================================================================
# Elasticidade a longo prazo do preco
#==============================================================================
eps2       <-  yfitAR1$coef[3]/(1-yfitAR1$coef[2])
print(eps2)
lnprecot 
 -0,4114 
In [7]:
#==============================================================================
# Elasticidade a longo prazo da renda
#==============================================================================
eps3       <-  yfitAR1$coef[4]/(1-yfitAR1$coef[2])
print(eps3)
lnR.Popt 
  0,9705 
In [8]:
#==============================================================================
# O vetor de derivadas avaliado nas estimativas e' dado por
#==============================================================================
aux1       <- 1 / (1-yfitAR1$coef[2])
aux2       <- yfitAR1$coef[3]*(aux1^2)
aux3       <- yfitAR1$coef[4]*(aux1^2)
d.eps2     <- c(0,aux2,aux1,0,0,0)
d.eps3     <- c(0,aux3,0,aux1,0,0)
s2         <- sum(yfitAR1$res^2) / yfitAR1$df.res
W          <- cbind(rep(1,n-1),lnG.Poptm1,lnprecot,lnR.Popt,lnPCNt,lnPCUt)
WtW        <- t(W)%*%W
WtWinv     <- solve(WtW)
se.eps2    <- sqrt(d.eps2%*%(s2*WtWinv)%*%d.eps2)
se.eps3    <- sqrt(d.eps3%*%(s2*WtWinv)%*%d.eps3)
print(c(se.eps2))  # erro assintotico estimado   
print(c(se.eps3))  # erro assintotico estimado

# Verificar o intervalo de confianca no fim deste arquivo para epsilon_2 e 
# epsilon_3 baseados em teoria assintotica (normalidade).
[1] 0,1523
[1] 0,1624
In [9]:
#==============================================================================
# O modelo de regressao SEM a variavel defasada (dados completos)
#==============================================================================
yfit     <- lm(lnG.Pop ~ lnpreco + lnR.Pop + lnPCN + lnPCU)
In [10]:
#==============================================================================
# Resumo da estimacao
#==============================================================================
print(summary(yfit))
Call:
lm(formula = lnG.Pop ~ lnpreco + lnR.Pop + lnPCN + lnPCU)

Residuals:
    Min      1Q  Median      3Q     Max 
-0,1383 -0,0284  0,0130  0,0337  0,0743 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -21,2111     0,7532  -28,16   <2e-16 ***
lnpreco      -0,0212     0,0438   -0,48    0,630    
lnR.Pop       1,0959     0,0777   14,10   <2e-16 ***
lnPCN        -0,3736     0,1571   -2,38    0,021 *  
lnPCU         0,0200     0,1033    0,19    0,847    
---
Signif. codes:  0 ‘***’ 0,001 ‘**’ 0,01 ‘*’ 0,05 ‘.’ 0,1 ‘ ’ 1

Residual standard error: 0,0507 on 47 degrees of freedom
Multiple R-squared:  0,958,	Adjusted R-squared:  0,955 
F-statistic:  271 on 4 and 47 DF,  p-value: <2e-16

In [11]:
#==============================================================================
# Analise de variancia
#==============================================================================
print(anova(yfit))
Analysis of Variance Table

Response: lnG.Pop
          Df Sum Sq Mean Sq F value  Pr(>F)    
lnpreco    1  2,103   2,103  817,58 < 2e-16 ***
lnR.Pop    1  0,621   0,621  241,51 < 2e-16 ***
lnPCN      1  0,064   0,064   24,86 8,8e-06 ***
lnPCU      1  0,000   0,000    0,04    0,85    
Residuals 47  0,121   0,003                    
---
Signif. codes:  0 ‘***’ 0,001 ‘**’ 0,01 ‘*’ 0,05 ‘.’ 0,1 ‘ ’ 1
In [12]:
#==============================================================================
# Teste de hipoteses e intervalos de confianca
#==============================================================================
talpha   <- qt(0.95,47)
tk       <- (1.09587 - 1)/0.07771
print((tk > talpha))   # Nao rejeita H_0   
print(confint(yfit))   # Intervalo de confianca
[1] FALSE
               2,5 %    97,5 %
(Intercept) -22,7264 -19,69580
lnpreco      -0,1093   0,06685
lnR.Pop       0,9395   1,25220
lnPCN        -0,6896  -0,05762
lnPCU        -0,1878   0,22786
In [13]:
#==============================================================================
# Utilizando teoria assintotica podemos facilmente calcular o intervalo de 
# confianca de 95% para a elasticidade de longo prazo:
#==============================================================================
# (a) do consumo em relacao ao preco
#==============================================================================
print(c(eps2-qnorm(0.975)*se.eps2,eps2+qnorm(0.975)*se.eps2))
[1] -0,7099 -0,1129
In [14]:
#==============================================================================
# (b) do consumo em relacao a renda per capita 
#==============================================================================
print(c(eps3-qnorm(0.975)*se.eps3,eps3+qnorm(0.975)*se.eps3))
[1] 0,6523 1,2888