Private Quadtree Decomposition for Spatial Histograms

Real location data is brutally skewed: a downtown census tract may hold a hundred thousand check-ins while the rural cell beside it holds three, and a uniform differentially private grid spends the same noise budget on both — burning privacy on empty space while under-resolving the dense core. This page, a deep-dive inside the differentially private spatial aggregation guide within the Differential Privacy for Geospatial Data section, builds the alternative: an adaptive quadtree that subdivides only where a noisy count says density is worth resolving, splits the ε\varepsilon budget across tree levels, and reuses each node’s released count as its own split decision so no budget leaks to structure. Where the sibling guide on building differentially private heatmaps fixes the cell grid before any data arrives, everything here is data-dependent — the tree shape itself is a private release, and every noisy count it touches must be debited against the privacy budget.

Parameter Configuration and Calibration

A private quadtree is governed by four knobs, and each one trades noise against spatial resolution rather than raw performance. The defaults below assume EPSG:4326 decimal-degree input and event-level sensitivity Δf=1\Delta f = 1 (one point per subject inside the extent); if a subject can contribute multiple points, scale every Laplace scale by that per-user cap.

Parameter Symbol Default Rationale
Total budget εtotal\varepsilon_{\text{total}} 0.7 Cumulative ε\varepsilon debited along the deepest root-to-leaf path — the composition ceiling, not a per-node value.
Max depth HH 3 Bounds the finest cell; at WGS84 a depth-3 quarter of a 1° extent is ~13.9 km, still above the HIPAA 20,000-population geographic floor for most metros.
Split threshold τ\tau ≈ 3·b/ε_level Subdivide only when a node’s noisy count clears τ\tau; sized above the Laplace noise floor so pure noise rarely triggers a spurious split.
Budget ratio ρ\rho 4^{1/3} ≈ 1.587 Geometric allocation across levels; the cube-root-of-fanout ratio minimises range-query error for a quadtree (fanout 4).

The critical design choice is how εtotal\varepsilon_{\text{total}} is split across levels. Every node emits one noisy count, and along any root-to-leaf path the counts compose sequentially, so the per-level budgets must sum to at most εtotal\varepsilon_{\text{total}}:

∑i=0H−1εi≤εtotal,εi=εtotal⋅ρ i∑j=0H−1ρ j.\sum_{i=0}^{H-1} \varepsilon_i \le \varepsilon_{\text{total}}, \qquad \varepsilon_i = \varepsilon_{\text{total}} \cdot \frac{\rho^{\,i}}{\sum_{j=0}^{H-1}\rho^{\,j}}.

Geometric allocation with ρ>1\rho>1 deliberately gives deeper levels more budget. Deep cells hold small counts where a fixed noise magnitude does the most relative damage, so they need the tighter (larger-ε\varepsilon) scale; the coarse root can tolerate more noise because its count is large. For H=3H=3 and εtotal=0.7\varepsilon_{\text{total}}=0.7 this yields roughly ε1≈0.14\varepsilon_1\approx0.14, ε2≈0.22\varepsilon_2\approx0.22, ε3≈0.35\varepsilon_3\approx0.35. Uniform allocation (εi=εtotal/H\varepsilon_i=\varepsilon_{\text{total}}/H) is the safe fallback when your dominant query is a single point lookup rather than a range aggregate.

The split threshold τ\tau is where adaptivity pays for itself. A uniform grid materialises 4H4^H cells unconditionally; the quadtree only recurses where the noisy count exceeds τ\tau, so sparse rural quadrants stop at depth 0 or 1 and never spend the deeper-level budget at all. Because the thresholded statistic is the very count we already released for that node, testing it is free — the split decision consumes no additional budget beyond the count itself.

Reference Implementation

The class below builds the tree recursively. Each QuadNode carries the Laplace-noised count that both is its released statistic and drives its split test, together with the exact ε\varepsilon spent, so the budget audit in the next section can walk the tree and prove the composition bound. Comments flag every line where privacy budget is consumed.

python
from __future__ import annotations

from dataclasses import dataclass, field
from typing import List, Tuple

import numpy as np

Bounds = Tuple[float, float, float, float]  # (min_lon, min_lat, max_lon, max_lat)


