Objetivo

El objetivo de este laboratorio será introducir a los estudiantes a los siguientes elementos clave que son utilizados en el análisis exploratorio de datos espaciales (ESDA) y modelos econométricos. Estos elementos son:

Exploración preliminar de datos espaciales.
Matrices de peso espacial

Paquetes que se utilizarán

Antes de comenzar, debemos incluir los paquetes estadísticos necesarios para llevar a cabo el análisis exploratorio espacial, la construcción de matrices de peso espacial y medidas de autocorrelación. Preliminarmente, estos paquetes son:

(Estos paquetes se incluyen de acuerdo al supuesto de que ya están instalados en R. En el caso de que no estén instalados, usar el comando install.package(“nombre del paquete”). También, es posible instalar los paquetes en el menú Herramientas (tools) –> instalar paquetes)

library(spdep)
library(spData)
library(RColorBrewer)
library(leaflet)
library(dplyr)
library(ggplot2)
library(tmap)
library(tmaptools)
library(sf)
library(spDataLarge)

Información Espacial

Para agregar información espacial de tipo Shapes (en este caso, es un polígono) al ambiente global de trabajo debemos usar el siguiente código:

nepal.shp<-st_read('nepal.shp')

## Reading layer `nepal' from data source 
##   `/Users/yasnacortesgarriga/Library/CloudStorage/Dropbox/cursos/Cursos 2025/Econometría Espacial/Laboratorios_Actualizados/Contraste-autocorrelacion-espacial/nepal.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 75 features and 61 fields
## Geometry type: POLYGON
## Dimension:     XY
## Bounding box:  xmin: 80.06014 ymin: 26.34752 xmax: 88.20401 ymax: 30.44702
## Geodetic CRS:  WGS 84

class(nepal.shp)

## [1] "sf"         "data.frame"

Análisis Exploratorio Preliminar de los Datos.

Podemos usar otra interfaz para ver el mapa, a través del paquete leaflet

leaflet(nepal.shp) %>%
  addPolygons(stroke = FALSE, fillOpacity = 0.5, smoothFactor = 0.5) %>%
  addTiles() #adds a map title, the default is OpenStreetMap

Mapa de quantiles

tm_shape(nepal.shp) +
  tm_polygons(
    fill = "pcinc",
    fill.scale = tm_scale_intervals(style = "quantile"),
    fill.legend = tm_legend(title = "Ingresos Per Capita")
  ) +
  tm_layout(legend.outside = TRUE)

Usando el paquete leaflet, podemos usar un mapa de quantiles

qpal<-colorQuantile("OrRd", nepal.shp$pcinc, n=8) 
leaflet(nepal.shp) %>%
  addPolygons(stroke = FALSE, fillOpacity = .8, smoothFactor = 0.2, color = ~qpal(pcinc)
  ) %>%
  addTiles()

tmap_mode("view") #Para ver el mismo mapa de manera interactiva#

## ℹ tmap mode set to "view".

tm_shape(nepal.shp) +
  tm_polygons(
    fill = "pcinc",
    fill.scale = tm_scale_intervals(style = "quantile"),
    fill.legend = tm_legend(title = "Ingresos Per Capita")
  ) +
  tm_layout(legend.outside = TRUE)

Medidas de Autocorrelación Espacial Global

El Diagrama de Moran: El Camino Largo

El diagrama de Moran (Moran scatterplot) es una forma rápida para medir la estructura espacial de las observaciones de una variable bajo estudio. Este diagrama se construye centrando la variable bajo estudio de acuerdo a su media (aunque no es necesario) y el promedio de los valores de las observaciones vecinas, \(Wy\). En este caso, la matriz \(W\) está estandarizada (normalizada) por filas. Al centrar la variable de interés a su media, podemos delinear los cuadrantes del diagrama de acuerdo a \(y=0\) y \(Wy=0\).

