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

Spatial Leakage Explained: Why Random Cross-Validation Lies

Problem statement

The standard advice about machine learning on spatial data is: random cross-validation is optimistic because nearby points are correlated, so use spatial block cross-validation instead.

That advice is half right, and the missing half causes as much damage as the problem it solves.

Measured on a real modelling problem — predicting vegetation index from terrain over a 10 km window, 3,000 training cells, repeated over ten random samples:

                    mean R²    sd
random 5-fold CV     +0.555   0.019
spatial block CV     +0.213   0.106
true error over the
whole unsampled map  +0.557   0.014

Random CV was accurate to within 0.002. Block CV understated the model by 0.34 R² and was five times more variable. Following the standard advice here would have led you to reject a model that works.

Change the deployment scenario, though, and the picture inverts completely — see below. Neither scheme is "the honest one". The scheme has to match what you will actually do with the model.

Quick answer

Ask what the model will be asked to predict, and validate that:

you will predict… validate with why
unsampled locations inside the sampled area random k-fold matches deployment; block CV is pessimistic
a new region you have not sampled leave-one-region-out, held out for real both random and block CV are optimistic
a future time forward-chaining split by date space is not the leaking dimension 
from sklearn.model_selection import KFold, GroupKFold
from sklearn.cluster import KMeans

# interpolation within the sampled area
folds = KFold(n_splits=5, shuffle=True, random_state=0)

# extrapolation to unsampled ground
blocks = KMeans(n_clusters=5, random_state=0, n_init=10).fit_predict(coords)
folds = GroupKFold(n_splits=5)          # groups=blocks
Random cross-validation matching interpolation within a sampled area and block cross-validation matching extrapolation to a new region, with each being wrong for the other case.
The question is not which scheme is honest. It is which scheme resembles what you will do next.

Step-by-step solution

1. Recognise what leakage actually is

Leakage is information reaching the test set that will not be available at prediction time.

With spatially autocorrelated data and a random split, a test point often has a training point twenty metres away. The model can effectively memorise the neighbourhood.

That is only leakage if your real predictions will be far from any training point. If you will predict at locations surrounded by training data — filling gaps in a survey, densifying a monitoring network — then having a near neighbour is not leakage. It is the deployment condition.

2. Measure the interpolation case

Training on 3,000 randomly placed cells and predicting the remaining 256,776:

random 5-fold CV   +0.555 (sd 0.019 over 10 samples)
block CV           +0.213 (sd 0.106)
true map R²        +0.557 (sd 0.014)

Random CV reproduced the truth almost exactly. Block CV was pessimistic by 0.34 and unstable: across ten runs it ranged from +0.061 to +0.429, while the truth barely moved.

The instability has a mechanical cause. With five blocks, each fold trains on 80% of the area and predicts a large contiguous region, and how hard that is depends heavily on which region the clustering happened to carve out.

3. Measure the extrapolation case

Now train on the west half only and predict the east half — a genuinely unsampled region:

random 5-fold CV on the training half   +0.649
block CV on the training half           +0.385
actual R² on the eastern half           -0.771

Both lie, in the same direction. Random CV promised +0.649; block CV, the supposedly honest scheme, still promised +0.385. Reality was −0.771.

A negative R² means the model is worse than predicting the training mean everywhere. The terrain-to-vegetation relationship in the western, coastal half simply does not hold in the eastern, mountainous half — and no resampling of the western data could reveal that.

4. Draw the actual conclusion

Block CV is not a truth serum. It is a different estimator with a different bias:

  • For interpolation, it is strongly pessimistic (0.213 against 0.557) and noisy.
  • For extrapolation, it is less optimistic than random CV but still badly optimistic (0.385 against −0.771).

The only way to estimate transfer to a new region is to hold out a region of the kind you will actually predict, and accept that you cannot do this at all if you have only sampled one region.

5. Report the scheme with the number

R² = 0.55 is not a result. R² = 0.55 under random 5-fold CV, for prediction within the sampled extent is.