@dataclass
class QuadNode:
    """One cell of the private decomposition.

    `noisy_count` is BOTH the differentially private release for this cell
    and the statistic we threshold on to decide whether to subdivide — so
    the split decision is free. `eps_spent` records the budget this single
    noisy count consumed, for the path-composition audit downstream.
    """
    bounds: Bounds
    depth: int
    noisy_count: float
    eps_spent: float
    children: List["QuadNode"] = field(default_factory=list)

    @property
    def is_leaf(self) -> bool:
        return not self.children


def geometric_budget(eps_total: float, max_depth: int, ratio: float) -> np.ndarray:
    """Split eps_total across levels 0..max_depth-1 so deeper levels get more.

    Returns per-level epsilons summing to eps_total; ratio = fanout**(1/3)
    (~1.587 for a quadtree) is near-optimal for range-query error.
    """
    weights = ratio ** np.arange(max_depth, dtype=float)
    return eps_total * weights / weights.sum()  # sum == eps_total exactly


class PrivateQuadtree:
    """Adaptive, differentially private quadtree over 2-D points.

    Budget is spent ONLY on the per-node Laplace counts. Along any
    root-to-leaf path the counts compose sequentially, and the geometric
    per-level split guarantees the path total never exceeds eps_total.
    """

    def __init__(
        self,
        eps_total: float,
        max_depth: int = 3,
        split_threshold: float | None = None,
        ratio: float = 4.0 ** (1.0 / 3.0),
        sensitivity: float = 1.0,
        rng: np.random.Generator | None = None,
    ) -> None:
        if not (eps_total > 0 and max_depth >= 1):
            raise ValueError("require eps_total > 0 and max_depth >= 1")
        self.eps_total = eps_total
        self.max_depth = max_depth
        self.sensitivity = sensitivity          # Delta f: per-subject count cap
        self.level_eps = geometric_budget(eps_total, max_depth, ratio)
        # Default threshold sits above the noise floor of the finest level
        # so pure Laplace noise rarely forces a spurious subdivision.
        finest_scale = self.sensitivity / self.level_eps[-1]
        self.split_threshold = (
            split_threshold if split_threshold is not None else 3.0 * finest_scale
        )
        self.rng = rng if rng is not None else np.random.default_rng()

    def _noisy_count(self, true_count: int, level: int) -> Tuple[float, float]:
        """Laplace mechanism on a count. THIS is where budget is spent."""
        eps = float(self.level_eps[level])                 # budget for this level
        scale = self.sensitivity / eps                     # b = Delta f / eps
        noise = self.rng.laplace(loc=0.0, scale=scale)     # calibrated DP noise
        return true_count + noise, eps                     # (released, eps_spent)

    def build(self, points: np.ndarray, bounds: Bounds, depth: int = 0) -> QuadNode:
        """Recursively decompose `points` (shape (N, 2): lon, lat) in `bounds`."""
        min_x, min_y, max_x, max_y = bounds
        inside = (
            (points[:, 0] >= min_x) & (points[:, 0] < max_x)
            & (points[:, 1] >= min_y) & (points[:, 1] < max_y)
        )
        cell_points = points[inside]
        true_count = int(cell_points.shape[0])             # sensitive — never released raw

        noisy, eps = self._noisy_count(true_count, depth)  # <-- budget debit
        node = QuadNode(bounds=bounds, depth=depth, noisy_count=noisy, eps_spent=eps)

        # Data-dependent stop: recurse ONLY where the NOISY count justifies it.
        # We threshold the already-released count, so this test costs 0 budget.
        if depth + 1 < self.max_depth and noisy > self.split_threshold:
            mid_x = 0.5 * (min_x + max_x)
            mid_y = 0.5 * (min_y + max_y)
            quadrants: List[Bounds] = [
                (min_x, min_y, mid_x, mid_y),  # SW
                (mid_x, min_y, max_x, mid_y),  # SE
                (min_x, mid_y, mid_x, max_y),  # NW
                (mid_x, mid_y, max_x, max_y),  # NE
            ]
            node.children = [self.build(cell_points, q, depth + 1) for q in quadrants]
        return node

    def enforce_consistency(self, node: QuadNode) -> None:
        """Post-process (top-down): make each parent equal the sum of children.

        Constrained inference (Hay et al. 2010) — a pure post-processing step
        that spends NO budget yet reduces variance ("boosting") by exploiting
        the redundancy between a parent's count and its childrens' sum. We fix
        each node's children to its already-settled value, THEN recurse, so
        the whole tree ends internally consistent at every level.
        """
        if node.is_leaf:
            return
        child_sum = sum(c.noisy_count for c in node.children)
        residual = node.noisy_count - child_sum            # inconsistency to spread
        for child in node.children:
            child.noisy_count += residual / len(node.children)
        for child in node.children:                        # recurse AFTER fixing
            self.enforce_consistency(child)

