7 votos

Econometría espacial: cálculo de los residuos

Consideremos el modelo de econometría espacial denominado SAR por James P LeSage:

$y = \rho W y + X \beta+ \epsilon$

Utilizo el R paquete spdep y la caja de herramientas econométricas de LeSage de www.spatial-econometrics.com . Sin embargo, al comparar la forma en que estos dos calculan los residuos hay diferencias. En R de Roger Bivand los calcula mediante

$\dot{y}=y - \rho W y$

y luego regresión lineal $\dot{y}$ en el $X$ matriz. Obteniendo los residuos de este ajuste. ( lm(y-rWy~X-1 ). Sin embargo, LeSage calcula los residuos de la manera (que en un principio pensé que sería la correcta):

$y-(\rho W y)^{-1}Xb$

Así que aquí están mis preguntas:

  1. ¿Por qué los dos autores calculan los residuos de forma diferente?
  2. ¿Hay alguna objeción contra una u otra forma de computarlos?

Gracias por su amable ayuda.

7voto

Spencer Kormos Puntos 3082

Los residuos que calcula cada uno de ellos son diferentes. He aquí por qué:

El modelo es el siguiente:

$y = \rho Wy + xb + e$ con $e \sim n(0,1)$

Ahora bien, si jugamos con ello obtenemos:

$y = (I - \rho W)^{-1}(xb + e)$

Ahora lo que el Prof. LeSage hace es:

$y - (I-\hat{\rho} \cdot W)^{-1} \cdot x\hat{b} = (I-\rho W)^{-1}\cdot e$

Así que lo que se obtiene es el residuo con la correlación automática.

Por otra parte, transformando y:

$y - \hat{\rho}\cdot W \cdot y = xb + e$

Estimación, $xb$ y calcular los residuos, lo que Bivand está haciendo es darte $e$ en lugar de $(I-\rho W)^{-1}\cdot e$

La elección dependerá de su aplicación.

He aquí un código que lo demuestra. No estoy usando SPDEP directamente porque no estoy seguro de cómo crear mapas aleatorios... Pero eso está bien el código es bastante simple de todos modos:

#------------------ GENERATE SAMPLE DATA
rm(list=ls())   #clean
require(igraph) #random graphs
require(AER)    #get ivreg ...

n<-700   #700 locations
p=0.2
g <- erdos.renyi.game(n=n, p.or.m=p, type="gnp", directed=F, loops=F)
graph.density(g) 

w <- get.adjacency(g) #get an adjacency matrix
w <- w/rowSums(w)     #row standardize because of eigen vectors and eigen values
sum(rowSums(w)==0) 

rho <- 0.5
intercept   <- rep(1,n)
rvariable   <- rnorm(n)
y <- solve(diag(n) - rho*w) %*% ( 2*intercept + 3*rvariable + rnorm(n))

Una vez generados los datos según un modelo SAR LAG, los estimaremos mediante 2SLS (como ya te dije que podíamos hacer).

#------------------ GENERATE INSTRUMENTS
#get some instrumental variables 
z0 <- w%*%rvariable 
z1 <- w%*%w%*%rvariable

#check to see if there is a minimum of correlation
cor(z0, w%*%y)
cor(z1, w%*%y)

Los instrumentos funcionan porque rvariable es exógena. Por tanto, mientras w es exógena ¡tenemos un juego!

#------------------ NOW ONTO ESTIMATION

#The wrong way ...
summary(out<-lm(y ~ rvariable)) 
confint(out)

#The not so bad, but still very wrong way
summary(out<-lm(y ~ w%*%y + rvariable)) 
confint(out)

#ok now this should do it  ... not perfect beacuse 2sls is not efficient. 
#I am doing it this way because i did not want to generate random maps...
#Plus random graphs are easily available !

summary(out<-ivreg( y ~ w%*%y + rvariable, instruments=~ z0 + z1 + rvariable )) 
confint(out)

Ahora pasemos a lo que realmente importa, el cálculo de los residuos:

#residuals LeSage way
y_hat0    <- solve(diag(n) - coef(out)[2]*w ) %*% ( coef(out)[1]*intercept + coef(out)[3]*rvariable )
u_hat0    <- y - y_hat0

#residuals BiVand way
y_tilda   <- y - coef(out)[2]*w%*%y
summary(out_biv   <- lm( y_tilda ~ rvariable ))
#ok they are not the same due to rounding error ...
coef(out)[3] == coef(out_biv)[2]; round(coef(out)[3],5) == round(coef(out_biv)[2],5)

u_hat1 <- residuals(out_biv)
u_hat1 <- solve(diag(n) - coef(out)[2]*w)%*%u_hat1

#If we give Bivand some taste of autocorrelation it is the same as LeSage ...
round( u_hat0 - u_hat1, 5)

¡Al final debería ver la diferencia de residuos == 0 !

Hay que tener en cuenta que, dependiendo de la estructura del $W$ el efecto podría no ser identificable, por lo que la estrategia de utilizar el generador de gráficos aleatorios podría ser falsa en algunas ocasiones.

De todos modos espero que esto realmente resuelto su pregunta

0 votos

Entiendo tu razonamiento pero por otro lado la matriz $(I-\rho W)$ debe invertirse para la estimación de todos modos, ¿no? lagsarlm() función implementada.

0 votos

No, eso no es correcto (I-pW) no necesita invertirse. Podría ser el caso de que es en ese código, no sé (nunca mirado), pero este modelo puede ser estimado por "2SLS" como en Prucha sin necesidad de invertir (I-pW) ya que se instrumenta Wy con Wx ... por máxima verosimilitud en cuyo caso se maximiza la función de verosimilitud concentrada lnL() = + ln|In W| (n/2) ln(S()) sobre "p" o mediante MCMC como en LeSage donde evitan la inversión directa solvin (InW) = X a través de la descomposición cholsky ...

0 votos

Gracias por su extensa respuesta y su ayuda. si me permite resumirlo, la conclusión es la siguiente: el método bivand es más sencillo desde el punto de vista informático y, por tanto, preferible. bueno, es una buena pista, ya que los resultados de mi aplicación difieren bastante, lo que me parece bastante extraño, pero probablemente se deba a un error de especificación (?). de todos modos, habría esperado un razonamiento estadístico ;)

i-Ciencias.com

I-Ciencias es una comunidad de estudiantes y amantes de la ciencia en la que puedes resolver tus problemas y dudas.
Puedes consultar las preguntas de otros usuarios, hacer tus propias preguntas o resolver las de los demás.

Powered by:

X