Random CV accurate for interpolation and optimistic for extrapolation, block CV pessimistic for interpolation and still optimistic for extrapolation.
Both schemes are biased. The direction and size of the bias depend on which task you are actually doing.

Code examples

Example 1 — score the model both ways and against the truth 

import numpy as np
from sklearn.cluster import KMeans
from sklearn.ensemble import RandomForestRegressor
from sklearn.metrics import r2_score
from sklearn.model_selection import GroupKFold, KFold


def compare_schemes(X, y, coords, X_holdout=None, y_holdout=None,
                    n_blocks=5, seed=0):
    """Random CV, block CV, and — if you have it — the truth."""
    def model():
        return RandomForestRegressor(n_estimators=200, min_samples_leaf=2,
                                     n_jobs=-1, random_state=seed)

    def cv(splitter):
        predicted = np.empty(len(y))
        for train, test in splitter:
            predicted[test] = model().fit(X[train], y[train]).predict(X[test])
        return r2_score(y, predicted)

    random_r2 = cv(KFold(n_splits=5, shuffle=True, random_state=seed).split(X))

    blocks = KMeans(n_clusters=n_blocks, random_state=seed,
                    n_init=10).fit_predict(coords)
    block_r2 = cv(GroupKFold(n_splits=n_blocks).split(X, y, blocks))

    print(f"  random 5-fold CV   R2 {random_r2:+.3f}")
    print(f"  {n_blocks}-block spatial CV R2 {block_r2:+.3f}")

    result = {"random": random_r2, "block": block_r2}
    if X_holdout is not None:
        true_r2 = r2_score(y_holdout, model().fit(X, y).predict(X_holdout))
        print(f"  true held-out      R2 {true_r2:+.3f}")
        result["true"] = true_r2
    return result

Where you can construct a genuine hold-out — a dense reference dataset, a pilot area, a later survey — do it once. It calibrates how much to trust each CV scheme for your variable and landscape, and that calibration transfers to future projects far better than any rule of thumb.

Example 2 — check how stable your block CV is 

import numpy as np
from sklearn.cluster import KMeans


def block_cv_stability(X, y, coords, n_repeats=10, n_blocks=5):
    """Block CV depends on which blocks the clustering happened to make."""
    scores = []
    for seed in range(n_repeats):
        blocks = KMeans(n_clusters=n_blocks, random_state=seed,
                        n_init=10).fit_predict(coords)
        scores.append(cv_with_groups(X, y, blocks))

    scores = np.array(scores)
    print(f"  block CV over {n_repeats} block layouts: "
          f"mean {scores.mean():+.3f}, sd {scores.std():.3f}, "
          f"range {scores.min():+.3f} to {scores.max():+.3f}")
    if scores.std() > 0.05:
        print("  ! this estimate is unstable — report the spread, not one number")
    return scores
  block CV over 10 block layouts: mean +0.213, sd 0.106, range +0.061 to +0.429
  ! this estimate is unstable — report the spread, not one number

A single block-CV number carries a standard deviation of 0.106 here. Quoting it to three decimal places, as everyone does, is false precision by a wide margin.

Example 3 — the only honest extrapolation test 

import numpy as np
from sklearn.metrics import r2_score


def region_transfer_test(X, y, coords, split_axis=0):
    """Train on one half of the study area, predict the other. No resampling."""
    threshold = np.median(coords[:, split_axis])
    left = coords[:, split_axis] < threshold

    model = RandomForestRegressor(n_estimators=200, min_samples_leaf=2,
                                  n_jobs=-1, random_state=0)
    forward = r2_score(y[~left], model.fit(X[left], y[left]).predict(X[~left]))
    backward = r2_score(y[left], model.fit(X[~left], y[~left]).predict(X[left]))

    print(f"  train left  -> predict right: R2 {forward:+.3f}")
    print(f"  train right -> predict left : R2 {backward:+.3f}")
    if min(forward, backward) < 0:
        print("  ! the model does not transfer — it is worse than the mean "
              "in at least one direction")
    return {"forward": forward, "backward": backward}

