共计 3274 个字符，预计需要花费 9 分钟才能阅读完成。
技术背景
LISA（Local Indicators of Spatial Association）图是空间自相关分析的重要工具，能够帮助我们识别数据中的空间聚类模式（高值聚集或低值聚集）和异常值。在实际应用中，LISA 分析常用于：
- 公共卫生领域识别疾病高发区
- 城市规划中分析房价空间分布
- 环境科学中检测污染热点区域
数据准备
在开始分析前，需要确保数据满足以下要求：
- 空间数据应为面要素（Polygon）或点要素（Point）
- 属性表中包含待分析的数值型字段
- 数据采用适当的投影坐标系（推荐使用等面积投影）
- 检查并处理缺失值
推荐预处理步骤：
- 使用 ArcMap 的 ”Project” 工具统一坐标系
- 使用 ”Calculate Geometry” 添加面积 / 长度字段（如需标准化）
- 使用 ”Select by Attributes” 处理异常值
核心实现
1. 空间权重矩阵构建
空间权重矩阵定义了要素之间的空间关系，常见构建方法：
- 邻接法（Contiguity）
- 距离法（Distance Band）
- K 最近邻（K-Nearest Neighbors）
推荐使用 Queen 邻接法处理面数据：
import arcpy
import pysal as ps
# 构建空间权重矩阵
def create_weights(in_features):
w = ps.weights.Queen.from_shapefile(in_features)
return w2. 局部 Moran’s I 计算
关键参数说明：
- 输入字段：待分析的数值变量
- 空间权重：上一步构建的矩阵
- 排列次数：建议 999 次（蒙特卡罗模拟）
完整计算代码：
def calculate_lisa(in_features, field_name, weights):
# 读取属性值
y = np.array([row[0] for row in arcpy.da.TableToNumPyArray(in_features, [field_name])])
# 计算 LISA
lisa = ps.explore.moran.Moran_Local(y, weights, permutations=999)
# 将结果写入属性表
arcpy.AddField_management(in_features, "LISA_Score", "DOUBLE")
with arcpy.da.UpdateCursor(in_features, ["LISA_Score"]) as cursor:
for i, row in enumerate(cursor):
row[0] = lisa.Is[i]
cursor.updateRow(row)
return lisa3. 显著性评估
Moran’s I 结果的解读：
- 正 I 值：相似值聚集（高 - 高或低 - 低）
- 负 I 值：相异值聚集（高 - 低或低 - 高）
- p 值：统计显著性（通常取 0.05）
分类实现代码：
def classify_lisa(lisa, in_features):
# 创建分类字段
arcpy.AddField_management(in_features, "LISA_Cluster", "TEXT", field_length=10)
# 定义分类标准
sig = lisa.p_sim <= 0.05
quad = lisa.q
# 更新分类结果
with arcpy.da.UpdateCursor(in_features, ["LISA_Cluster"]) as cursor:
for i, row in enumerate(cursor):
if not sig[i]:
row[0] = "不显著"
elif quad[i] == 1:
row[0] = "高 - 高"
elif quad[i] == 2:
row[0] = "低 - 低"
elif quad[i] == 3:
row[0] = "低 - 高"
elif quad[i] == 4:
row[0] = "高 - 低"
cursor.updateRow(row)可视化技巧
在 ArcMap 中优化显示效果的实用方法：
- 使用分类渲染（Categories -> Unique Values）展示聚类类型
- 设置透明度（30-50%）避免重叠区域遮挡
- 添加显著性过滤（Definition Query: p_value <= 0.05）
- 使用暖色（红）表示高值，冷色（蓝）表示低值
- 添加空间参考比例尺（重要！）
避坑指南
1. 投影系统问题
常见错误：
- 未投影数据导致距离计算失真
- 不同图层坐标系不一致
解决方案：
- 统一使用等面积投影（如 Albers）
- 运行前检查坐标系：
sr = arcpy.Describe(in_features).spatialReference
if not sr.projectionName: # 未投影
raise ValueError("数据必须使用投影坐标系")2. 边界效应
边缘要素因邻居较少可能导致分析偏差，推荐：
- 使用缓冲区分析扩展研究区域
- 在结果解释时注明边界要素
- 考虑使用自适应带宽（Adaptive Kernel）
3. 性能优化
大数据量处理建议：
- 使用空间索引（arcpy.AddSpatialIndex_management）
- 分块处理（建议单次不超过 10,000 个要素）
- 关闭不必要的 ArcMap 窗口
- 考虑使用 64 位背景地理处理
实际应用建议
将 LISA 结果用于决策时应注意：
- 统计显著性≠实际意义，需结合领域知识判断
- 考虑时间维度变化（可进行时空聚类分析）
- 注意尺度效应（MAUP 问题）
- 结果需与其他空间分析方法交叉验证
完整的端到端实现脚本示例：
import arcpy
import numpy as np
import pysal as ps
# 参数设置
input_shp = r"C:\data\study_area.shp"
field_name = "COVID_RATE"
output_shp = r"C:\output\lisa_result.shp"
# 主流程
try:
# 1. 准备数据
arcpy.CopyFeatures_management(input_shp, output_shp)
# 2. 构建空间权重
w = ps.weights.Queen.from_shapefile(output_shp)
# 3. 计算 LISA
y = np.array([row[0] for row in arcpy.da.TableToNumPyArray(output_shp, [field_name])])
lisa = ps.explore.moran.Moran_Local(y, w, permutations=999)
# 4. 保存结果
arcpy.AddField_management(output_shp, "LISA_I", "DOUBLE")
arcpy.AddField_management(output_shp, "LISA_P", "DOUBLE")
arcpy.AddField_management(output_shp, "LISA_Type", "TEXT", 10)
with arcpy.da.UpdateCursor(output_shp,
["LISA_I", "LISA_P", "LISA_Type"]) as cursor:
for i, row in enumerate(cursor):
row[0] = lisa.Is[i]
row[1] = lisa.p_sim[i]
if lisa.p_sim[i] > 0.05:
row[2] = "NS"
else:
row[2] = ["HH", "LL", "LH", "HL"][lisa.q[i]-1]
cursor.updateRow(row)
print("LISA 分析完成，结果已保存至：", output_shp)
except Exception as e:
print("错误：", str(e))
arcpy.AddError(str(e))结语
通过上述方法，我们可以在 ArcMap 环境中高效完成 LISA 空间聚类分析。实际应用中，建议先在小规模测试数据上验证流程，再处理完整数据集。记得定期保存中间结果，并使用 Python 的 logging 模块记录处理日志。
LISA 分析结果应与领域知识结合解读，避免过度依赖统计显著性。当需要更复杂的空间模式分析时，可考虑进一步探索 Getis-Ord Gi* 统计量或时空聚类方法。
正文完
