计算莫兰指数和Geary’s C 空间自相关程度

卷积核类型

常见的卷积核为Rook,Bishop,Queen，如上图所示。

Molan’s I

Geary’s C

代码实现为

# 利用空间统计量Moran和Geary计算遥感数据的自相关程度
import numpy as np
import pandas as pddef getMoranV(path,t=0,method="Moran"):''':param path: 图像路径:param t: 卷积核类型:param method: 空间自相关方法:return:'''# 读取数据wb1 = pd.read_excel(path, index_col=0)data = wb1.valuesn, m = wb1.shape# 定义卷积核# queen 卷积dir2 = [[1, 0], [-1, 0], [0, 1], [0, -1], [1, 1], [1, -1], [-1, 1], [-1, -1]]# Rook 卷积dir1 = [[1, 0], [0, 1], [-1, 0], [0, -1]]# Bishop 卷积dir3=[[1,1],[1,-1],[-1,1],[-1,-1]]if t==0:dir=dir1elif t==1:dir=dir2else:dir=dir3# 开辟O(n^4)权重矩阵# 这玩意很耗内存，所以可以考虑用局部计算代替全局计算a = [[0]*(n*m) for _ in range((n*m))]for i in range(n):for j in range(m):# 二维的点转到一维loc = i * m + j  # 一维上的坐标# 找到他在二维上的临界点for detX, detY in dir:x = i + detXy = j + detY# 转化为一维# 边界处理if 0 <= x < m and 0 <= y < m:a[loc][x * m + y] = 1  # 表示临接def moranIndex(data,weight):# 莫兰指数计算s_0=np.sum(np.sum(weight))n,m=len(data),len(data[0])x_hat=np.mean(data)up_sum=0down_sum=0for i in range(n):for j in range(m):# 下一轮for x in range(n):for y in range(m):if not (val:=weight[i*m+j][x*m+y]):continueup_sum+=val*(data[i][j]-x_hat)*(data[x][y]-x_hat)down_sum+=(data[i][j]-x_hat)**2return (n*m/s_0)*(up_sum/down_sum)def Geary_C(data,weight):# Geary's C计算s_0=np.sum(np.sum(weight))*2n, m = len(data), len(data[0])x_hat = np.mean(data)up_sum = 0down_sum = 0for i in range(n):for j in range(m):# 下一轮for x in range(n):for y in range(m):if not (val := weight[i * m + j][x * m + y]):continueup_sum += val * (data[i][j] - data[x][y])**2down_sum += (data[i][j] - x_hat) ** 2return ((n*m-1)/(s_0))*(up_sum/down_sum)if (M:=method.upper())=="MORAN":return moranIndex(data,a)elif M=="GEARY":return Geary_C(data,a)def getHelp():print("Method: {'Moran','Geary'}")print("t: {0:Rook,1:Queen,2:bishop}")if __name__ == '__main__':getHelp()path = r"data.xlsx"print("墨兰指数为: %s"%getMoranV(path,0))print("Geary's C为: %s"%getMoranV(path,0,"Geary"))

输出结果为

Method: {'Moran','Geary'}
t: {0:Rook,1:Queen,2:bishop}
墨兰指数为: 0.5579039646993681
Geary's C为: 0.4355131948652942

优化

我们发现，这样子需要开辟n4n^4n4大小的权重矩阵，且需要O(4n2+n4)O(4n^2+n^4)O(4n2+n4)也就是O(n4)O(n^4)O(n4)的时间复杂度，那能不能对其做一些剪枝，降低复杂度呢？

我们可以动态开辟空间权重矩阵，从而将空间复杂度降低到O(M∗n2)O(M*n^2)O(M∗n2)，其中，MMM表示卷积核的有效大小。(非空洞大小)

这样子也能将时间复杂度降低为O(M∗n2)O(M*n^2)O(M∗n2) ，也就是O(n2)O(n^2)O(n2)级别。

