Matriz de covarianza no escalar

Si la matriz\(\color{red}{E[\epsilon*\epsilon^T]=\Omega}\), es decir que se ha detectado la presencia de Heterocedasticidad y/o de Autocorrelación.

Entonces la Matriz de varianza covarianza de los parámetros del modelo estimado será:

\(\color{purple}{Var[\beta]=\left(X^TX\right)^{-1}X^T\Omega X\left(X^TX\right)^{-1}}\)

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("https://stats.idre.ucla.edu/stat/data/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

\(\color{red}{\text{Hay evidencia de heterocedasticidad ya que}}\) \(p_{value}<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

\(\color{green}{\text{No hay evidencia de Autocorrelación de 2º orden ya que }}\) \(p_{value}>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.942
##  Alternative hypothesis: rho != 0

\(\color{green}{\text{No hay evidencia de Autocorrelación de 1º orden ya que }}\) \(p_{value}>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")