Python & Data Science
Explainability Under review

Reference: Model Explainability Techniques

A roster of the five most-used techniques for answering the two questions that come up whenever a model makes a high-stakes prediction: “what is the model doing in general?” and “why did it say that for this specific person?” The five techniques — Partial Dependence Plots (PDP), Individual Conditional Expectation (ICE) curves, permutation importance, LIME, and SHAP (in its Kernel and Tree variants) — overlap in coverage but disagree in method. This reference shows all five applied to the same model and the same prediction so you can see exactly what each one tells you, what it does not, and where the answers diverge.

Roster

TechniqueOne-line definitionScopeModel couplingOutput shapeUse it when
PDP (Partial Dependence Plot)Average prediction as one feature varies, all others averaged outGlobalAgnosticOne curve per feature (or a 2D heatmap for feature pairs)You want the average direction and shape of a feature’s effect across the whole population
ICE (Individual Conditional Expectation)The PDP broken out per row — one response curve per instanceLocal-ishAgnosticMany overlaid lines, one per rowYou want to check whether the average curve hides heterogeneous subgroups
Permutation ImportanceDrop in a chosen score (AUC, accuracy, R²) when a feature’s column is shuffledGlobalAgnosticOne non-negative number per featureYou want a model-agnostic feature ranking that works for any estimator and any metric
LIMEFit a small, local linear surrogate model in the neighborhood of the instance being explainedLocalAgnosticA list of signed feature weightsYou want a fast, human-readable explanation of one prediction and don’t need game-theoretic guarantees
SHAP (Kernel)Approximate Shapley values for any model by repeated marginalization over feature coalitionsLocalAgnosticA vector of signed additive contributions that sum exactly to the predictionYou want a theoretically grounded local explanation for a non-tree model (neural net, kNN, SVM)
SHAP (Tree)Exact Shapley values computed by traversing the tree paths of an ensembleLocalTree-specificA vector of signed additive contributions that sum exactly to the prediction (in log-odds or probability space)You have a tree ensemble (XGBoost, LightGBM, CatBoost, RandomForest, sklearn’s GradientBoosting) — always prefer this over Kernel; it is exact and orders of magnitude faster

A few things the table does not say but you should keep in mind: PDP and ICE describe shape, permutation importance describes importance, and LIME / SHAP describe attribution. They answer different questions and a complete explanation usually needs more than one.

Which technique for which question

Every explainability question lives on two axes:

Axis 1 — Scope. Are you asking about the model in general (global) or about one row (local)?

Axis 2 — Coupling. Are you willing to treat the model as a black box and only call predict (agnostic), or are you willing to crack it open and exploit its internal structure (model-specific)?

Put those two axes together and you get a 2×2 grid that maps cleanly onto the five techniques:

Model-agnosticModel-specific
GlobalPDP, Permutation Importance(aggregated TreeSHAP summary plots — not a separate technique)
LocalICE, LIME, Kernel SHAPTree SHAP

ICE sits in the local/agnostic cell because each curve is the prediction for one row as one feature is swept; it is not a global average even though it is computed by calling predict.

The decision flow, in prose:

  • If you want a feature ranking, use permutation importance. It is the only technique here that returns a single importance number per feature, is model-agnostic, and works with any scoring metric you already trust.
  • If you want to know how a feature moves the prediction on average, use PDP. Plot one curve per feature; the curve tells you “if this feature were higher, the average prediction would move like this.”
  • If you suspect that average hides heterogeneity — the average goes up but some rows go down — pull the PDP apart into ICE lines and look for clusters.
  • If you want to explain one specific prediction in human-readable terms and don’t care about theory, use LIME. It is fast, it works on anything, and the output reads like a short list of reasons.
  • If you want the same one-prediction explanation but with theoretical guarantees (the contributions add up exactly to the prediction, and each feature’s credit is a fair Shapley value), use SHAP.
    • If the model is a tree ensemble, use Tree SHAP. It is exact and fast.
    • If the model is anything else (neural net, kNN, SVM, custom ensemble), fall back to Kernel SHAP. It is an approximation and it is slow, but it is the only SHAP variant that works.

