*資料前處理
library(sf);library(spdep);library(units)
setwd("D:/1092SA/Data")
TW=st_read("Popn_TWN2.shp",options="ENCODING=BIG5")
pop=TW$A0A14_CNT+TW$A15A64_CNT+TW$A65UP_CNT
area=st_area(TW)%>%set_units("km^2")%>%drop_units
dens=pop/area
定義鄰近,建立鄰近關係
TW.c=poly2nb(TW)
TW.cw=nb2listw(TW.c, zero.policy=T)
coords=st_centroid(TW)%>%st_coordinates
TW.k=knn2nb(knearneigh(coords, k=2))
TW.kw=nb2listw(TW.k)
TW.d=dnearneigh(coords, d1=0, d2=20000)
TW.dw=nb2listw(TW.d, zero.policy=T)
一、計算以下統計量與繪製圖表,說明其參數設定,並解釋其意義。
Moran’s I Statistic
M=moran.test(dens, listw=TW.cw, zero.policy=T); M
##
## Moran I test under randomisation
##
## data: dens
## weights: TW.cw n reduced by no-neighbour observations
##
##
## Moran I statistic standard deviate = 21.508, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.703816518 -0.002808989 0.001079383
Monte-Carlo Simulation
mc=moran.mc(dens,listw=TW.cw,nsim=999,zero.policy=T)
hist(mc$res, freq=TRUE, breaks=20, xlab="Simulated Moran's I",main="Monte-Carlo simulation")
abline(v=M$estimate[1], col="red")
Moran Scatter Plot
moran.plot (dens, TW.cw, xlab="dens", ylab="SpatialLag dens",zero.policy=T)
Moran Spatial Correlogram
cor=sp.correlogram(TW.c, dens, order=10, method="I", style="W",zero.policy=T)
print(cor); plot(cor)
## Spatial correlogram for dens
## method: Moran's I
## estimate expectation variance standard deviate Pr(I) two sided
## 1 (357) 0.70381652 -0.00280899 0.00107938 21.5081 < 2.2e-16
## 2 (357) 0.37701617 -0.00280899 0.00061233 15.3494 < 2.2e-16
## 3 (353) 0.14626861 -0.00284091 0.00039596 7.4934 6.71e-14
## 4 (349) 0.02460139 -0.00287356 0.00025198 1.7308 0.0834825
## 5 (349) -0.00634159 -0.00287356 0.00020052 -0.2449 0.8065285
## 6 (349) -0.04681396 -0.00287356 0.00016801 -3.3900 0.0006990
## 7 (349) -0.04513285 -0.00287356 0.00014538 -3.5048 0.0004569
## 8 (349) -0.01006903 -0.00287356 0.00013443 -0.6206 0.5348668
## 9 (349) -0.03484390 -0.00287356 0.00014026 -2.6995 0.0069441
## 10 (344) -0.12162522 -0.00291545 0.00016661 -9.1968 < 2.2e-16
##
## 1 (357) ***
## 2 (357) ***
## 3 (353) ***
## 4 (349) .
## 5 (349)
## 6 (349) ***
## 7 (349) ***
## 8 (349)
## 9 (349) **
## 10 (344) ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
General G
G=globalG.test(dens, listw=TW.cw,zero.policy=T); G
##
## Getis-Ord global G statistic
##
## data: dens
## weights: TW.cw
##
## standard deviate = 20.78, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Global G statistic Expectation Variance
## 1.098029e-02 2.808989e-03 1.546298e-07
二、利用以下三種不同的空間鄰近定義,計算Moran’s I coefficient,比較其數值的差異,並討論可能的原因。
Contiguity:相接相鄰
(Mc=moran.test(dens, listw=TW.cw, zero.policy=T))
##
## Moran I test under randomisation
##
## data: dens
## weights: TW.cw n reduced by no-neighbour observations
##
##
## Moran I statistic standard deviate = 21.508, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.703816518 -0.002808989 0.001079383
K-nearest Neighbors (KNN):最近的前K個
Mk=moran.test(dens, listw=TW.kw, zero.policy=T);Mk
##
## Moran I test under randomisation
##
## data: dens
## weights: TW.kw
##
## Moran I statistic standard deviate = 19.452, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.904261534 -0.002724796 0.002174078
Distance-based:距離在閾值內
Md=moran.test(dens, listw=TW.dw, zero.policy=T);Md
##
## Moran I test under randomisation
##
## data: dens
## weights: TW.dw n reduced by no-neighbour observations
##
##
## Moran I statistic standard deviate = 14.61, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.3677143051 -0.0027777778 0.0006431049