Modelo Estimado

Ejemplo tomado de Statistical Methods for Social Sciences, Third Edition by Alan Agresti and Barbara Finlay (Prentice Hall, 1997), capítulo 9.

Variables en el dataset: Identificación del estado (sid), Nombre del estado (state), Crímenes violentos por cada 100,000 habitantes (crime), Asesinatos por cada millón de habitantes (murder), Porcentaje de población que vive en áreas metropolitanas (pctmetro), Porcentaje de población blanca (pctwhite), Porcentaje de población con un mínimo de segundaria como educación (pcths), Porcentaje de población por debajo de la línea de pobreza (poverty), Porcentaje de la población que es padre soltero (single).

Se estimará el el número de crimenes violentos, en función de la población por debajo de la línea de pobreza y el Porcentaje de la población que es padre soltero. (single): 

options(scipen = 999999)
library(foreign)
datos_regresion <- read.dta("crime.dta")
modelo_crime<-lm(crime~poverty+single,data=datos_regresion)
print(modelo_crime)
## 
## Call:
## lm(formula = crime ~ poverty + single, data = datos_regresion)
## 
## Coefficients:
## (Intercept)      poverty       single  
##   -1368.189        6.787      166.373


Pruebas de Heterocedasticidad y Autocorrelación.

Prueba de White (prueba de Breusch Pagan)

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.938
##  Alternative hypothesis: rho != 0

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

Estimación Robusta (uso del estimador HAC)

Se usará las librerias “lmtest” y “sandwich”, al respecto puede leer el siguiente paper: Econometric Computing with HC and HAC Covariance Matrix Estimators.


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 un estimador HAC)

options(scipen = 99999)
library(lmtest)
library(sandwich)
#Corregido
#HC0 Corrige Sólo Heterocedasticidad, use HC1 para corregir también Autocorrelación de Primer Orden
estimacion_omega<-vcovHC(modelo_crime,type = "HC0") 

coeftest(modelo_crime,vcov. = estimacion_omega)
## 
## t test of coefficients:
## 
##               Estimate Std. Error t value      Pr(>|t|)    
## (Intercept) -1368.1887   276.4111 -4.9498 0.00000956181 ***
## poverty         6.7874    10.6010  0.6403        0.5251    
## single        166.3727    25.4510  6.5370 0.00000003774 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


Estimador HAC de NeweyWest

En presencia de autocorrelación, es preferible usar el Estimador HAC de NeweyWest, indicando el orden de correlación detectada.

Si se hubiera detectado la correlación de segundo orden se tendría que corregir así:

Nota:(si la correlación es de orden 1, es más simple usar el estimador HC1) 

library(lmtest)
library(sandwich)

#Corregido:

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


Estimación Robusta

options(scipen = 999999)
library(robustbase)
library(stargazer)
modelo_crime_robust<-lmrob(crime~poverty+single,data=datos_regresion)
# print(summary(modelo_crime_robust))
stargazer(modelo_crime,modelo_crime_robust,type = "html",title = "Comparativa")
Comparativa
Dependent variable:
crime
OLS MM-type
linear
(1) (2)
poverty 6.787 11.466
(8.989) (9.263)
single 166.373*** 176.569***
(19.423) (23.223)
Constant -1,368.189*** -1,539.640***
(187.205) (235.765)
Observations 51 51
R2 0.707 0.795
Adjusted R2 0.695 0.787
Residual Std. Error (df = 48) 243.610 191.864
F Statistic 57.964*** (df = 2; 48)
Note: p<0.1; p<0.05; p<0.01