R로 공간적 자기상관 나타내기

공간적 자기상관(spatial autocorrelation)은 어떤 지역의 특성이 주변 지역과 얼마나 유사한지를 나타내는 개념이다. 토블러의 지리학 제1법칙(모든 것은 먼 것보다 가까운 것에 영향을 더 받는다)에 따라 공간 데이터는 통계학의 제1법칙인 데이터 관측값의 독립성을 위배한다. 따라서 공간데이터를 분석할 때는 나타나는 전역적/국지적 공간적 자기상관을 파악하여야 한다. 공간적 자기상관은 전체 데이터 집합에 대해서 살펴보느냐, 아니면 특정 위치에서 주변과의 관계를 분석하느냐에 따라 전역적(Global) 공간적 자기상관과 국지적(Local) 공간적 자기상관으로 나뉜다. 국지적 공간적 자기상관 분석은 LISA(Local Indicator of Spatial Association) 클러스터라고도 불린다.

전역적 공간적 자기상관

전역적 공간적 자기상관은 전체 데이터 집합이 얼마나 서로 연관을 가지고 분포해있는가를 나타낸 것이다. 즉, 주어진 데이터의 전체 공간에서 값들이 얼마나 유사하거나 상이한지를 평가한 것이다. 이를 통해 특정 특성(예: 기온, 집값 등)이 공간적으로 군집되어 있는지, 아니면 무작위로 분포되어 있는지를 알 수 있다. 전역적 공간적 자기상관을 측정하기 위한 대표적인 통계지표는 전역적 모란지수(Global Moran's I)이다. 전역적 모란지수의 값이 양수라면, 유사한 값들이 가까운 지역에 모여 있다는 것을 의미한다. 예를 들어, 집값이 최고로 높은 아파트 주변 집값도 높은 경우를 들 수 있다. 값이 음수라면, 서로 다른 값들이 가까운 지역에 모여 있다는 것을 의미한다. 아래 사진의 체스판같은 배열이 그것이다. 값이 0에 가깝다면 좌측 하단에서 볼 수 있듯 무작위에 가깝게 분포되어 있음을 뜻한다.

출처: 위키피디아

전역적 공간적 자기상관 분석은 전체 연구 영역을 요약하는 유일한 통계를 생성한다.  다시 말해, 전역적 공간적 자기상관 분석은 전체 영역의 동질성을 전제한다는 것이다. 이 전제가 들어맞지 않는다면 통계가 공간에 따라 달라야 하므로 통계가 하나만 있는 것은 의미가 없다. 그렇기 때문에 지리학적 연구에서 전역적 공간적 자기상관은 항상 국지적(지역적) 공간적 자기상관과 함께 분석되는 경향이 있다.

국지적 공간적 자기상관 분석 (LISA 클러스터)

특정 지역의 주변 지역과의 상관을 분석하는 척도로는 Global Moran's I에서 파생된 국지적 모란지수(Local Moran's I)와 게티스-오드 지아이스타(Getis-Ord's Gi*)가 있다. Local Moran's I는 특정 지점에서 주변 지점들과의 관계를 분석하여 공간적 자기상관을 평가한다. 즉, 특정 지점의 값이 주변 지점들의 값과 얼마나 유사한지를 나타낸다. 즉, '패턴 찾기'에 특화된 것이다. Moran's I에서 파생된 Getis-Ord's Gi*는 특정 지점 주변의 값들이 전반적으로 높거나 낮은지를 평가한다. 즉, 특정 지점을 중심으로 주변 지점들의 값이 집합적으로 클러스터링 되어 있는지(핫스팟) 또는 분산되어 있는지(콜드스팟)를 나타낸다. '덩어리 찾기'에 특화된 셈이다. 

 

R 코드

전국 행정동 shp: https://www.vworld.kr/dtmk/dtmk_ntads_s002.do?dsId=30017

브이월드

활용한 데이터: 저소득 독거노인 비율

저소득독거노인.csv

0.01MB

 

데이터 전처리

library(sf)
library(tmap)
library(tidyverse)
korea<-st_read('C:/path/BND_ADM_DONG_PG.shp', options="ENCODING=EUC-KR")
seoul <- korea %>% filter(substr(korea$ADM_CD,1,2)=="11") %>% arrange(ADM_NM)
seoul$ADM_NM <- gsub("·", "", seoul$ADM_NM)
old_alone<-read.csv("C:path/저소득독거노인.csv")
old_alone <- old_alone %>% arrange(ADM_NM)
seoul_oa<-merge(seoul, old_alone, all.x=TRUE)
dim(seoul_oa)
seoul_oa$저소득독거노인<-as.numeric(seoul_oa$저소득독거노인)
seoul_oa[is.na(seoul_oa)] <- 0
tm_shape(seoul_oa)+tm_polygons(col="저소득독거노인", style="quantile", n=8, palette="Greens")+tm_legend(outside=TRUE)

Global Moran's I

library(spdep)
library(sp)

#폴리곤에서 이웃 행렬 추출
nb <- poly2nb(seoul_oa, queen=TRUE)

#이웃 가중치 설정
#각 행의 가중치가 그 행의 합으로 나눠져 합이 1이 되도록 표준화 (W)
nb_w <- nb2listw(nb, style="W", zero.policy = TRUE)    #0.46
#인접 관계만을 나타내며 가중치는 0 또는 1로 설정 (B)
nb_b <- nb2listw(nb, style="B", zero.policy = TRUE)   #0.44

#Global Moran's I 계산
moran.test(seoul_oa$저소득독거노인, nb_w, randomisation = TRUE)
moran.test(seoul_oa$저소득독거노인, nb_b, randomisation = TRUE)

 

Local Moran's I

# LISA 클러스터 데이터 생성
MC <- moran.mc(seoul_oa$저소득독거노인, nb_w, nsim=999)
oa.lmi <- localmoran(seoul_oa$저소득독거노인, nb_w, na.action = na.exclude)
local_moran.map <- cbind(seoul_oa, oa.lmi)

# 클러스터 할당
quadrant <- vector(mode="numeric", length = nrow(oa.lmi))
m.qualification <- seoul_oa$저소득독거노인 - mean(seoul_oa$저소득독거노인)
m.local <- oa.lmi[,1] - mean(oa.lmi[,1])
signif <- 0.05
quadrant[m.qualification > 0 & m.local > 0] <- 4 
quadrant[m.qualification < 0 & m.local < 0] <- 1     
quadrant[m.qualification < 0 & m.local > 0] <- 2
quadrant[m.qualification > 0 & m.local < 0] <- 3
quadrant[oa.lmi[,5] > signif] <- 0

# 클러스터를 tmap용 데이터프레임에 추가
local_moran.map$LISA <- factor(quadrant, levels = 0:4, labels = c("Non significant", "Low-Low", "Low-High", "High-Low", "High-High"))

# 색상 팔레트 정의
palette <- c("white", "blue", rgb(0, 0, 1, alpha = 0.4), rgb(1, 0, 0, alpha = 0.4), "red")

# tmap으로 시각화
tm_shape(local_moran.map) +
  tm_polygons("LISA", palette = palette, title = "LISA Cluster Map") +
  tm_layout(legend.title.size = 0.8, legend.text.size = 0.6)

Getis-Ord's Gi*

library(RColorBrewer)
library(sp)

#Local G 계산
local_g <- localG(seoul_oa$저소득독거노인, nb_w)

#벡터데이터에 벡터 형식으로 붙이기
seoul_oa <- cbind(seoul_oa, as.vector(local_g))

#열 번호 확인하고 열 이름 변경
head(seoul_oa)
names(seoul_oa)[6] <- c("Gi")

#색상 지정
MoranColours<-rev(brewer.pal(7, "RdBu"))
breaks1<-c(-3.00,-2.58,-1.96,-1.65,1.65,1.96,2.58,15.00)

#지도 출력
tm_shape(seoul_oa) + tm_borders() +
  tm_polygons(col="Gi", breaks=breaks1, palette = MoranColours, title="Z-score (Gi*)")