Measuring Spatial Query Error Under Differential Privacy

An analyst asks “how many devices were in this polygon between 08:00 and 09:00?” and gets a number back from a differentially private store. The number is wrong — that is the deal — but how wrong, and wrong in which direction, are questions the pipeline must be able to answer before anyone builds a decision on top of it. This page is the measurement recipe for range and polygon queries specifically, sitting under privacy–utility trade-off measurement within Core Fundamentals & Architecture for Spatial Privacy. The mechanism itself — grid, noise, post-processing — is covered in differentially private spatial aggregation; here we only measure what it produced.

The distinguishing feature of a spatial range query is that its error does not come from one noisy number. A polygon covering $m$ cells sums $m$ independent noise draws, so error grows as $\sqrt{m}$ while the count itself grows as $m$ — which means large queries are relatively more accurate and small ones are relatively worse, the opposite of most people’s intuition about “big queries leak more”.

Parameter Configuration and Calibration

Three settings determine what error you will measure, and all three must be recorded next to any number you report.

• epsilon_per_release — the budget covering the whole grid, not one cell. Under parallel composition, disjoint cells each get the full $\varepsilon$, because a device contributes to exactly one cell per release. Splitting $\varepsilon$ across cells is a common and expensive mistake: it multiplies the noise scale by the cell count for no privacy gain. The corollary is that per-cell noise scale is $b = \Delta/\varepsilon$ regardless of grid size.
• report_floor — the count below which a cell is not worth reporting. Relative error is meaningless on a cell holding four people; conditioning the metric on this floor is what makes the numbers comparable between a dense downtown grid and a rural one. Set it from the smallest count the downstream decision actually distinguishes, typically 25–50 for mobility work.
• query_corpus — the set of shapes you measure against. Error is a property of the query, so a single “mean error” figure is only defined relative to a corpus. Build the corpus from real queries pulled from the access log, bucketed by area, rather than from synthetic circles: real polygons are irregular and straddle cell boundaries, which is where the extra error lives.

The expected absolute error of a $m$-cell polygon under the Laplace mechanism is $\mathbb{E}|e| = \sqrt{m}\,\sigma$ with $\sigma = \sqrt{2}\,\Delta/\varepsilon$, so a 400-cell neighbourhood at $\varepsilon = 1$ carries about 28 counts of expected error, while a 4-cell intersection carries about 2.8. Relative to counts of, say, 12,000 and 40 respectively, that is 0.2% and 7% — the same mechanism, two orders of magnitude apart in usefulness.

Reference Implementation

The function below measures a query corpus against a held-out truth grid and returns the error distribution rather than a single mean, because the mean of a heavy-tailed error distribution is the number least useful for deciding whether to trust an answer.

from __future__ import annotations

from dataclasses import dataclass
from typing import Sequence

import numpy as np


@dataclass(frozen=True)
class QueryError:
    """Error distribution for one bucket of the query corpus."""

    cells: int
    n_queries: int
    mean_abs_error: float
    p95_abs_error: float
    mean_rel_error: float
    p95_rel_error: float
    sign_bias: float          # mean signed error; non-zero means the pipeline is biased


def measure_query_error(
    truth: np.ndarray,
    queries: Sequence[np.ndarray],
    *,
    epsilon: float,
    sensitivity: float = 1.0,
    draws: int = 200,
    report_floor: float = 25.0,
    clamp_negative: bool = True,
    seed: int = 11,
) -> QueryError:
    """Measure answer error for one bucket of same-sized polygon queries.

    `truth` is the unprotected per-cell count vector; each entry of `queries` is an
    index array naming the cells a polygon covers. Noise is drawn per cell per draw,
    exactly as the store would, so summation over the polygon accumulates the same
    sqrt(m) growth the analyst experiences.
    """
    rng = np.random.default_rng(seed)
    scale = sensitivity / epsilon

    abs_errors: list[float] = []
    rel_errors: list[float] = []
    signed: list[float] = []

    for _ in range(draws):
        noisy = truth + rng.laplace(0.0, scale, size=truth.shape)
        if clamp_negative:
            # Post-processing: free in privacy terms, but it makes the estimator
            # one-sided, which is exactly what `sign_bias` below is here to expose.
            noisy = np.maximum(noisy, 0.0)
        for idx in queries:
            true_answer = float(truth[idx].sum())
            noisy_answer = float(noisy[idx].sum())
            err = noisy_answer - true_answer
            abs_errors.append(abs(err))
            signed.append(err)
            if true_answer >= report_floor:
                rel_errors.append(abs(err) / true_answer)

    rel = np.asarray(rel_errors) if rel_errors else np.asarray([np.nan])
    return QueryError(
        cells=int(np.mean([len(q) for q in queries])),
        n_queries=len(queries),
        mean_abs_error=float(np.mean(abs_errors)),
        p95_abs_error=float(np.percentile(abs_errors, 95)),
        mean_rel_error=float(np.mean(rel)),
        p95_rel_error=float(np.percentile(rel, 95)),
        sign_bias=float(np.mean(signed)),
    )


