Explainable AI: SHAP, LIME, and Feature Importance in Practice

Explainable AI (XAI) helps a bank answer a critical question after a model makes a decision: why did this customer receive this credit decision? This article uses a loan-default scenario to compare three methods, with tested Python for each. SHAP uses Shapley values to split one prediction into additive feature contributions, LIME fits a simple local surrogate model around one prediction, and permutation importance measures how much shuffling a feature degrades model performance across a population.

Banking problem statement

A retail bank trains a model to predict whether an applicant will default within 12 months. Credit officers, risk teams, regulators, and customers need more than a probability: why an application was declined or flagged for review, which attributes drove the score, whether the model relies on unstable, correlated, biased, or inappropriate variables, and how explanations are produced consistently in production. A useful explanation is faithful to the model, stable enough for its purpose, and understandable to a non-technical analyst. These properties must be tested; no explanation method guarantees them by itself.

Running example: loan default

The model is a gradient-boosting classifier trained on synthetic data where default = 1 means the customer defaulted. Because the data is synthetic, we know the ground truth — including two features, age and num_open_accounts, that have no effect at all — and can check whether each method recovers it.

# Pin versions: SHAP output shapes and defaults have changed across releases
pip install "scikit-learn>=1.4" "shap>=0.46" lime pandas numpy matplotlib

Shared data and model

Two problems in the earlier version of this example changed every result downstream. First, income and loan amount entered the synthetic default signal in raw dollars (-0.00003 × income), which pushed the default rate down to about 3% and made income dominate everything. Second, the model was judged with classification_report(y_test, model.predict(X_test)), which uses a 0.5 cut-off: at that rate it flagged no defaulters at all (recall 0.00 on 39 defaults) while reporting 97% accuracy. It also overfitted, giving the two no-effect features more SHAP weight than credit utilization.

The version below puts income and loan size on a log scale (default rate about 12%), stops boosting early, adds monotonic constraints so risk can only rise with utilization or debt-to-income and only fall with income — a common requirement for credit models because it keeps explanations directionally sensible — and evaluates at the bank’s review threshold with ROC-AUC, PR-AUC, and the Brier score. In testing it reached ROC-AUC 0.76 and PR-AUC 0.32, and flagged 26% of defaulters at the 0.30 threshold. Save it as data_model.py; the later sections import from it.

# Save as data_model.py: the later sections import from it.
import numpy as np
import pandas as pd
from sklearn.ensemble import HistGradientBoostingClassifier
from sklearn.metrics import (average_precision_score, brier_score_loss,
                             classification_report, roc_auc_score)
from sklearn.model_selection import train_test_split

RANDOM_STATE = 42
REVIEW_THRESHOLD = 0.30  # set by credit policy, not predict()'s 0.5 default
rng = np.random.default_rng(RANDOM_STATE)

n = 20_000
df = pd.DataFrame({
    "age": rng.integers(21, 70, n),
    "annual_income_k": rng.lognormal(4.1, 0.45, n).round(1),  # in $000s
    "loan_amount_k": rng.lognormal(2.3, 0.55, n).round(1),    # in $000s
    "credit_utilization": rng.beta(2, 5, n).round(3),
    "months_since_delinquency": rng.integers(0, 121, n),
    "num_open_accounts": rng.integers(1, 15, n),
    "employment_years": rng.integers(0, 41, n),
    "debt_to_income": np.clip(rng.beta(2, 6, n) * 1.5, 0, 1).round(3),
    "prior_defaults": rng.binomial(1, 0.12, n),
})

# Known ground truth. Income and loan size enter on a log scale so no
# single dollar-denominated column swamps the signal. age and
# num_open_accounts have NO effect: a sound explanation should say so.
logit = (
    -2.6
    + 2.5 * df["credit_utilization"]
    + 2.0 * df["debt_to_income"]
    + 1.4 * df["prior_defaults"]
    - 0.012 * df["months_since_delinquency"]
    - 0.6 * np.log(df["annual_income_k"] / 60)
    + 0.5 * np.log(df["loan_amount_k"] / 10)
    - 0.03 * df["employment_years"]
)
df["default"] = rng.binomial(1, 1 / (1 + np.exp(-logit)))

