机器学习+城市规划第十二期:利用SHAP来分析GWR的回归过程
引言
随着城市化进程的不断推进,城市规划需要依赖大量的数据分析和建模来实现更加精细化、精准化的决策。地理加权回归(GWR)作为空间数据分析中的一种重要方法,已经被广泛应用于城市规划领域。GWR通过考虑空间异质性,能够揭示出每个地理位置的局部回归关系,进而帮助我们更好地理解城市内部的空间结构和特征。
在实际应用中,如何理解模型的预测结果和每个特征对模型输出的贡献是至关重要的。SHAP(Shapley Additive Explanations)值提供了一种直观的方式来解释模型的预测结果,并可以帮助我们理解模型如何根据每个特征做出决策。在本篇博客中,我们将结合GWR模型和SHAP值,探索如何分析GWR的回归过程,揭示空间数据中的局部特征对预测结果的贡献。
1. 依赖库与数据处理
首先,我们需要导入一些常用的库,包括pandas
、geopandas
、numpy
等,以及GWR和SHAP相关的工具。接下来,我们加载并清洗数据,并构建一个GeoDataFrame对象,用于后续的地理空间分析。
# ===== 导入依赖库 =====
import pandas as pd
import geopandas as gpd
from shapely.geometry import Point
import numpy as np
np.float = float
from mgwr.gwr import GWR
from mgwr.sel_bw import Sel_BW
import matplotlib
matplotlib.use('TkAgg')
import matplotlib.pyplot as plt
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus'] = False
from gwr_tools import save_gwr_summary
在数据预处理阶段,我们加载一个名为“shuju bike.csv
”的数据文件,并对其中的经纬度坐标、因变量(Y)和自变量(X)进行清洗。
# ===== 1. 读取数据 =====
df = pd.read_csv("shuju bike.csv", header=None)
df.columns = ["lon", "lat", "Y", "X"]
# 清洗数据
for col in ["lon", "lat", "Y", "X"]:
df[col] = pd.to_numeric(df[col], errors="coerce")
df.dropna(subset=["lon", "lat", "Y", "X"], inplace=True)
df.reset_index(drop=True, inplace=True)
在这个步骤中,我们将列名指定为lon
、lat
、Y
、X
,并确保数据中的无效值(如NaN
和Inf
)被清除。
接着,我们将经纬度坐标转换为空间数据,并将其转换为适合后续分析的投影坐标系(EPSG: 3857)。
# 构造空间数据
geometry = [Point(xy) for xy in zip(df["lon"], df["lat"])]
gdf = gpd.GeoDataFrame(df, geometry=geometry)
gdf.set_crs(epsg=4326, inplace=True)
gdf = gdf.to_crs(epsg=3857)
2. 数据建模与GWR拟合
为了分析GWR模型的回归过程,我们需要使用mgwr
库中的GWR
类来构建和拟合模型。在GWR中,模型的关键参数之一是带宽(Bandwidth),它控制了地理邻近样本的影响范围。我们首先通过Sel_BW
类来选择最优带宽,然后利用GWR模型进行拟合。
# ===== 2. 提取建模数据 =====
coords = np.array([(point.x, point.y) for point in gdf.geometry])
y = gdf["Y"].values.reshape((-1, 1))
X = gdf["X"].values.reshape((-1, 1))
# ===== 3. 带宽选择 + 模型拟合 =====
print("正在选择最优带宽,请稍候...")
selector = Sel_BW(coords, y, X)
bw = selector.search()
print(f"选定的最优带宽为:{bw}")
model = GWR(coords, y, X, bw)
results = model.fit()
print("\n======= GWR 模型结果汇总 =======")
print(results.summary())
在这里,Sel_BW
类的search
方法会自动选择最优带宽,GWR
类则用于根据选择的带宽拟合回归模型。模型拟合完成后,我们打印出模型的汇总结果,帮助我们更好地理解每个自变量的回归系数和模型的拟合效果。
3. 模型结果保存与可视化
GWR模型的结果中包含了每个地理位置的局部回归系数,这些系数揭示了空间数据的空间异质性。我们将这些结果保存为GeoJSON文件,并进行可视化展示。
# ===== 4. 保存结果到 GeoJSON,用于地图展示 =====
gdf["GWR_coef_X"] = results.params[:, 0]
gdf["GWR_resid"] = results.resid_response
gdf["local_R2"] = results.localR2
gdf.to_file("gwr_result.geojson", driver="GeoJSON")
print("\n📍 局部系数结果已保存为 'gwr_result.geojson'。")
这样,您可以将每个观测点的局部回归系数、残差以及局部R²值保存为GeoJSON文件,方便进一步的空间数据可视化和分析。
4. SHAP分析与可视化
虽然GWR模型可以提供局部回归系数,但我们还需要使用SHAP值来解释模型的预测过程。SHAP值可以帮助我们了解每个特征对模型输出的具体贡献。在这部分中,我们将利用SHAP来解释GWR模型。
# ===== 5. 保存模型统计摘要为 CSV 文件(调用封装的函数) =====
save_gwr_summary(results, model, filename_prefix="bike_gwr")
利用SHAP,我们可以绘制每个特征的影响力图,帮助我们更好地理解空间特征如何影响GWR模型的预测结果。以下代码将展示如何计算SHAP值并绘制SHAP值的总结图。
# ===== 8. SHAP 值计算 =====
# 创建SHAP解释器
def gwr_predict(X):
local_predictions = np.dot(X, results.params[:, 1:].T) + results.params[:, 0]
return local_predictions.flatten()
# 标准化X数据
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)
# 使用SHAP计算值
explainer = shap.KernelExplainer(gwr_predict, X_scaled)
shap_values = explainer.shap_values(X_scaled)
# ===== 9. SHAP 可视化 =====
shap.summary_plot(shap_values, X_scaled, feature_names=["X"])
在这段代码中,我们首先定义了一个自定义的gwr_predict
函数来计算GWR模型的预测值,并用shap.KernelExplainer
来计算SHAP值。最后,我们使用shap.summary_plot
来可视化每个特征对预测结果的贡献。
5. 结论
通过将GWR模型与SHAP结合使用,我们不仅可以了解空间回归模型中每个变量的回归系数,还能够通过SHAP值解释模型的预测过程。SHAP值的可视化帮助我们深入理解每个自变量对模型输出的影响,从而为城市规划决策提供更加直观的支持。
未来,随着地理信息系统(GIS)和机器学习技术的进一步发展,GWR与SHAP结合的分析方法将为空间数据分析和城市规划提供更多的可能性和价值。
注:本文代码基于mgwr
库和shap
库,分别用于地理加权回归建模和机器学习模型解释。如果您希望了解更多关于GWR和SHAP的内容,建议参考相关文档并进行深入学习。
原创声明:本教程由课题组内部教学使用,利用CSDN平台记录,不进行任何商业盈利。