def expected_abs_error(cells: int, epsilon: float, sensitivity: float = 1.0) -> float:
    """Analytic expectation: sqrt(m) * sqrt(2) * Delta / epsilon."""
    return float(np.sqrt(cells) * np.sqrt(2.0) * sensitivity / epsilon)

Validation Checkpoint

The measurement is only trustworthy if it reproduces the analytic expectation on data where the expectation is known. Run this before trusting any number the function produces on real data.

def _validate() -> None:
    rng = np.random.default_rng(3)
    # A dense grid, so clamping never fires and the estimator stays unbiased.
    truth = rng.integers(400, 900, size=2_000).astype(float)
    queries = [rng.choice(2_000, size=64, replace=False) for _ in range(40)]

    measured = measure_query_error(truth, queries, epsilon=1.0, draws=120)
    predicted = expected_abs_error(cells=64, epsilon=1.0)

    # 1. Measured error tracks the sqrt(m) prediction within sampling tolerance.
    assert abs(measured.mean_abs_error - predicted) / predicted < 0.15, (
        f"measured {measured.mean_abs_error:.1f} vs predicted {predicted:.1f}"
    )
    # 2. On dense data with no clamping effect the estimator is unbiased.
    assert abs(measured.sign_bias) < 0.1 * measured.mean_abs_error

    # 3. Halving epsilon must double the error, not change it by some other factor.
    tighter = measure_query_error(truth, queries, epsilon=0.5, draws=120)
    ratio = tighter.mean_abs_error / measured.mean_abs_error
    assert 1.8 < ratio < 2.2, f"error should scale as 1/epsilon, got {ratio:.2f}"

    # 4. Error grows as sqrt(m), so 4x the cells is 2x the error.
    big = [rng.choice(2_000, size=256, replace=False) for _ in range(40)]
    wide = measure_query_error(truth, big, epsilon=1.0, draws=120)
    assert 1.7 < wide.mean_abs_error / measured.mean_abs_error < 2.3

    # 5. On a sparse grid, clamping introduces a measurable upward bias.
    sparse = rng.integers(0, 3, size=2_000).astype(float)
    biased = measure_query_error(sparse, queries, epsilon=0.5, draws=120)
    assert biased.sign_bias > 0, "clamping negatives must bias sparse answers upward"

    print("query-error measurement: all assertions passed")


_validate()

Assertion 5 is the one worth keeping in CI permanently. It encodes the single most consequential property of a spatial DP store: on sparse geography the answers are systematically too high, and the size of that inflation depends on the query’s cell count and the local density, neither of which the analyst can see.

Incident Response and Edge Cases

• An analyst reports an impossible answer — a count larger than the population. Expected on small polygons: the noise tail is unbounded. Do not clip the answer to a “plausible” range, which would create a second, harder-to-model bias. Publish the error bar instead: return the answer alongside $\pm 1.96\sqrt{m}\sigma$ so an out-of-range value reads as low confidence rather than as a bug.
• Two overlapping queries return inconsistent totals. Also expected: independent draws over overlapping cell sets need not be additively consistent. If consistency matters more than freshness, answer from a single materialised noisy grid rather than re-noising per query — this makes answers consistent by construction and costs nothing extra in budget, since the grid was one release.
• Relative error looks fine but decisions are flipping. The queries in the corpus are larger than the queries in production. Re-bucket the corpus by area from the real access log; the failure is almost always concentrated in the smallest bucket.
• Measured error is much worse than $\sqrt{m}\sigma$. Sensitivity is understated. A device appearing in many cells within one release breaks the sensitivity-1 assumption, and no amount of budget compensates. Check that per-device contributions were clipped before binning, then re-measure.

Frequently Asked Questions

Related

• Privacy–Utility Trade-off Measurement — the parent method this recipe plugs into.
• Choosing Epsilon from a Utility Target — inverting these formulas to pick a budget.
• Differentially Private Spatial Aggregation — the mechanism producing the answers.
• Privacy Budget Management — accounting for the releases measured here.

Up one level: Privacy–Utility Trade-off Measurement · Section: Core Fundamentals & Architecture for Spatial Privacy.