How to Explain Machine Learning Predictions in Python: A Beginner's Guide to Feature Importance and SHAP (Delivery ML Series, Part 6)

The short answer

A prediction people can't question is a prediction people won't trust. In this guide you'll find out why the model predicted 63.8 minutes for a rainy Saturday biryani order: rain added 9.4 minutes, Old Town 7.6 and Saturday dinner time 5.3. Then you'll check whether the model treats riders fairly, and make the API explain itself.

Quick answer: To explain a machine learning model in Python: (1) for a global view, use the model's weights (linear models) or permutation importance (any model), which shuffles one feature at a time and measures how much worse predictions get; (2) to explain one prediction, use SHAP values, which split the prediction into a starting average plus one contribution per feature that adds up exactly to the prediction; (3) check fairness by comparing errors across groups; (4) remember that explanations describe the model, not the real world.

What you will build

By the end of this part you can answer three questions about the delivery model, the kind that managers, customers and regulators actually ask:

Then you'll add an optional ?explain=true to the Part 4 API, so the app can show customers why their delivery will take longer:

{
  "predicted_minutes": 63.8,
  "average_minutes": 40.5,
  "top_reasons": [
    {"feature": "weather", "minutes": 9.5},
    {"feature": "city_zone", "minutes": 7.6},
    {"feature": "order_time", "minutes": 5.6}
  ]
}

(The figures differ slightly from the 9.4 and 7.6 in Step 4 because the API uses a smaller sample of past orders, to stay fast. Step 7 and the theory break "how to judge an explanation" explain why.)

How this part connects to the series

