Validating polygon topology with GeoPandas

6 min read

Validating polygon topology with GeoPandas is the problem of proving — deterministically and at scale — that every polygon in a GeoDataFrame is well-formed under the OGC Simple Features model before it reaches a spatial join, a tiling engine, or a rasteriser. The specific tools are GeoPandas 0.14+, the Shapely 2.x vectorised geometry API, and the GEOS engine underneath both; the predicate at the centre of it all is GeoSeries.is_valid, paired with shapely.validation.explain_validity for diagnostics and shapely.make_valid for repair. This page sits beneath geometry validation patterns, the structural-integrity layer of the Spatial Test Pattern Design & Implementation discipline, and walks through the exact code an engineer runs to turn a noisy is_valid boolean into a versioned, auditable QA gate.

A boolean alone is not a test. A production check has to say which feature failed, why GEOS rejected it, whether the failure survives precision snapping, and whether it can be repaired without distorting the geometry beyond an agreed tolerance. The rest of this page builds that check from the predicate up.

Why polygon topology degrades in the first place

Polygon invalidity is almost never authored directly; it is manufactured by transformations that the engineer did not realise were lossy. Understanding the mechanism is what lets you choose the right tolerance instead of guessing.

• Floating-point overlay residue. union, intersection, and difference are computed in IEEE-754 doubles. When two near-collinear edges are clipped, the intersection point lands a few ULPs off the true line, producing a sliver — a polygon with non-zero but sub-millimetre area. GEOS still considers many slivers valid, but they break R-tree indexing and inflate vertex counts.
• Self-intersection and bowties. Aggressive generalisation (Douglas–Peucker) or hand digitisation can drag one edge across another, so the exterior ring crosses itself. GEOS reports this as Self-intersection and the geometry fails is_valid.
• Ring degeneracy on reprojection. A linear ring needs at least four coordinates (the first repeated as the last). Reprojecting dense WGS84 vertices into a coarse projected grid can collapse coincident points, dropping a ring below four vertices or leaving it unclosed.
• Interior rings escaping the shell. After an overlay, a hole can poke outside its exterior ring or two holes can overlap. GEOS flags Holes are nested or Interior is disconnected.
• Precision collapse during CRS transformation. Converting high-precision EPSG:4326 degrees to a metre-based CRS without an explicit grid introduces noise at the last representable digit, which is enough to flip is_valid non-deterministically between GEOS builds.

The remedy that ties these together is precision snapping: by quantising every coordinate to a fixed grid before validation, you collapse numerical noise into exact equality and make the check reproducible across machines.

Parameter and signature reference

The three Shapely 2.x callables below are the whole toolkit. Note that all of them are vectorised — apply them to a GeoSeries directly rather than looping in Python.

The grid_size you pass to set_precision must be expressed in the units of the active CRS, so derive it the same way you would when setting up spatial tolerance thresholds in assertions. A practical sliver-area floor is the square of the linear grid size:

where $g$ is the grid_size. To confirm a repair did not move the boundary too far, bound the one-sided Hausdorff distance between the original ring $P$ and the repaired ring $P'$:

with $\varepsilon$ set to one grid cell. A repair that exceeds $\varepsilon$ should be quarantined, not silently accepted.

Step-by-step implementation

Step 1 — Snap precision deterministically

Quantise coordinates to a grid before anything else. This is the single most effective step for reproducibility, because it removes the floating-point noise that makes is_valid machine-dependent.

import geopandas as gpd
from shapely import set_precision

# GeoPandas 0.14+, Shapely 2.x, GEOS 3.11+
def snap_grid(gdf: gpd.GeoDataFrame, grid_size: float = 0.001) -> gpd.GeoDataFrame:
    gdf = gdf.copy()
    # grid_size is in CRS units; 0.001 == 1 mm in a metre-based projected CRS
    gdf["geometry"] = set_precision(gdf.geometry.values, grid_size=grid_size)
    return gdf

Step 2 — Assess validity and capture the reason

Run the vectorised predicate, then attach the diagnostic string so a failure is self-explaining in CI logs.

from shapely.validation import explain_validity

def assess(gdf: gpd.GeoDataFrame) -> gpd.GeoDataFrame:
    gdf = gdf.copy()
    gdf["is_valid"] = gdf.geometry.is_valid
    # explain_validity returns "Valid Geometry" or e.g. "Self-intersection[12.0 4.5]"
    gdf["validity_reason"] = gdf.geometry.apply(explain_validity)
    return gdf

Step 3 — Isolate, repair, and bound the drift

Separate the invalid rows, repair them with make_valid, and measure how far the repair moved the boundary so you can accept or quarantine each one.