# 利用空间统计量Moran和Geary计算遥感数据的自相关程度
import numpy as np
import pandas as pddef getMoranV(path,val=1,t=0,method="Moran"):''':param path: 图像路径:param t: 卷积核类型:param method: 空间自相关方法:param val: 空间权重:return:'''# 读取数据wb1 = pd.read_excel(path, index_col=0)data = wb1.valuesn, m = wb1.shape# 定义卷积核# queen 卷积dir2 = [[1, 0], [-1, 0], [0, 1], [0, -1], [1, 1], [1, -1], [-1, 1], [-1, -1]]# Rook 卷积dir1 = [[1, 0], [0, 1], [-1, 0], [0, -1]]# Bishop 卷积dir3=[[1,1],[1,-1],[-1,1],[-1,-1]]if t==0:dir=dir1elif t==1:dir=dir2else:dir=dir3# 开辟O(n^4)权重矩阵# 这玩意很耗内存，所以可以考虑用局部计算代替全局计算a = [[] for _ in range(n*m)]a_weight=0for i in range(n):for j in range(m):# 找到他在二维上的临界点loc=i*m+jfor detX, detY in dir:x = i + detXy = j + detY# 转化为一维# 边界处理if 0 <= x < n and 0 <= y < m:a[loc].append([x,y,val])a_weight+=valdef moranIndex(data,weight):# 莫兰指数计算# s_0=np.sum(np.sum(weight))s_0=a_weightn,m=len(data),len(data[0])x_hat=np.mean(data)up_sum=0down_sum=0for i in range(n):for j in range(m):# 下一轮loc=i*m+jfor v in weight[loc]:up_sum+=v[2]*(data[i][j]-x_hat)*(data[v[0]][v[1]]-x_hat)down_sum+=(data[i][j]-x_hat)**2return (n*m/s_0)*(up_sum/down_sum)def Geary_C(data,weight):# Geary's C计算s_0=a_weight*2n, m = len(data), len(data[0])x_hat = np.mean(data)up_sum = 0down_sum = 0for i in range(n):for j in range(m):# 下一轮loc = i * m + jfor v in weight[loc]:up_sum += v[2] * (data[i][j] - data[v[0]][v[1]])**2down_sum += (data[i][j] - x_hat) ** 2return ((n*m-1)/(s_0))*(up_sum/down_sum)if (M:=method.upper())=="MORAN":return moranIndex(data,a)elif M=="GEARY":return Geary_C(data,a)def getHelp():print("Method: {'Moran','Geary'}")print("t: {0:Rook,1:Queen,2:bishop}")if __name__ == '__main__':getHelp()path = r"data.xlsx"print("墨兰指数为: %s"%getMoranV(path,val=1,t=0))print("Geary's C为: %s"%getMoranV(path,val=1,t=0,method="Geary"))

【空间分析】莫兰指数与Geary‘s Python实现相关推荐

空间分析 | 莫兰指数的计算
什么是莫兰指数? 根据百度百科的定义是"空间自相关系数的一种,其值分布在[-1,1],用于判别空间是否存在自相关." 简单的说就是判定一定范围内的空间实体相互之间是否存在相关关系, ...
全局莫兰指数_空间自相关 | 莫兰指数
空间自相关:是指一些变量在同一个分布区内的观测数据之间潜在的相互依赖性.其中,自相关中的"自"表示当你进行相关性观察统计量,是来源于不同对象的同一属性.Tobler(1970)曾指 ...
arcgis设置nodata值_新版白话空间统计（6）：在ArcGIS中实现莫兰指数计算
上一篇简单说了一下莫兰指数的计算原理和计算公式,如果是学生或者基础研究者,鼓励好好的学习一下手算或者编程计算,所谓的基础不牢,地动山摇--但是对于工程界和应用人士,特别是非基础学,重复造轮子是没有啥意 ...
空间分析：3-3.geoda计算莫兰指数
莫兰指数是一个地学统计概念,用来表示空间自相关性. 我们使用geoda来计算一下北京二手房的莫兰指数,看看它的空间自相关. 一.莫兰指数莫兰指数是最常用的空间自相关指标,最早由统计学家莫兰提出,所以 ...
Stata做空间杜宾模型、莫兰指数等操作
以下内容完全由本人在实际操作中搜集整理总结得到,很细致的介绍:从如何在stata中导入数据,怎么定义面板数据,再到如何做局部和全局空间相关性检验(莫兰指数)和空间杜宾模型等. 1.导入面板数据在ex ...
spdep | 除了莫兰指数，还有哪些指数可以衡量空间自相关性？
前面的推文中介绍了莫兰指数:spdep | 如何在R语言中计算空间自相关指数. 它的全局指数形式如下: 其中,表示的z得分(减去均值再除以标准差). 局部指数形式如下: 下面介绍另外两种空间自相关指数 ...
白话空间统计之：Moran's I（莫兰指数）
前两天聊了空间统计学里面的两个经典概念,今天来说说第一篇文章留下的大坑: Moran's I . 首先,Moran's I这个东西,官方叫做:莫兰指数,是澳大利亚统计学家帕特里克·阿尔弗雷德·皮尔斯· ...
空间统计：Moran's I（莫兰指数）
前两天聊了空间统计学里面的两个经典概念,今天来说说第一篇文章留下的大坑:Moran's I. 首先,Moran's I这个东西,官方叫做:莫兰指数,是澳大利亚统计学家帕特里克·阿尔弗雷德·皮尔斯·莫兰 ...
新版白话空间统计（6）：在ArcGIS中实现莫兰指数计算
CSDN的被爬虫专用声明:虾神原创,公众号\知乎:虾神说D 转发.转载和爬虫,请主动保留此声明. 上一篇简单说了一下莫兰指数的计算原理和计算公式,如果是学生或者基础研究者,鼓励好好的学习一下手算或者编 ...

【空间分析】莫兰指数与Geary‘s Python实现

计算莫兰指数和Geary’s C 空间自相关程度

【空间分析】莫兰指数与Geary‘s Python实现相关推荐

最新文章

热门文章