También, podemos dibujar el grado de ajuste de la relación entre la variable espacial bajo estudio y su rezago espacial a través de una línea de regresión. Con ello podemos ver con mayor claridad la fuerza de la relación entre la variable de interés y su rezago espacial.

Si las observaciones están distribuidas de manera aleatoria en el espacio, no se podría observar una relación particular entre \(y\) y \(Wy\). En este caso, la pendiente de la relación debería ser igual a cero.

En el caso contrario, el parámetro de regresión es distinto de cero. De esta forma, cada cuadrante se ajustaría a un tipo de asociación espacial específica.

Para construir manualmente este diagrama, podemos seguir el siguiente procedimiento:

Crear una matriz W estandarizada por filas (no estocástica por filas)

queen.v2<-nb2listw(poly2nb(nepal.shp, row.names = nepal.shp$id, queen = TRUE),style ="W",zero.policy = TRUE)
summary(queen.v2)

## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 75 
## Number of nonzero links: 366 
## Percentage nonzero weights: 6.506667 
## Average number of links: 4.88 
## Link number distribution:
## 
##  2  3  4  5  6  7  9 
##  3 15 15 15 15  9  3 
## 3 least connected regions:
## 31 42 43 with 2 links
## 3 most connected regions:
## 18 20 39 with 9 links
## 
## Weights style: W 
## Weights constants summary:
##    n   nn S0       S1       S2
## W 75 5625 75 32.85108 309.2524

queen.v2$weights[1]

## [[1]]
## [1] 0.1666667 0.1666667 0.1666667 0.1666667 0.1666667 0.1666667

Centramos la varible a estudiar a su media.

var<-scale(nepal.shp$pcinc, center = TRUE, scale = FALSE)
var

Calculamos el rezago espacial de la variable bajo estudio.

var.lag<-lag.listw(queen.v2, var)
var.lag

Para estimar la línea de regresión para visualizar el grado de ajuste, calculamos la siguiente regresión.

mor<-lm(var.lag ~ var)
summary(mor)

## 
## Call:
## lm(formula = var.lag ~ var)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -335.12  -84.95  -36.02   91.53  505.59 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  4.08087   18.09491   0.226    0.822    
## var          0.38460    0.06644   5.789 1.66e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 156.7 on 73 degrees of freedom
## Multiple R-squared:  0.3146, Adjusted R-squared:  0.3052 
## F-statistic: 33.51 on 1 and 73 DF,  p-value: 1.657e-07

Finalmente, podemos graficar el diagrama de moran de la siguiente forma:

plot(var.lag ~ var,  pch=20, asp=1, las=1)
abline(mor, lwd=2)
abline(h=mean(var), lt=2)
abline(v=mean(var.lag), lt=2)

Donde la pendiente de la curva de regresión nos indica la relación líneal entre ambas variables.

El Índice de Moran

El índice de moran, o estadístico de moran mide la relación lineal entre la variable de interés y sus valores espacialmente rezagados. Es decir,

\[ corr(y, Wy)= \dfrac{Cov(y, Wy)}{\sqrt{Var(y)Var(Wy)}} \] De acuerdo a esta ecuación, usamos la expresión del rezago espacial. También, asumimos que bajo el supuesto de estacionaridad (el mismo de series de tiempo estacionarias):

\[ Var(y)=Var(Wy) \] Por tanto, este coeficiente de correlación es igual a:

\[ corr(y, Wy)= \dfrac{Cov(y, Wy)}{Var(y)} \] En términos empíricos, el índice de moran puede estimarse de la siguiente forma:

\[ I= \dfrac{\sum_{i=1}^{n} \sum_{j=1}^n w_{ij}z_{i}z_{j}}{\sum_{i=1}^{n}z_{i}^{2}} \] donde \(n\) es el número de observaciones, \(z_{i}\) es el valor de la variable para la observación \(i\), \(z_{j}\) es el valor de la variable para la observación \(j\). Tanto \(z_{i}\) y \(z_{j}\) están centradas en la media. Finalmente, \(w_{ij}\) es el elemento \(ij\) de la matriz de peso espacial estandarizada por filas.