features = [c for c in df.columns if c != "default"]
X, y = df[features], df["default"]
X_train, X_test, y_train, y_test = train_test_split(
    X, y, test_size=0.25, stratify=y, random_state=RANDOM_STATE)

model = HistGradientBoostingClassifier(
    max_iter=500, learning_rate=0.05, max_depth=4,
    early_stopping=True, n_iter_no_change=20,  # stop before fitting noise
    # Risk may only rise with utilization, DTI, etc. (+1) or fall (-1)
    monotonic_cst={"credit_utilization": 1, "debt_to_income": 1,
                   "prior_defaults": 1, "loan_amount_k": 1,
                   "months_since_delinquency": -1, "annual_income_k": -1,
                   "employment_years": -1},
    random_state=RANDOM_STATE,
)
model.fit(X_train, y_train)

probs = model.predict_proba(X_test)[:, 1]
print(f"Default rate: {y.mean():.1%}  trees: {model.n_iter_}")
print(f"ROC-AUC: {roc_auc_score(y_test, probs):.3f}  "
      f"PR-AUC: {average_precision_score(y_test, probs):.3f}  "
      f"Brier: {brier_score_loss(y_test, probs):.4f}")
flagged = (probs >= REVIEW_THRESHOLD).astype(int)
# The riskiest test applicant: explained by every method below
APPLICANT_ID = X_test.index[np.argmax(probs)]
print(classification_report(y_test, flagged, zero_division=0,
                            target_names=["repaid", "default"]))

The rows and labels are synthetic and say nothing about real credit risk. A real bank needs governed feature and label definitions, point-in-time correctness, privacy controls, calibration checks (the threshold assumes the probabilities mean what they say), and monitoring.

1. SHAP

Problem statement

The bank needs customer-level explanations. If an application is declined, the credit officer needs the main reasons — for example high credit utilization, high debt-to-income ratio, or a prior default.

Conceptual explanation

SHAP assigns each feature a contribution to the difference between the model’s output for one applicant and a baseline. For applicant i with M features, f(xᵢ) = φ₀ + Σⱼ φᵢⱼ, where f is the model output on the scale the explainer uses, φ₀ is the base value on that same scale, and φᵢⱼ is feature j’s contribution. With default settings, TreeExplainer on this scikit-learn model works on the raw log-odds scale: positive values push toward default, negative values away from it, and they are not percentage-point changes in probability. A probability-scale explanation is possible but must be requested explicitly, and it then needs a background dataset — which is where the choice of a representative background matters. The default tree-path-dependent algorithm uses the training distribution stored in the trees instead.

Runnable Python code

The code computes the explanation once and reuses it, asserts that base value plus contributions reproduces the model’s raw score, and handles both 2-D and 3-D SHAP outputs. The older pattern if isinstance(shap_values, list) no longer covers multi-output models: current SHAP returns a 3-D array (for a random forest, shape (n, features, 2)), which that check passes straight through to code expecting 2-D input.

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import shap

from data_model import (APPLICANT_ID, RANDOM_STATE, X_test, X_train,
                        features, model)

# Default TreeExplainer: path-dependent, no background data needed,
# output on the model's raw (log-odds) scale.
explainer = shap.TreeExplainer(model)
X_explain = X_test.sample(500, random_state=RANDOM_STATE)

explanation = explainer(X_explain)  # compute once; reuse for all outputs
if explanation.values.ndim == 3:     # some models return one slice per class
    explanation = explanation[:, :, 1]
values = explanation.values

# Additivity check: base value + contributions == model raw score
raw_score = model.decision_function(X_explain)
assert np.allclose(explanation.base_values + values.sum(axis=1),
                   raw_score, atol=1e-6)

global_shap = pd.Series(np.abs(values).mean(axis=0), index=features)
print(global_shap.sort_values(ascending=False).round(3))

# Local explanation for one applicant, in log-odds
applicant = X_test.loc[[APPLICANT_ID]]
local_expl = explainer(applicant)
local = pd.DataFrame({"value": applicant.iloc[0],
                      "shap_log_odds": local_expl.values[0]})
print(local.sort_values("shap_log_odds", key=abs, ascending=False))

# Probability-scale explanation needs a background sample and the
# interventional algorithm; contributions then sum to the probability.
background = shap.sample(X_train, 100, random_state=RANDOM_STATE)
prob_explainer = shap.TreeExplainer(
    model, data=background, feature_perturbation="interventional",
    model_output="probability")