This is Part 6 of the Delivery ML series. It explains the pipeline you built in Part 3 and served in Part 4, and compares it with the gradient boosting model from Part 3, Step 6. (Part 5, on tuning, isn't required: this part works directly from Parts 3 and 4.) In Part 8, these explanations become the raw material for a generative AI model that writes friendly delay messages to customers.

Prerequisites

You need Part 3's data and feature files, Part 4's API folder for the last step, and one new library: SHAP.

From earlier in the series

Knowledge

Software

Hardware, time and cost

• Computer: any laptop or Google Colab (Steps 1–6 work in Colab; Step 7 needs your own computer, as in Part 4).
• Time: about 75 minutes.
• Cost: free.
• Patience: one optional calculation in Step 5 takes about 30 seconds. Everything else runs in a second or two.

Checklist before you start

• I have delivery_data.py and delivery_features.py from Part 3
• I can rebuild Part 3's linear and gradient boosting pipelines (Step 1 shows the code)
• I have Part 4's part4-api folder (for Step 7)
• I can install a new library with pip

Theory break: why explain a model, and what kind of explanation?

Explanations help you trust, debug and improve a model, and in many industries they're required by law. There are two kinds: global ("what does the model care about in general?") and local ("why this particular prediction?").

Why explain at all?

A model that's accurate on the test set can still be right for the wrong reasons. Explanations are how you find out.

Global vs local explanations

The two can disagree, and that's fine. Distance matters most overall, but for our rainy Old Town biryani the biggest push comes from the rain. A feature that rarely matters can matter hugely for one order.

Interpretable models vs black boxes

• Interpretable models can be read directly. Part 2's linear regression literally is a recipe: "−4.74 + 3.63 minutes per km + …". A small decision tree is a readable flowchart.
• Black-box models are too complex to read: a random forest of 200 trees, gradient boosting with 100 trees, a neural network with millions of weights.

Black boxes are often a little more accurate, so we need tools that work on any model. That's what permutation importance (Step 3) and SHAP (Steps 4–5) do: they treat the model as a black box, change its inputs, and watch how the outputs move.

An honest warning before we start

Every technique in this guide explains what the model does, not how the world works. If the model says rain adds 9.4 minutes, that's what the model believes. It matches reality here only because our data has a known recipe. With real data, an explanation can reveal a model's mistakes just as easily as the world's truths, and that's one of the most useful things it can do.

Step 1: Set up and rebuild Part 3's two models

Install SHAP, copy Part 3's two helper files, and rebuild the linear and gradient boosting pipelines. They predict 63.8 and 66.7 minutes for the messy biryani order, the two answers we're about to explain.

Install SHAP and create the folder

cd first-ml-model
source venv/bin/activate          # Windows: venv\Scripts\activate
pip install shap

mkdir part6-explain
cd part6-explain
cp ../part3-clean-data/delivery_data.py .
cp ../part3-clean-data/delivery_features.py .

In Colab, run !pip install shap and create the two files with %%writefile, as in Part 3.

Rebuild the pipelines

This is Part 3's code, unchanged: the same data, split, preprocessing and models. The linear pipeline it builds is identical to the delivery_pipeline.joblib you saved in Part 3. Save it as explain_delivery.py, or run it in a notebook cell:

import numpy as np
import pandas as pd
from sklearn.compose import ColumnTransformer
from sklearn.ensemble import GradientBoostingRegressor
from sklearn.impute import SimpleImputer
from sklearn.linear_model import LinearRegression
from sklearn.model_selection import train_test_split
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import FunctionTransformer, OneHotEncoder, StandardScaler

from delivery_data import make_messy_orders
from delivery_features import clean_orders

orders = make_messy_orders().drop_duplicates()
X = orders.drop(columns="delivery_min")
y = orders["delivery_min"]
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

numeric_features = ["distance_km", "prep_time_min", "rider_experience_months", "is_raining",
                    "is_peak_hour", "is_weekend", "distance_x_rain", "distance_x_peak",
                    "distance_x_rain_x_peak"]
categorical_features = ["city_zone", "cuisine"]


def make_preprocess():
    return ColumnTransformer([
        ("num", make_pipeline(SimpleImputer(strategy="median"), StandardScaler()), numeric_features),
        ("cat", make_pipeline(SimpleImputer(strategy="most_frequent"),
                              OneHotEncoder(handle_unknown="ignore")), categorical_features),
    ])


linear = make_pipeline(FunctionTransformer(clean_orders), make_preprocess(), LinearRegression())
linear.fit(X_train, y_train)

boosted = make_pipeline(FunctionTransformer(clean_orders), make_preprocess(),
                        GradientBoostingRegressor(random_state=42))
boosted.fit(X_train, y_train)

new_order = pd.DataFrame([{
    "order_time": "2026-10-03 20:15", "city_zone": "  old TOWN", "cuisine": "Biryani",
    "distance_km": 6000, "prep_time_min": None, "weather": "RAINY", "rider_experience_months": 6,
}])
print("Linear:", linear.predict(new_order).round(1))
print("Boosted:", boosted.predict(new_order).round(1))

Output:

Linear: [63.8]
Boosted: [66.7]

Why rebuild a **make_preprocess()** function? Each pipeline needs its own preprocessing object, because fitting one pipeline fills in the learned medians and categories. A function that builds a fresh one each time avoids accidentally sharing it between the two models.

Two models, two different answers for the same order. By the end of this guide you'll know exactly where the 2.9-minute difference comes from.

Step 2: The global view from the linear model's weights

A linear model explains itself: its weights show what it cares about. Ranked by size, distance comes first, then Old Town and prep time. But the weights also hide things, and rain seems to count for nothing.

Read the weights

names = linear.named_steps["columntransformer"].get_feature_names_out()
weights = linear.named_steps["linearregression"].coef_
weights_table = pd.Series(weights, index=names).sort_values(key=abs, ascending=False)
print(weights_table.round(2).head(8))

Output:

num__distance_km                10.06
cat__city_zone_Old Town          6.65
num__prep_time_min               5.75
cat__city_zone_Suburbs          -4.61
num__distance_x_rain             4.51
num__distance_x_peak             3.51
num__rider_experience_months    -1.55
num__is_weekend                  1.28

How to read these numbers

These aren't "minutes per km" as in Part 2, because Part 3's pipeline scales the number columns first (Part 3's theory break). So each number weight means minutes per one typical spread of that feature:

• **distance_km** +10.06: distance has a spread (standard deviation) of 3.3 km, so 10.06 ÷ 3.3 ≈ 3 minutes per km. That's Part 2's recipe again.
• **city_zone_Old Town** +6.65 and **city_zone_Suburbs** −4.61: one-hot columns are 0 or 1, so these compare zones with each other: Old Town orders are estimated about 11 minutes slower than Suburbs orders with everything else equal.
• **rider_experience_months** −1.55: more experience means slightly faster deliveries.

Scaling is what makes these weights comparable: without it, a weight per kilometre and a weight per month of experience couldn't be ranked against each other.

Where weights mislead

Look further down the full list and you'll find num__is_raining with a weight of −0.00. Does the model think rain doesn't matter? No. Rain's effect lives in the interaction columns (distance_x_rain +4.51, distance_x_rain_x_peak +1.24), because rain slows every kilometre (Part 2, Step 6). The flat is_raining column has nothing left to explain.

This is the weakness of reading weights:

• One real-world factor can be spread over several columns (rain over three; zone over four one-hot columns), so no single weight tells you how much "rain" matters.
• Weights only exist for linear models. The gradient boosting model has no weights to read at all.

We need a method that asks about the original features (weather, zone, distance), works for any model, and measures importance in minutes of error. That's permutation importance.

Step 3: Permutation importance, for any model

Shuffle one feature's values so they no longer match their orders, and see how much worse the predictions get. Shuffling distance makes the average error jump from 4.49 to 15.40 minutes; shuffling cuisine changes almost nothing.

The idea in one experiment

If a feature matters, scrambling it should hurt. Try it by hand with distance:

from sklearn.metrics import mean_absolute_error

print("Normal MAE:", round(mean_absolute_error(y_test, linear.predict(X_test)), 2))

shuffled = X_test.copy()
shuffled["distance_km"] = np.random.default_rng(0).permutation(shuffled["distance_km"].values)
print("Distance shuffled:", round(mean_absolute_error(y_test, linear.predict(shuffled)), 2))
Normal MAE: 4.49
Distance shuffled: 15.4

Shuffling keeps the same distances (same average, same spread), but gives each order someone else's. The model's error more than triples. That jump, about 10.9 minutes, is distance's importance. Notice that we shuffle the raw column, so the pipeline's cleaning, interactions and encoding all see the scrambled distance, exactly as they would in real use.

Do it properly for every feature

scikit-learn repeats this for every column, several times each, and averages the results:

from sklearn.inspection import permutation_importance

result = permutation_importance(linear, X_test, y_test, n_repeats=10, random_state=42,
                                scoring="neg_mean_absolute_error")
importance = pd.Series(result.importances_mean, index=X_test.columns).sort_values(ascending=False)
print(importance.round(2))

Output:

distance_km                10.83
prep_time_min               3.25
weather                     2.50
order_time                  1.72
city_zone                   1.47
rider_experience_months     0.29
cuisine                     0.02

Use test data, not training data. Importance on the training set tells you what the model memorised; on the test set it tells you what actually helps it predict new orders.

What this reveals that the weights didn't

• Weather is the third most important feature (2.50 minutes), even though the is_raining weight was 0.00. Permutation importance asks about the original weather column, so all of rain's effect, including through the interactions, is counted together.
• Order time matters (1.72), though the model never sees the timestamp directly: it's the peak-hour and weekend features hidden inside it.
• Cuisine is almost useless (0.02). That's not because cuisine doesn't affect delivery times. It's because prep time already carries cuisine's effect: biryani is slow because it takes long to cook. When prep time is known, cuisine adds nothing new.

Same result for the black box

Run the same code with boosted instead of linear. Gradient boosting has no weights to read, but permutation importance works just the same:

distance_km                10.72
prep_time_min               3.05
weather                     2.48
order_time                  1.71
city_zone                   1.51
rider_experience_months     0.31
cuisine                     0.08

The two models agree on the ranking and nearly on the sizes. When very different models agree on what matters, you can be more confident it reflects the data, not one model's quirks.

The limit of a global view: this tells you distance matters most on average. It can't tell you why one particular order got 63.8 minutes. That needs a local explanation.

Step 4: Explain one prediction: why 63.8 minutes?

SHAP splits a single prediction into a starting point, the average order (40.8 minutes), plus one push per feature. For the rainy Saturday biryani: rain +9.4, Old Town +7.6, Saturday dinner time +5.3. The pushes add up exactly to 63.8.

Theory break: Shapley values in one story

Three friends run a lemonade stand together and make ₹900. How should they split it fairly, when the stand only works with all three together? In the 1950s the economist Lloyd Shapley found a fair answer: imagine the friends joining in every possible order, note how much each one adds when they join, and average those additions. His idea later earned him a share of a Nobel Prize.

SHAP (SHapley Additive exPlanations) applies this to predictions:

To find out what "rain" adds, SHAP compares predictions with the order's real weather against predictions with weather values taken from a sample of ordinary past orders, the background. It does this for every combination of the other features and averages the results.

Two properties make SHAP values especially easy to trust:

• They add up. Average prediction + all the SHAP values = this order's prediction, exactly.
• They're in the target's units: minutes here, so anyone can read them.

A small helper, because our pipeline takes raw text

SHAP's built-in tools expect every input column to be a number. Our pipeline deliberately takes raw orders, typos and text included, so we give SHAP a small helper that swaps whole columns between the order and the background orders. Save this as explain_tools.py:

import pandas as pd
import shap


def make_explainer(pipeline, background):
    """A SHAP explainer that works on raw, messy orders (text columns included)."""
    cols = list(background.columns)

    def predict_raw(rows):
        return pipeline.predict(pd.DataFrame(rows, columns=cols))

    def masker(mask, x):
        # Keep the explained order's value where mask is True; use background orders elsewhere
        out = background.copy()
        for col, keep in zip(cols, mask):
            if keep:
                out[col] = x[cols.index(col)]
        return out.values

    masker.shape = (len(background), len(cols))
    return shap.Explainer(predict_raw, masker, algorithm="exact", feature_names=cols)

Explain the order

from explain_tools import make_explainer

background = X_train.sample(300, random_state=0).reset_index(drop=True)
explainer = make_explainer(linear, background)
explanation = explainer(new_order.values.astype(object))

print("Average prediction:", round(float(explanation.base_values[0]), 2))
for name, value in zip(X.columns, explanation.values[0]):
    print(f"{name:25s} {value:+.2f}")
print("Sum:", round(float(explanation.base_values[0] + explanation.values[0].sum()), 2))

Output:

Average prediction: 40.77
order_time                +5.27
city_zone                 +7.63
cuisine                   +1.01
distance_km               -2.26
prep_time_min             -0.11
weather                   +9.42
rider_experience_months   +2.04
Sum: 63.76

• **background**: 300 random training orders that define "an ordinary order". The base value of 40.77 minutes is simply the model's average prediction for them.
• **.astype(object)**: lets one array hold both text and numbers.
• **explanation.values[0]**: the 7 SHAP values for our one order.

The sum is 63.76, which rounds to the 63.8 the model predicted.

Read it as a waterfall

Start at the top with an ordinary order, 40.8 minutes, and follow each feature's push down to this order's 63.8:

• Rain (+9.4): the single biggest reason. Rain slows every kilometre.
• Old Town (+7.6): its narrow streets, the effect Part 3 found the model handles less well.
• Saturday 20:15 (+5.3): peak hour and weekend.
• Distance (−2.3), the surprise: 6 km sounds ordinary, so why does it reduce the time? Because the average background order is about 6.4 km, and in rain at rush hour every kilometre is expensive. Being 0.4 km shorter than average saves about 2 minutes under these conditions. SHAP values are always relative to the background, never absolute.
• Rider of 6 months (+2.0): a fairly new rider, so a little slower than average.
• Biryani (+1.0): a small push.
• Prep time missing (−0.1): almost nothing. The missing value was filled with the median, 17 minutes, which is close to a typical order's prep time, so it can't push the prediction anywhere. That's Part 3's imputation lesson made visible: a filled-in value carries no information, so the model effectively ignored the fact that biryani takes long to cook.

Step 5: Open the black box: SHAP for gradient boosting

The same helper explains the gradient boosting model, which has no weights at all. It predicts 66.7 minutes for the same order, 2.9 more than the linear model, and SHAP shows exactly where that gap comes from.

One order, two models

Only the pipeline changes. The background stays the same, so both explanations start from the same "ordinary order":

lin_explanation = make_explainer(linear, background)(new_order.values.astype(object))
gb_explanation = make_explainer(boosted, background)(new_order.values.astype(object))

print("Average prediction (boosted):", round(float(gb_explanation.base_values[0]), 2))
print(f"{'feature':25s} {'linear':>7s} {'boosted':>8s}")
for name, a, b in zip(X.columns, lin_explanation.values[0], gb_explanation.values[0]):
    print(f"{name:25s} {a:+7.2f} {b:+8.2f}")
print("Sum boosted:", round(float(gb_explanation.base_values[0] + gb_explanation.values[0].sum()), 2))

Output:

Average prediction (boosted): 40.76
feature                    linear  boosted
order_time                  +5.27    +5.57
city_zone                   +7.63    +7.60
cuisine                     +1.01    +2.72
distance_km                 -2.26    -3.79
prep_time_min               -0.11    +1.16
weather                     +9.42   +11.02
rider_experience_months     +2.04    +1.70
Sum boosted: 66.73

• **zip(X.columns, a, b)** walks the feature names and both models' SHAP values side by side, one row per feature.
• **{a:+7.2f}** prints a number with its sign (+ or -), 2 decimals, 7 characters wide, so the columns line up.

Where the 2.9 minutes come from

Both models start at almost the same average (40.77 vs 40.76) and agree on the big picture: rain, Old Town and Saturday night are the top three reasons in both. They disagree in four places:

The prep time row is the interesting one. In the linear model, a value of 17 minutes always has the same effect, so a typical value pushes nothing. In a tree model, the effect of one feature depends on the others: 17 minutes of prep on a rainy Saturday night can count differently from 17 minutes on a quiet Tuesday. That's an interaction, and it's exactly what makes trees more flexible, and harder to read without SHAP.

From one order to the whole model

Average the size of the SHAP values (ignoring the sign) over many orders and you get a global view, much like permutation importance in Step 3:

import time

start = time.perf_counter()
bg100 = X_train.sample(100, random_state=0).reset_index(drop=True)
sample = X_test.sample(100, random_state=1)
many = make_explainer(boosted, bg100)(sample.values.astype(object))
print(f"Explained 100 orders in {time.perf_counter() - start:.1f} s")

mean_abs = np.abs(many.values).mean(axis=0)
for name, v in sorted(zip(X.columns, mean_abs), key=lambda t: -t[1]):
    print(f"{name:25s} {v:5.2f}")

Output (your timing will differ):

Explained 100 orders in 27.3 s
distance_km               11.14
prep_time_min              4.38
weather                    3.62
order_time                 2.89
city_zone                  2.03
rider_experience_months    1.07
cuisine                    0.39

The ranking is the same as the boosted model's permutation importance in Step 3: distance, prep time, weather, order time, zone, rider, cuisine. Two different methods, one answer. That agreement is your best evidence that the explanations are reliable.

Notice how different the global and local pictures are. Across all orders, distance matters most (11.1 minutes on average). For our one rainy biryani, distance reduced the time a little and rain dominated. Both are true. That's why you need both views.

Step 6: Use explanations to check fairness

The model adds about 5 minutes to an order with a brand-new rider compared with a 4-year veteran. That's a fair estimate for customers, but it would be unfair if the same number were used to judge or assign riders.

Explanations are not only for curiosity. They let you ask, "Is the model treating a group of people in a way we're comfortable with?" Here the people are riders, and the feature is rider_experience_months.

Ask "what if?"

Keep everything about the rainy biryani the same and change only the rider:

for months in [0, 6, 24, 48]:
    what_if = new_order.assign(rider_experience_months=months)
    print(months, round(float(linear.predict(what_if)[0]), 1))

Output:

0 64.4
6 63.8
24 61.8
48 59.1

.assign(...) returns a copy of the order with one column replaced. The model says each extra year of experience saves the customer roughly 1.3 minutes on this order.

Check the errors by group

An explanation tells you what the model believes. Whether the belief is right shows up in the errors. Group the test orders by rider experience:

results = clean_orders(X_test).assign(actual=y_test.values, predicted=linear.predict(X_test))
results["error"] = results["predicted"] - results["actual"]
results["experience"] = pd.cut(results["rider_experience_months"], bins=[-1, 6, 12, 24, 1000],
                               labels=["0-6 months", "7-12 months", "13-24 months", "25+ months"])
summary = results.groupby("experience", observed=True).agg(
    orders=("error", "size"), actual=("actual", "mean"), predicted=("predicted", "mean"),
    mean_error=("error", "mean"), mae=("error", lambda e: e.abs().mean())).round(2)
print(summary)

Output:

              orders  actual  predicted  mean_error   mae
experience                                               
0-6 months        88   41.64      40.57       -1.07  3.91
7-12 months       65   40.99      41.93        0.94  4.31
13-24 months     149   43.38      43.65        0.27  4.37
25+ months       271   37.69      37.57       -0.12  4.70

The bands cover 573 of the 600 test orders. The other 27 have a missing experience value, which pd.cut leaves out of every band.

What this tells us

• Accuracy is similar for every group. The MAE stays between 3.9 and 4.7 minutes. No group is badly served.
• The newest riders get slightly optimistic estimates (−1.07): their orders really take about a minute longer than predicted. Small, but worth watching.
• Is it fair? It depends on how the number is used:To tell a customer when food will arrive: using experience is reasonable. New riders really are a little slower, and an honest estimate helps the customer.
• To rate riders, or decide who gets orders: this is risky. If the app sends fewer orders to new riders because the model predicts they're slower, they never get the experience to improve. That's a feedback loop: the model's prediction creates the data that confirms it.

The same table by city_zone repeats Part 3's finding: Old Town has the highest MAE (5.75). An explanation won't fix that, but it tells you where to collect better data.

Rule of thumb: a feature can be accurate and inappropriate. Decide what the prediction will be used for before deciding which features are fair to use.

Step 7: Add a "why" to the Part 4 API

Add **?explain=true** to the request and the API returns the three biggest reasons alongside the prediction. Without it, the API behaves exactly as in Part 4, so existing apps keep working.

Save a background file

The API needs its "ordinary orders" on disk. In the Part 4 project folder, run once:

X_train.sample(100, random_state=0).to_csv("background_orders.csv", index=False)

Copy explain_tools.py into the same folder. We use 100 background orders, not 300, to keep each explanation fast. That's why the API's average (40.5) differs slightly from Step 4's (40.8). The theory break below shows what happens if you go too small.

The updated `app.py`

The new or changed lines are marked with # NEW:

"""Delivery ML Series, Part 6: the Part 4 API, now able to explain its predictions."""
import joblib
import pandas as pd
from fastapi import FastAPI, HTTPException
from pydantic import BaseModel, Field

from delivery_features import clean_orders
from explain_tools import make_explainer  # NEW

MODEL_VERSION = "part3-linear-v1"

pipeline = joblib.load("delivery_pipeline.joblib")
background = pd.read_csv("background_orders.csv")  # NEW
explainer = make_explainer(pipeline, background)  # NEW

app = FastAPI(title="Delivery Time API", version="1.1")


class Order(BaseModel):
    """One raw order, exactly as the delivery app sends it."""
    order_time: str = Field(examples=["2026-10-03 20:15"])
    city_zone: str = Field(examples=["  old TOWN"])
    cuisine: str = Field(examples=["Biryani"])
    distance_km: float = Field(gt=0, examples=[6000])
    prep_time_min: float | None = Field(default=None, ge=0, examples=[None])
    weather: str = Field(examples=["RAINY"])
    rider_experience_months: float | None = Field(default=None, ge=0, examples=[6])


class Reason(BaseModel):  # NEW
    feature: str
    minutes: float


class Prediction(BaseModel):
    predicted_minutes: float
    model_version: str
    average_minutes: float | None = None  # NEW
    top_reasons: list[Reason] | None = None  # NEW


@app.get("/health")
def health():
    return {"status": "ok", "model_version": MODEL_VERSION}


@app.post("/predict", response_model=Prediction)
def predict(order: Order, explain: bool = False):  # NEW: explain
    row = pd.DataFrame([order.model_dump()])

    try:
        cleaned = clean_orders(row)
    except (ValueError, TypeError):
        raise HTTPException(status_code=422, detail="order_time must look like 2026-10-03 20:15")
    if cleaned["city_zone"].isna().any():
        raise HTTPException(status_code=422, detail=f"Unknown city_zone: {order.city_zone!r}")
    if cleaned["is_raining"].isna().any():
        raise HTTPException(status_code=422, detail=f"Unknown weather: {order.weather!r}")

    minutes = float(pipeline.predict(row)[0])
    result = Prediction(predicted_minutes=round(minutes, 1), model_version=MODEL_VERSION)

    if explain:  # NEW: the whole block below
        shap_values = explainer(row[background.columns].values.astype(object))
        contributions = dict(zip(background.columns, shap_values.values[0]))
        top = sorted(contributions.items(), key=lambda item: abs(item[1]), reverse=True)[:3]
        result.average_minutes = round(float(shap_values.base_values[0]), 1)
        result.top_reasons = [Reason(feature=name, minutes=round(float(value), 1)) for name, value in top]

    return result

The validation checks run before the explanation, so a bad order is still rejected with a clear 422 and never reaches SHAP.

Try it

Start the server as in Part 4 (uvicorn app:app --reload), then send the messy order to /predict?explain=true. The response:

{
  "predicted_minutes": 63.8,
  "model_version": "part3-linear-v1",
  "average_minutes": 40.5,
  "top_reasons": [
    {"feature": "weather", "minutes": 9.5},
    {"feature": "city_zone", "minutes": 7.6},
    {"feature": "order_time", "minutes": 5.6}
  ]
}

The delivery app can turn that into a friendly message: "About 64 minutes. It's raining, Old Town is busy, and it's Saturday dinner time."

Test it

Add these to test_app.py:

def test_predict_without_explanation_is_unchanged():
    body = client.post("/predict", json=MESSY_ORDER).json()
    assert body["predicted_minutes"] == 63.8
    assert body["top_reasons"] is None


def test_explanation_lists_three_reasons():
    body = client.post("/predict?explain=true", json=MESSY_ORDER).json()
    assert len(body["top_reasons"]) == 3
    assert body["top_reasons"][0]["feature"] == "weather"


def test_unknown_zone_is_still_rejected():
    response = client.post("/predict?explain=true", json={**MESSY_ORDER, "city_zone": "Old-Town"})
    assert response.status_code == 422

MESSY_ORDER is the rainy biryani as a dictionary, the same one used in Part 4's tests. Run pytest: Part 4's six tests plus these three give 9 passed. The first test matters most: it proves that adding explanations didn't change any existing answer.

The cost of "why"

On our laptop, a normal prediction took about 20 ms and an explained one about 250–280 ms: more than ten times slower, because SHAP calls the model 128 × 100 times. That's why explain is off by default. Ask for it when a person will actually read the answer (a customer tapping "Why so long?", a support agent), not on every request.

Theory break: how to judge an explanation

An explanation is itself a result, and results can be wrong. Before you show one to anyone, check five things: does it add up, is it stable, do different methods agree, are features sharing credit, and are you reading it as cause and effect when it isn't?

In Parts 1–3 we learned never to trust a model without measuring it. The same goes for explanations.

1. Does it add up? (faithfulness)

SHAP values must add up to the prediction. If they don't, something is wrong: the wrong model, the wrong columns, a mismatched background. One line catches it:

total = explanation.base_values[0] + explanation.values[0].sum()
assert abs(total - linear.predict(new_order)[0]) < 0.01, "SHAP values must add up to the prediction"

This is Part 3's assertion habit applied to explanations.

2. Is it stable?

SHAP compares the order with a sample of background orders, so a different sample gives slightly different values. We explained the rainy biryani five times with five different random backgrounds of each size:

for size in [50, 100, 300, 1000]:
    values = np.array([
        make_explainer(linear, X_train.sample(size, random_state=seed).reset_index(drop=True))(
            new_order.values.astype(object)).values[0]
        for seed in range(5)
    ])
    print(size, "weather", values[:, 5].min().round(1), values[:, 5].max().round(1),
          "distance", values[:, 3].min().round(1), values[:, 3].max().round(1))

(values[:, 5] is the weather column of all five runs, values[:, 3] the distance column.) The results, with timings from our laptop:

Three lessons:

• Big reasons are stable; small ones are noisy. Rain is +8 to +10 whatever you do. Distance's push is so small for this order that with 50 background orders it even flips sign. So read Step 4's "distance −2.3" as "distance barely mattered for this order", not as a precise number.
• Size is a trade-off between speed and stability. The API uses 100 because it only reports the top 3 reasons, and those were the same in every run. For a report that quotes exact values, use 300 or more.
• Never rely on a background of 50 or fewer. With 50, one run in five even changed the top-3 list.

3. Do different methods agree?

If SHAP said rain mattered most overall while permutation importance said distance, you'd have a problem. Ours agree: Step 3 (shuffling) and Step 5 (average SHAP size) rank the seven features in the same order. Agreement between methods that work differently is strong evidence.

4. Are features sharing credit?

When two features carry the same information, the model can lean on either one, and the explanation splits the credit between them, sometimes unevenly. Imagine adding "estimated driving minutes" next to distance_km: the two would share distance's importance, and each would look less important than distance really is.

We avoided this by explaining the raw columns. Our pipeline creates distance_x_rain and friends inside the cleaning step, so all of the distance credit lands on distance_km and all of the rain credit on weather. If you add overlapping features, check their combined importance, not each one alone.

5. Is it the model, or the world?

The explanation says the model adds 2 minutes for a rider of 6 months. It does not prove that inexperience causes slow delivery. Maybe new riders are given the hardest shifts, or the worst zones. SHAP describes what the model learned from patterns in the data. Proving cause and effect needs a controlled experiment, such as an A/B test, which is beyond this series.

The complete code

One script, **explain_delivery.py**, rebuilds both models and prints every explanation from this tutorial: permutation importance for both models, SHAP for the rainy biryani with an add-up check, and the fairness table. It runs in about 5 seconds.

Put it in the same folder as delivery_data.py and delivery_features.py (from Part 3) and explain_tools.py (from Step 4):

"""Delivery ML Series, Part 6: explain the delivery-time models."""
import numpy as np
import pandas as pd
from sklearn.compose import ColumnTransformer
from sklearn.ensemble import GradientBoostingRegressor
from sklearn.impute import SimpleImputer
from sklearn.inspection import permutation_importance
from sklearn.linear_model import LinearRegression
from sklearn.model_selection import train_test_split
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import FunctionTransformer, OneHotEncoder, StandardScaler

from delivery_data import make_messy_orders
from delivery_features import clean_orders
from explain_tools import make_explainer

# 1. Rebuild Part 3's data and both models
orders = make_messy_orders().drop_duplicates()
X = orders.drop(columns="delivery_min")
y = orders["delivery_min"]
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

numeric_features = ["distance_km", "prep_time_min", "rider_experience_months", "is_raining",
                    "is_peak_hour", "is_weekend", "distance_x_rain", "distance_x_peak",
                    "distance_x_rain_x_peak"]
categorical_features = ["city_zone", "cuisine"]


def make_preprocess():
    return ColumnTransformer([
        ("num", make_pipeline(SimpleImputer(strategy="median"), StandardScaler()), numeric_features),
        ("cat", make_pipeline(SimpleImputer(strategy="most_frequent"),
                              OneHotEncoder(handle_unknown="ignore")), categorical_features),
    ])


linear = make_pipeline(FunctionTransformer(clean_orders), make_preprocess(), LinearRegression())
linear.fit(X_train, y_train)
boosted = make_pipeline(FunctionTransformer(clean_orders), make_preprocess(),
                        GradientBoostingRegressor(random_state=42))
boosted.fit(X_train, y_train)

new_order = pd.DataFrame([{
    "order_time": "2026-10-03 20:15", "city_zone": "  old TOWN", "cuisine": "Biryani",
    "distance_km": 6000, "prep_time_min": None, "weather": "RAINY", "rider_experience_months": 6,
}])

# 2. Global: permutation importance (works for any model)
for name, model in [("linear", linear), ("boosted", boosted)]:
    result = permutation_importance(model, X_test, y_test, n_repeats=10, random_state=42,
                                    scoring="neg_mean_absolute_error")
    importance = pd.Series(result.importances_mean, index=X_test.columns).sort_values(ascending=False)
    print(f"\nPermutation importance ({name}, minutes of extra error):")
    print(importance.round(2).to_string())

# 3. Local: SHAP for one order, both models
background = X_train.sample(300, random_state=0).reset_index(drop=True)
for name, model in [("linear", linear), ("boosted", boosted)]:
    explanation = make_explainer(model, background)(new_order.values.astype(object))
    total = explanation.base_values[0] + explanation.values[0].sum()
    assert abs(total - model.predict(new_order)[0]) < 0.01, "SHAP values must add up to the prediction"
    print(f"\nWhy {total:.1f} minutes? ({name}, average order = {explanation.base_values[0]:.1f})")
    for feature, value in sorted(zip(X.columns, explanation.values[0]), key=lambda t: -abs(t[1])):
        print(f"  {feature:25s} {value:+.2f}")

# 4. Fairness: errors by rider experience
results = clean_orders(X_test).assign(actual=y_test.values, predicted=linear.predict(X_test))
results["error"] = results["predicted"] - results["actual"]
results["experience"] = pd.cut(results["rider_experience_months"], bins=[-1, 6, 12, 24, 1000],
                               labels=["0-6 months", "7-12 months", "13-24 months", "25+ months"])
print("\nErrors by rider experience (linear):")
print(results.groupby("experience", observed=True)["error"]
      .agg(orders="size", mean_error="mean", mae=lambda e: e.abs().mean()).round(2))

Run python explain_delivery.py. The numbers match Steps 3 to 6; the SHAP section prints each model's reasons sorted from biggest to smallest. (SHAP also shows a progress bar; you can ignore it.)

Common mistakes

Practice project: explain the electricity model

Explain the home-electricity model from the earlier practice projects: find which inputs matter most, explain why the May home in "Bombay" gets 636 kWh, and check whether the model serves homes with solar panels as well as homes without.

In Parts 1–4 you predicted monthly electricity use (kWh) for homes in Mumbai, Bengaluru, Chennai and Pune, cleaned the messy smart-meter export in Part 3, and served the model as an API in Part 4. Now explain it.

What you start with

In a new folder, copy:

• energy_data_messy.py (Part 3 practice) and energy_features.py with your clean_homes function (Part 4 practice);
• explain_tools.py from Step 4 of this tutorial. It works unchanged: nothing in it is specific to deliveries.

Rebuild the model exactly as in the Part 4 practice: drop duplicates, divide Wh readings (≥ 20,000) by 1,000, train_test_split(test_size=0.2, random_state=42), then the cleaning + ColumnTransformer + LinearRegression pipeline. Your test set should have 400 homes.

The home to explain:

home = pd.DataFrame([{"reading_date": "2026-05-20", "city": "Bombay", "home_size_sqm": 120,
                      "residents": 4, "solar_panels": "TRUE", "works_from_home": 1,
                      "appliance_age_years": None, "ac_hours_per_day": 8}])

Your tasks

• 1. Global importance (Step 3). Run permutation_importance on the test set with scoring="neg_mean_absolute_error", n_repeats=10, random_state=42. Which input matters most? Which two inputs matter almost not at all?
• 2. Explain one home (Step 4). Use a background of 300 training homes (random_state=0). What's the base value? What are the two biggest pushes, and in which direction? Check that everything adds up to the prediction.
• 3. Think about it (theory break). The reading is from May, which clean_homes marks as summer, yet reading_date gets almost no credit. Why? Hint: what else in the data already tells the model it's summer?
• 4. Check stability. Repeat task 2 with random_state 0 to 4. Do the top two reasons stay the same?
• 5. Fairness (Step 6). Group the test errors by city and by has_solar (after clean_homes). Is any group served noticeably worse?
• Stretch: add ?explain=true to your Part 4 energy_app.py and return the top 3 reasons, with a test that the prediction without explain is unchanged.

Expected results

Your numbers should match these to the first decimal (timings will differ):

Check your thinking

• Task 3: in this data, homes run their AC much longer in summer. ac_hours_per_day already carries the summer information, so once the model knows the AC hours, the date adds nothing. That's theory-break point 4, features sharing credit, in action. City is near zero for the same reason: hot cities show up as more AC hours.
• Task 5: solar panels cut a home's bill by a percentage (35% in this data), but a linear model can only subtract a fixed number of kWh. So it under-credits solar for big homes and over-credits it for small ones. The explanation points you to the fix: an interaction such as has_solar × home_size_sqm, or a model that learns interactions on its own (Part 7).

Frequently asked questions

What is the difference between feature importance and SHAP?

Feature importance (such as permutation importance) is global: it says how much a feature matters across all predictions. SHAP is local: it splits one prediction into one contribution per feature. Averaging the size of SHAP values over many predictions gives a global view too.

Do I need SHAP if my model is linear?

Not always. A linear model's scaled weights already show global importance. But SHAP still helps when your pipeline cleans and encodes the data: it explains the original inputs ("weather", "city_zone") in minutes, instead of 19 encoded columns.

Why does SHAP need background data?

A SHAP value answers, "How different is this prediction because of this feature, compared with an ordinary input?" The background defines "ordinary". Use 100 or more rows from the training data; 300 or more if you quote exact values.

Is SHAP slow?

The exact method is fine for a handful of features: about 0.2–0.3 seconds per explanation here. With dozens of features, SHAP's approximate methods (such as its permutation or tree explainers) are much faster.

Can explanations prove the model is fair?

No. They show what the model relies on and help you find groups that are served worse. Whether that's acceptable depends on how the prediction is used, which is a human decision.

Does a big SHAP value mean the feature causes the outcome?

No. It means the model uses that feature. Cause and effect needs an experiment.

Glossary

Series roadmap

Next up, Part 7: the gradient boosting model got a lot of attention here: it leaned harder on rain and biryani, and found interactions the linear model can't. Part 7 takes boosting further with XGBoost and LightGBM, and you'll use this part's SHAP tools to see what the stronger models learn. In Part 8, the API's top_reasons become the input to a generative AI model that writes the customer message for you.

Conclusion

You can now answer the question every model eventually faces: why?

• Weights explain a linear model, once the features are scaled.
• Permutation importance explains any model globally, by breaking one feature at a time.
• SHAP explains a single prediction as a fair split that adds up exactly: 40.8 minutes for an average order, plus 9.4 for rain, 7.6 for Old Town and 5.3 for Saturday dinner time, and so on, to 63.8.
• Explanations need checking like any result: do they add up, are they stable, do methods agree?
• Explanations reveal fairness questions but don't answer them: a feature can be accurate and still inappropriate for some uses.

Your API can now say not just "about 64 minutes" but "because it's raining, Old Town is busy and it's Saturday night". That's the difference between a number people tolerate and one they trust.