Python空间分析| 01 利用Python计算全局莫兰指数(Global Moran's I)

酱肉包- 2021-04-07 20:18:30
Python spatial


全局空间自相关

空间自相关(spatial autocorrelation)是指一些变量在同一个分布区内的观测数据之间潜在的相互依赖性。Tobler(1970)曾指出“地理学第一定律:任何东西与别的东西之间都是相关的,但近处的东西比远处的东西相关性更强”

全局莫兰指数(Global Moran's I)是最常用的空间自相关指数,用来反映全局的空间相关性,其计算公式为:

I = N W i j w i j ( x i x ˉ ) ( x j x ˉ ) i ( x i x ˉ ) 2 I = \frac N W \frac {\sum_i \sum_j w_{ij} (x_i-\bar x) (x_j-\bar x)} {\sum_i (x_i-\bar x)^2}

式中, N N 表示空间单元的数量, i i j j 表示空间单元的位置索引, x x 为分析的变量, x ˉ \bar x 表示分析变量的均值, w i j w_{ij} 表示空间权重矩阵, W W 表示空间权重矩阵的和

Moran's I大于0时,表示数据呈现空间正相关,其值越大空间相关性越明显;Moran's I小于0时,表示数据呈现空间负相关,其值越小空间差异越大;Moran's I为0时,空间呈随机性

解读莫兰指数的时候,需要有P值和Z得分来判定,P值小于0.05(通过95%置信度检验),且Z得分超过临界值1.65(拒绝零假设设定的阈值)

工具介绍

  • PySAL是一个基于Python进行探索性空间数据分析的开放源码库,是由亚利桑那州立大学GeoDa Center for Geospatial Analysis and Computation赞助的社区项目

  • libpysal:是Python空间分析库核心库,它提供了四个模块,构成了PySAL系列的空间分析工具:

    • 空间权重: libpysal.weights

    • 输入和输出:libpysal.io

    • 计算几何学:libpysal.cg

    • 内置示例数据集 libpysal.examples

  • esda 是一个用于空间数据探索性分析的开放源码 Python 库,它是 PySAL (Python Spatial Analysis Library)的一个子包,包括全局和局部空间自相关分析的方法。

import esda
import pandas as pd
import geopandas as gpd
from geopandas import GeoDataFrame
import libpysal as lps
import numpy as np
import matplotlib.pyplot as plt
from shapely.geometry import Point
import contextily as ctx
from pylab import figure, scatter, show
%matplotlib inline
root_dir="/root/learning/pysal/实验1/"
gdf = gpd.read_file(root_dir+'data/econ90040212.shp') # 读取数据

数据概况

gdf.columns.values #字段名
array(['CODE', 'COUNT', 'SUM_AREA', 'FIRST_ANAM', 'OID_', 'CODE_1',
'DATAFLAG', 'TOTPOP', 'TOTPOP_10K', 'RURPOP_10K', 'TOWNPOP_10',
'AGRPRODUCT', 'AGRLBR_10K', 'AGRSTOTGDP', 'FSTGDPRATE',
'SCNDGDPRAT', 'THRDGDPRAT', 'Province', 'geometry'], dtype=object)
gdf.head(1) # 查看一条记录

计算全局空间自相关指数Moran's I

  • 第一产业占GDP比重FSTGDPRATE为变量
ax=gdf.plot(figsize=(8,8),column="FSTGDPRATE")
ax.set_axis_off()

可以看出有一定的空间聚集性,但聚集性程度的大小需要通过Moran's I的计算进一步衡量

计算空间权重矩阵

df = gdf
wq = lps.weights.Queen.from_dataframe(df)# 使用Quuen式邻接矩阵
wq.transform = 'r' # 标准化矩阵
('WARNING: ', 98, ' is an island (no neighbors)')
('WARNING: ', 108, ' is an island (no neighbors)')
('WARNING: ', 137, ' is an island (no neighbors)')
('WARNING: ', 138, ' is an island (no neighbors)')
('WARNING: ', 139, ' is an island (no neighbors)')
/usr/local/lib/python3.6/dist-packages/libpysal/weights/weights.py:172: UserWarning: The weights matrix is not fully connected:
There are 6 disconnected components.
There are 5 islands with ids: 98, 108, 137, 138, 139.
warnings.warn(message)

可以看到有5个单元没有邻接对象

接下来将计算得到的权重矩阵绘制出来进行观察

centroids = gdf.geometry.centroid # 计算多边形几何中心
fig = figure(figsize=(8,8))
plt.plot(centroids.x, centroids.y,'.')
for k,neighs in wq.neighbors.items():
#print(k,neighs)
origin = centroids[k]
for neigh in neighs:
segment = centroids[[k,neigh]]
plt.plot(segment.x, segment.y, '-')
plt.title('Queen Neighbor Graph')
plt.axis('off')
plt.show()

wr = lps.weights.Rook.from_dataframe(df) # 使用Rook式邻接矩阵
# wr.transform = 'r' # 标准化矩阵
fig = figure(figsize=(8,8))
plt.plot(centroids.x, centroids.y,'.')
for k,neighs in wr.neighbors.items():
#print(k,neighs)
origin = centroids[k]
for neigh in neighs:
segment = centroids[[k,neigh]]
plt.plot(segment.x, segment.y, '-')
plt.title('Rook Neighbor Graph')
plt.axis('off')
plt.show()

计算Moran'I

y=gdf["FSTGDPRATE"]
mi = esda.moran.Moran(y, wq)
print("Moran's I 值为:",mi.I)
print("随机分布假设下Z检验值为:",mi.z_rand)
print("随机分布假设下Z检验的P值为:",mi.p_rand)
print("正态分布假设下Z检验值为:",mi.z_norm)
print("正态分布假设下Z检验的P值为:",mi.p_norm)
Moran's I 值为: 0.7133088262106428
随机分布假设下Z检验值为: 16.7858372659958
随机分布假设下Z检验的P值为: 0.0
正态分布假设下Z检验值为: 16.826039087303325
正态分布假设下Z检验的P值为: 0.0

p<0.001,说明能够在99.9%的置信度下证明GDP中第一产业占比存在空间聚集性,空间相关性显著

绘制Moran's I散点图

Moran散点图分为四个象限,分别表示四种局部空间相关关系:

  • 第一象限(HH):高值区被高值包围
  • 第二象限(LH):低值区被高值包围
  • 第三象限(LL):低值区被低值包围
  • 第四象限(HL):高值区被低值包围
from splot.esda import plot_moran
plot_moran(mi, zstandard=True, figsize=(10,4))
plt.show()

根据概率分布图,可以看出莫兰指数的值远大于大于随机正态分布的期望,这表明空间现象并不是随机分布的,而是具有一定相关性

根据散点图结果,可以看出全局Moran’s I 指数为0.708,并且大于0.5,因此其空间相关性较高

参考链接

版权声明
本文为[酱肉包-]所创,转载请带上原文链接,感谢
https://my.oschina.net/jiangroubao/blog/5011462

  1. Python web menu project takes another step forward to learn about the built-in user authentication system from the application layer
  2. Python classic interview questions (with answers)!
  3. 【Python从零到壹】Python的循环结构详解
  4. 【Python从零到壹】Python列表详解
  5. 【Python从零到壹】Python的字典详解
  6. 【Python从零到壹】Python的字符串详解
  7. 【Python从零到壹】Python基础之函数的应用
  8. 【Python从零到壹】用Python实现植物大战僵尸里的面向对象
  9. Detailed explanation of Python loop structure
  10. Detailed explanation of Python list
  11. Detailed explanation of Python dictionary
  12. Detailed explanation of Python string
  13. [Python from zero to one] the application of Python basic functions
  14. [Python from zero to one] using Python to realize object-oriented in plant vs. zombie
  15. 用 Python 实现微信版飞机大战
  16. 用 Python 实现***帝国中的数字雨落既视感
  17. 想知道未来孩子长相?Python人脸融合告诉你
  18. 我用 Python 做了一个全球疫情数据大屏
  19. Using Python to realize wechat aircraft war
  20. Using Python to realize the visual sense of digital rain in the Empire of the Communist Party of China
  21. Want to know what kids will look like in the future? Python face fusion tells you
  22. I made a big screen of global epidemic data with Python
  23. python你TM太皮了——区区30行代码就能记录键盘的一举一动
  24. Python you TM too skinny - just 30 lines of code can record every move of the keyboard
  25. python的装饰器概念学习基础基础版
  26. Python decorator concept learning basic edition
  27. SQL配合Python-Flask的中转注入
  28. python3使用kivy生成安卓程序
  29. 不到 150 行代码写一个 Python 版的贪吃蛇
  30. Transfer injection of SQL and python flash
  31. Using Kivy to generate Android program in Python 3
  32. Less than 150 lines of code to write a python version of the snake
  33. Python面向对象练习题
  34. Python数据分析入门(八):Pandas统计计算和描述
  35. Python面向对象练习题
  36. Python object oriented exercises
  37. Introduction to Python data analysis (8): Pandas statistical calculation and description
  38. Python object oriented exercises
  39. WEB4-通过python获得flag
  40. python-web5
  41. Pandas-二进制操作
  42. python入门教程14-01 (python语法入门之python内存泄露)
  43. Web4 - get flag through Python
  44. python-web5
  45. Pandas binary operation
  46. python入门教程13-06 (python语法入门之视图、触发器、事务、存储过程、函数)
  47. python入门教程13-07 (python语法入门之ORM框架SQLAlchemy)
  48. python入门教程13-08 (python语法入门之python索引原理与慢查询优化)
  49. 定投指数到底能不能赚钱?Python 来告诉你答案
  50. Python入门学习之:10分钟1500访问量
  51. Getting started with Python 14-01
  52. 用 Python 画哆啦 A 梦
  53. Python 图表利器 pyecharts
  54. 用 Python 抓取公号文章保存成 HTML
  55. Introduction to Python 13-06 (view, trigger, transaction, stored procedure, function of introduction to Python syntax)
  56. Getting started with Python 13-07 (ORM framework Sqlalchemy for getting started with Python syntax)
  57. Introduction to Python 13-08
  58. Can fixed investment index make money? Python will tell you the answer
  59. Introduction to Python: 1500 visits in 10 minutes
  60. 用 Python 获取股市交易数据