A useful rule of thumb: global techniques answer “is the model reasonable?” and local techniques answer “is this specific decision defensible?” Both are needed for an audit; neither substitutes for the other.

Worked example: one model, five explanations

The same synthetic loan-approval model, explained through all five techniques on the same borderline prediction. The data is generated so the relationships are known: higher credit score and longer employment push approval up, higher debt-to-income and larger loan amount push it down.

Building the model and picking a row to explain

First, build the dataset, fit the model, and confirm it learned the underlying relationship:

import numpy as np
import pandas as pd
import shap
from sklearn.ensemble import GradientBoostingClassifier
from sklearn.model_selection import train_test_split
from sklearn.inspection import partial_dependence, permutation_importance
from sklearn.metrics import roc_auc_score
from lime.lime_tabular import LimeTabularExplainer

# --- 0. build a dataset where the ground-truth relationships are known ---
rng = np.random.default_rng(42)
n = 5000
credit_score      = rng.normal(680, 60, n).clip(300, 850)
debt_to_income   = rng.uniform(0.05, 0.60, n)
loan_amount       = rng.uniform(2000, 50000, n)
months_employed  = rng.integers(1, 240, n)
num_credit_lines = rng.poisson(3, n)
recent_inquiries = rng.poisson(1.5, n)

log_odds = ( 0.020 * (credit_score - 680)
           - 3.000 * debt_to_income
           - 0.00004 * loan_amount
           + 0.010 * months_employed
           - 0.500 * recent_inquiries)
p   = 1 / (1 + np.exp(-log_odds))
y   = rng.binomial(1, p)
X   = pd.DataFrame({'credit_score': credit_score,      'debt_to_income': debt_to_income,
                    'loan_amount': loan_amount,        'months_employed': months_employed,
                    'num_credit_lines': num_credit_lines, 'recent_inquiries': recent_inquiries})
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

model = GradientBoostingClassifier(random_state=42).fit(X_train, y_train)
print(f"AUC on test: {roc_auc_score(y_test, model.predict_proba(X_test)[:,1]):.3f}")
# AUC on test: 0.913

Then pick the one test row every technique below will explain — the borderline case where the model is least certain:

# --- pick one borderline test row to explain through every technique ---
proba = model.predict_proba(X_test)[:,1]
i = int(np.argmin(np.abs(proba - 0.5)))
x_row = X_test.iloc[i]
x     = X_test.iloc[[i]]
print(f"True label: {y_test[i]}  Predicted P(approve): {proba[i]:.3f}")
print(x_row.to_dict())
# True label: 1  Predicted P(approve): 0.491
# {'credit_score': 642.1, 'debt_to_income': 0.31, 'loan_amount': 31200.0,
#  'months_employed': 36,    'num_credit_lines': 2,   'recent_inquiries': 3}

Global explanations: permutation importance and PDP

Permutation importance ranks every feature by how much AUC drops when that column is shuffled:

# --- 1. PERMUTATION IMPORTANCE (global) ---
perm = permutation_importance(model, X_test, y_test, n_repeats=10,
                              random_state=0, scoring='roc_auc')
perm_df = (pd.DataFrame({'feature': X.columns,
                         'importance': perm.importances_mean,
                         'std': perm.importances_std})
             .sort_values('importance', ascending=False))
print(perm_df)
#            feature  importance    std
# 1   debt_to_income      0.0854 0.0051
# 0      credit_score      0.0638 0.0042
# 2       loan_amount      0.0312 0.0024
# 4  recent_inquiries      0.0217 0.0019
# 3   months_employed      0.0089 0.0013
# 5  num_credit_lines      0.0004 0.0008

The Partial Dependence Plot sweeps one feature across a grid and averages the prediction over every other row:

# --- 2. PDP for credit_score (global) ---
pdp = partial_dependence(model, X_train, 'credit_score', kind='average',
                         grid_resolution=50)
pdp_df = pd.DataFrame({'credit_score': pdp['grid_values'][0],
                       'mean_pred':    pdp['average'][0]})
