Location Intelligence Retail Site Selection Automation

Fixing sliver polygons in spatial join operations

This page solves one precise failure in the catchment enrichment pipeline: the microscopic gap and overlap polygons — slivers — that appear when retail trade areas are overlaid against census geometry, and the deterministic GeoPandas routine that detects and removes them before they corrupt demographic aggregation.

Slivers emerge during overlay operations, coordinate reference system transformations, or when merging datasets with differing vertex densities and floating-point precision limits. When drive-time isochrones, custom buffer zones, or municipal boundaries are intersected with census block groups and tract geometry, sub-metre coordinate mismatches generate polygons often smaller than 500 m². Left in place, these artifacts silently double-count households, misallocate median income, and distort revenue forecasts. Resolving them is a hard data-integrity gate for any team running Performing Point-in-Polygon Joins for Store Catchments inside an automated demographic data integration workflow.

Anatomy of a sliver polygon Two adjacent boundaries with non-coincident vertices produce two artifact polygons: a thin gap sliver where the edges pull apart, and an overlap sliver where they cross. After snapping vertices and a morphological close, the shared edge becomes coincident and the slivers vanish. Where slivers come from — and what removes them Before · mismatched vertices catchment census block gap sliver overlap sliver (double-counts) snap + close After · coincident shared edge catchment census block one coincident edge — no artifacts

Prerequisites

Before running the remediation routine below, confirm the following are in place:

  • Python packages: geopandas>=0.14, shapely>=2.0 (the 2.x make_valid and vectorized buffer API are assumed), and pyproj for CRS handling. Shapely 2.x is required because the area and buffer operations here are vectorized over the GeoSeries.
  • A projected CRS on every input layer. Area thresholds are meaningless in degrees. Reproject both the catchment layer and the reference layer to a metric CRS — a regional UTM zone (e.g. EPSG:32617 for the eastern US) or EPSG:5070 (NAD83 / Conus Albers) for national footprints. CRS standardization is covered in the source schema design at Setting Up PostGIS for Retail Analytics.
  • Valid input geometry, or a tolerance for repairing it. Self-intersections and bad ring orientation precede most sliver formation, so the routine repairs topology with make_valid() before measuring area.
  • The two layers being joined. Typically a trade area polygon layer (isochrones or buffers) and a census layer such as block groups synced via Syncing US Census ACS Data via API.

Configuration and execution parameters

Sliver remediation is configuration-driven, not a manual cartographic edit. The following parameters control detection sensitivity and the remediation aggressiveness. Tune them per market density rather than hard-coding a single value.

Parameter Type Valid range Retail default Purpose
crs EPSG code any projected CRS EPSG:5070 Metric CRS for accurate area in m²
sliver_area_threshold_m2 float 250 – 5000 1000.0 Area below which a polygon is a sliver
snap_tolerance_m float 0.5 – 2.0 1.0 Vertex snapping band between layers
buffer_close_dist_m float 0.25 – 1.0 0.5 Negative-then-positive buffer distance
max_allowed_sliver_pct float 0.0 – 0.01 0.005 Fraction of total area allowed as slivers
validation_mode enum strict / warn / bypass strict Pipeline behaviour when the gate trips

Calibrate sliver_area_threshold_m2 to market density: in dense urban markets 1000 m² isolates true artifacts without consuming real blocks, while rural and exurban zones tolerate 5000 m². Set snap_tolerance_m to the positional accuracy of the source data — 0.5–2.0 m for municipal parcel data or Census TIGER/Line files. A YAML representation of these settings keeps the stage reproducible across planning cycles:

yaml
spatial_join_config:
  crs: "EPSG:5070"
  sliver_area_threshold_m2: 1000.0
  snap_tolerance_m: 1.0
  buffer_close_dist_m: 0.5
  max_allowed_sliver_pct: 0.005  # 0.5% of total catchment area
  validation_mode: "strict"      # strict | warn | bypass

Detection and remediation implementation

The routine has three stages — detect, snap-and-close, and gate — wired so the same config dict drives all of them. Each function asserts the CRS is projected before it touches an area calculation, which is the single most common cause of silent failure.

python
import geopandas as gpd
import logging
from shapely.ops import snap
from shapely.validation import make_valid

logger = logging.getLogger(__name__)


def _assert_projected(gdf: gpd.GeoDataFrame) -> None:
    """Fail fast if area math would run in degrees."""
    if gdf.crs is None or not gdf.crs.is_projected:
        raise ValueError("CRS must be projected to calculate area in square metres.")


def detect_slivers(gdf: gpd.GeoDataFrame, threshold_m2: float = 1000.0) -> gpd.GeoDataFrame:
    """Return rows whose polygon area falls below the configured threshold."""
    _assert_projected(gdf)

    # Repair topology BEFORE measuring area: self-intersections inflate area errors.
    gdf = gdf.copy()
    gdf["geometry"] = gdf.geometry.apply(lambda g: make_valid(g) if g else g)

    area = gdf.geometry.area
    mask = area < threshold_m2
    if mask.any():
        logger.warning(
            "Detected %d sliver polygons below %.0f m2; max sliver area %.2f m2",
            int(mask.sum()), threshold_m2, float(area[mask].max()),
        )
    return gdf[mask]


