期刊图片复现|Python绘制SHAP重要性蜂巢图+条形图+相对贡献饼图组图

# =========================================================================================# ====================================== 1. 库的导入 =========================================# =========================================================================================import numpy as npimport pandas as pdimport matplotlib.pyplot as pltimport matplotlib.gridspec as gridspecimport matplotlibimport xgboost as xgbimport shap

# =========================================================================================# ====================================== 2.颜色库 =========================================# =========================================================================================COLOR_SCHEMES = {    1: (['#e41a1c', '#377eb8', '#4daf4a'], LinearSegmentedColormap.from_list('c1', ['#008080', '#FFFFFF', '#FF7F50'])),}

# =========================================================================================# ======================================3.蜂群图辅助函数=======================================# =========================================================================================def simple_beeswarm(x_values, nbins=40, width=0.1):    np.random.seed(42)    hist_range = (np.min(x_values), np.max(x_values))  #数据的最小值和最大值范围    if hist_range[0] == hist_range[1]:  # 如果最大值等于最小值        hist_range = (hist_range[0] - 0.1, hist_range[1] + 0.1)  #手动扩展范围        current_width = (counts[i] / max_count) * width  # 根据当前箱子的密度计算抖动宽度        ys = np.linspace(-current_width, current_width, len(idxs))  # 在宽度范围内生成均匀分布的Y值        np.random.shuffle(ys)  # 打乱Y值顺序        y_values[idxs] = ys  # 将计算好的Y值赋给对应的数据点    return y_values  # 返回计算好的Y轴抖动坐标

# =========================================================================================# ====================================== 4.3D饼图 =========================================# =========================================================================================def draw_3d_pie_chart_on_ax(ax, sizes, colors, labels, tilt_angle=0.4, depth_layers=80,total_thickness=0.5):    step = total_thickness / depth_layers  #计算每一层的厚度    # 循环按层绘制实现堆叠的3D效果    for i in range(depth_layers):    #绘制表层饼图    wedges, texts, autotexts = ax.pie(        sizes, #数据        colors=colors, #颜色        startangle=90, #角度        radius=1.2, #半径        autopct='%1.2f%%', #文本        pctdistance=0.6, #文本距圆心的相对距离        shadow=False #关闭默认的阴影效果    )    #文本样式    for autotext in autotexts:        autotext.set_fontsize(11)  #字体大小        autotext.set_fontweight('bold')  #加粗        autotext.set_color('black')  #颜色        autotext.set_zorder(depth_layers + 3)  #层    ax.set_aspect(tilt_angle)  #角度    #加图例    ax.legend(wedges,#图例对象              labels,  #标签              loc="center left",  #位置              bbox_to_anchor=(0.92, 0.5),  #坐标              frameon=False,  #边框              fontsize=14,  #字体大小              handletextpad=0.5)  #间距

# =========================================================================================# ======================================5.shap重要性散点组图=========================================# =========================================================================================def plot_custom_shap_dual_axis(ax_bar, ax_scatter, shap_values, X_df, scheme_colors, cmap,max_features=15):    mean_abs_shap = np.abs(shap_values.values).mean(axis=0)  #每个特征SHAP绝对值的平均数    order = np.argsort(mean_abs_shap)[::-1]  #排序    top_indices = order[:max_features]  #前15个最重要的特征    other_indices = order[max_features:]  #其它特征    #拼接名字    top_features = [X_df.columns[i] for i in top_indices] + ["other features"]    top_features.reverse()  #倒序    y_pos = np.arange(len(top_features))  #每个特征的y坐标    total_mean = sum(mean_abs_shap)  #所有特征的总重要性，用来求百分比    bar_color_pos = matplotlib.colors.to_rgba(scheme_colors[0], alpha=0.25)  #正向色    bar_color_neg = matplotlib.colors.to_rgba(scheme_colors[1], alpha=0.25)  #负向色    bar_color_neutral = matplotlib.colors.to_rgba('grey', alpha=0.15)  #其它特征颜色        if feature == "other features":            bar_color = bar_color_neutral  #使用灰色        else:            bar_color = bar_color_pos if corr >= 0 else bar_color_neg  #大于0用正向色，小于0用负向色        #百分比        ax_bar.text(max(plot_mean_abs) * 1.1,  #x                    idx,  #y                    f"{pct:.2f}%",  #文本                    va='center',  #垂直                    ha='left',  #水平                    fontsize=14)  #大小    ax_bar.set_yticks(y_pos)  #左侧y轴刻度位置    ax_bar.set_yticklabels(top_features, fontsize=14)  #y轴刻度标注特征文本    ax_bar.set_xlabel('mean(|SHAP value|)', fontsize=16, fontweight='bold')  #x轴标题    #边框    ax_scatter.spines['bottom'].set_visible(False)    ax_scatter.spines['right'].set_visible(False)