prob_expl = prob_explainer(applicant)
print("P(default):", model.predict_proba(applicant)[0, 1].round(3),
      "=", (prob_expl.base_values[0] + prob_expl.values[0].sum()).round(3))

# Plots for review packs: show=False so they also work outside notebooks
shap.plots.beeswarm(explanation, show=False)
plt.savefig("shap_beeswarm.png", dpi=150, bbox_inches="tight")
plt.close("all")
shap.plots.waterfall(local_expl[0], show=False)
plt.savefig("shap_waterfall.png", dpi=150, bbox_inches="tight")
plt.close("all")

On the test data, the two no-effect features had the smallest mean |SHAP| (about 0.04, against 0.20–0.35 for real drivers) and the probability-scale contributions summed to the model’s probability. shap.plots.beeswarm replaces the legacy summary_plot.

Real-world use case and reason codes

The riskiest applicant in the test set has a predicted default probability of 0.92, above the review threshold. Its largest log-odds contributions come from a prior default, high debt-to-income, and short employment. Raw SHAP values are not something to show a customer; a reviewed mapping turns them into reasons, citing only factors that pushed risk up, largest first. Features that must never appear as a reason — age here — cause the case to be escalated rather than silently skipped, and a feature with no reviewed wording is an error, not a blank.

import numpy as np

# Reviewed by credit policy and compliance. None = must never be cited;
# a decision driven by it is escalated rather than explained away.
REASON_TEXT = {
    "credit_utilization": "Proportion of revolving balances to limits "
                          "is too high",
    "debt_to_income": "Debt obligations are high relative to income",
    "prior_defaults": "Previous default on a credit obligation",
    "months_since_delinquency": "Recent delinquency on an account",
    "annual_income_k": "Income insufficient for amount of credit requested",
    "loan_amount_k": "Amount of credit requested is high",
    "employment_years": "Length of employment",
    "num_open_accounts": "Number of open accounts",
    "age": None,
}

class NeedsManualReview(Exception):
    pass

def reason_codes(feature_names, shap_row, max_reasons=4,
                 min_contribution=0.05):
    """Map one applicant's log-odds SHAP values to adverse-action reasons.

    Only features that pushed risk UP are cited, largest first.
    """
    missing = set(feature_names) - REASON_TEXT.keys()
    if missing:
        raise KeyError(f"no reviewed reason text for: {sorted(missing)}")
    reasons = []
    for j in np.argsort(shap_row)[::-1]:
        if shap_row[j] < min_contribution or len(reasons) == max_reasons:
            break
        text = REASON_TEXT[feature_names[j]]
        if text is None:
            raise NeedsManualReview(
                f"{feature_names[j]} is a principal driver of this decision")
        reasons.append(text)
    return reasons

# Usage with section 2: reason_codes(features, local_expl.values[0])

For this applicant the mapper returned: previous default; debt obligations high relative to income; length of employment; recent delinquency. A credit officer can then say the application was routed for review mainly because of a previous default and high debt relative to income.

Common mistakes and production notes

Wrong output or scale. Confirm the explanation refers to the default class and document whether it is log-odds or probability; never present log-odds values as probability changes.

Causal reading. SHAP explains the model, not the real-world cause of default.

Correlated features. Attribution spreads across related variables (income, loan size, and debt-to-income in real data), so individual values can look small even when the group matters.

Global ≠ local. Mean |SHAP| ranks features across a population; it does not describe any single applicant.

Governance. Explain with the same feature pipeline used for scoring; store model version, feature schema, explainer settings, and SHAP version with each explanation; monitor drift in both inputs and SHAP importance, across time and segments; require human review for high-impact decisions such as declines.

Protected attributes and proxies. Do not use protected attributes in scoring unless law and policy explicitly allow it. Where permitted, keep access-controlled data for fairness testing, and test whether other features act as proxies — excluding a column does not remove proxy effects.

2. LIME

Problem statement

The bank uses a complex or proprietary scoring service where exact SHAP values are hard or expensive to obtain, and still needs a model-agnostic explanation for a particular application.

Conceptual explanation

