模型可解释性

# 对比可解释模型和黑箱模型的决策透明度
from sklearn.linear_model import LogisticRegression
from sklearn.ensemble import RandomForestClassifier
import pandas as pd

# 线性模型：可以直接查看权重
lr = LogisticRegression()
lr.fit(X_train, y_train)
for feat, coef in zip(feature_names, lr.coef_[0]):
    print(f"  {feat}: {coef:+.4f}")

# 随机森林：无法直接读取规则
rf = RandomForestClassifier(n_estimators=100)
rf.fit(X_train, y_train)
print(f"  树数量: {rf.n_estimators}, 最大深度: {rf.max_depth_}")
print("  无法直接理解单棵树的决策规则")

# 评估模型可解释性需求的决策框架
def assess_explainability_need(domain, impact_level, stakeholder):
    """评估场景对可解释性的需求等级"""
    risk_matrix = {
        "medical": {"high": "必须", "medium": "必须", "low": "推荐"},
        "finance": {"high": "必须", "medium": "推荐", "low": "可选"},
        "marketing": {"high": "推荐", "medium": "可选", "low": "可选"},
        "judicial": {"high": "必须", "medium": "必须", "low": "必须"},
    }
    need = risk_matrix.get(domain, {}).get(impact_level, "可选")
    print(f"领域: {domain}, 影响级别: {impact_level}")
    print(f"利益相关者: {stakeholder}")
    print(f"可解释性需求: {need}")
    return need

assess_explainability_need("medical", "high", "患者和医生")

# 全局解释：特征重要性（所有样本的平均贡献）
from sklearn.ensemble import GradientBoostingClassifier
import numpy as np

model = GradientBoostingClassifier(n_estimators=100)
model.fit(X_train, y_train)

# 全局特征重要性
importances = model.feature_importances_
for feat, imp in sorted(zip(feature_names, importances), key=lambda x: -x[1]):
    print(f"  {feat}: {imp:.4f}")

# 这告诉我们哪些特征"总体上"最重要
# 但无法解释对单个样本的预测

# 局部解释：为什么这个特定样本被预测为正例？
def explain_single_prediction(model, x_instance, feature_names):
    """为单个样本生成本地解释（简化版）"""
    base_pred = model.predict_proba(x_instance.reshape(1, -1))[0][1]
    explanations = []
    for i in range(x_instance.shape[0]):
        x_perturbed = x_instance.copy()
        x_perturbed[i] = np.median(X_train[:, i])
        perturbed_pred = model.predict_proba(x_perturbed.reshape(1, -1))[0][1]
        contribution = base_pred - perturbed_pred
        explanations.append((feature_names[i], contribution))
    explanations.sort(key=lambda x: -abs(x[1]))
    for feat, contrib in explanations[:5]:
        direction = "+" if contrib > 0 else "-"
        print(f"  {direction} {feat}: |{contrib:.4f}|")
    return explanations

from sklearn.inspection import permutation_importance
from sklearn.ensemble import RandomForestClassifier
import matplotlib.pyplot as plt

# 训练模型
model = RandomForestClassifier(n_estimators=200, random_state=42)
model.fit(X_train, y_train)

# 计算置换重要性（在测试集上）
result = permutation_importance(
    model, X_test, y_test,
    n_repeats=30,
    random_state=42,
    scoring="accuracy"
)

# 排序并可视化
sorted_idx = result.importances_mean.argsort()
plt.boxplot(result.importances[sorted_idx].T,
            vert=False, labels=feature_names[sorted_idx])
plt.title("Permutation Importances (test set)")
plt.show()

for feat, imp, std in sorted(
    zip(feature_names, result.importances_mean, result.importances_std),
    key=lambda x: -x[1]
)[:5]:
    print(f"  {feat}: {imp:.4f} (+/- {std:.4f})")

# 置换重要性 vs 内置重要性：对比分析
from sklearn.inspection import permutation_importance

def compare_importance_methods(model, X_train, y_train, X_test, y_test, names):
    """对比内置特征重要性和置换重要性"""
    if hasattr(model, "feature_importances_"):
        built_in = model.feature_importances_
    else:
        built_in = np.zeros(len(names))
    perm_result = permutation_importance(
        model, X_test, y_test, n_repeats=20, random_state=42
    )
    perm_imp = perm_result.importances_mean
    comparison = pd.DataFrame({
        "内置重要性": built_in,
        "置换重要性": perm_imp,
        "差异": np.abs(built_in - perm_imp)
    }, index=names).sort_values("置换重要性", ascending=False)
    print(comparison.to_string())
    return comparison

compare_importance_methods(rf, X_train, y_train, X_test, y_test, feature_names)

import lime
import lime.lime_tabular
from sklearn.ensemble import RandomForestClassifier

# 训练黑箱模型
model = RandomForestClassifier(n_estimators=200, random_state=42)
model.fit(X_train, y_train)

# 创建 LIME 解释器
explainer = lime.lime_tabular.LimeTabularExplainer(
    X_train,
    feature_names=feature_names,
    class_names=["拒绝", "批准"],
    mode="classification",
    discretize_continuous=True
)

# 解释单个样本
idx = 42
explanation = explainer.explain_instance(
    X_test[idx],
    model.predict_proba,
    num_features=5,
    num_samples=5000
)

print(f"样本 {idx} 的预测: {model.predict_proba([X_test[idx]])}")
for feat, contrib in explanation.as_list():
    print(f"  {feat}: {contrib:+.4f}")

