它的全局指数形式如下：

其中，表示的z得分（减去均值再除以标准差）。

局部指数形式如下：

下面介绍另外两种空间自相关指数。

1 Geary's C指数

全局Geary's C指数的计算公式如下：

402 Payment Required

  • Geary's C指数的形式和莫兰指数很相似，区别主要在于交叉乘积项不同：Moran's I的交叉乘积项为，Geary's C的交叉乘积项为；
  • Geary's C指数的范围在0-2之间，其中0 < C<1表示正相关，1<C<2表示负相关。

对应的计算函数为geary.test()：

geary.test(x, listw, randomisation = TRUE,
           zero.policy = NULL,
           alternative = "greater",
           spChk = NULL, adjust.n = TRUE)

局部Geary's C指数的计算公式如下：

对应的计算函数为localC()：

localC(x, listw, ..., zero.policy=NULL)
  • 需要注意的是，使用localC()函数需要升级spdep工具包至最新版（version 1.2-2）,而该版本需要在R 4.0版本上才能安装，因此还需要更新R的版本。

示例数据：

library(sf)
library(spdep)
nc <- st_read(system.file("shape/nc.shp", package = "sf"))
nb <- poly2nb(nc)
w <- nb2listw(nb)

全局C指数计算：

geary.test(nc$PERIMETER, w)
##  Geary C test under randomisation
## 
## data:  nc$PERIMETER 
## weights: w 
## 
## Geary C statistic standard deviate = 4.8111, p-value = 7.504e-07
## alternative hypothesis: Expectation greater than statistic
## sample estimates:
## Geary C statistic       Expectation          Variance 
##       0.647379018       1.000000000       0.005371826

局部C指数的计算：

localC(nc$PERIMETER, w)
##   [1]  0.51871787  1.10439458  0.41156194  5.31073565  0.97425908  0.95669222
##   [7]  3.19510123  0.28153922  1.21040904  0.15333337  0.08101971  0.15217679
##  [13]  0.71528016  0.12772676  0.95567198  1.17663802  0.26277716  1.65514050
##  [19]  0.60799549  0.19324480  0.08778058  0.52086725  0.82627745  0.67515755
##  [25]  0.18556866  0.15568586  0.62861626  0.83662284  0.63068066  1.42171090
##  [31]  0.48092260  0.18685971  0.63813298  0.44709774  0.73871203  1.66111416
##  [37]  1.11051027  0.35764638  0.80139349  0.79071385  1.80516436  0.37353909
##  [43]  0.40289037  5.48715784  3.63503582  0.41605715  0.24066795  1.87750650
##  [49]  0.60811979  0.24171111  2.51226979  0.25992709  0.93348509  0.43351900
##  [55]  0.45803593  2.83362107  7.68988661  2.74305490  1.08559576  3.19789972
##  [61]  0.55841615  0.59652820  0.73013432  0.13412889  0.41134715  0.61242858
##  [67]  1.03866518  0.67725251  0.67055591  0.11752121  0.22153419  0.52679953
##  [73]  2.43721783  1.73913298  0.54755929  0.26162366  1.59008583  0.84028967
##  [79]  0.61778575 12.05212920  0.83458465  0.44389302  1.63895880  0.27868399
##  [85]  0.10929799  1.50588998  5.03877814  0.26278862  0.75834165  1.47374432
##  [91]  4.60894025  1.87370371  1.46066858  1.19153837  2.82518333  0.14154167
##  [97]  0.99003582  0.38531353  5.11559152  1.69555793

2 Getis's G指数

全局Getis's G（Getis-Ord G）指数的计算公式如下：

402 Payment Required

对应的计算函数为globalG.test()：

globalG.test(x, listw, zero.policy = NULL,
             alternative = "greater", spChk = NULL,
             adjust.n = TRUE, B1correct = TRUE,
             adjust.x = TRUE, Arc_all_x = FALSE)
  • 当G指数的p值高于一定的显著性水平时，可以认为属性数据在空间上的分布是随机的；
  • 当G指数的p值低于一定的显著性水平时：
  • 若G大于其期望，表明属性数据在空间上存在聚集关系；
  • 若G小于其期望，表明属性数据在空间上存在离散关系；

具体可参见下方链接：

https://desktop.arcgis.com/zh-cn/arcmap/10.3/tools/spatial-statistics-toolbox/h-how-high-low-clustering-getis-ord-general-g-spat.htm

局部Getis's G指数的计算公式如下：

  • 当与可以相等时，记为指数。

对应的计算函数为localG()：

localG(x, listw, zero.policy = NULL,
       spChk = NULL, return_internals = FALSE,
       GeoDa = FALSE, alternative = "two.sided")

全局G指数的计算：

globalG.test(nc$PERIMETER, w)
##  Getis-Ord global G statistic
## 
## data:  nc$PERIMETER 
## weights: w 
## 
## standard deviate = 5.6461, p-value = 8.207e-09
## alternative hypothesis: greater
## sample estimates:
## Global G statistic        Expectation           Variance 
##       1.071096e-02       1.010101e-02       1.167074e-08
  • 可以看出，p < 0.05，G指数大于其期望。

局部G指数的计算：