En términos matriciales, esta formula puede expresarse de la siguiente forma: \[ I=\dfrac{z^{T}Wz}{z^{T}z} \] Donde \(W\) es la matriz de peso espacial estandarizada por filas.

Tal como puede notar, este índice es igual a los estimados que se obtienen a través de OLS. Por tanto, la pendiente de la regresión entre el rezago espacial de la variable bajo estudio (y) y la variable bajo estudio (x) es igual al índice de moran.

El índice de moran se encuentra en el intervalo entre \([-1,1]\), De esta forma, podemos notar tres situaciones:

Cuando la variable z está distribuida aleatoriamente, el índice de moran se encuentra cercano a 0.
Cuando el valor observado del índice de moran es positivo, indica la existencia de autocorrelación positiva.
Cundo el valor observado del índice moran es negativo, indica la existencia de autocorrelación negativa.

Contraste de Hipótesis para el Índice de Moran

El contraste de hipótesis para determinar la significancia del índice de moran se lleva a cabo teniendo en cuenta lo siguiente:

\[ H_{0}: \mbox{Las observaciones se distruyen aleatoriamente. No hay autocorrelación espacial} \\ H_{1}: \mbox{Existe autocorrelación espacial} \] Para llevar a cabo el contraste de hipótesis, podemos usar dos formas:

Bajo Aleatorización.
Bajo Normalidad.

Contraste de Hipótesis bajo Aleatorización

Este método no asume una distribución a priori de la variable bajo estudio. Por tanto, el estadístico de moran se estima usando los datos es comparado con la distribución de las variables que surge de una reorganización de los datos de forma aleatoria. Si la hipótesis nula es cierta, todas las posibles combinaciones de datos son equiprobales (todos los valores son igualmente probables).

El procedimiento para realizar este constraste de hipótesis es el siguiente:

Creamos la distribución de índices de moran, el cual corresponde a cada uno de los valores permutados en el área de análisis. En otras palabras, por cada permutación de valores, obtenemos un índice de moran por cada permtuación. Esto lo podemos crear de la siguiente forma:

n <- 599   # Define the number of simulations
I.r <- vector(length=n)  # Create an empty vector

for (i in 1:n){
  # Randomly shuffle Values
  x <- sample(var, replace=FALSE)
  # Compute a new set of lagged values
  x.lag <- lag.listw(queen.v2, x)
  # Compute the regression slope and store its value
  M.r    <- lm(x.lag ~ x)
  I.r[i] <- coef(M.r)[2]
}

I.r