# LIME 解释稳定性分析
def lime_stability_check(explainer, model, x_instance, n_runs=10):
    """检查 LIME 解释的稳定性（多次运行）"""
    all_explanations = []
    for i in range(n_runs):
        exp = explainer.explain_instance(
            x_instance, model.predict_proba,
            num_features=5, num_samples=5000
        )
        all_explanations.append(dict(exp.as_list()))
    all_feats = set()
    for exp in all_explanations:
        all_feats.update(exp.keys())
    for feat in all_feats:
        contributions = [exp.get(feat, 0) for exp in all_explanations]
        mean_c = np.mean(contributions)
        std_c = np.std(contributions)
        stability = "稳定" if std_c < 0.05 else "不稳定"
        print(f"  {feat}: mean={mean_c:+.4f}, std={std_c:.4f} [{stability}]")
    return all_explanations

lime_stability_check(explainer, model, X_test[42])

import shap

# 使用 TreeSHAP（专为树模型优化，精确且快速）
model = RandomForestClassifier(n_estimators=200, random_state=42)
model.fit(X_train, y_train)

explainer = shap.TreeExplainer(model)
shap_values = explainer.shap_values(X_test)

# 全局摘要图（所有样本所有特征）
shap.summary_plot(shap_values, X_test, feature_names=feature_names)

# 单样本 Force Plot（瀑布图）
shap.force_plot(
    explainer.expected_value[1],
    shap_values[1][0],
    X_test[0],
    feature_names=feature_names
)

# 依赖图：特征值 vs SHAP 值
shap.dependence_plot("income", shap_values[1], X_test,
                     feature_names=feature_names)

# KernelSHAP：适用于任何模型的 SHAP 近似
# 当没有专用 explainer 时使用（如 SVM、神经网络）

import shap

# 创建 KernelSHAP 解释器
explainer = shap.KernelExplainer(model.predict_proba, X_train[:100])
# background data 使用训练集子集

# 计算 SHAP 值
shap_values = explainer.shap_values(X_test[:5], nsamples=1000)

# 解释单个样本
idx = 0
shap.force_plot(
    explainer.expected_value[1],
    shap_values[1][idx],
    X_test[idx],
    feature_names=feature_names,
    matplotlib=True
)

# shap_values[1] 对应正类的 SHAP 值矩阵
print(f"SHAP values shape: {shap_values[1].shape}")

from sklearn.inspection import PartialDependenceDisplay
import matplotlib.pyplot as plt

# PDP 图：展示单个特征的边际效应
fig, ax = plt.subplots(figsize=(8, 5))
PartialDependenceDisplay.from_estimator(
    model, X_test, features=["income"],
    kind="average",
    ax=ax,
    feature_names=feature_names
)
ax.set_title("Partial Dependence: Income")
plt.show()

# PDP + ICE 组合图
fig, ax = plt.subplots(figsize=(8, 5))
PartialDependenceDisplay.from_estimator(
    model, X_test, features=["income"],
    kind="both",
    ax=ax,
    feature_names=feature_names,
    subsample=50
)
plt.show()

# 二维 PDP：展示两个特征的交互效应
from sklearn.inspection import PartialDependenceDisplay

features = [("income", "debt_ratio"), ("age", "credit_history")]
fig, ax = plt.subplots(1, 2, figsize=(14, 5))

for i, feat_pair in enumerate(features):
    PartialDependenceDisplay.from_estimator(
        model, X_test,
        features=feat_pair,
        kind="average",
        ax=ax[i],
        feature_names=feature_names,
        grid_resolution=50
    )
    ax[i].set_title(f"2D PDP: {feat_pair[0]} x {feat_pair[1]}")

plt.tight_layout()
plt.show()

import shap
import lime.lime_tabular
import xgboost as xgb

# 1. 训练 XGBoost 模型
model = xgb.XGBClassifier(n_estimators=200, max_depth=5, random_state=42)
model.fit(X_train, y_train)

# 2. SHAP 全局解释
explainer = shap.TreeExplainer(model)
shap_values = explainer.shap_values(X_test)

# Summary plot（蜜蜂图）
shap.summary_plot(shap_values, X_test, feature_names=feature_names,
                   plot_type="dot", show=False)

# SHAP 特征重要性排序
shap.summary_plot(shap_values, X_test, feature_names=feature_names,
                   plot_type="bar", show=False)

# 3. LIME 局部解释
lime_explainer = lime.lime_tabular.LimeTabularExplainer(
    X_train, feature_names=feature_names,
    class_names=["拒绝", "批准"], mode="classification"
)
rejected_idx = np.where(model.predict(X_test) == 0)[0][0]
exp = lime_explainer.explain_instance(
    X_test[rejected_idx], model.predict_proba,
    num_features=6, num_samples=5000
)
exp.show_in_notebook()

# SHAP 高级可视化：决策图和依赖图

# Decision plot：展示从基础值到预测的完整路径
shap.decision_plot(
    explainer.expected_value[1],
    shap_values[1][:10],
    X_test[:10],
    feature_names=feature_names,
    link="logit"
)

# Dependence plot：收入对 SHAP 值的影响
shap.dependence_plot(
    "income", shap_values[1], X_test,
    feature_names=feature_names,
    interaction_index="debt_ratio",
    show=False
)

# 生成解释报告
def generate_explanation_report(model, X_sample, idx):
    shap_val = shap_values[1][idx]
    base_val = explainer.expected_value[1]
    pred = model.predict_proba([X_sample[idx]])[0][1]
    report = f"预测概率: {pred:.2%}\n"
    report += f"基础值: {base_val:.4f}\n"
    report += "--- 特征贡献 (SHAP) ---\n"
    for feat, val in sorted(zip(feature_names, shap_val), key=lambda x: -abs(x[1]))[:5]:
        direction = "+" if val > 0 else "-"
        report += f"  {direction} {feat}: {val:+.4f}\n"
    return report

print(generate_explanation_report(model, X_test, 0))