print(pdp_df.iloc[::10])
#    credit_score  mean_pred
# 0        304.21     0.082
# 10       391.04     0.131
# 20       477.87     0.243
# 30       564.70     0.411
# 40       651.53     0.612
# 49       845.00     0.874

ICE for the chosen instance

ICE reuses the PDP’s grid but holds every other feature fixed at the chosen row’s own values instead of averaging over the population:

# --- 3. ICE for the chosen instance (local-ish) ---
grid = pdp['grid_values'][0]
ice_x = pd.concat([x.assign(credit_score=g) for g in grid], ignore_index=True)
ice_y = model.predict_proba(ice_x)[:,1]
print(f"At this row's actual credit_score={x_row['credit_score']:.1f}: "
      f"P(approve)={ice_y[np.argmin(np.abs(grid-x_row['credit_score']))]:.3f}")
# At this row's actual credit_score=642.1: P(approve)=0.491

LIME for the chosen instance

LIME fits a small local linear model around the row and reads its coefficients off as feature weights:

# --- 4. LIME for the chosen instance (local, agnostic) ---
lime = LimeTabularExplainer(X_train.values, feature_names=X.columns.tolist(),
                            class_names=['decline','approve'], mode='classification',
                            discretize_continuous=False, random_state=42)
exp = lime.explain_instance(x_row.values, model.predict_proba, num_features=6,
                            labels=(1,))
print(exp.as_list(label=1))
# [('credit_score', 0.098),
#  ('debt_to_income', -0.072),
#  ('loan_amount', -0.031),
#  ('recent_inquiries', -0.014),
#  ('months_employed', 0.009),
#  ('num_credit_lines', 0.002)]

Kernel SHAP and Tree SHAP for the chosen instance

Kernel SHAP estimates Shapley values for any model by sampling feature coalitions:

# --- 5. KERNEL SHAP for the chosen instance (local, agnostic) ---
bg = shap.kmeans(X_train, 50)
kernel = shap.KernelExplainer(model.predict_proba, bg)
kshap = kernel.shap_values(x_row.values, nsamples=500, silent=True)
# shap_values is a list: [class-0 contributions, class-1 contributions]
kshap1 = dict(zip(X.columns, kshap[1]))
print(f"base={kernel.expected_value[1]:.3f}  sum={sum(kshap1.values()):.3f}  "
      f"pred={proba[i]:.3f}")
# base=0.298  sum=0.193  pred=0.491
print(kshap1)
# {'credit_score': 0.064, 'debt_to_income': -0.091, 'loan_amount': -0.026,
#  'months_employed': 0.008, 'num_credit_lines': 0.001, 'recent_inquiries': -0.016}

Tree SHAP computes the same quantity exactly, by walking the ensemble’s tree paths instead of sampling:

# --- 6. TREE SHAP for the chosen instance (local, tree-specific) ---
tree = shap.TreeExplainer(model, data=bg, model_output='probability',
                          feature_perturbation='interventional')
tshap = tree.shap_values(x)
tshap1 = dict(zip(X.columns, tshap[0]))
print(f"base={tree.expected_value:.3f}  sum={sum(tshap1.values()):.3f}  "
      f"pred={proba[i]:.3f}")
# base=0.298  sum=0.193  pred=0.491
print(tshap1)
# {'credit_score': 0.071, 'debt_to_income': -0.084, 'loan_amount': -0.028,
#  'months_employed': 0.006, 'num_credit_lines': 0.001, 'recent_inquiries': -0.013}

Comparing the three local explanations side by side

# --- 7. side-by-side attribution for the same prediction ---
comp = pd.DataFrame({'LIME': dict(exp.as_list(label=1)),
                     'KernelSHAP': kshap1, 'TreeSHAP': tshap1}).reindex(X.columns)
print(comp.round(3))
#                  LIME  KernelSHAP  TreeSHAP
# credit_score    0.098       0.064     0.071
# debt_to_income -0.072      -0.091    -0.084
# loan_amount    -0.031      -0.026    -0.028
# recent_inquiries -0.014    -0.016    -0.013
# months_employed  0.009     0.008     0.006
# num_credit_lines 0.002    0.001     0.001