##   [1] -8.340192e-03  1.330821e-01  8.832500e-03 -6.501069e-02  8.439835e-02
##   [6]  6.450232e-02 -5.668037e-02  7.952456e-03 -1.769770e-02 -1.474528e-01
##  [11] -4.793226e-03  6.763761e-02  6.809256e-02  1.217361e-03 -3.592993e-02
##  [16] -4.276188e-02  6.115488e-02 -1.511978e-02 -2.160567e-02  6.682959e-03
##  [21]  2.677613e-02  1.735097e-01  1.455513e-02  1.347758e-02  9.658942e-02
##  [26]  7.431402e-02 -7.595409e-02 -1.561698e-01  5.680614e-02  3.008875e-02
##  [31] -8.242900e-02 -6.469847e-03  1.141841e-01 -4.934563e-02 -3.613062e-02
##  [36]  4.779920e-02  6.398806e-02 -1.739808e-02 -2.672493e-02 -1.621914e-02
##  [41] -7.062662e-02  1.272042e-01 -9.747532e-02  6.069758e-02  3.050740e-01
##  [46] -6.006333e-02 -1.093974e-01 -6.680527e-02 -3.608464e-02 -9.941121e-02
##  [51]  1.009505e-02  2.945911e-02 -4.548175e-02 -5.816494e-02 -1.020466e-01
##  [56]  3.055358e-02 -3.559254e-02 -1.561695e-02  1.392127e-01  2.829109e-02
##  [61]  9.880937e-02  1.066533e-01 -1.059908e-01  5.043292e-02 -4.335678e-02
##  [66] -1.036883e-02 -4.470548e-02 -2.673567e-02 -2.462302e-04  1.727647e-02
##  [71] -1.440027e-01  6.956376e-02  1.556005e-01 -1.228742e-02 -1.601196e-01
##  [76]  1.596055e-01 -1.563878e-01  3.091732e-02 -4.461004e-02  2.016488e-03
##  [81] -9.103700e-02 -1.060429e-02  1.043951e-01 -5.003996e-02  3.121360e-03
##  [86] -2.937789e-02 -1.561300e-02 -1.485136e-01 -6.991715e-02 -1.148072e-01
##  [91]  1.261958e-02 -5.353616e-02  7.817533e-02 -9.279197e-03 -9.040529e-03
##  [96] -2.377714e-02  4.904242e-02 -3.764237e-02 -1.291518e-02  3.261059e-02
## [101] -1.247958e-01 -2.471935e-02 -2.709647e-02 -2.943383e-02 -8.626278e-02
## [106] -8.857281e-02  8.640954e-02 -8.195811e-02 -2.573202e-02  2.139832e-02
## [111] -5.736144e-02  5.189769e-02 -1.127321e-01 -4.803059e-02 -6.798467e-02
## [116]  6.694077e-02  5.769448e-02  1.042800e-01  1.715750e-02  6.193873e-02
## [121] -7.763408e-02 -4.453424e-03 -3.364657e-02 -4.746064e-02 -1.962035e-02
## [126] -3.862657e-02 -5.338941e-02  5.750314e-02 -1.252533e-01  4.818850e-02
## [131]  4.421915e-03  6.503006e-02 -2.901811e-03 -4.196471e-03 -6.538116e-02
## [136]  3.094437e-02  1.091638e-02  7.736717e-05  1.459476e-01  3.142901e-03
## [141] -1.330623e-01 -6.659576e-02  1.920305e-02 -6.574538e-02 -8.789999e-02
## [146] -4.396219e-02 -1.343744e-01 -2.881669e-02  4.052035e-02 -1.273579e-01
## [151] -3.677886e-03 -4.311875e-03 -7.356278e-03  1.117896e-01 -2.248800e-02
## [156] -2.521799e-02 -8.925021e-02 -8.924060e-03  1.691132e-02  8.373076e-02
## [161] -1.544583e-01 -4.564824e-02 -4.724852e-02  8.403253e-02 -8.786686e-02
## [166] -2.590173e-02  4.463862e-02 -1.793311e-02 -3.935217e-02 -2.541662e-02
## [171] -5.324673e-02  3.976512e-03 -4.942764e-02 -4.257906e-02 -2.830354e-02
## [176] -2.447668e-02 -4.699989e-02  2.952475e-01  8.051685e-03  2.522102e-03
## [181] -2.920894e-02  6.264095e-02 -7.937690e-03 -3.476528e-03 -3.796452e-02
## [186]  1.333747e-01 -3.376025e-02 -1.917688e-02 -1.406977e-01 -1.263193e-02
## [191] -6.113945e-02 -2.702837e-02 -5.761922e-02 -5.697697e-02 -6.496509e-02
## [196] -4.397517e-02  6.622655e-02 -6.669631e-02 -7.763956e-02 -1.589225e-02
## [201] -1.315202e-02 -8.040290e-02 -7.669548e-02  4.854059e-02  1.247652e-01
## [206]  9.992006e-03  3.192556e-02  2.457564e-02 -5.210365e-02  1.278891e-01
## [211] -3.398664e-02 -1.061196e-01 -2.444244e-02 -6.796547e-02 -2.151634e-02
## [216] -3.708663e-02 -1.416326e-01  2.503321e-02 -6.728566e-02  1.873535e-02
## [221]  3.406698e-02  3.695849e-02  1.745487e-02 -8.252041e-02  2.133396e-01
## [226] -1.233386e-01 -2.087973e-02 -9.688620e-02 -6.510703e-02  1.615308e-02
## [231]  9.286614e-02 -3.230468e-02  3.488794e-02  2.286747e-02 -1.131148e-01
## [236]  4.191891e-02 -8.175945e-03 -4.024841e-03 -8.983300e-02 -3.416878e-02
## [241] -1.047712e-01 -4.442618e-02 -6.316159e-02  3.242933e-03 -3.560826e-02
## [246]  7.727019e-02 -1.204637e-01  2.326197e-02 -9.658363e-02  1.535842e-02
## [251]  3.604395e-02  3.045816e-02  2.248586e-02 -5.653324e-02 -6.613046e-02
## [256] -2.490770e-02  6.032074e-02 -6.729517e-02  3.303361e-02 -4.226532e-02
## [261]  2.616784e-02 -4.960997e-02 -4.900038e-02 -4.310185e-02 -1.620579e-02
## [266] -4.520856e-02  4.687165e-02  7.524188e-02 -1.167679e-01 -5.127239e-02
## [271] -1.194731e-01 -1.209580e-01 -1.478169e-01  2.015635e-01  5.369709e-02
## [276]  1.002652e-01 -1.106540e-01 -2.230553e-02 -8.888246e-02 -6.049492e-02
## [281]  5.251158e-02  2.487862e-02 -4.684730e-02  4.794752e-03  1.160511e-02
## [286] -5.112922e-03 -6.072746e-02 -2.848011e-02  1.569602e-02 -1.187042e-01
## [291] -1.291685e-01 -1.355189e-01  9.328871e-02 -2.228717e-02 -1.028688e-01
## [296] -1.084479e-01 -7.984789e-02  1.054419e-01 -3.212838e-02 -8.227190e-02
## [301] -3.600915e-02 -1.131898e-01 -3.005863e-02  4.429999e-02 -1.073487e-01
## [306] -1.711109e-03 -1.083208e-01  1.987310e-02 -8.137491e-02 -5.029034e-02
## [311] -2.105058e-02 -5.532209e-03 -1.426358e-02 -6.347569e-02  1.739361e-02
## [316] -1.506694e-02 -6.157727e-02 -3.145035e-03  3.122829e-02  1.957961e-02
## [321]  1.374602e-03  8.731298e-02  5.991715e-02 -5.238851e-02  5.696860e-02
## [326] -5.840703e-02  1.440251e-01  1.098005e-01 -2.006708e-02 -6.271965e-02
## [331] -2.176123e-02 -1.294549e-01  5.350835e-02 -2.677381e-02  5.815963e-02
## [336]  1.016967e-01 -3.665471e-02 -6.072883e-03  3.507454e-02 -4.807645e-02
## [341]  7.249030e-02  1.244376e-03 -7.973046e-02  4.151375e-02  6.095331e-02
## [346] -5.032533e-02  8.152772e-02 -1.157908e-01  1.902394e-02 -9.788985e-02
## [351] -7.564343e-02  3.293105e-02 -1.558930e-01 -7.788126e-02 -7.787293e-02
## [356] -5.833866e-02  1.511883e-01 -1.055460e-01 -6.211755e-02  1.126260e-01
## [361] -6.680426e-02  6.602034e-03 -4.887781e-02  1.506564e-02  4.905040e-02
## [366] -9.372278e-02  7.592186e-02  1.826911e-02 -7.447885e-02 -1.368633e-02
## [371] -9.105265e-02 -1.146723e-01 -1.886440e-03 -1.051122e-01  5.698269e-02
## [376]  1.966768e-01 -5.279624e-02 -1.293598e-01 -4.444588e-02 -3.504356e-02
## [381] -9.115284e-02  4.915321e-02 -1.372534e-01  7.960812e-02 -7.794920e-02
## [386]  2.836629e-01  1.731620e-02 -4.975433e-02 -2.819442e-02 -3.010343e-02
## [391] -6.799611e-02 -7.745316e-02 -1.996549e-02  9.113079e-03 -7.571472e-03
## [396] -8.208781e-02 -3.144090e-02 -7.762582e-02  1.199527e-04 -8.688759e-02
## [401] -3.356633e-02 -3.606445e-02 -9.387513e-03 -6.345869e-02  7.145345e-02
## [406] -7.542977e-02  2.602310e-03  8.122820e-03 -2.076616e-02 -1.929204e-02
## [411] -1.926539e-02  1.794378e-02 -6.716412e-02 -7.675700e-02  4.720371e-02
## [416]  8.074247e-04  2.779253e-02  9.849419e-02 -1.273272e-01 -9.255546e-02
## [421] -8.628877e-02 -1.730873e-03  9.434429e-02 -6.652431e-02  5.753046e-02
## [426] -1.460745e-02 -1.665561e-02  8.943846e-02 -1.135874e-01 -7.890156e-02
## [431]  1.089510e-02 -1.547413e-01 -9.657457e-02 -9.546605e-02 -8.027396e-02
## [436] -7.468846e-02 -1.530554e-01  2.934203e-02 -1.028623e-01 -3.414044e-02
## [441] -1.877781e-02  2.813206e-02 -6.351701e-02  1.271455e-01  4.694267e-02
## [446]  1.104715e-02 -2.400187e-02  2.004578e-02 -9.725553e-02 -1.173979e-01
## [451] -7.150087e-02  9.370887e-03  8.274586e-04  3.180107e-02  1.453732e-01
## [456] -1.289184e-01 -6.707001e-02 -5.198793e-02  6.632200e-02 -3.920222e-03
## [461] -3.988758e-02  2.603739e-01 -7.661773e-02 -5.837442e-03 -3.947543e-02
## [466] -2.610042e-02 -6.683134e-02  1.429350e-01 -2.316128e-02 -4.784011e-02
## [471] -7.515251e-02 -4.745493e-02 -8.623205e-02  1.098234e-01 -7.074904e-02
## [476]  5.422995e-02 -6.654059e-02 -3.102956e-02 -2.424173e-02 -6.241076e-02
## [481]  2.643169e-02 -8.106580e-02 -4.842585e-02 -1.112177e-01 -7.942664e-02
## [486] -9.501238e-02 -6.511093e-02 -1.336261e-02 -1.633922e-02 -1.789513e-02
## [491]  1.845137e-02 -6.037589e-02 -4.420083e-03  2.365326e-04 -4.691711e-02
## [496] -1.986820e-02  1.331669e-02 -8.322215e-02  3.641444e-02  7.558804e-02
## [501] -6.143432e-02 -7.224380e-02 -5.186587e-02 -1.503623e-02  2.292414e-02
## [506]  2.001900e-01  1.255460e-01  1.292104e-02 -1.047118e-02  2.551340e-02
## [511] -1.294576e-02  6.876941e-02  2.274656e-02 -3.014442e-02  1.014953e-01
## [516] -1.857049e-03 -5.898802e-02  8.951726e-02 -1.179642e-01 -7.704737e-03
## [521]  2.064706e-01 -3.951684e-02 -5.325141e-02  2.664404e-02 -2.541099e-02
## [526]  2.976877e-02 -1.658799e-03 -1.275565e-01 -1.221048e-01 -3.655916e-02
## [531] -4.065182e-02  8.247514e-02 -1.588230e-03  1.838113e-02 -2.405630e-02
## [536] -5.857242e-02 -8.727239e-02 -1.278594e-01  2.477759e-02  1.435742e-01
## [541] -3.953203e-02 -8.560301e-02 -8.678251e-02 -3.504055e-02 -3.909464e-02
## [546]  1.880046e-02  9.670053e-02 -4.849597e-04 -4.521335e-02 -2.616649e-02
## [551] -5.043546e-03 -5.725210e-02  2.319659e-02 -6.333702e-02 -3.476342e-02
## [556]  1.099041e-02 -4.472605e-03 -1.399970e-01  2.212885e-02 -3.125752e-02
## [561]  4.697148e-03  1.221257e-01  2.119686e-02  2.076914e-02  2.309848e-02
## [566] -5.990519e-02  5.711605e-02 -6.115319e-02 -9.838410e-02  6.810007e-02
## [571] -2.595382e-02  1.131288e-01  2.566021e-02 -1.373899e-02 -8.066162e-02
## [576]  9.010974e-03 -1.731256e-01  2.130047e-01 -1.217246e-01 -4.075293e-03
## [581]  6.275944e-02 -3.861183e-02 -2.822627e-02  1.277578e-01 -5.768726e-02
## [586]  5.492956e-02  1.458666e-02 -1.173202e-01 -3.039034e-02 -1.003251e-01
## [591] -5.951389e-02  4.621716e-02 -6.165212e-02 -3.216849e-02  4.228613e-02
## [596]  5.893672e-02 -5.312468e-02  6.126082e-03  2.851082e-02