The top-down equal-share spread in enforce_consistency is the simplest estimator that leaves every level internally consistent for zero budget. Hay’s variance-optimal version does more: a bottom-up weighted-mean pass first (weighting each node by its inverse noise variance, which differs per level here because each εi\varepsilon_i differs), then a top-down distribution — squeezing more boosting out of the same counts. The equal-share pass trades a little of that variance reduction for a much shorter implementation while never manufacturing a consistency violation.

Validation Checkpoint

Two invariants matter and neither can be eyeballed: the composition bound (no path may overspend εtotal\varepsilon_{\text{total}}) and the adaptivity property (dense regions must out-resolve sparse ones). The harness below asserts both on a deliberately skewed dataset — a tight urban cluster beside a thin rural scatter.

python
def _validate() -> None:
    rng = np.random.default_rng(42)
    # Dense cluster (12k pts) tight in the NE quadrant; sparse scatter (15 pts)
    # in the SW quadrant. An explicit threshold keeps the contrast well above
    # the noise floor so the sparse quadrant provably stops early.
    dense = rng.normal(loc=[0.78, 0.78], scale=0.015, size=(12_000, 2))
    sparse = rng.uniform(low=0.0, high=0.4, size=(15, 2))
    points = np.vstack([dense, sparse]).clip(0.0, 1.0)

    qt = PrivateQuadtree(eps_total=0.7, max_depth=5, split_threshold=80.0, rng=rng)
    root = qt.build(points, bounds=(0.0, 0.0, 1.0, 1.0))

    # 1. Composition: eps summed along EVERY root-to-leaf path <= eps_total.
    def max_path_eps(node: QuadNode, acc: float) -> float:
        acc += node.eps_spent
        if node.is_leaf:
            return acc
        return max(max_path_eps(c, acc) for c in node.children)

    assert max_path_eps(root, 0.0) <= qt.eps_total + 1e-9, "path overspends budget"

    # 2. Adaptivity: the dense NE subtree subdivides deeper than sparse cells.
    def max_depth_reached(node: QuadNode) -> int:
        return node.depth if node.is_leaf else max(map(max_depth_reached, node.children))

    ne_child = root.children[3]                    # NE quadrant holds the cluster
    sw_child = root.children[0]                     # SW quadrant is near-empty
    assert max_depth_reached(ne_child) > max_depth_reached(sw_child), \
        "dense region must resolve deeper than sparse region"

    # 3. Consistency: after boosting, each parent equals the sum of its children.
    qt.enforce_consistency(root)

    def check_consistent(node: QuadNode) -> None:
        if node.is_leaf:
            return
        assert abs(node.noisy_count - sum(c.noisy_count for c in node.children)) < 1e-6
        for c in node.children:
            check_consistent(c)

    check_consistent(root)
    print("quadtree invariants hold: bounded budget, adaptive depth, consistent counts")


if __name__ == "__main__":
    _validate()

Invariant 1 is the whole privacy argument: because geometric_budget sums to exactly eps_total and the deepest path visits one node per level, the worst-case path spends precisely eps_total and no configuration of the data can exceed it. Invariant 2 confirms the tree earned its complexity — if a uniform grid would have resolved every cell equally, the adaptive machinery bought nothing.

