Python 空间统计分析实战3种方法实现地理加权回归(GWR)与OLS结果对比地理加权回归Geographically Weighted Regression, GWR作为空间统计分析的重要工具近年来在计量地理学、城市规划、环境科学等领域获得广泛应用。与传统的普通最小二乘回归OLS相比GWR能够捕捉空间异质性为局部关系建模提供更精细的视角。本文将深入探讨三种Python实现GWR的方法并与OLS进行系统对比帮助研究人员选择适合的分析工具。1. 空间统计分析基础与环境准备空间统计分析的核心在于考虑地理要素的空间依赖性和异质性。地理学第一定律指出任何事物都与其他事物相关但邻近的事物比遥远的事物更相关。这为空间统计分析提供了理论基础。1.1 必要库安装与数据准备首先确保安装以下Python库pip install geopandas pysal mgwr statsmodels matplotlib我们使用波士顿房价数据集作为示例数据该数据集包含506个社区的各项特征import geopandas as gpd from libpysal.examples import load_example from mgwr.sel_bw import Sel_BW from mgwr.gwr import GWR import numpy as np import statsmodels.api as sm # 加载数据 boston load_example(Boston) gdf gpd.read_file(boston.get_path(boston.shp)) y gdf[MEDV].values.reshape((-1,1)) X gdf[[CRIM, RM, LSTAT]].values coords list(zip(gdf[X], gdf[Y]))1.2 空间权重矩阵构建空间权重矩阵是GWR分析的关键常用的构建方法包括K最近邻法每个点选择固定数量的邻近点距离阈值法设定距离阈值范围内的点为邻居反距离法权重随距离增加而衰减from libpysal.weights import KNN, DistanceBand # K最近邻法 knn KNN.from_array(coords, k10) # 距离阈值法 dist_band DistanceBand.from_array(coords, threshold1000)2. 三种GWR实现方法对比2.1 MGWR库实现MGWR是专门为地理加权回归设计的Python库提供了完整的GWR分析流程# 选择最优带宽 sel_bw Sel_BW(coords, y, X).search() print(f最优带宽{sel_bw}) # 构建GWR模型 gwr_model GWR(coords, y, X, sel_bw) results gwr_model.fit() # 结果分析 print(results.summary()) print(f调整R²{results.adj_rsq}) print(fAICc{results.aicc})MGWR库的优势在于自动带宽选择完整的结果输出支持多尺度GWR2.2 PySAL实现PySAL是空间分析的标准库也提供了GWR实现from spreg import GWR # 转换为PySAL格式 X sm.add_constant(X) gwr_pysal GWR(coords, y, X, bwsel_bw) results_pysal gwr_pysal.fit() # 获取局部参数估计 params results_pysal.betasPySAL版本更适合与其他空间分析工具链集成。2.3 手动实现核心算法理解GWR核心算法有助于深入掌握其原理def gwr_manual(y, X, coords, bw): n len(y) params np.zeros((n, X.shape[1])) for i in range(n): # 计算高斯权重 dists np.array([np.sqrt((coords[i][0]-c[0])**2 (coords[i][1]-c[1])**2) for c in coords]) weights np.exp(-0.5*(dists/bw)**2) # 加权最小二乘 W np.diag(weights) XtWX X.T W X XtWy X.T W y params[i] np.linalg.inv(XtWX) XtWy return params # 使用示例 manual_params gwr_manual(y, X, coords, sel_bw)3. OLS与GWR模型对比3.1 全局OLS模型实现# 全局OLS回归 ols_model sm.OLS(y, X).fit() print(ols_model.summary())3.2 模型性能对比指标我们使用以下指标对比模型性能指标说明计算公式R²模型解释方差比例1 - RSS/TSS调整R²考虑变量数的R²调整1 - (1-R²)(n-1)/(n-p-1)AICc修正的Akaike信息准则AIC 2k(k1)/(n-k-1)残差空间自相关残差的Morans I指数衡量残差的空间聚集程度3.3 实际对比结果下表展示了一个典型的对比结果模型R²调整R²AICc残差Morans I计算时间(s)OLS0.6730.6713024.50.32**0.02GWR0.8210.8152897.30.084.56注意**表示p0.01Morans I检验显著4. 空间异质性分析与可视化4.1 参数空间分布import matplotlib.pyplot as plt # 绘制RM变量的GWR系数空间分布 fig, ax plt.subplots(figsize(10,8)) gdf[gwr_rm] results.params[:,2] gdf.plot(columngwr_rm, legendTrue, schemequantiles, k5, cmapcoolwarm, axax) plt.title(RM变量的GWR系数空间分布) plt.show()4.2 带宽选择影响带宽是GWR的关键参数影响模型的空间尺度bws np.linspace(500, 5000, 10) aics [] for bw in bws: model GWR(coords, y, X, bw).fit() aics.append(model.aicc) plt.plot(bws, aics) plt.xlabel(带宽) plt.ylabel(AICc) plt.title(带宽选择与模型拟合)4.3 局部R²分析gdf[local_r2] results.localR2 fig, ax plt.subplots(figsize(10,8)) gdf.plot(columnlocal_r2, legendTrue, schemequantiles, k5, cmapviridis, axax) plt.title(局部R²空间分布) plt.show()5. 实际应用建议与注意事项5.1 方法选择指南OLS适用场景空间异质性不显著样本量较小(100)需要快速初步分析GWR适用场景明显的空间异质性样本量较大(100)需要探究局部关系5.2 常见问题解决方案带宽选择不稳定尝试多种选择方法(CV、AICc)考虑使用自适应带宽多重共线性问题计算局部VIF指标移除高度相关变量计算效率优化使用并行计算考虑空间子集分析5.3 高级技巧半参数GWR结合参数和非参数部分时空GWR加入时间维度多尺度GWR(MGWR)不同变量使用不同带宽from mgwr.sel_bw import Sel_BW from mgwr.gwr import GWR # 多尺度GWR sel_bw Sel_BW(coords, y, X, multiTrue).search() mgwr_model GWR(coords, y, X, sel_bw, multiTrue) results mgwr_model.fit()在实际项目中我发现GWR对空间权重矩阵的选择非常敏感。一次分析中使用KNN权重和距离阈值权重得到了截然不同的结果这提示我们需要通过敏感性分析来验证结论的稳健性。同时当处理大规模数据时MGWR库的内存优化版本或分布式计算框架可能更为适合。