Matriz de ponderação espacial com
spdep
André Luiz Ferreira e Silva
2023-09-20
1 Introdução
Em estatística espacial, uma questão relevante consiste em definir o
arranjo espacial que melhor acomoda os efeitos da autocorrelação
espacial. Busca-se identificar o arranjo que proporciona uma estimativa
consistente e paciomoniosa do grau de autocorrelação espacial. Para esse
fim, deve-se escolher uma matriz de ponderação espacial que melhor
represente as interações do fenômeno em estudo. Assim, esta nota visa
demonstrar de modo simples, como construir diferentes arranjos
espaciais. Para tanto, diferentes critérios de ponderação espacial são
abordados por meio do pacote spdep.
2 Formalização
Usualmente, a matriz de ponderação espacial W se baseia em critérios de proximidade geográfica, seja de contiguidade ou distância entre pontos. Desse modo, quatro diferentes tipos de matrizes são avaliadas.
2.1 Contiguidade de primeira-ordem
A primeira matriz se baseia no critério Queen de contiguidade de primeira-ordem \[ w_{ij} = \left\{ \begin{array}{cc} 1 & \mbox{se i e j são contíguos} \\ 0 & \mbox{outros casos} \\ \end{array} \right. \] Trata-se de uma matriz binária, não necessariamente simétrica.
2.2 k-vizinhos mais próximos
A segunda matriz é definida de acordo com o critério dos k-vizinhos mais próximos. \[ w_{ij}(k) = \left\{ \begin{array}{ccc} 1 & \mbox{se} & d_{ij} \leq d_i(k) \\ 0 & \mbox{se} & d_{ij} > d_i(k) \\ \end{array} \right. \] Em que, \(d_i(k)\) é a distância de corte para uma localidade específica. É a menor distância entre \(i,j\), tal que o arranjo suporte exatamente k-vizinhos.
2.3 Distância limite de corte
A terceira matriz também é binária, mas toma a distância máxima de corte como critério de vizinhança de primeira ordem. \[ w_{ij} = \left\{ \begin{array}{ccc} 1 & \mbox{se} & 0 \leq d_{ij} \leq d_{max} \\ 0 & \mbox{se} & d_{ij} > d_{max} \\ \end{array} \right. \] A distância de corte \(d_{max}\) deve ser escolhida, tal que cada localidade tenha pelo menos um vizinho, evitando ilhas.
2.4 Função inversa da distância
Esta matriz assemelha-se a anterior, porém, seus pesos são determinados por uma função inversa da distância. \[ w_{ij} = \left\{ \begin{array}{ccc} f(d_{ij}) & \mbox{se} & \leq d_{ij} \leq d_{max} \\ 0 & \mbox{se} & d_{ij} > d_{max} \\ \end{array} \right. \] Em que, \(f(d_{ij}) = d_{ij}^{-\alpha}\) é uma função inversa da distância, tal que \(\alpha=1\) é um parâmetro de alisamento.
3 Desenvolvimento
3.1 Carregando os pacotes e lendo o polígono
library(rgdal)
library(spdep)O pacote spdep fornece as funções necessárias à
construção das matrizes de ponderação espacial. O pacote
rgdal auxilia na leitura e visualização do polígono.
O arquivo Mun_Amazonia_Legal_2022.shp é um polígono com
a distribuição geográfica dos 772 municípios da Amazônia Legal.
bmap <- readOGR(dsn = "Mun_Amazonia_Legal_2022.shp", layer = "Mun_Amazonia_Legal_2022")## OGR data source with driver: ESRI Shapefile
## Source: "/home/andre/Mun_Amazonia_Legal_2022.shp", layer: "Mun_Amazonia_Legal_2022"
## with 772 features
## It has 9 fields
## Integer64 fields read as strings: CD_UFhead(bmap@data)## CD_MUN NM_REGIAO CD_UF NM_UF SIGLA_UF NM_MUN AREA_TOT
## 0 1100015 Norte 11 Rondônia RO Alta Floresta D'Oeste 7067.127
## 1 1100023 Norte 11 Rondônia RO Ariquemes 4426.571
## 2 1100031 Norte 11 Rondônia RO Cabixi 1314.352
## 3 1100049 Norte 11 Rondônia RO Cacoal 3793.000
## 4 1100056 Norte 11 Rondônia RO Cerejeiras 2783.300
## 5 1100064 Norte 11 Rondônia RO Colorado do Oeste 1451.060
## AREA_INT PORC_INT
## 0 7067.127 100
## 1 4426.571 100
## 2 1314.352 100
## 3 3793.000 100
## 4 2783.300 100
## 5 1451.060 1003.2 Matriz de contiguidade tipo Queen
A função poly2nb gera uma lista de contiguidade para
cada unidade espacial do objeto bmap, classe
"sp". A opção queen = TRUE considera que a
vizinhança deve ter pelo menos um ponto de fronteira comum. O critério
de contiguidade Rook é atribuído com
queen = FALSE.
# Lista de contiguidade baseada no critério Queen
lnb <- poly2nb(pl = bmap, queen = TRUE)
str(lnb[1:10])## List of 10
## $ : int [1:7] 14 20 21 23 26 29 46
## $ : int [1:9] 11 13 19 27 28 30 40 49 51
## $ : int [1:3] 6 43 665
## $ : int [1:7] 9 16 18 20 33 38 744
## $ : int [1:3] 6 7 43
## $ : int [1:7] 3 5 7 22 34 43 665
## $ : int [1:4] 5 6 34 43
## $ : int [1:3] 10 46 47
## $ : int [1:5] 4 16 22 690 744
## $ : int [1:6] 8 23 24 31 36 47summary(lnb)## Neighbour list object:
## Number of regions: 772
## Number of nonzero links: 4382
## Percentage nonzero weights: 0.7352546
## Average number of links: 5.676166
## Link number distribution:
##
## 1 2 3 4 5 6 7 8 9 10 11 12 13 14
## 1 20 62 137 179 155 92 55 40 17 6 3 4 1
## 1 least connected region:
## 638 with 1 link
## 1 most connected region:
## 758 with 14 linksExaminado pontualmente o polígono, observa-se que Ariquemes em Rondônia, faz fronteira com outros nove municípios, os quais podem ser identificados da seguinte forma
lnb[[2]]## [1] 11 13 19 27 28 30 40 49 51bmap$CD_MUN[c(11,13,19,27,30,40,49,51)]## [1] "1100114" "1100130" "1100262" "1100403" "1100601" "1101401" "1101609"
## [8] "1101757"bmap$NM_MUN[c(11,13,19,27,30,40,49,51)]## [1] "Jaru" "Machadinho D'Oeste" "Rio Crespo"
## [4] "Alto Paraíso" "Cacaulândia" "Monte Negro"
## [7] "Theobroma" "Vale do Anari"A função summary(lnb) retorna um resumo da matriz. Entre
as 595984 ligações possíveis, existem 4382 (0,74%) não nulas, tal que a
frequência média se aproxima de 6 ligações. Baseado no critério
Queen de contiguidade, o gráfico abaixo oferece um panorama
geral das ligações entre municípios.
plot(bmap, border = 'grey', main = "Figura 1: Arranjo espacial Queen contiguidade")
plot(lnb, coordinates(bmap), add = TRUE, col = 'red', pch = 1, cex = .5)
Usando a função nb2listw, inserindo no argumento a lista
lnb, cria-se, então, uma matriz de contiguidade de primeira
ordem na forma de lista. Usando a função nb2mat, então é
possível convertê-la em uma matriz numérica. A opção
style = "W" indica que a matriz é normalizada na linha.1
wq <- nb2listw(neighbours = lnb, style = "W")
print(wq)## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 772
## Number of nonzero links: 4382
## Percentage nonzero weights: 0.7352546
## Average number of links: 5.676166
##
## Weights style: W
## Weights constants summary:
## n nn S0 S1 S2
## W 772 595984 772 291.0116 3209.968wqm <- nb2mat(lnb)
wqm[1:6, 1:6]## [,1] [,2] [,3] [,4] [,5] [,6]
## 0 0 0 0.0000000 0 0.0000000 0.0000000
## 1 0 0 0.0000000 0 0.0000000 0.0000000
## 2 0 0 0.0000000 0 0.0000000 0.3333333
## 3 0 0 0.0000000 0 0.0000000 0.0000000
## 4 0 0 0.0000000 0 0.0000000 0.3333333
## 5 0 0 0.1428571 0 0.1428571 0.00000003.3 Matriz k-vizinhos mais próximos.
Esta matriz adota a distância entre centróides como critério de
vizinhança. Usando a função coordinates do pacote
rgdal, o primeiro passo é determinar pontos de coordenadas
para os 772 municípios da Amazônia Legal. Essa informação é utilizada na
criação de uma lista de objetos (municípios), formada pelos k=2
vizinhos mais próximos. A função knearneigh retorna uma
matriz com a identificação dos objetos pertencentes ao conjunto
k=2.
coords <- coordinates(bmap)
lk2nb <- knn2nb(knearneigh(coords, k=2))
wk2nb <- nb2listw(lk2nb,zero.policy = TRUE)
print(wk2nb)## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 772
## Number of nonzero links: 1544
## Percentage nonzero weights: 0.2590674
## Average number of links: 2
## Non-symmetric neighbours list
##
## Weights style: W
## Weights constants summary:
## n nn S0 S1 S2
## W 772 595984 772 648 3304str(wk2nb$neighbours[1:10])## List of 10
## $ : int [1:2] 21 26
## $ : int [1:2] 19 30
## $ : int [1:2] 5 6
## $ : int [1:2] 9 38
## $ : int [1:2] 7 43
## $ : int [1:2] 3 7
## $ : int [1:2] 5 34
## $ : int [1:2] 46 47
## $ : int [1:2] 4 16
## $ : int [1:2] 8 24Ao fixar k=2, esta matriz pondera exatamente os dois vizinhos mais próximos. Ressalta-se que, no escopo da econometria espacial, a determinação do parâmetro k é uma questão que merece cuidado. Na falta de um bom critério de seleção, deve-se escolher o número k que proporciona o melhor ajuste aos resíduos do modelo espacial.2
plot(bmap, border = 'grey', main = "Figura 2: Arranjo espacial k=2 vizinhos mais próximos")
plot(lk2nb, coords, add = TRUE, col = 'red', pch = 1, cex = .5)
3.4 Matriz distância de corte
A matriz distância de corte é definida com a função
dnearneigh, que tem recebe no seu argumento: uma matriz de
pontos de coordenadas dos objetos (municípios), a distância inferior e
superior de corte. Neste caso, a distância máxima de corte é extraída de
um vetor de distância entre objetos. Assim, o primeiro passo é gerar uma
lista de vizinhança de objetos.
# Lista de vizinhança
lknb <- knn2nb(knearneigh(coords))
# Vetor de distância entre objetos
dnb <- unlist(nbdists(lknb,coords))
summary(dnb)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.0570 0.2060 0.3305 0.4116 0.5272 2.0968dmax <- max(dnb)
# Lista de vizinhança
dlnb <- dnearneigh(coords, 0, dmax)
summary(dlnb)## Neighbour list object:
## Number of regions: 772
## Number of nonzero links: 39786
## Percentage nonzero weights: 6.675683
## Average number of links: 51.53627
## Link number distribution:
##
## 1 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21
## 4 5 9 12 16 9 14 10 10 11 9 14 10 7 12 12 12 13 12 9
## 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41
## 5 13 6 14 8 7 4 8 12 5 5 8 13 10 14 9 6 11 11 7
## 42 43 44 45 46 47 48 49 50 51 52 53 55 57 58 59 60 61 62 63
## 7 1 2 3 5 1 2 2 3 2 4 4 1 3 2 6 7 4 4 8
## 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83
## 6 7 6 9 8 9 6 9 12 7 8 3 10 6 7 7 6 6 9 3
## 84 85 86 87 88 89 90 91 92 93 95 96 97 98 99 100 101 102 103 104
## 11 7 9 7 7 5 2 4 5 4 3 1 6 3 4 10 4 5 3 4
## 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124
## 1 3 6 4 6 5 3 4 1 1 3 2 3 2 2 3 2 3 5 2
## 125 126 128
## 4 1 1
## 4 least connected regions:
## 82 124 126 159 with 1 link
## 1 most connected region:
## 569 with 128 links# Converte em matriz
wdc <- nb2listw(dlnb)
print(wdc)## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 772
## Number of nonzero links: 39786
## Percentage nonzero weights: 6.675683
## Average number of links: 51.53627
##
## Weights style: W
## Weights constants summary:
## n nn S0 S1 S2
## W 772 595984 772 69.33981 3114.637Observe que as ligações variam entre objetos (no limite, existem 4
municípios com uma ligação e um município com 128 ligações), para uma
distância de corte equivalente a dmax. Com resultado, o
grau de conectividade da matriz distância de corte, estimado em 6,67%, é
superior às duas primeiras.
plot(bmap, border = 'grey', main = "Figura 3: Arranjo espacial da matriz distância de corte")
plot(dlnb, coords, add = TRUE, col = 'red', pch = 1, cex = .5)
3.5 Matriz função inversa da distância
Este arranjo espacial assemelha-se ao anterior, porém, difere em um
aspecto: os pesos são definidos como uma função inversa da distância.
Assim, o primeiro argumento da função nb2listwdist é
preenchido com a lista de vizinhança da matriz anterior. A opção
type="idw" define a função inversa da distância, ao passo
que style = "W" indica que a matriz é normalização na
linha. 3
wid <- nb2listwdist(neighbours = dlnb, x = bmap, type="idw", style = "W")
str(wid$weights[1:10])## List of 10
## $ : num [1:35] 0.0181 0.0229 0.0277 0.0187 0.0274 ...
## $ : num [1:33] 0.0143 0.04 0.0229 0.0246 0.0163 ...
## $ : num [1:20] 0.0298 0.0846 0.1719 0.0799 0.0354 ...
## $ : num [1:39] 0.0173 0.0136 0.0129 0.0158 0.0485 ...
## $ : num [1:23] 0.0348 0.0647 0.0226 0.0617 0.1334 ...
## $ : num [1:22] 0.0271 0.1511 0.0247 0.071 0.0854 ...
## $ : num [1:24] 0.032 0.0566 0.0243 0.1236 0.0689 ...
## $ : num [1:15] 0.0553 0.1022 0.0541 0.0831 0.0532 ...
## $ : num [1:34] 0.0196 0.0676 0.0189 0.0199 0.0226 ...
## $ : num [1:16] 0.0456 0.0976 0.0458 0.0602 0.1015 ...4 Referências
ALMEIDA, E. Econometria Espacial Aplicada. Campinas, SP: Editora Alínea, 2012.
GOLGHER, A. B. Introdução à Econometria Espacial. Jundiaí, SP: Paco Editorial, 2015.
KOPCZEWSKA, K. Applied spatial statistics and econometrics: data analysis in R. London and New York: Routledge, 2021.