LIME generates perturbed versions of one applicant, scores them with the original model, and fits a weighted linear model to those scores, weighting samples by how close they are to the applicant. The result is a local approximation: “for this applicant, these features most influenced the decision.”

A common description — that LIME “samples around the instance” — is not what the lime package does by default. With discretize_continuous=True, each numeric feature is binned by training quartiles and new values are drawn from the training distribution; without discretization, values are drawn from a normal distribution centred on the training mean unless sample_around_instance=True. The proximity weighting is what makes the fit local. The package itself has had no release since version 0.2.0.1 in June 2020, so pin it and test it like any unmaintained dependency.

Runnable Python code

Compared with the earlier version, this declares prior_defaults as categorical (the article advised it but the code did not), uses more samples, and actually measures the two properties the checklist asks for: local fidelity and stability across seeds.

import numpy as np
import pandas as pd
from lime.lime_tabular import LimeTabularExplainer

from data_model import (APPLICANT_ID, RANDOM_STATE, X_test, X_train,
                        features, model)

CATEGORICAL = ["prior_defaults"]
cat_idx = [features.index(c) for c in CATEGORICAL]

def predict_proba_np(X_array):
    # Restore column names and order: the model was fitted on a DataFrame
    return model.predict_proba(pd.DataFrame(X_array, columns=features))

def make_explainer(seed):
    return LimeTabularExplainer(
        training_data=X_train.to_numpy(dtype=float),
        feature_names=features,
        class_names=["repaid", "default"],
        categorical_features=cat_idx,
        categorical_names={cat_idx[0]: ["no prior default",
                                        "prior default"]},
        mode="classification",
        discretize_continuous=True,
        random_state=seed,
    )

row = X_test.loc[APPLICANT_ID].to_numpy(dtype=float)
exp = make_explainer(RANDOM_STATE).explain_instance(
    row, predict_proba_np, labels=(1,), num_features=5, num_samples=5_000)

for rule, weight in exp.as_list(label=1):
    print(f"{rule:<40} {weight:+.4f}")

# Local fidelity: does the surrogate reproduce the model here?
model_p = predict_proba_np(row[None, :])[0, 1]
print(f"model P(default) {model_p:.3f}  surrogate {exp.local_pred[0]:.3f}  "
      f"surrogate R^2 {exp.score:.2f}")

# Stability: the same applicant under different seeds
def top_features(seed, k=5):
    e = make_explainer(seed).explain_instance(
        row, predict_proba_np, labels=(1,), num_features=k,
        num_samples=5_000)
    return {features[j] for j, _ in e.as_map()[1]}

runs = [top_features(seed) for seed in range(5)]
overlap = np.mean([len(a & b) / len(a | b)
                   for n, a in enumerate(runs) for b in runs[n + 1:]])
print(f"mean Jaccard overlap of top-5 features across seeds: {overlap:.2f}")

For the same high-risk applicant, LIME ranked prior default, debt-to-income, and credit utilization highest, and the top five features were identical across five seeds. Fidelity was the weak point: the surrogate predicted 0.51 where the model predicted 0.92 (surrogate R² 0.57). That gap is exactly why LIME output needs a fidelity check before anyone relies on it. Note also that one explainer object advances its random state on every call: in testing, two calls on the same explainer gave different weights, while a fresh explainer with the same seed reproduced the first call. A “fixed seed” only gives reproducible output if the explainer is created per explanation.

Real-world use case

A fraud team uses a third-party transaction-risk API that returns only a score. Wrapping it with LIME answers “which observable patterns raised this score?” — but each explanation makes thousands of scoring calls (5,000 here), so check the API’s rate limits, cost, and data-sharing terms first, and batch the calls. LIME weights are on the probability scale of the explained class and are relative to the applicant’s own bins, so they are not comparable across applicants.

Common mistakes and production notes

Treating local weights as global importance, or comparing them across applicants as if standardized.

Too few samples. Unstable explanations usually mean num_samples is too low; measure variance across seeds rather than guessing.

Input mismatch. LIME passes NumPy arrays; restore column names, order, and types before calling a model fitted on a DataFrame, and use the same preprocessing as scoring.

Undeclared categoricals. Binary and categorical features must be declared, or they are perturbed as if continuous.

Caching. Cache an explanation only when inputs, model version, and explainer configuration are unchanged.