Graficamos la distribución obtenida y comparamos con el valor estimado con la relación espacial determinada a través de la matriz W.

hist(I.r, main=NULL, xlab="Moran's I", las=1, xlim =c(-0.2,0.5))
abline(v=coef(mor)[2], col="red")

Calculamos el p-valor.

Obtenemos el número de valores del índice de moran simulados que son más grandes que nuestro valor observado (un lado):

N.greater <- sum(coef(mor)[2] > I.r)
p <- min(N.greater + 1, n + 1 - N.greater) / (n + 1)
p

## [1] 0.001666667

En este caso, rechazamos la hipótesis nula en favor de la alternativa. Por tanto, existe autocorrelación espacial.

Contraste de Hipótesis bajo Normalidad

Para calcular esta hipótesis bajo normalidad, usamos

moran.test(nepal.shp$pcinc, listw = queen.v2,randomisation = FALSE)

## 
##  Moran I test under normality
## 
## data:  nepal.shp$pcinc  
## weights: queen.v2    
## 
## Moran I statistic standard deviate = 5.3883, p-value = 3.556e-08
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       0.384599553      -0.013513514       0.005458871

El Diagrama de Moran: El Camino Corto

plot_mor<-moran.plot(nepal.shp$pcinc, listw = queen.v2, xlab="Income PC", ylab = "Lag Income PC", pch=19, labels = FALSE)

