Splitting Chips into Spatially Disjoint Folds

To split chips so validation is honest, derive a block id from each chip’s map coordinates and hand that to a group-aware splitter:

import numpy as np
from sklearn.model_selection import GroupKFold

BLOCK_M = 20_000                                   # 20 km blocks

block = (np.floor(chip_x / BLOCK_M).astype(int) * 100_000
         + np.floor(chip_y / BLOCK_M).astype(int))  # one id per block

for train_idx, val_idx in GroupKFold(n_splits=5).split(X, y, groups=block):
    ...

The whole point is that GroupKFold never puts the same group on both sides of a split, so a block — and everything spatially close within it — moves as a unit. This page belongs to building training datasets from satellite imagery, part of Raster Machine Learning & Model Inference.

Why Random Splits Fail on Raster Data

Spatial data violates the independence assumption that random splitting relies on. Two pixels 30 m apart are nearly the same observation: same field, same soil, same weather, same atmospheric state on the same day. A random split puts one in training and the other in validation, and the model scores well by recognising a neighbourhood it has already seen.

Reported accuracy against block size With a random split the reported accuracy is about 0.94. As the block size grows from one kilometre to twenty the reported accuracy falls steeply, then flattens near 0.79 beyond about fifteen kilometres. The plateau is the honest estimate; the distance at which it flattens is the autocorrelation range of the target. The score falls until the blocks exceed the correlation range 0.95 0.75 random split, 0.94 plateau at 0.79 — the honest number 0 km 15 km 30 km Run this curve once per project; the plateau tells you the block size to use everywhere after.

The curve above is worth producing for any new target variable. Sweep the block size, plot the cross-validated score, and use the size at which the score stops falling. It costs a handful of training runs and settles an argument that otherwise runs for the life of the project.

Environment & Setup

Package Version pin Used for
numpy >=1.23 Block id arithmetic
pandas >=2.0 The chip index table that carries coordinates and folds
scikit-learn >=1.3 GroupKFold and StratifiedGroupKFold
geopandas >=0.14 Region polygons for leave-one-region-out designs 
pip install "numpy>=1.23" "pandas>=2.0" "scikit-learn>=1.3" "geopandas>=0.14"

Complete Working Example

import numpy as np
import pandas as pd
from sklearn.model_selection import StratifiedGroupKFold


def assign_spatial_folds(index: pd.DataFrame, *, block_m: float = 20_000,
                         n_folds: int = 5, buffer_blocks: int = 0,
                         seed: int = 0) -> pd.DataFrame:
    """Add block_id and fold columns to a chip index table.

    index must carry `x` and `y`: the chip centre in the projected CRS.
    Rows with fold == -1 fall on a buffer and should be dropped.
    """
    bx = np.floor(index["x"].to_numpy() / block_m).astype(int)
    by = np.floor(index["y"].to_numpy() / block_m).astype(int)
    index = index.assign(block_id=bx * 1_000_003 + by)   # collision-free enough

    # Assign whole blocks to folds, shuffled so folds are not geographic stripes
    rng = np.random.default_rng(seed)
    blocks = index["block_id"].unique()
    rng.shuffle(blocks)
    fold_of = {b: i % n_folds for i, b in enumerate(blocks)}
    index = index.assign(fold=index["block_id"].map(fold_of))

    if buffer_blocks:
        # Drop chips in a ring of blocks whose neighbours belong to another fold
        keep = np.ones(len(index), dtype=bool)
        coords = {(x, y): fold_of[x * 1_000_003 + y] for x, y in zip(bx, by)}
        for i, (x, y) in enumerate(zip(bx, by)):
            here = coords[(x, y)]
            for dx in range(-buffer_blocks, buffer_blocks + 1):
                for dy in range(-buffer_blocks, buffer_blocks + 1):
                    other = coords.get((x + dx, y + dy))
                    if other is not None and other != here:
                        keep[i] = False
        index.loc[~keep, "fold"] = -1

    return index