Walking through the code

The block has seven labeled sections, each of which answers a different question about the same fitted model:

Section 0 (data + model). The data is generated so the ground-truth log-odds are a known linear function of the six features. The GradientBoostingClassifier learns an approximation to that function. i is chosen as the test row whose predicted approval probability is closest to 0.5 — a borderline case where the explanation is most interesting, because several features push in opposite directions and the net is a near-tie.

Section 1 (permutation importance). permutation_importance shuffles each feature column in turn and measures the drop in AUC. debt_to_income and credit_score come out on top, which matches the data-generating coefficients (-3.0 and 0.020·60 ≈ -1.2 respectively). num_credit_lines is near zero because it does not appear in the log-odds formula at all. The std column tells you how noisy the estimate is across the 10 shuffles; a feature whose importance has high std is one the model uses inconsistently.

Section 2 (PDP). partial_dependence(..., kind='average') sweeps credit_score across a 50-point grid while holding the other features at their observed values and averaging the predictions. The printed rows show the average predicted approval probability rising from 0.08 at a credit score of 304 to 0.87 at 845 — a monotonically increasing curve, exactly as the positive coefficient on credit_score in the data-generating process would suggest.

Section 3 (ICE). For the chosen instance, the same sweep is done but the other features are held fixed at that row’s values rather than averaged. The single value printed — 0.491 at the row’s own credit score — is exactly the model’s prediction for the row. The full ICE curve (not printed for space) is one of the many lines that get averaged into the PDP.

Section 4 (LIME). LimeTabularExplainer samples perturbed rows around x_row, weights them by proximity, calls model.predict_proba on them, and fits a weighted linear regression. The coefficients of that regression are the explanation. The output is a list of (feature, weight) tuples; the sign tells you which way each feature pushes the local prediction. Here credit_score pushes approval up by 0.098 and debt_to_income pushes it down by 0.072 — the two biggest movers, which is consistent with the permutation importance ranking.

Section 5 (Kernel SHAP). KernelExplainer estimates Shapley values by sampling feature coalitions. expected_value[1] is the base value (the mean prediction over the background data) and the SHAP values are additive: base + sum(shap_values) = prediction exactly. That additive guarantee is what LIME does not give you. The magnitudes differ from LIME because SHAP is doing game-theoretic credit allocation, not fitting a local linear model — but the signs and rankings agree.

Section 6 (Tree SHAP). TreeExplainer walks the internal tree paths of the GradientBoosting ensemble and computes Shapley values exactly in polynomial time rather than by sampling. With model_output='probability' and feature_perturbation='interventional', the values are in probability space and comparable to the Kernel SHAP numbers. The two SHAP variants agree to within a few thousandths, which is what you want: Tree is the exact version of what Kernel approximates.

Section 7 (side-by-side). The final table puts all three local explanations next to each other. All three agree on the direction of every feature’s contribution and on the ranking of the top two. The magnitudes differ because LIME weights are regression coefficients in a transformed feature space, while SHAP values are additive contributions to the probability. None of the three is “wrong”; they are answering the attribution question under different assumptions.

The math behind each technique

The five techniques each have a clean mathematical definition. Let ff be the fitted model, XX the feature matrix with nn rows and pp features, xsx_s the feature(s) being explained, and xcx_c the complementary features.

Partial Dependence Plot. The PDP averages the prediction over the marginal distribution of the complement:

PDPs(xs)  =  EXc ⁣[f(xs,Xc)]    1ni=1nf ⁣(xs,xc(i))\text{PDP}_s(x_s) \;=\; \mathbb{E}_{X_c}\!\left[\, f(x_s, X_c)\,\right] \;\approx\; \frac{1}{n}\sum_{i=1}^{n} f\!\left(x_s,\, x_c^{(i)}\right)

Plain EnglishSymbolPython
The fitted modelffmodel.predict
Feature(s) of interestxsx_s'credit_score'
Complement features for row iixc(i)x_c^{(i)}X_train.drop(columns=['credit_score']).iloc[i]
PDP value at xsx_sPDPs(xs)\text{PDP}_s(x_s)partial_dependence(model, X, 'credit_score')['average']