Adaptive quadtree subdivision with a geometric per-level epsilon split On the left, a square map extent is recursively subdivided. The whole square is the depth-0 root (epsilon one). One line-pair splits it into four quadrants at depth 1; three of those quadrants are sparse (few dots) and stop, while the dense north-east quadrant (many dots) splits again into four at depth 2, and its densest sub-quadrant splits a third time at depth 3. On the right, a panel lists the geometric budget split: level 1 root epsilon one approximately 0.14, level 2 epsilon two approximately 0.22, level 3 epsilon three approximately 0.35, deeper levels receiving more budget, and the deepest root-to-leaf path summing to the total budget of 0.7. A note states that sparse quadrants never spend the deeper-level budget. Adaptive quadtree — subdivide where density warrants, split ε across levels depth 0 · root · ε₁ sparse · stop sparse · stop sparse · stop dense → depth 3 Geometric budget split (ρ = 4¹⁄³ ≈ 1.59) deeper levels resolve small counts → get more ε depth 0 · ε₁ ≈ 0.14 depth 1 · ε₂ ≈ 0.22 depth 2 · ε₃ ≈ 0.35 deepest path spends ε₁+ε₂+ε₃ = ε total 0.7 sparse quadrants never spend the deeper-level budget

A uniform grid would spend all three levels of budget on every cell, including the empty ones.

Incident Response and Edge Cases

Adaptive decompositions fail in quiet, data-dependent ways that a uniform grid never exhibits — the tree shape leaks, or the noise manufactures structure. The high-frequency production failures:

  • Threshold set below the noise floor. If τ\tau is smaller than a few multiples of the finest Laplace scale Δf/εH−1\Delta f / \varepsilon_{H-1}, pure noise routinely pushes an empty cell over the bar and the tree grows phantom dense regions in genuinely empty rural space. Detection: count leaves at max depth whose consistency-corrected count is near zero. Remediation: raise τ\tau to ≥3Δf/εH−1\ge 3\Delta f/\varepsilon_{H-1}, or switch to uniform allocation so the finest level keeps a usable scale.
  • Budget silently overspent by an extra count. Adding any statistic outside build — a total-population count, a re-scan for a dashboard — composes on top of the path budget and voids the εtotal\varepsilon_{\text{total}} guarantee. Detection: the max_path_eps audit plus a global ledger that debits every Laplace draw. Remediation: route all counts through the same accountant described in privacy budget management; never issue an off-book count.
  • Sensitivity miscalibrated for multi-point subjects. The default Δf=1\Delta f = 1 assumes one point per subject inside the extent. A trajectory dataset where a subject contributes dozens of pings has sensitivity equal to that per-user cap, and using 11 under-noises every node. Detection: assert a per-subject point cap at ingestion. Remediation: set sensitivity to the enforced cap and clip contributions before building.
  • Negative or non-monotone counts after boosting. Laplace noise and residual spreading can drive a leaf count below zero or make a parent smaller than a child. Detection: scan for negative corrected counts. Remediation: clamp with max(count, 0) as a post-processing step (it spends no budget) and prefer the variance-weighted consistency pass, which is far less prone to sign flips than equal-share spreading.

Frequently Asked Questions

Why does an adaptive quadtree beat a uniform grid on skewed data?

A uniform grid fixes its cells before seeing the data, so it spends the same noise budget on a dense downtown tile and an empty rural tile. On skewed distributions that wastes most of the budget on cells with near-zero signal while leaving the dense core under-resolved. The quadtree spends its finest-level budget only where a noisy count says density is worth resolving, so error concentrates where the data actually is and range queries over dense regions come back far more accurate for the same total epsilon.

Does the split decision consume extra privacy budget?

No. The statistic we threshold to decide whether to subdivide is the same Laplace-noised count we already released for that node. Because it is post-processing on an already-private value, testing it against the threshold spends zero additional budget. The only budget spent is the one noisy count per node, and those counts compose sequentially along each root-to-leaf path.

Why give deeper tree levels more epsilon instead of splitting evenly?

Deep cells hold small counts, where a fixed noise magnitude does the most relative damage; the coarse root holds a huge count that tolerates more noise. Allocating budget geometrically with ratio equal to the cube root of the fanout (about 1.587 for a quadtree) gives deeper levels the tighter, larger-epsilon noise scale and provably minimises error for range queries. Uniform allocation is the better default only when your dominant query is a single point lookup.

What does consistency enforcement (boosting) actually improve?

Each node is noised independently, so a parent's count and the sum of its children disagree even though they estimate related quantities. Constrained inference redistributes that disagreement to make every parent equal the sum of its children. It spends no budget because it is pure post-processing, yet it reduces variance by exploiting the redundancy between levels, producing more accurate and internally consistent counts for downstream range aggregation.

Up: Differentially Private Spatial Aggregation · Differential Privacy for Geospatial Data