🧭 Concept Spatial machine learning Spatial Machine Learning · step 12 of 14 Updated

Extrapolation in Space Explained: The Area of Applicability

Problem statement

A prediction map covers every pixel. Some of those pixels resemble the training data closely; others do not resemble it at all, and the model is guessing there with no way to say so.

Measured on a real problem, a model trained on the western half of a study area and applied to the eastern half:

random cross-validation on the training half   +0.649
block cross-validation on the training half    +0.385
actual R² on the eastern half                  -0.771

Both internal estimates were optimistic by more than a unit of R², and the true result was worse than predicting the training mean.

No cross-validation scheme fixes this, because the training data contains no information about the new area. What you can do is say where the model is extrapolating, which is what an area of applicability is.

Quick answer

Measure distance in scaled feature space from each prediction location to the training data, and threshold it against the distances within the training set:

import numpy as np
from scipy.spatial import cKDTree
from sklearn.preprocessing import StandardScaler


def applicability(X_train, X_new, quantile=0.95, k=5):
    scaler = StandardScaler().fit(X_train)
    train = scaler.transform(X_train)
    new = scaler.transform(X_new)

    tree = cKDTree(train)
    within, _ = tree.query(train, k=k + 1)        # k+1 because self is nearest
    threshold = np.quantile(within[:, 1:].mean(axis=1), quantile)

    outside, _ = tree.query(new, k=k)
    distance = outside.mean(axis=1)
    return distance <= threshold, distance, threshold

Ship the resulting mask with the prediction. It does not estimate the error outside; it says where the model is guessing.

Prediction locations inside and outside the training data's feature-space envelope, with the outside region flagged rather than assigned an error.
The mask is about feature space, not geography. A nearby pixel can be outside it and a distant one inside.

Step-by-step solution

1. Distinguish geographic distance from feature-space distance

A pixel a hundred kilometres away with the same slope, elevation and aspect is well inside the model's competence. A pixel next door on a different geology may be far outside it.

Applicability is about the covariates, not the coordinates. That is why the mask often looks nothing like a buffer around the training points.

2. Scale the features before measuring distance

Elevation in metres ranges over a thousand; northness ranges over two. Unscaled, elevation dominates the distance entirely.

StandardScaler puts each feature on a comparable footing. Weighting by feature importance is a refinement: features the model barely uses should contribute less to the distance.

3. Derive the threshold from the training data

The useful comparison is against how far apart the training points themselves are. If a new location is further from the training data than 95% of training points are from their own neighbours, it is outside the envelope the model learned.

That makes the threshold self-calibrating: a dense training set gives a tight envelope, a sparse one a generous one.

4. Check the univariate ranges too

Multivariate distance can hide a single feature far outside its range, because the other features compensate.

lo, hi = np.percentile(X_train[:, i], [1, 99])
outside = ((X_new[:, i] < lo) | (X_new[:, i] > hi)).mean()

A feature with 30% of its new values outside the training band is extrapolation regardless of what the distance says.

5. Ship the mask, not a corrected prediction

There is no defensible correction for extrapolation. The honest output is two rasters: the prediction, and a mask of where it is supported.

A nearby pixel outside the feature-space envelope and a distant pixel inside it, showing that applicability is not a buffer.
Applicability follows the covariates. Distance on the map is a poor proxy for it.

Code examples

Example 1 — an importance-weighted applicability mask 

import numpy as np
from scipy.spatial import cKDTree
from sklearn.preprocessing import StandardScaler


def applicability_mask(X_train, X_new, importances=None, quantile=0.95, k=5):
    """Feature-space distance, weighted by how much the model uses each feature."""
    scaler = StandardScaler().fit(X_train)
    train = scaler.transform(X_train)
    new = scaler.transform(X_new)

    if importances is not None:
        weights = np.asarray(importances, dtype=float)
        weights = weights / weights.sum() * len(weights)
        train = train * np.sqrt(weights)
        new = new * np.sqrt(weights)

    tree = cKDTree(train)
    within, _ = tree.query(train, k=min(k + 1, len(train)))
    reference = within[:, 1:].mean(axis=1)
    threshold = float(np.quantile(reference, quantile))

    distances, _ = tree.query(new, k=min(k, len(train)))
    mean_distance = distances.mean(axis=1)
    inside = mean_distance <= threshold

    print(f"  training neighbour distance: median {np.median(reference):.3f}, "
          f"p{quantile * 100:.0f} {threshold:.3f}")
    print(f"  deployment distance: median {np.median(mean_distance):.3f}, "
          f"p95 {np.quantile(mean_distance, 0.95):.3f}")
    print(f"  inside the area of applicability: {inside.mean():.1%}")
    return inside, mean_distance, threshold

Weighting by importance is the refinement that matters most. Without it, a feature the model barely uses contributes as much to the distance as the one it depends on, and the mask flags places the model would actually handle fine.

Example 2 — a per-feature extrapolation report 

import numpy as np