if __name__ == "__main__":
    idx = pd.read_parquet("chips_index.parquet")
    idx = assign_spatial_folds(idx, block_m=20_000, n_folds=5, buffer_blocks=1)
    usable = idx[idx.fold >= 0]
    print(idx.groupby("fold").size())
    print(f"{1 - len(usable) / len(idx):.1%} of chips dropped to the buffer")

    splitter = StratifiedGroupKFold(n_splits=5, shuffle=True, random_state=0)
    for train, val in splitter.split(usable, usable["class_id"], usable["block_id"]):
        print(len(train), "train /", len(val), "val")

Shuffling the block list before assigning folds matters more than it looks. Assigning folds in raster order produces folds that are geographic stripes, so fold 0 is the north of the study area and fold 4 the south. Cross-validation then measures north-to-south transfer specifically, which is a real question but not the one five-fold CV is supposed to answer.

StratifiedGroupKFold adds class balance on top of the grouping, which matters when a rare class clusters — without it a fold can end up with no examples of a class at all and the per-class metrics become undefined.

Variant Patterns

1. Leave-one-region-out

When the deployment question is “will this work in the next country”, the fold unit should be the region, not an arbitrary block.

Leave-one-region-out folds Five regions are drawn as adjacent panels. In each of five runs one region is held out entirely and the other four train the model, so every reported score answers the question of how the model behaves in a region it has never seen. Five runs, each holding out one whole region run 1 run 2 run 3 run 4 run 5 north central coast Orange is held out. The spread across the five scores is as informative as their mean. 
import geopandas as gpd
from sklearn.model_selection import LeaveOneGroupOut

regions = gpd.read_file("regions.gpkg").to_crs(chip_crs)
chips = gpd.sjoin(chip_points, regions[["region_id", "geometry"]], how="left")

for train, val in LeaveOneGroupOut().split(X, y, groups=chips["region_id"]):
    ...   # five scores; report the spread, not just the mean

2. Blocks by tile rather than by distance

When the imagery arrives as tiles, the tile is often the natural block: atmospheric state, acquisition date and processing baseline are constant within it and differ across it. Using the tile id as the group is one line and catches sensor-state leakage that distance-based blocks miss.

groups = index["scene_id"].to_numpy()     # e.g. S2A_36NYF_20260614

3. Nested folds for hyperparameter tuning

Selecting hyperparameters on the same fold that reports the score leaks just as surely as a random split does. Nest the search inside the outer split, with the same grouping at both levels.

Nested folds keep tuning out of the reported score The outer split holds one fifth of the blocks back as a test set that no tuning ever touches. Within the remaining four fifths an inner split searches hyperparameters. The score reported comes only from the outer test blocks, so hyperparameter choices cannot inflate it. Outer split reports; inner split tunes outer blocks used for fitting and tuning test inner search Both levels group by block id, so a block never appears on two sides of either split. Only the orange blocks contribute to the number you publish. 
from sklearn.model_selection import GridSearchCV, GroupKFold

outer, inner = GroupKFold(5), GroupKFold(4)
for train, test in outer.split(X, y, groups=block):
    search = GridSearchCV(model, grid, cv=inner)
    search.fit(X[train], y[train], groups=block[train])
    score = search.score(X[test], y[test])

Common Errors

ValueError: n_splits cannot be greater than the number of groups

There are fewer blocks than folds, which means the block size is larger than the study area. Reduce the block edge, or accept that this dataset supports only a hold-out rather than cross-validation.

Folds have wildly different sizes

GroupKFold balances by sample count, but blocks vary in how many chips they contain. If one block holds a third of the chips, no grouping scheme will produce even folds — the collection design is the problem, not the splitter.

A class vanishes from a validation fold

Use StratifiedGroupKFold, and if it still happens, the class exists in too few blocks to validate. Report it as unmeasured rather than quietly averaging over the folds where it appeared.

Frequently Asked Questions

Q: How big should a spatial block be? Larger than the distance over which your target is correlated, which for land cover is usually a few kilometres and for soil or climate variables can be tens. Estimate it from an empirical variogram of the target rather than guessing, then round up.

Q: Does blocking waste data? Buffering does, by a few per cent. Blocking itself does not — the same chips are used, only grouped differently. What blocking costs is the illusion of a high score; what it buys is a number that predicts performance on a new area.

Q: What if one block holds all of a rare class? Then that class cannot be validated, and the honest response is to say so rather than to fall back on a random split. Collect labels for the class in at least three separated locations; until then treat its accuracy as unmeasured.