Spatial Epi Public Health GIS

Regression Testing Cluster Detection Output

This guide, part of Testing Geospatial Analysis Code, builds the regression test for the output that matters most and changes most subtly: which areas a cluster statistic flags, and how strongly.

Problem Context & Constraints

A surveillance pipeline’s cluster output drives action — an investigation, a press statement, a resource shift — so an unintended change to it is the most consequential regression the pipeline can suffer. It is also awkward to test for three reasons.

Permutation inference is stochastic. Local Moran’s I, Getis-Ord statistics and scan statistics all assign significance by comparing the observed statistic to a distribution built from random permutations. Two runs with different permutation draws produce slightly different pseudo-p-values, and an area near the threshold can flip between significant and not.

A cluster is a set, not a number. Comparing two runs means comparing which areas were flagged, how they were grouped, and how strongly — a structured comparison, not an equality check.

Some change is legitimate. A bug fix, a corrected weights matrix, or a deliberate method improvement should change the output. The test’s job is to make every change visible and deliberate, not to forbid change.

The constraint that makes the test tractable is controlling the randomness: fix the permutation seed and count so that, for identical inputs and code, the output is bit-for-bit identical. Once that holds, any difference is caused by something that changed, and the test can report exactly what.

Threshold Flicker Without a Fixed Seed Two sets of five repeated runs of the same local cluster statistic on identical data. In the unseeded runs, most areas are consistently significant or consistently not, but two areas whose pseudo-p-values sit close to the threshold flip between flagged and unflagged from run to run, making any regression comparison noisy. In the seeded runs, every area gives the same result every time, so any difference in a later run can only have been caused by a change to data, code or dependencies. Fix the permutation seed and the flicker disappears Unseeded: five identical runs Seeded: five identical runs area 1 area 2 area 3 area 4 area 5 areas 2 and 4 flip between runs identical every time filled = flagged significant · outlined = not flagged

Prerequisites

  • python 3.11
  • libpysal 4.10, esda 2.5.1, geopandas 0.14.4, numpy 1.26.4, pandas 2.2.2
  • A synthetic regression fixture with a planted cluster, as in Building Synthetic Case Data for Pipeline Tests
  • A stored known-good output from a reviewed run
  • A fixed permutation count and seed, recorded with the stored output

Step-by-Step Solution

Run the statistic with a fixed seed, compare membership exactly, compare statistics with a tolerance, and produce a structured diff that says precisely what changed.

# Regression comparison for local cluster statistics.
# Pinned: libpysal==4.10, esda==2.5.1, geopandas==0.14.4, numpy==1.26.4, pandas==2.2.2
import logging

import numpy as np
import pandas as pd
import geopandas as gpd
from esda.moran import Moran_Local
from libpysal.weights import Queen

logging.basicConfig(level=logging.INFO, format="%(asctime)s %(levelname)s %(message)s")
log = logging.getLogger("tests.cluster_regression")

PERMUTATIONS = 9999
SEED = 20260918
ALPHA = 0.05


def run_local_moran(areas: gpd.GeoDataFrame, value: str) -> pd.DataFrame:
    """Deterministic local Moran's I: sorted input, fixed permutations, fixed seed."""
    areas = areas.sort_values("geoid").reset_index(drop=True)
    w = Queen.from_dataframe(areas, use_index=False, silence_warnings=True)
    w.transform = "r"
    lm = Moran_Local(areas[value].to_numpy(float), w,
                     permutations=PERMUTATIONS, seed=SEED)
    quad = np.where(lm.p_sim < ALPHA, lm.q, 0)          # 0 = not significant
    return pd.DataFrame({"geoid": areas["geoid"], "Is": lm.Is, "p_sim": lm.p_sim,
                         "quadrant": quad})


def diff_cluster_output(ref: pd.DataFrame, new: pd.DataFrame,
                        stat_tol: float = 1e-9) -> dict:
    """A structured diff: exactly which areas changed, and how."""
    m = ref.merge(new, on="geoid", how="outer", suffixes=("_ref", "_new"),
                  indicator=True)
    missing = m[m["_merge"] != "both"]["geoid"].tolist()
    both = m[m["_merge"] == "both"]
    changed_quad = both[both["quadrant_ref"] != both["quadrant_new"]]
    stat_moved = both[~np.isclose(both["Is_ref"], both["Is_new"], atol=stat_tol)]
    newly_sig = changed_quad[(changed_quad["quadrant_ref"] == 0)
                             & (changed_quad["quadrant_new"] != 0)]["geoid"].tolist()
    lost_sig = changed_quad[(changed_quad["quadrant_ref"] != 0)
                            & (changed_quad["quadrant_new"] == 0)]["geoid"].tolist()
    out = {"areas_missing": missing,
           "quadrant_changed": changed_quad["geoid"].tolist(),
           "newly_significant": newly_sig,
           "no_longer_significant": lost_sig,
           "statistic_moved": stat_moved["geoid"].tolist(),
           "max_stat_delta": float(np.abs(both["Is_new"] - both["Is_ref"]).max())
                             if len(both) else 0.0}
    out["identical"] = not (missing or out["quadrant_changed"] or out["statistic_moved"])
    (log.info if out["identical"] else log.error)(
        "cluster regression: identical=%s newly_sig=%d lost_sig=%d quad_changed=%d "
        "stat_moved=%d max_delta=%.3e",
        out["identical"], len(newly_sig), len(lost_sig), len(out["quadrant_changed"]),
        len(out["statistic_moved"]), out["max_stat_delta"])
    return out


