Aplicación de la librería stargazer para presentar modelos corregidos con estimadores HAC

1. Estimación del Modelo

options(scipen = 999999)
library (stargazer)
library(foreign)
datos_regresion <- read.dta("https://stats.idre.ucla.edu/stat/data/crime.dta")
modelo_crime<-lm(crime~poverty+single,data=datos_regresion)
stargazer(modelo_crime,type = "html",title = "Modelo Crimen")
Modelo Crimen
Dependent variable:
crime
poverty 6.787
(8.989)
single 166.373***
(19.423)
Constant -1,368.189***
(187.205)
Observations 51
R2 0.707
Adjusted R2 0.695
Residual Std. Error 243.610 (df = 48)
F Statistic 57.964*** (df = 2; 48)
Note: p<0.1; p<0.05; p<0.01

2. Pruebas de Heterocedasticidad y Autocorrelación

Prueba de Heterocedasticidad

Usando la libreria “lmtest”

library(lmtest)
white_test<-bptest(modelo_crime,~I(poverty^2)+I(single^2)+poverty*single,data = datos_regresion)
print(white_test)
studentized Breusch-Pagan test

data: modelo_crime BP = 10.73, df = 5, p-value = 0.057

Hay evidencia de heterocedasticidad ya que pvalue < 0.05

Prueba del Multiplicador de Lagrange (Breusch Godfrey)

Autocorrelación de 2º Orden:

library(lmtest)
prueba_LM<-bgtest(modelo_crime,order = 2)
print(prueba_LM)
Breusch-Godfrey test for serial correlation of order up to 2

data: modelo_crime LM test = 0.27165, df = 2, p-value = 0.873

No hay evidencia de Autocorrelación de 2º orden ya que pvalue > 0.05

Autocorrelación de 1º orden (prueba de Durbin Watson)

library(car)
durbinWatsonTest(model = modelo_crime)

lag Autocorrelation D-W Statistic p-value 1 -0.07014421 2.040007 0.904 Alternative hypothesis: rho != 0

No hay evidencia de Autocorrelación de 1º orden ya que pvalue > 0.05

3. Estimación Robusta (uso del estimador HAC)

Nota: No hay evidencia de autocorrelación, pero para la simulación del ejemplo se supondra que existe presencia de autocorrelación de primer orden, lo que se usa el Estimador HAC de NeweyWest.

Sin Corregir

options(scipen = 99999)
library(lmtest)
#Sin corregir:
coeftest(modelo_crime)
## 
## t test of coefficients:
## 
##               Estimate Std. Error t value         Pr(>|t|)    
## (Intercept) -1368.1887   187.2052 -7.3085 0.00000000247861 ***
## poverty         6.7874     8.9885  0.7551           0.4539    
## single        166.3727    19.4229  8.5658 0.00000000003117 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Corregido usando HAC

library(lmtest)
library(sandwich)

#Corregido:
# Para Autocorrelacion de Orden 1
estimacion_omega1<-vcovHC(modelo_crime,type = "HC1")
coeftest(modelo_crime,vcov. = estimacion_omega1) 
## 
## t test of coefficients:
## 
##               Estimate Std. Error t value      Pr(>|t|)    
## (Intercept) -1368.1887   284.9180 -4.8020 0.00001576624 ***
## poverty         6.7874    10.9273  0.6211        0.5374    
## single        166.3727    26.2343  6.3418 0.00000007519 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
robust.se<-sqrt(diag(estimacion_omega1))


#Para Autocorrelacion de Orden 2
estimacion_omega<-NeweyWest(modelo_crime,lag = 2)
coeftest(modelo_crime,vcov. = estimacion_omega)
## 
## t test of coefficients:
## 
##               Estimate Std. Error t value      Pr(>|t|)    
## (Intercept) -1368.1887   303.8466 -4.5029 0.00004279768 ***
## poverty         6.7874    10.5943  0.6407        0.5248    
## single        166.3727    25.9154  6.4198 0.00000005708 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

4. Estimación Robusta usando libreria Stargazer

options(scipen = 999999)
library(robustbase)
library(stargazer)
modelo_crime_corregido<-lmrob(crime~poverty+single,data=datos_regresion,)

#Para ver el modelo original y el modelo corregido en dos columnas
stargazer(modelo_crime,modelo_crime,se=list(NULL, robust.se),type = "html",title = "Modelo de Regresion Lineal", column.labels = c("Modelo Original", "Modelo Corregido HAC") )
Modelo de Regresion Lineal
Dependent variable:
crime
Modelo Original Modelo Corregido HAC
(1) (2)
poverty 6.787 6.787
(8.989) (10.927)
single 166.373*** 166.373***
(19.423) (26.234)
Constant -1,368.189*** -1,368.189***
(187.205) (284.918)
Observations 51 51
R2 0.707 0.707
Adjusted R2 0.695 0.695
Residual Std. Error (df = 48) 243.610 243.610
F Statistic (df = 2; 48) 57.964*** 57.964***
Note: p<0.1; p<0.05; p<0.01