El Índice de Moran bajo Aleatorización

moran.test(nepal.shp$pcinc, listw = queen.v2,  randomisation = TRUE)

## 
##  Moran I test under randomisation
## 
## data:  nepal.shp$pcinc  
## weights: queen.v2    
## 
## Moran I statistic standard deviate = 5.6986, p-value = 6.039e-09
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       0.384599553      -0.013513514       0.004880579

mc<-moran.mc(nepal.shp$pcinc, queen.v2, nsim=599)
plot(mc)

El Índice de Moran bajo Normalidad

moran.test(nepal.shp$pcinc, listw = queen.v2,randomisation = FALSE)

## 
##  Moran I test under normality
## 
## data:  nepal.shp$pcinc  
## weights: queen.v2    
## 
## Moran I statistic standard deviate = 5.3883, p-value = 3.556e-08
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       0.384599553      -0.013513514       0.005458871

El Índice de Moran Local

El diagrama de moran que presentamos anteriormente, no entrega indicaciones acerca de la significancia de los clústeres espaciales detectados. Por tanto, necesitamos calcular el estadístico LISA para cada observación. De esta forma, obtenemos un indicador de cuales son los clústeres espaciales de valores similares significativos.

El estadístico LISA es usado para la detección de clústeres espaciales, como también para la inestabilidad local, outliers significativos y regímenes espaciales. El estadístico LISA puede ser calculado de la siguiente forma:

