🎯 What You'll Learn
- Explain why interpretability matters: regulation (GDPR, ECOA), debugging, stakeholder trust, and fairness auditing
- Distinguish global interpretability (how the model behaves overall) from local interpretability (why this one prediction)
- Distinguish intrinsically interpretable models from post-hoc explanation of black-box models
- Understand Shapley values and the additive feature attribution property that SHAP is built on
- Choose the right SHAP explainer — TreeExplainer, KernelExplainer, or DeepExplainer — for a given model type
- Read and produce SHAP summary plots, force plots, dependence plots, and waterfall plots
- Use LIME to build local surrogate explanations and articulate the SHAP-vs-LIME tradeoffs in speed, theory, and stability
- Explain an individual credit-risk or fraud-model decision to a non-technical stakeholder and detect proxy-discrimination in a model
- Define and compute formal fairness metrics — Demographic Parity, Equalized Odds, Equal Opportunity — and understand why satisfying all of them at once is often mathematically impossible
You've spent Lessons 11–24 and 36–60 learning to build models that are accurate. This lesson is about a different, equally important question: can you explain why a model did what it did? A random forest or gradient-boosted tree with 500 estimators, or a transformer with hundreds of millions of parameters, is a function so complex that no human can trace through it by hand. Interpretability tools don't change the model — they build a second, much simpler model on top of it whose entire job is to answer "why?" in a way a human can audit. SHAP answers this with game theory: it treats each feature as a "player" contributing to a "payout" (the prediction) and fairly splits credit among them. LIME answers it more pragmatically: zoom in close enough to one prediction and even the wiggliest black box looks locally like a straight line — so just fit that line. Both turn opaque models into accountable ones.
1 Why Interpretability Matters
Through Lessons 11–24 you optimized for predictive performance — accuracy, AUC, F1, RMSE. In production, performance is necessary but not sufficient. A model that is 94% accurate but cannot explain a single one of its decisions is, in many regulated or high-stakes settings, not deployable at all. Four forces push interpretability from "nice to have" to "required":
1. Regulatory requirements
The EU's GDPR (Articles 13–15 and 22) grants individuals a qualified "right to explanation" for decisions made by automated systems that produce legal or similarly significant effects — loan denials, insurance pricing, automated hiring screens. In the US, the Equal Credit Opportunity Act (ECOA) and Regulation B require lenders to give applicants specific, accurate reasons for credit denial — "the model said no" is not an acceptable adverse action notice. Healthcare (FDA software-as-a-medical-device guidance) and insurance underwriting carry similar obligations. If you can't explain a prediction, you may not be legally allowed to act on it.
2. Debugging models
A model can achieve excellent validation accuracy while learning the wrong thing entirely — the classic example is an image classifier that "detects huskies" by learning to detect snow in the background. Interpretability tools surface this kind of shortcut learning before it ships. If a fraud model's top feature turns out to be an internal database timestamp rather than any genuine behavioral signal, that's a data leakage bug, and you will only catch it by asking the model "why?" on individual predictions, not by staring at an aggregate AUC score.
3. Stakeholder trust
A doctor will not act on a sepsis-risk score they cannot interrogate. A loan officer will not override their own judgment for a black-box number with no rationale. Interpretability is frequently the deciding factor in whether a model that performs well in offline evaluation is ever actually adopted in the field — trust is a prerequisite for the model creating any value at all.
4. Bias and fairness auditing
A model can be statistically accurate and still systematically disadvantage a protected group — for instance by relying heavily on a feature like zip code that correlates strongly with race (a classic proxy variable). Local and global explanation tools let you check, feature by feature, whether the model's reasoning is defensible, not just whether its aggregate metrics look acceptable. We'll build exactly this kind of audit in the Real-World Spotlight.
It's tempting to think "simpler model = more interpretable, but less accurate" and stop there. In practice, for tabular data, gradient-boosted trees (Lesson 24) are often both more accurate and, with SHAP's TreeExplainer, fully and exactly explainable — there is no accuracy tax to pay for interpretability on this data type. The real tradeoff is between models cheap to explain with exact methods (trees) and models that require approximate, sampling-based explanation (deep nets, ensembles of heterogeneous models, anything behind an API you don't control).
2 Global vs Local, Intrinsic vs Post-Hoc
Two independent axes organize the whole field of interpretability. Knowing where a technique sits on each tells you immediately what question it can and can't answer.
Global vs Local
- Global interpretability describes how the model behaves on average, across the whole input distribution: "Across all applicants, debt-to-income ratio is the single most influential feature, and higher values push predictions toward default." This is what you need for model documentation, regulatory model risk reviews, and sanity-checking that the model learned sensible relationships overall.
- Local interpretability describes why the model produced this specific output for this specific input: "Applicant #4471 was denied primarily because of a 41% debt-to-income ratio and two recent missed payments, partially offset by 12 years of credit history." This is what you need for an adverse-action notice, a single clinical decision, or debugging one wrong prediction a user reported.
Critically, a technique that's good at one is not automatically good at the other. Permutation feature importance (which you may have seen used for feature selection) is inherently global — it tells you nothing about any individual row. SHAP and LIME, as we'll see, are fundamentally local methods that can be aggregated upward into global summaries — which is exactly why they have become the dominant tools in this space.
Intrinsic vs Post-Hoc
- Intrinsically interpretable models are simple enough that the model itself is the explanation. A fitted linear regression's coefficients (Lesson 11) tell you exactly how much each unit of a feature moves the prediction, holding others constant. A shallow decision tree (Lesson 13) can be read top-to-bottom as a sequence of human-readable rules. Their feature-importance is built in, not bolted on, and requires no extra machinery.
- Post-hoc explanation is required for everything else: random forests, XGBoost/LightGBM ensembles, and any neural network. These models have no transparent internal structure a human can read directly, so a second, separate explanation technique — SHAP or LIME — is run after training to approximate or exactly decompose what the trained black box is doing.
The two axes that organize every interpretability technique in this lesson. SHAP and LIME both live natively in the bottom-right (local, post-hoc) cell — but SHAP's local values, uniquely among post-hoc methods, aggregate cleanly upward into the top-right (global, post-hoc) cell, which is exactly why summary plots exist and why SHAP has become the default choice.
import numpy as np
import pandas as pd
from sklearn.linear_model import LogisticRegression
from sklearn.tree import DecisionTreeClassifier, export_text
from sklearn.ensemble import RandomForestClassifier
from sklearn.datasets import make_classification
X, y = make_classification(n_samples=2000, n_features=6, n_informative=4,
random_state=42)
feature_names = ['debt_to_income', 'credit_score', 'income', 'years_employed',
'num_open_accounts', 'recent_inquiries']
X = pd.DataFrame(X, columns=feature_names)
# ── Intrinsically interpretable: logistic regression coefficients ──
logreg = LogisticRegression(max_iter=1000).fit(X, y)
print("Logistic regression — the model IS the explanation:")
for name, coef in sorted(zip(feature_names, logreg.coef_[0]),
key=lambda t: abs(t[1]), reverse=True):
direction = "increases" if coef > 0 else "decreases"
print(f" {name:20s} coef={coef:+.3f} ({direction} predicted default log-odds)")
# ── Intrinsically interpretable: shallow decision tree as readable rules ──
tree = DecisionTreeClassifier(max_depth=3, random_state=42).fit(X, y)
print("\nDecision tree — readable as if/else rules:")
print(export_text(tree, feature_names=feature_names, max_depth=2))
# ── Black box requiring post-hoc explanation ──
rf = RandomForestClassifier(n_estimators=500, max_depth=12, random_state=42).fit(X, y)
print(f"\nRandom forest: {rf.n_estimators} trees, no single readable rule set.")
print("Built-in .feature_importances_ is GLOBAL ONLY and biased toward high-cardinality")
print("features — it cannot tell you why any one applicant was scored as it was.")
print("This is exactly the gap SHAP and LIME are built to fill →")You met feature_importances_ on random forests and boosted trees back in Lessons 13 and 24. That's a global, gain- or split-count-based statistic with two well-documented weaknesses: it's biased toward high-cardinality and continuous features, and it cannot tell you anything about an individual prediction. SHAP values, covered next, are computed per-row, are consistent (a feature that contributes more can never get a lower attribution), and sum exactly to the model's output — properties the built-in importances do not have.
3 SHAP and the Shapley Value
SHAP (SHapley Additive exPlanations, Lundberg & Lee, 2017) is built on the Shapley value, a concept from 1950s cooperative game theory. The original question Lloyd Shapley answered: if a group of players cooperate to produce some total payout, how do you fairly split the credit among them, accounting for the fact that players contribute different amounts depending on which other players are already "in the room"?
From Game Theory to Feature Attribution
SHAP reframes a model prediction as a cooperative game: the "players" are the input features, and the "payout" is the difference between the model's prediction for this instance and the model's average prediction over the training data (the baseline). The Shapley value for feature i is its average marginal contribution to the prediction, computed across every possible order in which features could be "added" to the model:
# Conceptual definition of the Shapley value for feature i (not how SHAP
# actually computes it for trees — see TreeExplainer below for the real algorithm)
from itertools import permutations
import numpy as np
def shapley_value_conceptual(model_predict, x, baseline, feature_idx, all_features):
"""
Illustrates the DEFINITION of a Shapley value: average the marginal
contribution of `feature_idx` across all orderings of features being
'switched on' from the baseline to the real value. Exponential in the
number of features — this is why real implementations use smarter
algorithms (TreeExplainer's polynomial-time exact algorithm, or
KernelExplainer's weighted sampling).
"""
contributions = []
other_features = [f for f in all_features if f != feature_idx]
for perm in permutations(other_features):
# Build up the "coalition" in this ordering, then add feature_idx
coalition_before = np.array(baseline, dtype=float)
for f in perm:
coalition_before[f] = x[f]
pred_before = model_predict(coalition_before.reshape(1, -1))[0]
coalition_after = coalition_before.copy()
coalition_after[feature_idx] = x[feature_idx]
pred_after = model_predict(coalition_after.reshape(1, -1))[0]
contributions.append(pred_after - pred_before)
return np.mean(contributions)
# In practice: 10 features -> 9! = 362,880 orderings PER feature. SHAP's
# TreeExplainer exploits tree structure to compute the exact answer in
# polynomial time instead — that's the breakthrough that made SHAP practical.The Additive Feature Attribution Property
The reason SHAP values are interpretable — and not just "some numbers" — is that they satisfy a guarantee called additive feature attribution: the SHAP values for a single prediction always sum exactly to the gap between that prediction and the baseline (expected) prediction:
# f(x) = base_value + sum of SHAP values for every feature
#
# prediction = E[f(X)] + φ_1 + φ_2 + ... + φ_n
# (avg model (debt_to_ (credit_
# output) income score
# contrib.) contrib.)
#
# Example for one applicant:
# 0.81 (predicted default prob.) = 0.30 (baseline) + 0.28 (debt_to_income)
# + (-0.05) (credit_score) + 0.31 (recent_inquiries) + ...
print("SHAP guarantee: base_value + sum(shap_values) == model.predict(x)")
print("This is what makes a SHAP force plot or waterfall plot additive and exact —")
print("every bar you see literally adds up to the final prediction.")SHAP additionally guarantees two other properties that ad-hoc attribution schemes don't: consistency (if a model changes so that a feature's marginal contribution increases or stays the same regardless of other features, its attributed importance cannot decrease) and missingness (a feature that has no effect gets a SHAP value of exactly zero). These are the theoretical guarantees people mean when they say "SHAP values" instead of just "feature importance."
The conceptual definition above requires evaluating the model on 2n feature subsets (every possible coalition) — for 20 features that's over a million model evaluations, for 100 features it's astronomically infeasible. This is precisely why SHAP isn't one algorithm but a family of explainers, each exploiting structure in a specific model class to make the computation tractable. Picking the right explainer for your model type, covered next, is the single most important practical decision when using SHAP.
4 Choosing an Explainer: Tree, Kernel, and Deep
The shap library ships several explainer classes. They all return SHAP values with the same additive guarantees, but differ enormously in speed and in whether the result is exact or approximate.
| Explainer | Model types | Exact? | Speed |
|---|---|---|---|
TreeExplainer |
XGBoost, LightGBM, CatBoost, sklearn trees/forests | Yes — polynomial-time exact algorithm | Very fast (milliseconds–seconds) |
DeepExplainer |
Keras / PyTorch neural networks | Approximate (DeepLIFT-style backprop) | Fast (one backward pass per sample) |
GradientExplainer |
Differentiable models (NNs) | Approximate (expected gradients) | Fast |
KernelExplainer |
Any model — fully model-agnostic | Approximate (weighted linear regression) | Slow (thousands of model calls per row) |
LinearExplainer |
Linear / logistic regression | Yes — closed form | Instant |
A practical explainer-selection flow: check the model class first, since that determines whether SHAP can be exact and fast (TreeExplainer), approximate but still efficient (DeepExplainer), or fully general but slow (KernelExplainer). Always try to fall into the TreeExplainer or LinearExplainer branch before reaching for the model-agnostic fallback.
TreeExplainer: Fast and Exact for Tree Models
Because XGBoost, LightGBM, and random forests are by far the most common production models for tabular data (Lessons 13, 18, 24), TreeExplainer is the workhorse of practical SHAP usage. It uses the tree structure itself — which features are split on, and where — to compute exact Shapley values without any sampling, in time proportional to the number of trees times leaf depth, not 2n.
import shap
import xgboost as xgb
import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split
# Simulated credit risk dataset
np.random.seed(7)
n = 5000
df = pd.DataFrame({
'debt_to_income': np.random.beta(2, 5, n) * 100,
'credit_score': np.random.normal(680, 80, n).clip(300, 850),
'income': np.random.lognormal(10.8, 0.5, n),
'years_employed': np.random.exponential(5, n).clip(0, 40),
'num_open_accounts': np.random.poisson(4, n),
'recent_inquiries': np.random.poisson(1.2, n),
'months_since_late_payment': np.random.exponential(24, n).clip(0, 120),
})
# Default risk increases with DTI and inquiries, decreases with credit score
logit = (0.05 * df.debt_to_income - 0.015 * (df.credit_score - 680)
+ 0.4 * df.recent_inquiries - 0.01 * df.months_since_late_payment - 1.0)
prob_default = 1 / (1 + np.exp(-logit))
y = np.random.binomial(1, prob_default)
X_train, X_test, y_train, y_test = train_test_split(df, y, test_size=0.2, random_state=42)
model = xgb.XGBClassifier(n_estimators=300, max_depth=4, learning_rate=0.05,
eval_metric='logloss', random_state=42)
model.fit(X_train, y_train)
# TreeExplainer: exact, fast SHAP values for the whole test set
explainer = shap.TreeExplainer(model)
shap_values = explainer(X_test) # shap.Explanation object
print(f"shap_values.values.shape: {shap_values.values.shape}") # (n_test, n_features)
print(f"Base value (expected model output): {shap_values.base_values[0]:.4f}")
print(f"\nVerifying additivity for row 0:")
row0_sum = shap_values.base_values[0] + shap_values.values[0].sum()
raw_pred = model.predict(X_test.iloc[[0]], output_margin=True)[0]
print(f" base_value + sum(shap_values) = {row0_sum:.4f}")
print(f" model raw margin output = {raw_pred:.4f} (should match)")KernelExplainer: Model-Agnostic, but Slow
KernelExplainer makes no assumption about model internals — it treats the model purely as a black-box function it can call. It estimates Shapley values by sampling many random feature coalitions, masking out the "absent" features (typically by replacing them with values from a background dataset), and fitting a weighted linear regression to the resulting predictions. This generality is its strength and its weakness: it works on literally any predict function (an ensemble of heterogeneous models, a model behind a REST API, a scikit-learn pipeline with custom preprocessing) but needs hundreds to thousands of model evaluations per explained row.
import shap
from sklearn.svm import SVC
from sklearn.preprocessing import StandardScaler
from sklearn.pipeline import Pipeline
# A model type TreeExplainer cannot handle: an RBF-kernel SVM
svm_pipeline = Pipeline([
('scale', StandardScaler()),
('svm', SVC(kernel='rbf', probability=True, random_state=42)),
]).fit(X_train, y_train)
# KernelExplainer needs a small representative "background" sample —
# this stands in for the missing features during coalition sampling.
background = shap.sample(X_train, 100, random_state=42)
kernel_explainer = shap.KernelExplainer(
svm_pipeline.predict_proba, background
)
# Explain just 20 test rows — KernelExplainer is too slow for the full set
sample_to_explain = X_test.iloc[:20]
kernel_shap_values = kernel_explainer.shap_values(sample_to_explain, nsamples=200)
print(f"KernelExplainer evaluated the model ~{200 * 20:,} times for 20 rows.")
print(f"TreeExplainer would have computed all {len(X_test):,} rows exactly,")
print(f"with zero sampling noise, in a fraction of the time.")DeepExplainer: For Neural Networks
For the deep learning models from Lessons 37–56, DeepExplainer implements an approximation based on DeepLIFT, propagating attributions backward through the network's layers rather than sampling coalitions — much faster than KernelExplainer for high-dimensional inputs like images or embeddings, at the cost of being an approximation rather than an exact Shapley value.
import shap
import torch
import torch.nn as nn
# A small tabular neural net (PyTorch) — same credit risk task
class CreditNet(nn.Module):
def __init__(self, n_features):
super().__init__()
self.net = nn.Sequential(
nn.Linear(n_features, 32), nn.ReLU(),
nn.Linear(32, 16), nn.ReLU(),
nn.Linear(16, 1),
)
def forward(self, x):
return self.net(x)
torch_model = CreditNet(n_features=X_train.shape[1])
torch_model.eval()
X_train_t = torch.tensor(X_train.values, dtype=torch.float32)
X_test_t = torch.tensor(X_test.values, dtype=torch.float32)
# Background: a sample of training data to integrate over
background_t = X_train_t[torch.randperm(len(X_train_t))[:100]]
deep_explainer = shap.DeepExplainer(torch_model, background_t)
deep_shap_values = deep_explainer.shap_values(X_test_t[:50])
print(f"DeepExplainer SHAP values shape: {deep_shap_values.shape}")
print("Used for tabular NNs, CNN image classifiers, and embedding-based models —")
print("trades exactness for tractability on models with millions of parameters.")Both KernelExplainer and DeepExplainer require a background dataset that defines the "default" or "absent" value of each feature. A SHAP value is always relative to this baseline — change the background sample and the attributed values shift, sometimes substantially. The convention is to use a representative random sample (50–200 rows) of the training distribution, not a single row of zeros or means, which can put the baseline in an unrealistic part of feature space and distort attributions.
Quick Check
5 Reading SHAP Plots: Summary, Force, Dependence, Waterfall
Raw SHAP values are a matrix of numbers; the SHAP library's real value is the standard set of plots built on top of that matrix, each answering a different question.
Summary Plot — Global Importance, with Local Detail
The summary plot is usually the first thing you look at after fitting an explainer. Each row is a feature, sorted by overall importance (mean absolute SHAP value); each point is one row of test data, colored by that feature's raw value (red = high, blue = low) and positioned left-right by its SHAP value (impact on the prediction). It compresses a global view and thousands of local explanations into a single chart.
import shap
import matplotlib.pyplot as plt
shap_values = explainer(X_test) # from Section 4, TreeExplainer on XGBoost
shap.summary_plot(shap_values, X_test, show=False)
plt.tight_layout()
plt.savefig("shap_summary.png", dpi=150)
plt.close()
# A bar-style global summary (mean |SHAP value| per feature) for executives
shap.summary_plot(shap_values, X_test, plot_type="bar", show=False)
plt.tight_layout()
plt.savefig("shap_summary_bar.png", dpi=150)
print("Summary plot reading guide:")
print(" - Vertical position: feature, ranked by overall impact")
print(" - Horizontal position: SHAP value (pushes prediction up or down)")
print(" - Color: the feature's actual value for that row (red=high, blue=low)")
print(" e.g. 'debt_to_income' shows red dots clustered on the right ->")
print(" high DTI consistently pushes default probability UP")Illustrative SHAP summary (beeswarm) plot for the credit-risk model — not computed from a real trained model. Feature order and each row's spread are hand-tuned to match the mean |SHAP value| ranking reported by generate_explainability_report() in Section 7 (recent_inquiries ≈0.41, debt_to_income ≈0.39, credit_score ≈0.30, months_since_late_payment ≈0.18, years_employed ≈0.09). Each dot is one applicant, colored by that applicant's raw feature value (blue = low, amber = high). Notice credit_score and years_employed run the opposite direction from recent_inquiries and debt_to_income: high credit_score and long tenure push risk down (negative SHAP), while frequent inquiries and high debt push it up (positive SHAP) — exactly the signs from the logistic regression coefficients in Section 2.
Force Plot — One Prediction, Visually
A force plot shows a single row's explanation as a tug-of-war: features pushing the prediction higher than the base value in red, features pushing it lower in blue, with bar width proportional to magnitude. It's the most stakeholder-friendly SHAP visualization — many production dashboards embed a force plot directly next to a model's decision.
import shap
shap.initjs() # enables the interactive force plot in a Jupyter notebook
row_idx = 0
shap.force_plot(
explainer.expected_value,
shap_values.values[row_idx],
X_test.iloc[row_idx],
matplotlib=True, # static image instead of interactive JS widget
show=False,
).savefig("force_plot_applicant_0.png", dpi=150, bbox_inches='tight')
print(f"Applicant {row_idx}: base rate = {explainer.expected_value:.3f}")
top_contributors = sorted(
zip(X_test.columns, shap_values.values[row_idx]),
key=lambda t: abs(t[1]), reverse=True
)[:3]
for feat, val in top_contributors:
direction = "increased" if val > 0 else "decreased"
print(f" {feat:25s} {direction} risk by {abs(val):.3f}")Dependence Plot — How One Feature Drives Predictions
A dependence plot scatters one feature's raw value (x-axis) against its SHAP value (y-axis) across every row, revealing the shape of the relationship the model learned — linear, threshold-like, or U-shaped — and automatically colors points by a second feature SHAP detects the strongest interaction with.
import shap
import matplotlib.pyplot as plt
shap.dependence_plot(
"debt_to_income", shap_values.values, X_test,
interaction_index="auto", # auto-picks the feature with strongest interaction
show=False,
)
plt.tight_layout()
plt.savefig("dependence_dti.png", dpi=150)
print("Dependence plot for debt_to_income reveals:")
print(" - Roughly monotonic: higher DTI -> higher SHAP value (more predicted risk)")
print(" - A visible 'knee' around DTI=40%: the model learned a soft threshold")
print(" effect, not a smooth linear one -- info a linear model could never show")Waterfall Plot — The Modern Single-Prediction View
The waterfall plot is the newer, typically preferred alternative to the force plot for explaining one row: it shows the same additive decomposition as a vertical staircase from the base value up (or down) to the final prediction, ordered by magnitude, which many people find easier to read than the force plot's horizontal bars.
import shap
import matplotlib.pyplot as plt
shap_explanation = explainer(X_test)
shap.plots.waterfall(shap_explanation[0], show=False)
plt.tight_layout()
plt.savefig("waterfall_applicant_0.png", dpi=150)
plt.close()
print("Waterfall plot for applicant 0:")
print(f" E[f(X)] = {shap_explanation.base_values[0]:.3f} (start)")
for feat, val in sorted(zip(X_test.columns, shap_explanation.values[0]),
key=lambda t: abs(t[1]), reverse=True):
sign = "+" if val >= 0 else "-"
print(f" {sign} {feat:25s} {val:+.3f}")
print(f" f(x) = {shap_explanation.base_values[0] + shap_explanation.values[0].sum():.3f} (final prediction)")A full waterfall built from the applicant example in Section 3's additive-property callout: 0.81 = 0.30 (baseline) + 0.28 (debt_to_income) − 0.05 (credit_score) + 0.31 (recent_inquiries) + …. Two small remaining contributions (years_employed, income) are added here, sized only to close the gap to the stated final 0.81 — numbers chosen to match the lesson's own text exactly, not output from a real trained model. Red bars push the predicted probability up, blue bars pull it down, and the violet bars are the running totals (base rate and final prediction).
In practice: ship the bar summary plot to executives and model-risk committees who want "what matters most, overall" in one glance; use the full summary (beeswarm) plot in model documentation and validation reports because it shows both importance and direction; embed a waterfall or force plot directly in customer-facing or analyst-facing tools that need to justify one specific decision; reach for dependence plots only during model debugging and feature engineering, when you need to understand a relationship's shape, not communicate it externally.
6 LIME: Local Surrogate Models
LIME (Local Interpretable Model-agnostic Explanations, Ribeiro, Singh & Guestrin, 2016) predates SHAP and takes a more pragmatic, less theoretically grounded approach to the same problem: explain one prediction by approximating the model's behavior in the immediate neighborhood of that one input.
How LIME Works
- Perturb: generate many synthetic samples by randomly perturbing the instance being explained (for tabular data: sample nearby feature values; for text: randomly remove words; for images: turn superpixels on/off)
- Predict: run the black-box model on every perturbed sample to get its predictions
- Weight: weight each perturbed sample by its proximity to the original instance — closer perturbations matter more
- Fit a surrogate: fit a simple, inherently interpretable model (almost always weighted linear regression, sometimes a shallow tree) to this weighted, perturbed dataset
- Explain: the surrogate model's coefficients, which are only ever claimed to be valid in this local neighborhood, become the explanation for the original prediction
import lime
import lime.lime_tabular
import numpy as np
# LIME needs to know the training distribution to generate realistic perturbations
explainer_lime = lime.lime_tabular.LimeTabularExplainer(
training_data=X_train.values,
feature_names=list(X_train.columns),
class_names=['repaid', 'default'],
mode='classification',
discretize_continuous=True, # bins continuous features for readability
random_state=42,
)
row_idx = 0
instance = X_test.iloc[row_idx].values
explanation = explainer_lime.explain_instance(
data_row=instance,
predict_fn=model.predict_proba, # any function: array -> class probabilities
num_features=7,
num_samples=5000, # perturbed samples used to fit the surrogate
)
print(f"LIME explanation for applicant {row_idx} (local surrogate model):")
for feature_rule, weight in explanation.as_list():
direction = "increases" if weight > 0 else "decreases"
print(f" {feature_rule:35s} {direction} default probability by {abs(weight):.3f}")
# Surrogate model's own fit quality -- LOW R^2 means the local linear
# approximation is unreliable in this region of feature space
print(f"\nLocal surrogate R^2 (fit quality): {explanation.score:.3f}")SHAP vs LIME: The Practical Tradeoffs
| Dimension | SHAP | LIME |
|---|---|---|
| Theoretical guarantees | Strong — uniquely satisfies additivity, consistency, missingness (Shapley axioms) | None — explanation quality depends on the local fit being good |
| Speed (tree models) | Very fast and exact with TreeExplainer | Moderate — thousands of perturbed predictions per explanation, every time |
| Stability | Deterministic for TreeExplainer; same input always gives same output | Stochastic — re-running on the same row can give noticeably different weights |
| Model coverage | Needs a matching explainer type for full speed; KernelExplainer covers the rest, slowly | Always model-agnostic, same algorithm regardless of model type |
| Global aggregation | Natural — sum/average local values into summary and dependence plots | Awkward — each explanation is an independent local fit, not designed to aggregate |
| Best fit | Default choice for tabular tree models and when guarantees matter (regulatory) | Quick sanity checks, text/image explanations, or when SHAP is too slow to run at all |
Because LIME refits a brand-new linear surrogate from freshly sampled random perturbations every time you call explain_instance, running it twice on the exact same row with the exact same model can surface different "top reasons" — sometimes flipping the sign of a feature's contribution. This is a serious problem for regulatory use cases where an adverse-action notice must be reproducible. Mitigate it by fixing the random seed, increasing num_samples substantially (5,000+), and treating the local surrogate's R² (explanation.score) as a confidence check — a low R² means the linear surrogate isn't a trustworthy local approximation at all, regardless of what coefficients it spits out.
7 Putting It Together: An Explainability Workflow
In a real MLOps pipeline (which you'll formalize with MLflow in the next lesson), interpretability isn't a one-off notebook cell — it's a repeatable step that runs alongside evaluation, every time a model is trained or retrained.
import shap
import pandas as pd
import numpy as np
import json
def generate_explainability_report(model, X_test, y_test, feature_names,
top_k_global=5, sample_local=3):
"""
Produces a structured interpretability report suitable for attaching
to a model card or a model-risk-review document.
"""
explainer = shap.TreeExplainer(model)
shap_exp = explainer(X_test)
# Global summary: mean absolute SHAP value per feature
mean_abs_shap = np.abs(shap_exp.values).mean(axis=0)
global_ranking = sorted(
zip(feature_names, mean_abs_shap), key=lambda t: t[1], reverse=True
)
report = {
"global_top_features": [
{"feature": f, "mean_abs_shap": round(float(v), 4)}
for f, v in global_ranking[:top_k_global]
],
"local_examples": [],
}
# A handful of representative local explanations
rng = np.random.RandomState(42)
sample_idx = rng.choice(len(X_test), size=sample_local, replace=False)
for idx in sample_idx:
contributions = sorted(
zip(feature_names, shap_exp.values[idx]),
key=lambda t: abs(t[1]), reverse=True
)[:3]
report["local_examples"].append({
"row_index": int(idx),
"predicted_prob": float(model.predict_proba(X_test.iloc[[idx]])[0, 1]),
"base_value": float(shap_exp.base_values[idx]),
"top_reasons": [
{"feature": f, "shap_value": round(float(v), 4)} for f, v in contributions
],
})
return report
report = generate_explainability_report(
model, X_test, y_test, list(X_test.columns)
)
print(json.dumps(report, indent=2))This report is the kind of artifact you attach to a model card, feed into automated regression tests ("did the top global feature change unexpectedly after retraining?"), and surface in a customer-facing explanation API. Notice it deliberately separates the global section (model-level summary, useful for model risk review) from the local section (per-decision detail, useful for individual adverse-action notices) — exactly the distinction from Section 2.
A model can keep the same accuracy after retraining while quietly shifting which features it relies on — for example because a correlated feature became more predictive due to a data pipeline change. Tracking the global SHAP feature ranking over time, the same way you'd track AUC or RMSE, catches this kind of silent behavioral drift that pure performance metrics miss entirely.
8 Fairness Metrics: Making "Biased" Measurable
Section 1 flagged proxy variables and Section 7's workflow (and the Real-World Spotlight below) audits whether SHAP attributions differ across groups — both are ways of asking "is this model biased?" without ever defining what "biased" means as a number. Fairness metrics formalize the question, and — critically — different formalizations can directly disagree with each other, so knowing which one you're optimizing for is part of the job, not an afterthought.
- Demographic Parity (Statistical Parity): the model should approve (or predict the positive class for) the same proportion of each group, regardless of their true underlying qualification rates. P(ŷ=1 | group=A) ≈ P(ŷ=1 | group=B).
- Equalized Odds: the model should have the same True Positive Rate and the same False Positive Rate across groups — among people who truly qualify, an equal fraction of each group is correctly approved, and among people who truly don't qualify, an equal fraction of each group is correctly rejected.
- Equal Opportunity: a relaxation of Equalized Odds that only requires equal True Positive Rate across groups (equal odds for the qualified population), ignoring the false positive side.
import numpy as np
import pandas as pd
np.random.seed(42)
n = 2000
# Simulate a lending model's predictions across two demographic groups,
# with group B genuinely having a slightly lower TRUE qualification rate
group = np.random.choice(['A', 'B'], size=n)
true_qualified = np.where(group == 'A',
np.random.random(n) < 0.55,
np.random.random(n) < 0.45)
# A model that is well-calibrated to true qualification but not perfectly fair
pred_approved = np.where(true_qualified,
np.random.random(n) < 0.85, # 85% TPR-ish
np.random.random(n) < 0.15) # 15% FPR-ish
df = pd.DataFrame({'group': group, 'qualified': true_qualified, 'approved': pred_approved})
def demographic_parity(df):
return df.groupby('group')['approved'].mean()
def equalized_odds(df):
tpr = df[df['qualified']].groupby('group')['approved'].mean()
fpr = df[~df['qualified']].groupby('group')['approved'].mean()
return tpr, fpr
dp = demographic_parity(df)
tpr, fpr = equalized_odds(df)
print("Demographic Parity (approval rate per group):")
print(dp, f"\ngap: {abs(dp['A'] - dp['B']):.3f}\n")
print("Equalized Odds:")
print(f"TPR -- A: {tpr['A']:.3f}, B: {tpr['B']:.3f} (gap: {abs(tpr['A']-tpr['B']):.3f})")
print(f"FPR -- A: {fpr['A']:.3f}, B: {fpr['B']:.3f} (gap: {abs(fpr['A']-fpr['B']):.3f})")If the two groups have genuinely different true qualification rates in the underlying population (as simulated above), a model that satisfies Demographic Parity — approving both groups at the same rate — will necessarily violate Equalized Odds, and vice versa. This isn't a bug to engineer away; it's a mathematically proven impossibility result (Chouldechova, 2017; Kleinberg et al., 2016) whenever base rates genuinely differ between groups. Choosing which fairness definition to optimize for is a policy decision informed by the specific harm you're trying to prevent, not a technical detail — it should involve legal, ethics, and domain stakeholders, not just the modeling team.
SHAP connects directly back to these metrics: once you've measured a fairness gap with the formulas above, SHAP lets you ask why — pull the SHAP values for false positives/negatives in the disadvantaged group and check whether a proxy feature (zip code, a correlate of the sensitive attribute you deliberately excluded from training) is driving the disparity, exactly as the Real-World Spotlight below demonstrates.
Real-World Spotlight: Credit Denials and Hiring Bias
Case 1: Explaining a Credit Denial to a Loan Applicant
ECOA's Regulation B requires that an adverse action notice list the specific principal reasons for denial — generic statements like "your application did not meet our criteria" are not compliant. A SHAP waterfall decomposition for the denied applicant's row maps directly onto this requirement: take the top-magnitude negative contributors (the reasons pushing toward denial) and translate each into the applicant-facing language a compliance team pre-approves.
import shap
FEATURE_TO_PLAIN_LANGUAGE = {
"debt_to_income": "Your debt-to-income ratio is higher than our threshold",
"credit_score": "Your credit score is below the level we typically approve",
"recent_inquiries": "You have had several recent credit inquiries",
"months_since_late_payment": "You have a recent history of late payments",
"years_employed": "Your employment history is shorter than typical for approval",
"income": "Your reported income is a contributing factor",
"num_open_accounts": "The number of open credit accounts is a contributing factor",
}
def generate_adverse_action_notice(model, explainer, applicant_row, feature_names,
n_reasons=4, decision_threshold=0.5):
shap_exp = explainer(applicant_row)
pred_prob = model.predict_proba(applicant_row)[0, 1]
if pred_prob < decision_threshold:
return {"decision": "approved"}
# Only reasons that pushed TOWARD denial (positive SHAP value for default risk)
risk_factors = [
(f, v) for f, v in zip(feature_names, shap_exp.values[0]) if v > 0
]
risk_factors.sort(key=lambda t: t[1], reverse=True)
top_reasons = risk_factors[:n_reasons]
return {
"decision": "denied",
"predicted_default_probability": round(float(pred_prob), 3),
"principal_reasons": [
FEATURE_TO_PLAIN_LANGUAGE.get(f, f) for f, _ in top_reasons
],
}
notice = generate_adverse_action_notice(
model, explainer, X_test.iloc[[0]], list(X_test.columns)
)
print(json.dumps(notice, indent=2))Case 2: Detecting Proxy Discrimination in a Hiring Model
A hiring screening model may never see a "race" or "gender" column directly (often deliberately excluded) and still discriminate, if it relies heavily on a feature that correlates strongly with a protected attribute — zip code standing in for race, "years since graduation" standing in for age, or college name standing in for socioeconomic background. SHAP makes this auditable: compute SHAP values for the full applicant pool, then check whether the global importance or the sign/magnitude of a suspect feature differs systematically across protected groups, even though the model never used the protected attribute as an input.
import shap
import pandas as pd
import numpy as np
from scipy import stats
def audit_proxy_discrimination(model, X, protected_attribute, feature_names,
suspect_feature, alpha=0.05):
"""
Checks whether a feature's SHAP contribution differs significantly
across groups of a protected attribute NOT used as a model input.
A significant difference suggests the model is using `suspect_feature`
as a proxy for the protected attribute.
"""
explainer = shap.TreeExplainer(model)
shap_exp = explainer(X)
suspect_idx = feature_names.index(suspect_feature)
suspect_shap = shap_exp.values[:, suspect_idx]
groups = protected_attribute.unique()
group_shap = {g: suspect_shap[protected_attribute == g] for g in groups}
print(f"Mean SHAP contribution of '{suspect_feature}' by group:")
for g, vals in group_shap.items():
print(f" {g:15s} mean={vals.mean():+.4f} n={len(vals)}")
# Two-group case: independent t-test on the SHAP value distributions
if len(groups) == 2:
g1, g2 = groups
t_stat, p_value = stats.ttest_ind(group_shap[g1], group_shap[g2])
print(f"\nWelch's t-test on SHAP values: t={t_stat:.3f}, p={p_value:.4f}")
if p_value < alpha:
print(f"FLAG: '{suspect_feature}' contributes significantly differently "
f"across {protected_attribute.name} groups -- investigate as a "
f"potential proxy variable, even though the model never saw "
f"{protected_attribute.name} directly.")
else:
print("No significant disparity detected for this feature.")
# Simulated hiring dataset: 'zip_code_income_tier' as a stand-in for a
# feature that may correlate with race, audited against a protected attribute
# that was withheld from training but is available for fairness auditing
hiring_df = pd.DataFrame({
'years_experience': np.random.exponential(5, 2000),
'zip_code_income_tier': np.random.randint(1, 11, 2000),
'university_tier': np.random.randint(1, 5, 2000),
'interview_score': np.random.normal(70, 15, 2000),
})
protected_group = pd.Series(
np.random.choice(['Group A', 'Group B'], 2000, p=[0.6, 0.4]), name='demographic_group'
)
# Simulate the proxy effect: zip_code_income_tier correlates with the group
hiring_df.loc[protected_group == 'Group B', 'zip_code_income_tier'] -= 2
hire_label = (hiring_df.sum(axis=1) > hiring_df.sum(axis=1).median()).astype(int)
hiring_model = xgb.XGBClassifier(n_estimators=200, max_depth=4, random_state=42)
hiring_model.fit(hiring_df, hire_label) # demographic_group is NOT a feature
audit_proxy_discrimination(
hiring_model, hiring_df, protected_group,
list(hiring_df.columns), suspect_feature='zip_code_income_tier'
)This kind of audit is exactly why "we don't use protected attributes as features" is not, by itself, a defense against disparate impact — regulators and internal fairness reviews increasingly expect this SHAP-based proxy analysis as standard due diligence before a hiring or lending model goes into production.
✍️ Practice Exercises
- Train an
XGBClassifieron a tabular dataset of your choice (e.g. the UCI Adult/Census Income dataset, or a Kaggle credit dataset). Useshap.TreeExplainerto generate a summary plot and a waterfall plot for the single highest-confidence positive prediction and the single highest-confidence negative prediction. Write one sentence per plot describing what a non-technical reviewer should take away. - On the same model, explain the same 10 rows with both
shap.TreeExplainerandlime.lime_tabular.LimeTabularExplainer. For each row, compare the top-3 features each method reports. Where do they agree? Run LIME twice on the same row with different random seeds — how much does the top reason change, and does increasingnum_samplesfrom 500 to 5000 stabilize it? - Build a model where you intentionally engineer a "leaky" feature (one that's a near-copy of the target, e.g. derived using future information). Confirm that SHAP correctly identifies it as overwhelmingly dominant in the summary plot, then remove it, retrain, and confirm the global SHAP ranking changes sensibly. This simulates the debugging use case from Section 1.
- Using the proxy-discrimination audit pattern from the Real-World Spotlight, pick a public dataset that includes a sensitive attribute withheld from training (e.g. the COMPAS recidivism dataset or a hiring dataset with a `gender`/`race` column you exclude from features). Run the t-test-based audit against every remaining feature, not just one, and report which features show a statistically significant SHAP disparity across groups.
▶ Show Solution (Exercise 1 — Summary & Waterfall Plots for Extreme Predictions)
import shap
import xgboost as xgb
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from sklearn.model_selection import train_test_split
from sklearn.datasets import fetch_openml
# UCI Adult / Census Income dataset (predict income > $50K)
adult = fetch_openml(name="adult", version=2, as_frame=True)
X = adult.data.select_dtypes(include=[np.number]).fillna(0)
y = (adult.target == '>50K').astype(int)
X_train, X_test, y_train, y_test = train_test_split(
X, y, test_size=0.2, random_state=42, stratify=y
)
model = xgb.XGBClassifier(n_estimators=300, max_depth=5, learning_rate=0.05,
eval_metric='logloss', random_state=42)
model.fit(X_train, y_train)
explainer = shap.TreeExplainer(model)
shap_exp = explainer(X_test)
# Global summary plot
shap.summary_plot(shap_exp, X_test, show=False)
plt.tight_layout()
plt.savefig("ex1_summary.png", dpi=150)
plt.close()
# Find the most confident positive and negative predictions
probs = model.predict_proba(X_test)[:, 1]
most_confident_positive = int(np.argmax(probs))
most_confident_negative = int(np.argmin(probs))
print(f"Most confident '>50K' prediction: row {most_confident_positive}, "
f"P={probs[most_confident_positive]:.4f}")
print(f"Most confident '<=50K' prediction: row {most_confident_negative}, "
f"P={probs[most_confident_negative]:.4f}")
for label, idx in [("highest-confidence positive", most_confident_positive),
("highest-confidence negative", most_confident_negative)]:
shap.plots.waterfall(shap_exp[idx], show=False)
plt.tight_layout()
plt.savefig(f"ex1_waterfall_{idx}.png", dpi=150)
plt.close()
top = sorted(zip(X_test.columns, shap_exp.values[idx]),
key=lambda t: abs(t[1]), reverse=True)[:3]
print(f"\n{label.title()} (row {idx}) top reasons:")
for feat, val in top:
print(f" {feat:20s} SHAP={val:+.4f}")
# Takeaways:
# - Summary plot: shows which features matter most ACROSS the whole population,
# and whether high values push toward or away from '>50K' (color vs position).
# - Waterfall (positive case): typically shows education/hours-per-week/age
# pushing strongly UP from the base rate toward the '>50K' prediction.
# - Waterfall (negative case): shows the mirror image -- the same features
# at low values pushing DOWN, confirming the model is using them consistently
# in both directions, which is exactly what the consistency property guarantees.📚 Primary Source for This Lesson
Lundberg & Lee (2017) — "A Unified Approach to Interpreting Model Predictions"
The SHAP paper, unifying Shapley values with several earlier explanation methods under one additive attribution framework. For LIME, see Ribeiro, Singh & Guestrin (2016) "'Why Should I Trust You?': Explaining the Predictions of Any Classifier." For the fairness impossibility result in Section 8, see Chouldechova (2017) "Fair Prediction with Disparate Impact."