Explainability (SHAP)

from sklearn.ensemble import RandomForestClassifier, GradientBoostingRegressor
from sklearn.datasets import load_breast_cancer, fetch_california_housing
from sklearn.model_selection import train_test_split
from sklearn.inspection import permutation_importance, partial_dependence, PartialDependenceDisplay
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

# 1) Train a model
data = load_breast_cancer(as_frame=True)
X, y = data.data, data.target
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, stratify=y, random_state=42)

model = RandomForestClassifier(n_estimators=300, random_state=42, n_jobs=-1).fit(X_train, y_train)
print('test acc:', model.score(X_test, y_test))

# 2) Built-in feature importance (tree models only)
imp = pd.Series(model.feature_importances_, index=X.columns).sort_values(ascending=False)
imp.head(10).plot.barh(); plt.tight_layout()
# Caveat: trees prefer high-cardinality numeric features; can be misleading. Verify with permutation importance.

# 3) Permutation importance — model-agnostic
perm = permutation_importance(model, X_test, y_test, n_repeats=10, random_state=42, n_jobs=-1)
imp_perm = pd.Series(perm.importances_mean, index=X.columns).sort_values(ascending=False)
print(imp_perm.head(10))
# Permutation importance asks: 'how much does shuffling this feature hurt the score?'
# Slow but trustworthy. Repeat n_repeats=10+ for stability.

# 4) SHAP — Shapley values, per-row + per-feature
# pip install shap
import shap
explainer = shap.TreeExplainer(model)
shap_values = explainer.shap_values(X_test)
# For multi-class returns list[class_idx] of (n_samples, n_features)
if isinstance(shap_values, list): shap_values = shap_values[1]    # positive class

# Global summary
shap.summary_plot(shap_values, X_test, max_display=15)
# Bar version
shap.summary_plot(shap_values, X_test, plot_type='bar')

# Single prediction
shap.force_plot(explainer.expected_value[1] if isinstance(explainer.expected_value, np.ndarray) else explainer.expected_value,
                                shap_values[0], X_test.iloc[0], matplotlib=True)

# Decision plot
shap.decision_plot(explainer.expected_value, shap_values[:5], X_test.iloc[:5])

# 5) SHAP for non-tree models
# Use KernelExplainer (slow) or LinearExplainer / DeepExplainer per model type
expl2 = shap.KernelExplainer(model.predict_proba, shap.sample(X_train, 100))
sv = expl2.shap_values(X_test.iloc[:50], nsamples=200)

# 6) LIME — local surrogate model around one prediction
# pip install lime
from lime.lime_tabular import LimeTabularExplainer

expl = LimeTabularExplainer(
    X_train.values,
    feature_names=list(X.columns),
    class_names=['benign', 'malignant'],
    discretize_continuous=True,
)
row = X_test.iloc[0].values
exp = expl.explain_instance(row, model.predict_proba, num_features=10)
exp.show_in_notebook(show_table=True)
# Trains a small linear model on PERTURBATIONS around 'row' to explain that one prediction.

# 7) Partial Dependence Plots — average effect of a feature
features = [(0,), (1,), (0, 1)]                  # singles + interaction
PartialDependenceDisplay.from_estimator(model, X_test, features, kind='average')
plt.tight_layout()

# Individual Conditional Expectation
PartialDependenceDisplay.from_estimator(model, X_test, [0], kind='both')
# kind='both' overlays the average PDP with per-row ICE lines.

# 8) Counterfactual explanations — 'what would change the prediction?'
# pip install dice-ml
import dice_ml
d = dice_ml.Data(dataframe=pd.concat([X_train, y_train], axis=1), continuous_features=list(X.columns), outcome_name='target')
m = dice_ml.Model(model=model, backend='sklearn')
exp = dice_ml.Dice(d, m, method='random')
cf = exp.generate_counterfactuals(X_test.iloc[:1], total_CFs=3, desired_class='opposite')
cf.visualize_as_dataframe()

# 'What's the smallest change to this patient's features that flips the prediction?'
# Great for actionable explanations: 'reduce cholesterol from 240 to 200'.

# 9) Global vs local explanations
# • Global:   feature importance, summary SHAP, PDP — overall model behaviour
# • Local:    force plot, LIME, counterfactual — single prediction
# Stakeholders almost always want both: 'on average, X drives risk' + 'this patient's risk is driven by Y, Z'.

# 10) Reporting templates
# Build a one-page report per prediction:
#   • Predicted probability
#   • Top 5 contributing features (SHAP)
#   • Confidence (model calibration; prediction interval if regression)
#   • Counterfactual: smallest change that flips prediction
#   • Comparison to nearest neighbour in training

# 11) Audit-ready explanations
# • Save: model version + data snapshot + SHAP/LIME outputs + decision + reviewer
# • Store as Parquet/PDF with hashes for tamper evidence
# • Pair with a human review gate for high-stakes decisions

# 12) Calibration matters
# An explanation only helps if the underlying probabilities are calibrated.
from sklearn.calibration import CalibratedClassifierCV
cal = CalibratedClassifierCV(estimator=GradientBoostingClassifier(), cv=3, method='isotonic').fit(X_train, y_train)
# Reliability diagrams: sklearn.calibration.calibration_curve

# 13) Common bugs
# • Built-in importance on trees fooled by high-cardinality features → use permutation importance
# • SHAP TreeExplainer used on a non-tree model — wrong values; pick the right explainer
# • Explaining the wrong class — check shap_values[positive_idx]
# • LIME perturbations crossing feature constraints (negative ages) — provide categorical_features + bounds
# • PDP on a feature CORRELATED with another → misleading averages; use ALE plots instead
# • Sharing per-row explanations that leak PII — anonymise
# • Calling 'feature importance' a 'cause' — these are correlations within the model
# • Stakeholders demand 100% explainability — communicate trade-offs honestly
# • Slow KernelExplainer in prod → precompute, cache, or switch to a simpler surrogate

import shap
explainer = shap.TreeExplainer(clf)
values = explainer.shap_values(Xte[:100])
shap.summary_plot(values, Xte[:100], feature_names=feature_names)