from shapely import make_valid, hausdorff_distance

def repair(gdf: gpd.GeoDataFrame, eps: float = 0.001) -> tuple[gpd.GeoDataFrame, gpd.GeoDataFrame]:
    invalid = gdf[~gdf["is_valid"]].copy()
    if invalid.empty:
        return gdf, invalid

    original = invalid.geometry.copy()
    invalid["geometry"] = make_valid(invalid.geometry.values, method="linework")
    invalid["repair_valid"] = invalid.geometry.is_valid
    # one-sided Hausdorff drift in CRS units; reject repairs that move the boundary > eps
    invalid["drift"] = hausdorff_distance(original.values, invalid.geometry.values)
    invalid["accepted"] = invalid["repair_valid"] & (invalid["drift"] <= eps)
    return gdf, invalid

Step 4 — Route results and emit a structured audit trail

Valid features and accepted repairs proceed; anything else is quarantined with enough context for a human to trace it upstream.

import json, logging
logging.basicConfig(level=logging.INFO, format="%(levelname)s | %(message)s")

def report(invalid: gpd.GeoDataFrame) -> None:
    accepted = int(invalid["accepted"].sum()) if not invalid.empty else 0
    logging.info(json.dumps({
        "check": "polygon_topology",
        "invalid": int(len(invalid)),
        "repaired_accepted": accepted,
        "quarantined": int(len(invalid)) - accepted,
    }))

Verify the fix

A topology gate is only trustworthy if it is itself tested. The pytest below pins the contract: after snapping and repair, every accepted geometry is valid and within tolerance, and a known bowtie is actually caught.

import geopandas as gpd
from shapely import Polygon

def test_bowtie_is_caught_and_repaired():
    bowtie = Polygon([(0, 0), (1, 1), (1, 0), (0, 1), (0, 0)])  # self-intersecting
    gdf = gpd.GeoDataFrame(geometry=[bowtie], crs="EPSG:3857")

    gdf = assess(snap_grid(gdf, grid_size=0.001))
    assert not gdf["is_valid"].iloc[0]            # the gate sees the failure

    _, invalid = repair(gdf, eps=0.001)
    assert invalid["repair_valid"].all()          # make_valid produced valid output
    assert invalid["geometry"].is_valid.all()     # GEOS agrees independently

Run it as a pre-merge gate with pytest -q test_topology.py; a non-zero exit blocks the merge. For ad-hoc triage, the same predicate is one line of CLI: python -c "import geopandas as g; d=g.read_file('in.gpkg'); print((~d.is_valid).sum(), 'invalid')".

Failure modes and edge cases

1. Empty and null geometries. set_precision on an empty POLYGON EMPTY returns empty, and is_valid reports True for it — a vacuous pass. Filter gdf.geometry.is_empty and gdf.geometry.isna() explicitly before the gate, or empties slip into curated data and break downstream predicates.
2. Anti-meridian wrap. A polygon spanning ±180° longitude in EPSG:4326 is geometrically valid but renders as a planet-wide horizontal smear after any planar overlay. is_valid will not catch it; add a width check (bounds spanning > 180°) or split on the dateline before validation.
3. make_valid changes geometry type. Repairing a self-intersecting Polygon frequently yields a MultiPolygon or even a GeometryCollection mixing lines and polygons. Assert the output type (geom_type) and, if your schema requires polygons, extract polygonal parts with shapely.get_parts and drop stray linework.
4. Grid too coarse. If grid_size exceeds the smallest real feature dimension, set_precision collapses legitimate narrow polygons to invalid or empty geometries — a false failure. Tie $g$ to the data’s true minimum feature size, never to a round number.
5. GEOS version skew. is_valid results can differ between GEOS 3.10 and 3.11 on borderline boundary-touching cases. Pin the GEOS/Shapely build in the CI container so the gate is reproducible; an unpinned runner turns a deterministic check into a flaky one.

By snapping precision first, capturing explain_validity, bounding repair drift with Hausdorff distance, and pinning GEOS in CI, polygon topology validation becomes a deterministic, auditable gate rather than an opaque boolean — exactly the structural-integrity contract that the broader catalogue of geometry validation patterns is built on.

Related

• Geometry Validation Patterns — the parent reference for validity, ring, sliver, and precision checks.
• Topology Rule Enforcement — adjacency, containment, and exclusion rules that run after per-feature validity passes.
• Attribute & Metadata Checks — pairing geometry validity with CRS, dtype, and feature-id integrity.
• Cross-Format Parity Testing — confirming a valid polygon survives serialisation to GeoJSON and GeoPackage.
• Setting up spatial tolerance thresholds in assertions — deriving the grid_size and eps values used above from the CRS unit.