localG(nc$PERIMETER, w)
##   [1] -0.633137157  0.028820686 -0.922143434  2.858745245  0.821878839
##   [6]  0.752315078  0.484315290 -1.746351560  0.295300220 -0.308633464
##  [11] -1.137410882 -0.977756619 -0.947039263 -1.188634241 -0.485902087
##  [16]  0.585017963 -1.463186507 -1.626462213 -0.583522388 -1.742398091
##  [21] -1.589460318 -0.804721755 -0.291046287  0.222382978 -0.899412011
##  [26] -0.453853301 -0.144302278  0.922360278 -0.856613487  0.160474792
##  [31]  0.222384607 -1.025540410  0.443801229 -1.124325789 -0.134775806
##  [36]  1.888435511  0.099613050 -0.035151360 -1.369734697 -0.453476232
##  [41]  0.009338262 -0.607088241 -0.970638610  2.833257795  1.016537905
##  [46] -0.583292856  0.597369932 -0.535678840 -0.433061934 -0.656799428
##  [51]  1.448702286 -0.904375114 -0.638748412  0.390006877  0.292187543
##  [56]  3.834122380  1.914940477 -0.526895898 -0.438797941  0.932795565
##  [61] -0.525281971 -0.124121208  1.241330496 -0.513672733 -0.184781374
##  [66]  0.397259293 -0.315610851 -0.646323829  0.049768053  0.320677178
##  [71] -0.243188022 -0.461480372  0.784949321  1.104454153  0.056627839
##  [76] -0.300251763 -0.163243851 -0.438909775  1.203591657  4.076388139
##  [81] -1.262048100  1.059846089  2.501316984 -0.421288765  0.073772108
##  [86]  0.452753928  3.159267921  1.134225652 -0.542509368  0.008052234
##  [91]  2.526605364  0.412310916  2.297652120  0.388078153  2.116573237
##  [96]  2.117293317  1.370675716  1.727242961  1.215320044  0.615774058
## attr(,"gstari")
## [1] FALSE
## attr(,"call")
## localG(x = nc$PERIMETER, listw = w)
## attr(,"class")
## [1] "localG"

局部G*指数的计算：

w2 <- nb2listw(include.self(nb))
localG(nc$PERIMETER, w2)
##   [1] -0.79019323 -0.42834051 -0.88592128  3.81168957  1.22453322  0.65510155
##   [7]  0.29430688 -1.92745237  0.06306807 -0.50210428 -1.31833293 -0.94877982
##  [13] -0.88155610 -1.38873027 -1.02508008  0.92650914 -1.78793333 -1.33858746
##  [19] -0.87692256 -2.19951435 -1.91656475 -1.11989474 -0.60393089  0.12558069
##  [25] -0.98410262 -0.41907796 -0.34919372  1.24253873 -1.10156100 -0.18574699
##  [31]  0.18672027 -1.25785084  0.28307296 -1.17073995 -0.45548983  1.93667132
##  [37]  0.41828641 -0.10457755 -1.24442102 -0.80734811 -0.54205953 -0.48963630
##  [43] -0.85876649  2.32248337  0.52520544 -0.55234560  0.65925225 -0.19535013
##  [49] -0.68317542 -0.67748717  1.45230950 -1.01959494 -0.37347110  0.49471088
##  [55]  0.28626539  5.22400265  2.86060610  0.11529943 -0.91212124  0.21495321
##  [61] -0.41095189 -0.12978968  1.17794565 -0.61178820 -0.39307452  0.29319542
##  [67] -0.10891967 -0.41797437 -0.24598642  0.34491892 -0.38901149 -0.59258143
##  [73]  0.20014053  0.93716403 -0.20814782 -0.47066247 -0.89196929 -0.40138670
##  [79]  1.51421332  3.01680769 -1.09067911  1.11431248  2.55471076 -0.34133198
##  [85]  0.01921657  0.07529056  3.73220070  1.21101278 -0.34900458 -0.68048696
##  [91]  3.23034253 -0.07210971  2.17584866  0.63141404  3.01878276  2.30414963
##  [97]  1.62885003  2.17584866  0.18975070  0.89739474
## attr(,"gstari")
## [1] TRUE
## attr(,"call")
## localG(x = nc$PERIMETER, listw = w2)
## attr(,"class")
## [1] "localG"
  • include.self()函数的作用是在计算空间权重矩阵时包括的情况。

绘图示例：

library(ggplot2)
library(RColorBrewer)
nc$gstar <- localG(nc$PERIMETER, w2)
ggplot(nc) +
  geom_sf(aes(fill = gstar)) +
  scale_fill_stepsn(n.breaks = 5,
                    colors = brewer.pal(5, "Spectral")) +
  guides(fill = guide_colourbar(title = expression(G^{"*"}))) +
  theme_bw()



空间自相关 python 莫兰 空间自相关数据_slam

在以上绘图过程中，除了ggplot2绘图系统外，还涉及到之前介绍的两个知识点：

对连续型变量使用离散型配色，见推文ggplot2 | 如何对连续型变量使用离散型调色板进行配色

图例标题存在特殊格式，见推文grDevices | 如何在图形中使用数学表达式作为标注文本