Run both directions. A model that transfers one way and not the other tells you the two halves differ, which is more useful than a single number and is often the actual scientific finding.

Explanation

Why random CV was accurate here

The 3,000 training cells were spread randomly over the whole window, and the model was then asked to predict the other 256,776 cells in the same window.

A randomly held-out cell is, in expectation, exactly as far from the training set as a randomly chosen unpredicted cell. The CV task and the deployment task are the same task, so the CV estimate is unbiased — which is precisely what the measurement shows: +0.555 against +0.557.

Why block CV was pessimistic

Block CV asks a harder question than deployment does. Each fold trains on 80% of the area and predicts a large contiguous region with no nearby training data.

If you will never do that, the number is answering a question you did not ask. Here it understated the model by 0.34 R² — enough to make a working model look like a failure.

It is also noisy, for a structural reason: with five folds, each test set is one large region, so the score is an average of five hard, correlated problems rather than of many independent ones. Hence sd 0.106 against random CV's 0.019.

Why both schemes failed on the new region

The west-to-east test is the case everyone is worried about, and it is worth noting that block CV did not save it. It promised +0.385 where reality was −0.771.

The reason is that block CV resamples within the training region. If the relationship between terrain and vegetation is different in the east — different geology, different land use, different rainfall — nothing inside the western data contains that information. No cross-validation scheme can estimate an error it has no evidence about.

This is the practical limit of validation. Cross-validation measures how well a model generalises to data drawn from the same distribution as the training data. A new region is a different distribution, and the honest answer is "we do not know", supported by an area-of-applicability mask rather than a number.

Why a negative R² is worth stating plainly

R² = −0.771 means the model's squared error is 1.77 times the variance of the target in the new region — it would have been better to predict the training mean everywhere and not build a model at all.

Negative R² is common in genuine transfer tests and almost never reported, because most pipelines never run one. When you see it, the model has not merely degraded; it has learned a relationship that is actively wrong somewhere else.

Three deployment tasks and whether a near training point counts as leakage for each.
Nothing about the split changes. What changes is the question the score is standing in for.

Edge cases or notes

  • State the scheme with every score. A bare R² is not interpretable.
  • Block CV has a large variance. Repeat over several block layouts and report the spread.
  • The block size matters as much as the scheme. Blocks should be at least the autocorrelation range, or the folds still leak.
  • Block CV scores are not comparable across sample sizes, because the block geometry changes with the sample.
  • Autocorrelation is not the only leaking dimension. Time, sensor, operator and processing batch all leak.
  • Duplicate or near-duplicate rows — two samples from one field — leak under any scheme.
  • A negative R² is meaningful, not a bug. It means the model is worse than the mean.
  • Do not tune hyperparameters on the same folds you report. That leaks in the ordinary, non-spatial way.

FAQ

Is random cross-validation always wrong for spatial data?

No. When the model will predict at unsampled locations inside a well-sampled area, random CV matched the truth almost exactly here — +0.555 against +0.557.

When should I use spatial block cross-validation?

When you will predict in areas with no nearby training data. For interpolation within a sampled area it was pessimistic by 0.34 R² and five times noisier.

Does block cross-validation give the true error?

No. On a genuine new-region test it promised +0.385 where the actual result was −0.771. It is less optimistic than random CV, not honest.

How do I estimate error in a region I have not sampled?

Hold out a real region of the same kind and test on it. If you have only sampled one region, you cannot estimate it — use an area-of-applicability mask and say so.

Why is my block CV score so variable?

Because five folds means five large, correlated test regions, and the score depends on which regions the clustering produced. It had a standard deviation of 0.106 across ten layouts here.

What does a negative R² mean?

The model is worse than predicting the mean. It happens routinely in genuine transfer tests and is worth reporting rather than hiding.

How big should the blocks be?

At least the autocorrelation range of the target, or neighbouring folds still share information. Estimate it from a variogram rather than guessing.