Regulated communications. Use LIME for analyst diagnostics. Customer notices need validated, decision-specific reason codes that reflect the principal factors actually used; SHAP or LIME output alone does not satisfy adverse-action requirements. Test across segments, including new-to-credit and thin-file applicants.

3. Permutation feature importance

Problem statement

Before deployment, the model-risk team wants to know which features the model relies on across the whole validation population, to support feature governance, data-quality monitoring, and early detection of leakage or over-reliance on one variable.

Conceptual explanation

Permutation importance shuffles one feature at a time and measures how much a metric falls; scikit-learn defines it as the baseline score minus the score after permuting the column. It is a global, performance-based measure and does not explain individual decisions. Computing it on training data as well as test data is a cheap overfitting check: a feature that matters much more on training data than on held-out data has been memorized.

Runnable Python code

import matplotlib.pyplot as plt
import pandas as pd
from sklearn.inspection import permutation_importance

from data_model import (RANDOM_STATE, X_test, X_train, features, model,
                        y_test, y_train)

def importance(X, y):
    r = permutation_importance(model, X, y, scoring="roc_auc", n_repeats=20,
                               random_state=RANDOM_STATE, n_jobs=-1)
    return pd.DataFrame({"mean": r.importances_mean,
                         "std": r.importances_std}, index=features)

test_imp = importance(X_test, y_test).sort_values("mean")
train_imp = importance(X_train, y_train)
print(test_imp[::-1].round(4))

# Train >> test importance for a feature means the model memorised it
gap = (train_imp["mean"] - test_imp["mean"]).sort_values(ascending=False)
print("\nTrain-minus-test importance:\n", gap.round(4))

fig, ax = plt.subplots(figsize=(8, 5))
ax.barh(test_imp.index, test_imp["mean"], xerr=test_imp["std"])
ax.axvline(0, color="grey", linewidth=0.8)
ax.set_xlabel("Decrease in ROC-AUC when feature is permuted (test set)")
ax.set_title("Permutation feature importance: loan-default model")
fig.tight_layout()
fig.savefig("permutation_importance.png", dpi=150)
plt.close(fig)

Real-world use case

On the test set, prior defaults, debt-to-income, and credit utilization had the largest importance (each about 0.045–0.05 ROC-AUC), and the two no-effect features were near zero, matching the ground truth. The earlier version’s narrative — credit utilization most important — did not match what its own code produced: there, income dominated and credit utilization scored slightly below zero. The lesson for governance is to write findings from the run, not from expectations.

When a feature expected to matter shows near-zero importance, it should trigger investigation: is it populated correctly, missing for many applicants, transformed in a way that removed signal, or affected by leakage or an upstream defect? Importance tracked this way becomes a model-governance control, not just a chart.

Common mistakes and production notes

Impurity-based importances (feature_importances_ on many tree models) are biased toward high-cardinality and continuous features; prefer permutation importance on held-out data.

Correlated features. Permuting one of a correlated pair can understate its importance because the model falls back on the other; check correlated groups and consider permuting them together.

Single numbers. Report mean and standard deviation over repeats, on a representative, time-based validation set, and document metric, dataset, seed, and software versions.

Monitoring. Track importance after each retrain, alert when a critical feature collapses or an unexpected one dominates, and confirm production feature definitions match training. Do not compare importance across models with different feature sets or pipelines, and do not confuse “important to performance” with “important to one decision.”

Choosing the right method

End-to-end production architecture

A robust deployment separates scoring from explanation while tying them to the same inputs:

Feature service produces a versioned feature vector.

Model service returns a default probability and decision.

Explanation service computes SHAP values from the exact vector that was scored. Tree SHAP is fast enough to run synchronously for declines; LIME is slower and stochastic, so keep it for offline analysis.

Reason-code mapper converts attributions into reviewed, business-readable reasons and escalates disallowed drivers.

Audit store saves the input hash, model and explainer versions, score, explanation, threshold, and decision.

Monitoring service tracks score distribution, calibration, feature drift, explanation stability, and fairness metrics.

Human-review workflow presents explanations to underwriters for flagged cases.

Takeaway: For regulated lending, explanations support — not replace — sound credit policy, documented model validation, and appropriate human oversight. Test them against known ground truth before trusting them on real applicants.

Back to all blogs