\[ I_{i}=\dfrac{z_{i}}{m_{0}}\sum_{j=1}^{n}w_{ij}z_{j} \] Donde

\[ m_{0}=\sum_{i=1}^{n}\dfrac{z_{i}^2}{n} \] Donde \(z_{i}\) es el valor de la observación \(i\), \(z_{i}\) es el valor de la observación \(j\). Ambas variables están centradas de acuerdo a su media. Finalmente, \(w_{ij}\) es la matriz estandarizada por filas.

La suma del estadístico de LISA para cada observación es proporcional al índice de moran global. En el caso de la matriz W sea estandarizada por filas, la media del LISA es igual al índice de moran global.

Para calcular el LISA, podemos usar la siguiente función:

lmoran<-localmoran(nepal.shp$pcinc, listw = queen.v2, alternative = "two.sided", zero.policy= TRUE)
head(lmoran)

##            Ii          E.Ii      Var.Ii       Z.Ii Pr(z != E(Ii))
## 1  0.69054902 -0.0179842596 0.205639755  1.5624521   0.1181815215
## 2  2.88646029 -0.0494108305 1.142064186  2.7472092   0.0060104785
## 3 -0.09578242 -0.0001888826 0.001857056 -2.2182780   0.0265358823
## 4  0.01866158 -0.0003092250 0.004382862  0.2865543   0.7744536457
## 5  1.45557108 -0.0062888556 0.151952300  3.7501801   0.0001767076
## 6 -0.24486530 -0.0006823121 0.007939309 -2.7404637   0.0061352556

