Implementing Spatial k-Anonymity with GeoPandas

Spatially join every point to a zone, then repeatedly absorb any zone holding fewer than k points into its nearest neighbor until each published zone — and therefore each individual within it — is indistinguishable from at least k − 1 others.

Core Calculation permalink

k-anonymity grouping for location traces generalizes each record’s location to a zone shared by enough people. Let ZZ be a partition of the study area into zones and n(z)n(z) the count of records in zone zz. The dataset satisfies spatial k-anonymity when:

zZ:  n(z)korn(z)=0\forall z \in Z : \; n(z) \geq k \quad\text{or}\quad n(z) = 0

Every occupied zone must hold at least k records so that the zone identifier, the only location released, maps to at least k indistinguishable individuals. When a zone falls short we merge it with a neighbor zz', forming zzz \cup z' with count n(z)+n(z)n(z) + n(z'). Adaptive merging repeats the rule:

while z:0<n(z)<k    zzarg minzN(z)d(z,z)\text{while } \exists\, z : 0 < n(z) < k \;\Rightarrow\; z \leftarrow z \cup \operatorname*{arg\,min}_{z' \in \mathcal{N}(z)} d(z, z')

where N(z)\mathcal{N}(z) is the set of candidate neighbors and dd is the centroid distance in a metric CRS. Merging into the nearest neighbor minimizes the added spatial generalization — the area a record’s true location is blurred into — so utility loss is kept small.

Parameter table permalink

Parameter Meaning Typical range
k Minimum members per published zone 5–100
Starting zones Initial partition (census tracts, H3, admin units) 100 m–5 km scale
Neighbor rule Adjacency, then nearest-centroid fallback
Metric CRS Projection for distance/area UTM zone / national grid

Worked numeric example permalink

  • k = 25
  • Starting zones: 6 census blocks with counts [40, 30, 12, 8, 60, 5]

Blocks with 12, 8, and 5 members are below k. The 5-member block merges into its nearest neighbor, the 8-member block, giving 13 — still short. That combined zone merges again into the nearest remaining neighbor, the 12-member block, giving 25 — exactly k. The final zone counts are [40, 30, 25, 60], every one at or above 25, with three small blocks collapsed into a single generalized zone covering their union.

Python Implementation permalink

The recipe below projects to a metric CRS, joins points to zones, merges sub-k zones into their nearest neighbor iteratively, dissolves geometries, and reassigns points. Input is WGS84 (EPSG:4326); all distance and adjacency work happens in a projected metric CRS, and the output is reprojected back for publication.

import numpy as np
import geopandas as gpd
from typing import Optional


def spatial_k_anonymize(
    points: gpd.GeoDataFrame,
    zones: gpd.GeoDataFrame,
    k: int,
    metric_crs: str = "EPSG:32633",  # UTM zone 33N; pick the zone for YOUR data
) -> gpd.GeoDataFrame:
    """
    Generalize points to zones so that every occupied zone holds >= k points,
    merging any sub-k zone into its nearest neighbor.

    Privacy rationale: releasing only the (generalized) zone id means each
    record is hidden among at least k-1 others sharing that zone, satisfying
    spatial k-anonymity for a single-point release.

    Args:
        points: Point GeoDataFrame in EPSG:4326, one row per individual.
        zones:  Polygon GeoDataFrame in EPSG:4326, the starting partition.
        k:      Minimum members per published zone.
        metric_crs: Projected CRS for distance/area; use the UTM zone that
                    covers your study area so distances are in meters.

    Returns:
        Points GeoDataFrame (EPSG:4326) with a 'zone_id' column whose every
        occupied zone contains >= k members.
    """
    if k < 2:
        raise ValueError(f"k must be >= 2 to hide a record among others, got {k!r}")

    # --- Project to metric CRS: distances and adjacency need real meters ---
    pts = points.to_crs(metric_crs).copy()
    zn = zones.to_crs(metric_crs).reset_index(drop=True).copy()
    zn["zone_id"] = zn.index.astype("int64")

    # --- Assign each point to a starting zone via spatial join ---
    joined = gpd.sjoin(pts, zn[["zone_id", "geometry"]], how="inner", predicate="within")
    counts = joined.groupby("zone_id").size()
    zn["n"] = zn["zone_id"].map(counts).fillna(0).astype(int)

    # Drop empty zones; only occupied zones need to satisfy k.
    occupied = zn[zn["n"] > 0].copy()

    # --- Iteratively merge the smallest sub-k zone into its nearest neighbor ---
    while (occupied["n"] < k).any():
        occupied = occupied.reset_index(drop=True)
        small_idx = occupied["n"].idxmin()          # smallest count first
        small = occupied.loc[small_idx]

        others = occupied.drop(index=small_idx)
        # Prefer a touching (adjacent) neighbor; fall back to nearest centroid
        # so isolated zones (e.g. islands) still merge across gaps.
        touching = others[others.geometry.touches(small.geometry)]
        pool = touching if not touching.empty else others
        nearest_idx = pool.geometry.centroid.distance(small.geometry.centroid).idxmin()

        # Union geometry and pool the counts into the neighbor.
        occupied.loc[nearest_idx, "geometry"] = (
            occupied.loc[nearest_idx, "geometry"].union(small.geometry)
        )
        occupied.loc[nearest_idx, "n"] = occupied.loc[nearest_idx, "n"] + small["n"]
        occupied = occupied.drop(index=small_idx)

    # --- Reassign points to the final generalized zones ---
    generalized = occupied[["geometry"]].reset_index(drop=True)
    generalized["zone_id"] = generalized.index.astype("int64")
    result = gpd.sjoin(
        pts.drop(columns=[c for c in ("zone_id", "index_right") if c in pts.columns]),
        generalized, how="left", predicate="within",
    ).drop(columns="index_right")

    # Reproject the published points back to WGS84.
    return result.to_crs("EPSG:4326")

The decisions that matter for the privacy guarantee:

  • Smallest-first merging: absorbing the smallest zone each iteration keeps every merge minimal, so generalized zones grow only as much as k demands.
  • Adjacency with nearest-centroid fallback: touching zones are preferred so merged shapes stay contiguous, but the fallback guarantees termination even for isolated zones that touch nothing.
  • Empty zones excluded: a zone with zero records is not a disclosure risk and must not force a spurious merge, so only occupied zones enter the loop.

Verification permalink

def assert_k_anonymous(result: gpd.GeoDataFrame, k: int) -> dict:
    """Confirm every occupied generalized zone holds >= k members."""
    sizes = result.groupby("zone_id").size()
    min_group = int(sizes.min())
    assert min_group >= k, f"Zone below k: min group size {min_group} < k={k}"
    return {
        "min_group_size": min_group,
        "n_zones": int(len(sizes)),
        "records": int(len(result)),
        "k_satisfied": True,
    }

# Expected for k=25 on the worked example:
# {'min_group_size': 25, 'n_zones': 4, 'records': 155, 'k_satisfied': True}

The single assertion sizes.min() >= k is the whole guarantee. Track n_zones before and after to quantify how much spatial resolution the merging cost; a large drop signals that k is too high for the starting zone granularity.

Edge Cases & Adjustments permalink

  • Islands and disconnected zones: a sub-k zone that touches nothing still merges via the nearest-centroid fallback, even across water. If the combined group is still short, the loop keeps absorbing outward; consider suppressing rather than merging across implausibly long distances that would generalize a coastal record into an inland zone.
  • Sparse rural data: when whole regions sit below k, merging can produce a few enormous zones that destroy utility. Cap the merged zone area and suppress records that cannot reach k within the cap, or coarsen with grid aggregation and spatial binning strategies as the starting layer instead.
  • Boundary effects: points near a zone edge may be assigned to the “wrong” side of an administrative line. Where the true location matters less than the count guarantee, this is acceptable; where it does, pair generalization with spatial fuzzing and buffer-zone implementation so edge records are also displaced.
  • CRS: always merge in a metric CRS. Nearest-neighbor selection in EPSG:4326 uses degree distances that mean different ground distances at different latitudes, biasing which zones merge; the tuning of k itself is covered in calculating k-anonymity thresholds for mobile tracking.

FAQ permalink

Why project to a metric CRS before merging?

Nearest-neighbor selection and area computation both depend on distance. In EPSG:4326 a degree of longitude covers different ground distances at different latitudes, so the “nearest” zone by degree distance may not be nearest on the ground. Projecting to a UTM zone or national metric CRS makes distances and areas comparable, then you reproject the generalized output back to EPSG:4326 for publication.

What if a sub-k zone has no adjacent neighbor, like an island?

Fall back from adjacency to nearest-centroid distance so an isolated zone still merges with the closest zone across the gap. If the combined group is still below k, keep merging outward or suppress the zone. Never release a lone sub-k zone because it is isolated — isolation makes re-identification easier, not harder.

How is this different from grid aggregation?

Grid aggregation snaps points onto a fixed lattice and suppresses cells below the count. Spatial k-anonymity by merging keeps irregular administrative or organic zones and grows them adaptively only where the count is short, preserving meaningful boundaries and avoiding discarded sparse cells. You can combine them by using a grid as the starting zone layer and merging its sub-k cells.

Does spatial k-anonymity protect against trajectory linkage?

Not on its own. Per-point zone generalization protects a single snapshot, but a sequence of zone visits can be unique even when each zone holds k people. For moving individuals, combine zone generalization with temporal cloaking or trajectory k-anonymity so the sequence, not just each point, meets the threshold.


Related

Back to k-Anonymity Grouping for Location Traces