def test_cluster_output_unchanged(fixture_areas, stored_reference):
    new = run_local_moran(fixture_areas, "rate")
    ref = pd.read_parquet(stored_reference)
    d = diff_cluster_output(ref, new)
    assert d["identical"], (
        f"cluster output changed — newly significant {d['newly_significant']}, "
        f"no longer significant {d['no_longer_significant']}, "
        f"quadrant changed {d['quadrant_changed']}")


def test_planted_cluster_still_found(fixture_areas, manifest):
    """Independent of the stored run: the planted signal must still be detected."""
    out = run_local_moran(fixture_areas, "rate")
    hh = set(out.loc[out["quadrant"] == 1, "geoid"])
    assert hh & set(manifest["planted_area_ids"]), "planted high-high cluster lost"

The structured diff is the practical value of the test. “Output changed” is not actionable; “tract A017 is newly significant and A022 is no longer significant, and no statistic moved by more than 1e-12” tells the reviewer that a boundary threshold was crossed by something tiny, while “forty statistics moved by 0.3 and the high-high set changed completely” tells them the weights matrix changed. The same failure message would hide the difference.

The planted-cluster test runs alongside the stored-output comparison for a reason. When a legitimate change forces the stored output to be updated, the planted test is the independent check that the new output is still right — that the update did not ratify a regression.

What the Shape of the Diff Tells You Three diff signatures and their likely causes. When a few areas change significance but no statistic moves materially, the cause is a change in permutation count, seed or significance threshold. When most statistics shift slightly and a handful of areas change quadrant, the cause is a change in the input values, such as a denominator vintage. When statistics shift substantially and the cluster sets rearrange, the cause is a change in the spatial weights, such as a repaired boundary or an altered contiguity rule. Each signature narrows the investigation before any code is read. The diff's shape points at the layer that changed Few flips, no stat moves areas near the threshold crossed it cause: seed, permutation count, or alpha changed Many small stat moves every value shifted a little, a few quadrants changed cause: the input values — a denominator vintage Large moves, sets rearranged neighbourhood structure is different cause: the spatial weights — a repaired boundary None of these is automatically wrong; all of them must be seen and decided on

Validation & Edge Cases

Unsorted input. The weights matrix and the permutation draws depend on row order, so an input sorted differently gives different results even with a fixed seed. Sort on a stable identifier before building weights.

Permutation count changed for speed. Dropping from 9,999 to 999 permutations to speed up continuous integration changes the resolution of the pseudo-p-value and will flip borderline areas. Keep the regression test’s permutation count equal to production’s.

Library upgrades that change the permutation algorithm. Even with a fixed seed, a new library version may draw permutations differently. The regression test will catch it as a set of flips with no statistic movement — the first signature above — which is the cue to review and, if accepted, update the stored output.

Scan statistics rather than local indicators. A space-time scan returns cluster objects — a centre, a radius, a time window and a likelihood ratio — rather than per-area flags. Compare those objects by overlap: match each reference cluster to the new cluster with the greatest area overlap, then report matched clusters whose radius, window or rank changed, unmatched reference clusters that disappeared and new clusters that appeared. A per-area comparison on scan output misreports a small shift in a cluster’s centre as dozens of independent changes.

Tied p-values at the permutation resolution. With 9,999 permutations the smallest attainable pseudo-p-value is 0.0001, and several areas can share it. An ordering-dependent tie-break in a downstream ranking then varies between runs even when every statistic is identical. Break ties explicitly on the identifier wherever output is ranked.

Regression fixture outgrown by the method. A fixture built when the pipeline used one statistic may not exercise a second statistic added later — no space-time structure, for instance, when a temporal scan was introduced. Review the fixture whenever a method is added, and extend the generator rather than adding a second, unrelated fixture.

Multiple-testing correction applied inconsistently. If production applies false discovery rate control and the test does not, they disagree by design. The test must reproduce production’s full significance procedure, including any correction.

Stored reference with no provenance. A stored output whose seed, permutation count and library versions were not recorded cannot be regenerated, so a legitimate update cannot be verified. Store the reference with its full provenance.

ERROR tests.cluster_regression cluster regression: identical=False newly_sig=1 lost_sig=1 quad_changed=2 stat_moved=0 max_delta=0.000e+00

For scan statistics the comparison is between cluster objects rather than per-area flags, and matching them by overlap keeps a small shift from looking like wholesale change:

Match Scan Clusters by Overlap, Not by Area Flags Reference scan clusters drawn as dashed circles and new-run clusters as solid circles. The first reference cluster overlaps a new cluster whose centre has shifted slightly; they are matched and reported as one cluster with a changed radius. A second reference cluster has no overlapping new cluster and is reported as lost. A new cluster with no reference counterpart is reported as appeared. A per-area comparison would have reported the first case as dozens of independent flag changes. One shifted cluster is one change, not forty matched: radius changed reference only lost new only appeared dashed circles are the stored reference run; solid circles are the new run

Compliance Notes

Record the stored reference output with its provenance — fixture version, seed, permutation count, significance threshold, correction method and library versions — and every decision to update it, with the diff that prompted the update and the reviewer who accepted it.

Cluster output is what drives public health action, so a change to the method that produces it is a change to the evidence base for those actions. Where an accepted change would have altered past alerts — a corrected weights matrix, for instance — assess whether any previously issued alert would not have been raised under the corrected method, and record that assessment even when the answer is none. It is the question an oversight body will ask first.