nepal.shp$lmp <- lmoran[, 5] # p-values están en la columna 5

mean_moran<-mean(lmoran[, 1]) # comprobamos que la media del índice de moran es igual a la media del global
mean_moran

## [1] 0.3845996

mp <- moran.plot(as.vector(scale(nepal.shp$pcinc)), queen.v2)

nepal.shp$quadrant <- NA
# high-high
nepal.shp[(mp$x >= 0 & mp$wx >= 0) & (nepal.shp$lmp <= 0.05), "quadrant"]<- 1
# low-low
nepal.shp[(mp$x <= 0 & mp$wx <= 0) & (nepal.shp$lmp <= 0.05), "quadrant"]<- 2
# high-low
nepal.shp[(mp$x >= 0 & mp$wx <= 0) & (nepal.shp$lmp <= 0.05), "quadrant"]<- 3
# low-high
nepal.shp[(mp$x <= 0 & mp$wx >= 0) & (nepal.shp$lmp <= 0.05), "quadrant"]<- 4
# non-significant
nepal.shp[(nepal.shp$lmp > 0.05), "quadrant"] <- 5

tm_shape(nepal.shp) + tm_fill(col = "quadrant", title = "",
                        breaks = c(1, 2, 3, 4, 5, 6),
                        palette =  c("red", "blue", "lightpink", "skyblue2", "white"),
                        labels = c("High-High", "Low-Low", "High-Low",
                                   "Low-High", "Not significant")) +
  tm_legend(text.size = 1)  + tm_borders(alpha = 0.5) +
  tm_layout(frame = FALSE,  title = "Local Moran Income PC")  +
  tm_layout(legend.outside = TRUE)

Introducción a la Econometría Espacial - Contraste Autocorrelación Espacial

Yasna Cortés

2025-09-08