Python 空间解释性机器学习实战:基于 RF 与 GeoShapley 的全流程代码解析
- 2026-06-21 22:03:26
Python 空间解释性机器学习实战:基于 RF 与 GeoShapley 的全流程代码解析
绘图结果 
库的导入以及字体设置 
颜色库设置以及配色方案的选择与提取 
GeoShapley值空间分布图绘制函数 
COHESION特征的GeoShapley值空间分布图，揭示了该特征对模型预测结果实际影响的空间异质性。图中散点的大小代表了目标变量的数值高低，可以明显观察到东中部地区散点普遍较大，而西部较小；散点的颜色则直观量化了COHESION特征在当地的贡献方向与力度，其中红色系代表正向推高预测值的区域，而蓝绿色系代表产生负向拉低作用的区域。图中的两颗五角星代表了该特征发挥极值影响的坐标：最大正向贡献点位于东北地区，而最大负向贡献点则落在了西部青藏高原周边，这种说明了COHESION对预测目标的实际影响并非全局统一，而是高度依赖于其所处的具体地理环境。 
SVC值空间分布图绘制函数 
COHESION特征的空间变异系数SVC分布图，剥离了单个样本的绝对影响，通过平滑的地理加权系数揭示了该变量底层敏感度的空间渐变规律。图右侧的色条显示所有SVC系数值均位于0.012至0.025以上的正数区间内，这意味着COHESION特征与目标变量始终保持着正相关的边际效应。然而，这种正向驱动力在空间上呈现出极其显著的西强东弱梯度分布格局：深红色和浅红色的高值区高度集中在新疆、西藏等西部腹地，表明在这些区域，COHESION特征数值的增加就能带来强烈的目标变量响应；中东部及沿海广大地区，散点普遍呈现代表较低系数的蓝绿色，说明在这些区域，该特征对目标变量的敏感度和边际效能已经大幅减弱。 
蜂巢图绘制函数 
GeoShapley值的蜂群图，展示了各个特征对模型预测结果的整体影响程度，揭示了影响的正负方向，图中纵轴按特征重要性从上到下排序，横轴表示GeoShapley值的大小，即对预测结果的具体影响，而数据点的颜色从深蓝到深红代表了该特征原始数值从低到高的变化。可以明显看出，LUE是影响最大的核心特征，其高值显著正向推高预测结果，低值则产生负向拉低作用，紧随其后的ROAD和SLOPE也表现出较强的影响力；此外，该图不仅列出了基础特征，还清晰地展示了纯空间效应（GEO）以及各特征与地理空间的交互作用项（如LUE x GEO等），这种分离让我们能够直观地观察到空间位置属性是如何改变单个特征的原始影响力的，同时数据点的横向发散宽度进一步反映了这些效应在整个数据集所有样本中的波动范围和分布异质性。 
部分依赖图(PDP)绘制函数 
依赖图（PDP），它通过散点分布与红色平滑拟合曲线（GAM曲线）的结合，刻画了每一个特征的原始观测数值（x轴）与其对应的GeoShapley贡献值（y轴）之间复杂的非线性映射关系。通过观察这些曲线的走势，我们可以精准识别出各个变量的作用模式和临界阈值，例如LUE、POPE、COHESION和CR呈现出明显的正相关趋势，即特征原始值越大对模型的正向推力就越强；相比之下，GDPE、FRAC_CV、SLOPE和LJI则表现出负相关特征，随着数值升高其影响逐渐穿过零线转为负面拉低作用；特别值得注意的是类似ROAD和GDPE这样的特征，前者的曲线在特定数值（0.4至0.6之间）出现了急剧上升，后者的曲线则在0.5之后出现断崖式下跌，这揭示了这些变量在机器学习模型中存在着非常敏感的突变临界阈值。 
全局贡献度柱状图绘制函数 
全局特征贡献排序柱状图，它通过计算各个特征GeoShapley绝对值的平均值，直观地量化并对所有变量的整体重要性进行了从大到小的绝对排名。在这张图中，LUE、ROAD和SLOPE稳居总体贡献的前三位，该图通过不同颜色的堆叠，将每个特征的总贡献严格拆解为了非空间贡献（Non-Geo，蓝绿色柱体部分）和空间交互贡献（Geo，橙红色柱体部分）。图中名为Geo的单独条目作为一个整体展示了地理位置本身的重要性，而对于其他基础变量来说，柱子末端的橙色量化了空间异质性对该特征的附加影响程度，例如排名第一的LUE虽然绝大部分由其自身非空间属性决定，但仍包含了一定比例的空间交互作用，这种量化拆解方式使得我们能够精准评估地理环境因素在各个具体变量发挥作用时所扮演的权重角色。 
执行部分：设置输入输出路径和数据。读入矢量底图数据与建模用的Excel数据。划分训练集和测试集。使用 
期刊图片复现|Python绘制二维偏依赖PDP图 期刊复现|python绘制基于SHAP分析和GAM模型拟合的单特征依赖图 期刊图片复现|python绘制带有渐变颜色shap特征重要性组合图（条形图+蜂巢图） 期刊复现|用Python绘制SHAP特征重要性总览图、依赖图、双特征交互效应SHAP图，解锁XGBoost模型的终极奥秘 期刊图片复现|Python绘制shap重要性蜂巢图+单特征依赖图+交互效应强度气泡图+交互效应依赖图（回归+二分类+分类） 
公众号中的所有所有的免费代码都已经下架了，都并入到付费部分里了，付费合集代码和数据的购买通道已经开通，全部合集100元，后续将会持续更新，决定购买请后台私信我，注意只会分享练习数据和代码文件，不会提供答疑服务，代码文件中已经包含了每行代码的完整注释，购买前请确保真的需要！！!
代码绘制成果展示
论文：GeoShapley: A Game Theory Approach to Measuring Spatial Effects in Machine Learning Models
代码解释
第一部分
# =========================================================================================# ====================================== 1. 环境设置 =======================================# =========================================================================================import osimport numpy as npimport pandas as pdimport matplotlibimport matplotlib.pyplot as pltimport matplotlib.colors as mcolors
第二部分
# =========================================================================================# ======================================2.颜色库=======================================# =========================================================================================COLOR_SCHEMES = {1: ['#086c76', '#3dbbc7', '#ffd5cf', '#f19c83', '#cc2a2e'],}scheme_id = 5#设置要使用的配色方案colors = COLOR_SCHEMES[scheme_id] #提取配色方案
第三部分
注意：所有内容均为我的个人理解，可能存在错误或不足，使用时还要阅读原文
# =========================================================================================# ====================================== 3.GeoShapley值空间绘图函数===========================# =========================================================================================def plot_spatial_shap_china(df_plot, plot_values, china_boundary, feature_name, save_dir, plot_type):#创建画布fig, ax = plt.subplots(figsize=(12, 8), dpi=120)# 将shp绘制在图上china_boundary.plot(ax=ax, #坐标轴facecolor='none', #填充色edgecolor='black', #边界线颜色linewidth=1.0, #边界线宽度zorder=1) #层idx_max = np.argmax(plot_values) #数值最大位置idx_min = np.argmin(plot_values) #数值最小位置minx, miny, maxx, maxy = china_boundary.total_bounds #四至范围，最小经/纬度，最大经/纬度buffer = 2.0 #边缘的缓冲距离#标题ax.set_title(f'{feature_name} ({plot_type})', #文本fontsize=22, #大小pad=15, #标题与图的间距fontweight='bold') #加粗ax.grid(True, linestyle='-', color='#e0e0e0', alpha=0.7, zorder=0) #背景网格线#设置边框线for spine in ax.spines.values():spine.set_edgecolor('black')spine.set_linewidth(2.0)ax.tick_params(axis='both', which='major', length=8, width=2.0, labelsize=14) #主刻度样式color_labels = [f"{bounds[i]:.4f} ~ {bounds[i + 1]:.4f}" for i in range(5)]#颜色图例矩形color_legend_elements = [mpatches.Patch(color=colors[i], label=color_labels[i]) for i in range(5)]png_path = os.path.join(save_dir, f'spatial_{plot_type}_{feature_name}_{scheme_id}.png')pdf_path = os.path.join(save_dir, f'spatial_{plot_type}_{feature_name}_{scheme_id}.pdf')plt.savefig(png_path, dpi=300, bbox_inches='tight')plt.savefig(pdf_path, bbox_inches='tight')plt.close() #关闭
第四部分
注意：所有内容均为我的个人理解，可能存在错误或不足，使用时还要阅读原文
# =========================================================================================# ====================================== 4.SVC值空间绘图函数===========================# =========================================================================================def plot_spatial_svc_china(df_plot, plot_values, china_boundary, feature_name, save_dir, plot_type):#创建画布fig, ax = plt.subplots(figsize=(12, 8), dpi=120)# 将shp绘制在图上china_boundary.plot(ax=ax, #坐标轴facecolor='none', #填充色edgecolor='black', #边界线颜色linewidth=1.0, #边界线宽度zorder=1) #层minx, miny, maxx, maxy = china_boundary.total_bounds #四至ax.grid(True,linestyle='-',color='#e0e0e0',alpha=0.7,zorder=0)#边框线for spine in ax.spines.values():spine.set_edgecolor('black')spine.set_linewidth(2.0)#设置刻度样式ax.tick_params(axis='both', #双轴which='major', #主刻度length=8, #刻度线长度width=2.0, #刻度线粗细labelsize=14) #标签字号png_path = os.path.join(save_dir, f'spatial_{plot_type}_{feature_name}_{scheme_id}.png')pdf_path = os.path.join(save_dir, f'spatial_{plot_type}_{feature_name}_{scheme_id}.pdf')plt.savefig(png_path, dpi=300, bbox_inches='tight')plt.savefig(pdf_path, bbox_inches='tight')plt.close() #关闭
第五部分
注意：所有内容均为我的个人理解，可能存在错误或不足，使用时还要阅读原文
# =========================================================================================# ====================================== 5.GeoShapley蜂巢图绘图函数===========================# =========================================================================================def plot_custom_summary(geo_results, save_dir):ax.tick_params(axis='both',which='major',labelsize=14,width=1.8,length=6)#设置轴刻度标注大小for label in (ax.get_xticklabels() + ax.get_yticklabels()):label.set_fontweight('bold')ax.xaxis.label.set_size(16) #X轴标题大小ax.xaxis.label.set_weight('bold') #X轴标题加粗ax.yaxis.label.set_size(16) #Y轴标题大小ax.yaxis.label.set_weight('bold') #Y轴标题加粗#保存plt.savefig(os.path.join(save_dir, f"summary_plot_{scheme_id}.png"), dpi=300, bbox_inches='tight')plt.savefig(os.path.join(save_dir, f"summary_plot_{scheme_id}.pdf"),bbox_inches='tight')plt.close() #关闭
第六部分
注意：所有内容均为我的个人理解，可能存在错误或不足，使用时还要阅读原文
# =========================================================================================# ====================================== 6.部分依赖图(PDP)绘图函数===========================# =========================================================================================def plot_custom_pdp(geo_results, save_dir):#调用内置函数绘图geo_results.partial_dependence_plots(gam_curve=True, #平滑GAM曲线拟合edgecolors='white', #边颜色lw=0.8) #边线宽#画布fig = plt.gcf()fig.set_size_inches(15, 10) #调整子图尺寸# 手动控制所有子图离边框左边缘的留白比例fig.subplots_adjust(left=0.08, #左边缘right=0.88, #右边缘top=0.92, #上边缘bottom=0.08, #下边缘wspace=0.35, #水平间隙hspace=0.45) #垂直间隙ax.yaxis.label.set_weight('bold') #Y标题加粗path_col = ax.collections[0] #取出第一组点的渲染层offsets = path_col.get_offsets() #点的XY坐标target_vals = offsets[:, 1] #提取每个点的Y轴坐标数值norm = mcolors.Normalize(vmin=target_vals.min(), vmax=target_vals.max()) #映射[min, max]至[0,1]face_colors = cmap(norm(target_vals)) #使用规范化后的数值从之前定义的色带中提取出对应的一组颜色path_col.set_facecolors(face_colors) #覆盖默认颜色for label in cbar.ax.get_yticklabels():label.set_fontweight('bold')cbar.outline.set_linewidth(1.5)#保存plt.savefig(os.path.join(save_dir, f"partial_dependence_plots_{scheme_id}.png"), dpi=300, bbox_inches='tight')plt.savefig(os.path.join(save_dir, f"partial_dependence_plots_{scheme_id}.pdf"),bbox_inches='tight')plt.close() #关闭
第七部分
注意：所有内容均为我的个人理解，可能存在错误或不足，使用时还要阅读原文
# =========================================================================================# ====================================== 7.绘制特征全局贡献柱状图绘图函数===========================# =========================================================================================def plot_custom_contribution(geo_results, save_dir):for i, patch in enumerate(ax.patches):if i < half: #如果当前条柱属于前半部分，代表非空间解释部分patch.set_facecolor(colors[1]) #颜色else: # 如果属于后半部分，空间解释相关patch.set_facecolor(colors[-2])#创建图例non_geo_patch = mpatches.Patch(color=colors[1], label='Non-Geo')geo_patch = mpatches.Patch(color=colors[-2], label='Geo')#添加图例ax.legend(handles=[non_geo_patch, geo_patch], loc='lower right', fontsize=14, frameon=True)ax.xaxis.label.set_size(16) #X轴标题大小ax.xaxis.label.set_weight('bold') #x轴标题加粗ax.yaxis.label.set_size(16) #标题字体大小ax.yaxis.label.set_weight('bold') #标题加粗plt.tight_layout() #自动调整#保存plt.savefig(os.path.join(save_dir, f"contribution_bar_plot_{scheme_id}.png"), dpi=300, bbox_inches='tight')plt.savefig(os.path.join(save_dir, f"contribution_bar_plot_{scheme_id}.pdf"),bbox_inches='tight')plt.close() #关闭
第八部分
GridSearchCV寻找随机森林模型最佳超参数，进行精度评估。初始化GeoShapleyExplainer解释器并计算预测解释值。检验GeoShapley值的可加性。获取总体统计摘要、转换为标准的SHAP矩阵形式，并利用内置接口拟合GWR地理加权回归提取SVC系数矩阵，将计算出来的结果保存为本地的Excel文件。调用前面定义好的函数绘制分析图（蜂巢、PDP、贡献度）和它们对应的 GeoShapley 空间分布图和 SVC 空间分布图。# =========================================================================================# ====================================== 8.GeoShapley值空间绘图函数===========================# =========================================================================================if __name__ == "__main__":BASE_OUTPUT_DIR = r"haley的地理空间可解释性分析与可视化绘图" #保存路径GEOSHAPLEY_DIR = os.path.join(BASE_OUTPUT_DIR, "GeoShapley_Plots") #GeoShapley图路径SVC_DIR = os.path.join(BASE_OUTPUT_DIR, "SVC_Plots") #SVC图保存路径OTHER_DIR = os.path.join(BASE_OUTPUT_DIR, "Other_Analysis") #非空见图保存路径#不存在就创建os.makedirs(GEOSHAPLEY_DIR, exist_ok=True)os.makedirs(SVC_DIR, exist_ok=True)os.makedirs(OTHER_DIR, exist_ok=True)LOCAL_SHP_PATH = r"2023年省级.shp" #shp文件路径LOCAL_DATA_PATH = os.path.join(BASE_OUTPUT_DIR, "synthetic_china_data.xlsx") #模型数据路径china_boundary = gpd.read_file(LOCAL_SHP_PATH) # 读取shp数据df = pd.read_excel(LOCAL_DATA_PATH) # 读取模型数据features = ['LUE', 'GDPE', 'POPE', 'COHESION', 'FRAC_CV', 'IJI', 'ROAD', 'SLOPE', 'CR'] #特征geo_features = features + ['LON', 'LAT'] #补充经纬度X = df[geo_features] #提取特征数据y = df['Target_ER'] #提取目标数据#划分数据集X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)#设置网格参数param_grid = {'n_estimators': [50, 100,200,500],'max_depth': [6, 7,8,9],'min_samples_split': [2, 5]}#实例化网格搜索优化器对象grid_search = GridSearchCV(estimator=RandomForestRegressor(random_state=42),param_grid=param_grid,cv=3,scoring='r2',n_jobs=-1,verbose=2)
如何应用到你自己的数据
1.设置配色方案，颜色库部分：
scheme_id = 1#设置要使用的配色方案2.设置保存的路径，执行部分：
BASE_OUTPUT_DIR = r"基于RF与Geoshaley的地理空间可解释性分析与可视化绘图" #保存路径3.设置shp文件的保存地址，执行部分：
LOCAL_SHP_PATH = r"2023年省级.shp" #shp文件路径4.设置excel文件的保存地址，执行部分：
LOCAL_DATA_PATH = os.path.join(BASE_OUTPUT_DIR, "synthetic_china_data.xlsx") #模型数据路径5.设置常规特征，执行部分：
features = ['LUE', 'GDPE', 'POPE', 'COHESION', 'FRAC_CV', 'IJI', 'ROAD', 'SLOPE', 'CR'] #特征6.设置经纬度特征，执行部分：
geo_features = features + ['LON', 'LAT'] #补充经纬度7.设置目标变量，执行部分：
y = df['Target_ER'] #提取目标数据8.设置超参数，执行部分：
param_grid = {'n_estimators': [50, 100,200,500],'max_depth': [6, 7,8,9],'min_samples_split': [2, 5]}
9.设置是否进行批量绘图，执行部分：
plot_all = False推荐
获取方式
本文来自网友投稿或网络内容，如有侵犯您的权益请联系我们删除，联系邮箱：wyl860211@qq.com 。