# =========================================================================================# ======================================6.组图绘制函数=========================================# =========================================================================================def plot_advanced_forest_chart(analysis_results, scheme_id, selected_colors_tuple, grid_cols=3):    scheme_colors = selected_colors_tuple[0]  #饼图配色    cmap_name = selected_colors_tuple[1]  #蜂巢图配色    hr = [3, 1] * grid_rows  #每行高度比例，主图:饼图    #创建画布    gs = gridspec.GridSpec(grid_rows * 2,  #行                           grid_cols,  #列                           height_ratios=hr,  #高度比                           hspace=0.05,  #上下间距                           wspace=0.45)  #左右间距    titles = [f'({chr(97 + i)}) {target_name}' for i, target_name in enumerate(targets_keys)]  #标题    cmap = plt.get_cmap(cmap_name) if isinstance(cmap_name, str) else cmap_name  #获取对应的Colormap对象        ax_pie = fig.add_subplot(gs_bottom[0])  #饼图轴        #饼图绘制        draw_3d_pie_chart_on_ax(ax_pie,  #轴                                res['pie_sizes'],  #数据                                scheme_colors,  #颜色                                ['X_3D', 'X_ND', 'X_LM'],#图例标签                                tilt_angle=0.6,  #角度                                total_thickness=0.45,  #厚度                                depth_layers=60)  #层数        #饼图标题        ax_pie.set_title('Relative Contribution',  #文本                         x=1.38,#x                         y=-0.32,  #y                         fontsize=20,  #字体大小                         fontweight='bold')  #加粗        ax_txt = fig.add_subplot(gs_bottom[1])  #右侧评估结果部分        ax_txt.axis('off')  #清空        text_content = f"R² = {res['R2']:.3f}\n\nMSE = {res['MSE']:.3f}\n\nRMSE = {res['RMSE']:.3f}\n\nMAE = {res['MAE']:.3f}"  #拼接文本

# =========================================================================================# ======================================7.执行部分 =========================================# =========================================================================================if __name__ == '__main__':    #类别划分    features_3d = ['H_MN', 'H_SD', 'Volume', 'H_SUM']    features_nd = ['NDVI_MN', 'NDVI_SD', 'NDVI_SUM']    features_lm = ['NP', 'MESH', 'SHAPE', 'COHESION', 'FRAC', 'PD', 'AI', 'LSI', 'LPI', 'TE', 'PAFRAC', 'CA','DIVISION', 'PARA', 'PLAND', 'SPLIT', 'ED']    excel_data_path = r"data.xlsx"  #原始数据    excel_file = pd.ExcelFile(excel_data_path)  #读取    target_names = excel_file.sheet_names  #读取表名就是目标名    #遍历每个目标    for target_name in target_names:        df_sheet = pd.read_excel(excel_data_path, sheet_name=target_name)  #读取当前表        #提取特征数据和目标数据        if target_name in df_sheet.columns:            y = df_sheet[target_name]            X = df_sheet.drop(columns=[target_name])        else:            y = df_sheet.iloc[:, -1]            X = df_sheet.iloc[:, :-1]        #划分数据        X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)        model = xgb.XGBRegressor(random_state=42)  #初始化XGBoost回归模型        #初始化网格搜索        grid = GridSearchCV(model, param_grid, cv=3, scoring='neg_mean_squared_error')        grid.fit(X_train, y_train)  #拟合        best_model = grid.best_estimator_  #最佳模型    print("-" * 50)    target_columns = 3  #组合图每一行的子图数        plot_advanced_forest_chart(df_real, scheme_id, selected_colors_tuple, grid_cols=target_columns)

features_3d = ['H_MN', 'H_SD', 'Volume', 'H_SUM']features_nd = ['NDVI_MN', 'NDVI_SD', 'NDVI_SUM']features_lm = ['NP', 'MESH', 'SHAPE', 'COHESION', 'FRAC', 'PD', 'AI', 'LSI', 'LPI', 'TE', 'PAFRAC', 'CA','DIVISION', 'PARA', 'PLAND', 'SPLIT', 'ED']

2.设置原始数据的保存地址，执行部分：

excel_data_path = r"data.xlsx"  #原始数据

param_grid = {'max_depth': [3, 4],              'n_estimators': [50, 100],              'learning_rate': [0.05, 0.1]}

4.设置是否要进行批量绘图，执行部分：

plot_all = True

plt.savefig(fr"advanced_shap_scheme_{scheme_id}.png",            dpi=300, bbox_inches='tight')