Estimadores HAC

Carga de datos y estimación del modelo

# carga de datos:

load("C:/Users/hp/Downloads/Gabriela Alexandra Belloso Martínez - smoke.RData")

# Estimación y presentación del modelo:

options(scipen = 9999)
library(stargazer)

modelo <- lm(formula = cigs ~ cigpric + lcigpric + income + lincome + age + agesq + educ + white + restaurn, data = data) 

stargazer(modelo,title = "Modelo",type = "html")
Modelo
Dependent variable:
cigs
cigpric 2.002
(1.493)
lcigpric -115.273
(85.424)
income -0.00005
(0.0001)
lincome 1.404
(1.708)
age 0.778***
(0.161)
agesq -0.009***
(0.002)
educ -0.495***
(0.168)
white -0.531
(1.461)
restaurn -2.644**
(1.130)
Constant 340.804
(260.016)
Observations 807
R2 0.055
Adjusted R2 0.044
Residual Std. Error 13.413 (df = 797)
F Statistic 5.169*** (df = 9; 797)
Note: p<0.1; p<0.05; p<0.01

Pruebas de Heterocedasticidad y Autocorrelación.

(Verificar si la matriz de Varianza-Covarianza del modelo es escalar)

Prueba de Breusch Pagan (White): Heterocedasticidad

options(scipen = 9999)
library(lmtest)
prueba_white<-bptest_result <- bptest(modelo, 
                          ~ I(cigpric^2) + I(lcigpric^2) + I(income^2) + I(lincome^2) + 
                            I(age^2) + I(agesq^2) + I(educ^2) + I(white^2) + I(restaurn^2) + 
                            cigpric:lcigpric + cigpric:income + cigpric:lincome +
                            cigpric:age + cigpric:agesq + cigpric:educ + 
                            cigpric:white + cigpric:restaurn + lcigpric:income + 
                            lcigpric:age + lcigpric:educ + 
                            lcigpric:white + income:lincome + income:age + income:educ + 
                            income:white + income:restaurn + lincome:agesq +
                            lincome:educ + lincome:white + lincome:restaurn + 
                            age:agesq + age:educ + age:white + age:restaurn +
                            agesq:educ + agesq:white + agesq:restaurn +
                            educ:white + educ:restaurn + white:restaurn, 
                          data = data)
print(prueba_white)
## 
##  studentized Breusch-Pagan test
## 
## data:  modelo
## BP = 58.87, df = 40, p-value = 0.02748

Hay evidencia de heterocedasticidad ya que \(\text{p-value} < 0.05\)

Prueba de Breusch Godfrey (Autocorrelación)

Verificando autocorrelación de 1° orden

library(lmtest)
bgtest(modelo,order = 1)
## 
##  Breusch-Godfrey test for serial correlation of order up to 1
## 
## data:  modelo
## LM test = 0.069737, df = 1, p-value = 0.7917

No hay evidencia de autocorrelación de 1° orden ya que (\(\text{p-value} > 0.05\))

Verificando autocorrelación de 1° orden

library(lmtest)
bgtest(modelo,order = 2)
## 
##  Breusch-Godfrey test for serial correlation of order up to 2
## 
## data:  modelo
## LM test = 0.26889, df = 2, p-value = 0.8742

No hay evidencia de autocorrelación de 2° orden ya que (\(\text{p-value} > 0.05\))

Estimación Robusta: Obtener el estimador HAC apropiado

library(lmtest)
library(sandwich)

#Sin corregir:

sin_corregir <- coeftest(modelo)

#Corregido:

vcov_HAC<-vcovHC(modelo,type = "HC1")
corregido <- coeftest(modelo,vcov. = vcov_HAC)

presentación del modelo inicial y el corregido

stargazer(sin_corregir, corregido, type = "html", 
          title = "Comparación de modelos", 
          column.labels = c("Sin corregir", "Corregido"), 
          align = TRUE)
Comparación de modelos
Dependent variable:
Sin corregir Corregido
(1) (2)
cigpric 2.002 2.002
(1.493) (1.613)
lcigpric -115.273 -115.273
(85.424) (91.916)
income -0.00005 -0.00005
(0.0001) (0.0001)
lincome 1.404 1.404
(1.708) (1.237)
age 0.778*** 0.778***
(0.161) (0.138)
agesq -0.009*** -0.009***
(0.002) (0.001)
educ -0.495*** -0.495***
(0.168) (0.164)
white -0.531 -0.531
(1.461) (1.370)
restaurn -2.644** -2.644**
(1.130) (1.045)
Constant 340.804 340.804
(260.016) (280.307)
Note: p<0.1; p<0.05; p<0.01