def extrapolation_report(X_train, X_new, feature_names, band=(1, 99)):
    """Which features are being extrapolated, and by how much."""
    rows = []
    for i, name in enumerate(feature_names):
        train, new = X_train[:, i], X_new[:, i]
        lo, hi = np.percentile(train, band)
        below = float((new < lo).mean())
        above = float((new > hi).mean())
        span = hi - lo
        reach = float(max((lo - new.min()) / span if new.min() < lo else 0,
                          (new.max() - hi) / span if new.max() > hi else 0))

        rows.append({"feature": name, "below": below, "above": above,
                     "reach": reach})
        flag = "  <- extrapolating" if below + above > 0.1 else ""
        print(f"  {name:14} below {below:6.1%}  above {above:6.1%}  "
              f"reach {reach:5.2f} training ranges{flag}")

    worst = max(rows, key=lambda r: r["below"] + r["above"])
    print(f"  worst: {worst['feature']} "
          f"({worst['below'] + worst['above']:.1%} outside)")
    return rows

The reach column — how far beyond the training range the new data goes, in units of the training range — is the number that predicts how badly a model will behave. A reach of 0.1 is a mild extension; a reach of 2 means the model has never seen anything like it.

Example 3 — writing the prediction and the mask together 

import numpy as np
import rasterio


def write_prediction(prediction, inside, distance, profile, out_path):
    """Three bands: the prediction, the mask, and the distance behind it."""
    profile = dict(profile) | {
        "count": 3, "dtype": "float32", "nodata": np.nan,
        "compress": "deflate", "tiled": True,
    }
    with rasterio.open(out_path, "w", **profile) as dst:
        dst.write(prediction.astype("float32"), 1)
        dst.write(inside.astype("float32"), 2)
        dst.write(distance.astype("float32"), 3)
        dst.set_band_description(1, "prediction")
        dst.set_band_description(2, "1 = inside the area of applicability")
        dst.set_band_description(3, "feature-space distance to training data")
        dst.update_tags(
            applicability_note="Band 1 is only supported where band 2 is 1. "
                               "Outside that, the model is extrapolating and "
                               "the error is unknown.")

    print(f"  {inside.mean():.1%} of pixels are inside the area of "
          f"applicability")
    print(f"  wrote {out_path} with the mask as band 2")
    return out_path

Writing the note into the file's tags is worth the line. The mask travels with the prediction, but the interpretation of the mask usually does not, and a raster tag survives every copy.

Explanation

Why cross-validation cannot estimate extrapolation error

Cross-validation resamples the training data. Every fold is drawn from the same distribution, so the estimate describes generalisation within that distribution.

A new region is a different distribution. The training data contains no observations of it, so no resampling of that data can say how the model behaves there.

The measured case is unambiguous: block CV — the scheme specifically recommended for spatial data — promised +0.385 where the actual transfer gave −0.771. It was less optimistic than random CV and still wrong by more than a unit of R².

Why feature space, not geographic space

The model does not see coordinates unless you gave it coordinates. It sees slope, elevation, reflectance, distance to water.

So its competence extends to combinations of those it has seen, wherever they occur. A distant valley with familiar terrain is inside; an adjacent quarry with an unfamiliar combination is outside.

That is why an applicability mask often looks patchy rather than like a buffer, and why it is more useful than a distance-from-training-points map.

Why the threshold comes from the training data

An absolute distance threshold means nothing, because feature-space distances depend on the features, their scaling and their number.

Comparing against the within-training neighbour distances makes it relative and self-calibrating. If the training set is dense in feature space, the envelope is tight. If it is sparse, the envelope is generous, correctly reflecting that the model was fitted on scattered evidence anyway.

The 95th percentile is a convention. Reporting the mask at two or three thresholds is more honest than picking one.

Why tree models make this more important

A random forest predicts the mean of training targets in a leaf, so it cannot produce a value outside the training target range. Beyond the training feature space it reuses its most extreme splits and saturates.

The result is a smooth, in-range, plausible map that contains no information. It does not look wrong, which is exactly the problem.

Measured, a model transferred to a new region with coordinates included scored −0.095 — close to predicting the mean — because it had been broken into near-constancy. Without an applicability mask, that map would be indistinguishable from a good one.

An applicability threshold derived by scaling features, weighting by importance and taking the 95th percentile of within-training distances.
Comparing against the training set itself is what makes the envelope adapt to how densely it was sampled.

Edge cases or notes

  • Scale the features before measuring distance, or the largest-range feature dominates.
  • Weight by importance; features the model ignores should not drive the mask.
  • Check univariate ranges too — multivariate distance can hide one extreme feature.
  • Derive the threshold from within-training distances, not from an absolute number.
  • Report the mask at several thresholds rather than one.
  • The mask is not an error estimate. It says where the error is unknown.
  • Tree models saturate outside their range, producing plausible flat maps.
  • Write the interpretation into the file tags, not just the documentation.

FAQ

What is the area of applicability?

The part of a prediction map where the covariates resemble the training data closely enough that the model is interpolating rather than guessing.

How is it different from a buffer around the training points?

It is measured in feature space, not geographic space. A distant location with familiar terrain is inside; an adjacent location with an unfamiliar combination is outside.

Can I estimate the error outside the area?

No. The training data contains no information about it. The mask says where the error is unknown, which is the honest available statement.

Why not just use cross-validation?

Because it resamples the training distribution. Block CV promised +0.385 where a real transfer test gave −0.771.

What threshold should I use?

The 95th percentile of within-training neighbour distances is a reasonable convention. Reporting two or three thresholds is more informative than one.

Should I weight by feature importance?

Yes. Otherwise a feature the model barely uses contributes as much to the distance as the one it depends on.

What do I ship with the prediction?

The prediction, the mask, and the feature-space distance, with the interpretation written into the file's tags.