def apply_snapping(target: gpd.GeoDataFrame,
                   reference: gpd.GeoDataFrame,
                   tolerance: float = 1.0) -> gpd.GeoDataFrame:
    """Snap target vertices onto the reference layer to close sub-pixel gaps."""
    _assert_projected(target)
    snapped = [
        snap(t, r, tolerance)
        for t, r in zip(target.geometry, reference.geometry)
    ]
    return target.set_geometry(snapped, crs=target.crs)


def morphological_close(gdf: gpd.GeoDataFrame, buffer_dist: float = 0.5) -> gpd.GeoDataFrame:
    """Negative-then-positive buffer: collapses thin gaps without moving the boundary."""
    _assert_projected(gdf)
    closed = gdf.geometry.buffer(-buffer_dist).buffer(buffer_dist)
    # Small polygons can buffer to empty/invalid geometry — re-validate.
    closed = closed.apply(lambda g: make_valid(g) if g and not g.is_empty else g)
    return gdf.set_geometry(closed, crs=gdf.crs)


def remediate(catchments: gpd.GeoDataFrame,
              reference: gpd.GeoDataFrame,
              config: dict) -> gpd.GeoDataFrame:
    """Full deterministic pass: reproject, snap, close, then gate on residual slivers."""
    crs = config["crs"]
    catchments = catchments.to_crs(crs)
    reference = reference.to_crs(crs)

    total_area = catchments.geometry.area.sum()

    fixed = apply_snapping(catchments, reference, config["snap_tolerance_m"])
    fixed = morphological_close(fixed, config["buffer_close_dist_m"])

    residual = detect_slivers(fixed, config["sliver_area_threshold_m2"])
    residual_pct = residual.geometry.area.sum() / total_area if total_area else 0.0

    if residual_pct > config["max_allowed_sliver_pct"]:
        msg = (f"Residual sliver area {residual_pct:.4%} exceeds "
               f"{config['max_allowed_sliver_pct']:.4%} budget.")
        if config["validation_mode"] == "strict":
            raise RuntimeError(msg)
        logger.error(msg)

    return fixed

For point-in-polygon joins specifically, slivers are mitigated a second way: enforce a minimum intersection area when reconciling overlapping catchments. GeoPandas sjoin exposes the spatial predicate, but the area filter must run after the join so duplicate point attribution is removed deterministically rather than being silently merged.

Sliver remediation decision flow The catchment and reference layers are reprojected to a metric CRS, repaired with make_valid, snapped within tolerance, morphologically closed with a negative-then-positive buffer, then measured for slivers below the area threshold. The residual sliver percentage branches the pipeline: at or below the budget the layer passes to the join; above the budget strict mode raises while warn mode logs and continues. Deterministic sliver remediation pass Catchment layer isochrones · buffers Reference layer census block groups Reproject to metric CRS make_valid() repair topology Snap vertices tolerance band Morphological close −buffer then +buffer detect_slivers() area < sliver_area_threshold_m2 residual % vs budget? ≤ budget Pass to join clean geometry > budget strict → raise warn → log & continue

Failure modes and debugging

Symptom Root cause Fix
ValueError: CRS must be projected Layer still in EPSG:4326 (degrees) to_crs() to a metric CRS before area math
Area threshold removes real blocks Threshold too high for market density Lower to 1000 m² urban; inspect area histogram
make_valid returns a GeometryCollection Mixed point/line/polygon debris from a bad overlay Filter to polygonal parts with .geometry.buffer(0) or extract polygons before measuring
Buffer produces empty geometry Polygon narrower than 2 × buffer_close_dist_m Reduce buffer_close_dist_m; re-validate after closing
Slivers persist after snapping Layers have non-coincident vertices outside the tolerance band Raise snap_tolerance_m toward source positional accuracy
Centroid jumps after remediation Closing distance too aggressive Reduce buffer_close_dist_m; verify centroid stability (below)

A CRS mismatch is the dominant failure: if either layer is unprojected, areas come back in square degrees and every polygon reads as a sliver. The _assert_projected guard turns that from a silent corruption into an immediate exception. When debugging a stubborn overlay, log structured JSON containing pipeline_run_id, source_layer_version, geometry_hash, and timestamp so attribution drift across quarterly ACS refreshes can be traced to a specific layer version.

Verification

Confirm correct output against measurable ground-truth tolerances before the layer feeds any downstream model. Tie this into the broader checks documented in Validating Spatial Join Accuracy with Ground Truth.

  1. Area delta. Total catchment area must deviate less than 0.1% from the pre-remediation total: abs(after - before) / before < 0.001.
  2. Centroid stability. Euclidean distance between pre- and post-remediation centroids should stay under 50 m. Flag any polygon that shifts further — it indicates over-aggressive closing.
  3. Topology consistency. Run shapely.validation.explain_validity() on a random 5% sample and confirm zero residual self-intersections.
  4. Row-count parity. Sliver removal must not drop whole records: the output row count equals the input row count (slivers are merged away, not deleted as rows).
  5. Demographic continuity. Re-aggregate household counts and median income against the pre-join baseline; deviations should sit inside US Census Bureau margin-of-error guidance.

When these five checks pass, the remediated layer is safe to feed into the next stage — point assignment and the demographic weighting that follows in Weighting Demographic Variables for Target Audiences.

← Back to Performing Point-in-Polygon Joins for Store Catchments