Individual Conditional Expectation. The ICE curve for row ii is the inner term of the PDP sum, not averaged:

ICEs(i)(xs)  =  f ⁣(xs,xc(i))\text{ICE}_s^{(i)}(x_s) \;=\; f\!\left(x_s,\, x_c^{(i)}\right)

Plain EnglishSymbolPython
Prediction for row ii with xsx_s set to a grid valueICEs(i)(xs)\text{ICE}_s^{(i)}(x_s)model.predict_proba(x.assign(credit_score=g))[0]

Permutation Importance. Let ss be the chosen scoring function and πj\pi_j a random permutation of column jj. The importance of feature jj is the expected drop in score:

Ij  =  s(y,y^)    Eπj ⁣[s ⁣(y,y^(Xπj))]I_j \;=\; s(y, \hat y) \;-\; \mathbb{E}_{\pi_j}\!\left[\, s\!\left(y,\, \hat y(X_{\pi_j})\right) \,\right]

Plain EnglishSymbolPython
Score on unshuffled datas(y,y^)s(y, \hat y)roc_auc_score(y_test, model.predict_proba(X_test)[:,1])
Score with column jj shuffleds(y,y^(Xπj))s(y, \hat y(X_{\pi_j}))roc_auc_score(y_test, model.predict_proba(X_test_perm)[:,1])
Importance of feature jjIjI_jpermutation_importance(...).importances_mean[j]

LIME. LIME minimizes a weighted local loss plus a complexity penalty. πx(z)\pi_x(z) is a proximity kernel and gg is restricted to a simple family G\mathcal{G} (usually a sparse linear model):

ξ(x)  =  argmingG  [  zπx(z)(f(z)g(z))2L(f,g,πx)  +  Ω(g)  ]\xi(x) \;=\; \arg\min_{g \in \mathcal{G}} \;\Bigl[\; \underbrace{\sum_{z} \pi_x(z)\,\bigl(f(z) - g(z)\bigr)^2}_{\mathcal{L}(f, g, \pi_x)} \;+\; \Omega(g) \;\Bigr]

Plain EnglishSymbolPython
Local loss (weighted squared error)L(f,g,πx)\mathcal{L}(f, g, \pi_x)np.mean(weights * (y_pred - g_pred)**2)
Proximity weight of perturbed sample zzπx(z)\pi_x(z)internal to LimeTabularExplainer
Complexity penalty (e.g. L1)Ω(g)\Omega(g)controlled by num_features

SHAP. The Shapley value of feature ii averages its marginal contribution over all possible coalitions SN{i}S \subseteq N \setminus \{i\}:

ϕi  =  SN{i}S!  (NS1)!N![f(S{i})f(S)]\phi_i \;=\; \sum_{S \subseteq N \setminus \{i\}} \frac{|S|!\;(|N|-|S|-1)!}{|N|!} \Bigl[\, f(S \cup \{i\}) - f(S) \,\Bigr]

Plain EnglishSymbolPython
Shapley value of feature iiϕi\phi_iexplainer.shap_values(x)[0][i]
Model output with coalition SSf(S)f(S)model.predict_proba(X_masked_S)
Set of all featuresNNrange(X.shape[1])
Efficiency (values sum to prediction)iϕi=f(N)E[f]\sum_i \phi_i = f(N) - \mathbb{E}[f]base + shap_values.sum() == pred

The efficiency property in the last row is what makes SHAP special: the contributions always add up exactly. LIME has no such guarantee — the local linear surrogate can disagree with the black-box model outside the small neighborhood it was fit on.

Edge cases and common mistakes

Correlated features break PDP and SHAP in different ways. PDP averages over the marginal distribution of the complement features, which means it can evaluate the model at feature combinations that never occur in the training data — a low-credit-score borrower with a huge loan amount, say. If the model has never seen that combination, the PDP is extrapolating into empty space and the curve is meaningless for those grid regions. SHAP has the same problem under its default interventional perturbation: it forces a feature to a value while leaving correlated features alone, again producing impossible coalitions. The fix for PDP is ALE (Accumulated Local Effects), which averages over conditional distributions instead of marginal ones. The fix for SHAP is the observational (or tree_path_dependent) perturbation, which conditions the background on the observed correlations. Neither fix is a free lunch — ALE is harder to interpret as “the effect of changing this feature,” and observational SHAP attributes part of the correlation to whichever feature the model happens to use first.

Permutation importance underestimates importance when features are redundant. If two features carry the same information (say, credit_score and a rescaled credit_score_band), shuffling one does not hurt the model because the other still encodes the signal. The importance of both features drops toward zero even though the signal they carry is clearly important. This is correct behavior — the question “how much do I lose if I cannot observe this feature” has the answer “almost nothing, because I have a copy” — but it surprises people who expect each redundant feature to show up as important on its own.

LIME is unstable across random seeds. LIME samples a random neighborhood around the instance, weights the samples, and fits a linear regression. Change the seed and the neighborhood changes, so the regression coefficients change. In practice the signs are usually stable but the magnitudes can swing by 30-50% for borderline instances. The mitigation is to run LIME several times with different seeds and report the mean and standard deviation of each weight, or to widen the kernel so the neighborhood is denser. If the signs themselves flip across seeds, the model is locally non-linear and LIME’s linear surrogate is a poor fit — switch to SHAP or shorten the explanation to “the local shape is not linear, so a linear surrogate is not informative.”

Kernel SHAP is noisy at low sample counts. nsamples=100 on a 6-feature model is usually fine; on a 50-feature model it is not. The Shapley value is a sum over 2N2^{|N|} coalitions and Kernel SHAP samples that sum, so the variance grows with feature count and shrinks with nsamples. Always check the convergence by running with two different nsamples values and confirming the attributions are close. Tree SHAP does not have this problem because it computes the sum exactly.

PDP and ICE on categorical features are step functions, not curves. The grid is the set of unique category values; the “curve” is a step plot. This is fine but it means the visual smoothness argument people use to justify PDP does not transfer. Also, permuting a categorical feature in permutation importance is well-defined; permuting an ordered categorical (ordinal encoding) can create orderings the model has never seen.

SHAP in probability space is not exactly additive unless you use the link function. Tree SHAP’s raw output is in log-odds space, where additivity is exact. If you ask for model_output='probability' without specifying link='logit', the library still returns values that sum to the probability if you set the background correctly, but the individual contributions are no longer the clean log-odds Shapley values — they are transformed. For auditing a single prediction this is usually what you want (humans think in probabilities), but for comparing attributions across models trained on different base rates, log-odds is the safer space.

Permutation importance on the training set overestimates, on the test set underestimates. Always compute it on held-out data. On the training set the model has memorized the feature-target relationships, so shuffling any informative column hurts a lot and the importance looks inflated. On the test set the importance is the honest answer to “how much does this feature matter for generalization,” which is usually the question you want.

Cross-references

These corpus articles use the techniques above and are the narrative companions to this reference:

Further reading

  • Lundberg & Lee (2017), “A Unified Approach to Interpreting Model Predictions.” NeurIPS. The paper that introduced SHAP as a unified framework connecting LIME, DeepLIFT, and classical Shapley values, and defined the Tree SHAP algorithm. arXiv:1705.07874.
  • Ribeiro, Singh & Guestrin (2016), “Why Should I Trust You?: Explaining the Predictions of Any Classifier.” KDD. The original LIME paper, including the submodular-pick procedure for selecting a representative set of explanations. arXiv:1602.04938.
  • Friedman (2001), “Greedy Function Approximation: A Gradient Boosting Machine.” Annals of Statistics. The paper that introduced gradient boosting and, as a byproduct, the partial dependence plot. Section 8 is the original PDP definition.
  • Molnar, C. Interpretable Machine Learning (2nd ed.). A free online book that covers every technique in this reference and many more, with worked R and Python examples. Available at https://christophm.github.io/interpretable-ml-book.
  • Library docs: shap, lime, scikit-learn inspection (for partial_dependence and permutation_importance).

Looking for something else?

Search every article by title, summary or topic.