Setting Decimal Precision for Disease Coordinate Mapping

This guide solves one narrow operational problem: enforcing a fixed, defensible number of decimal places on case coordinates at ingestion, so that spatial resolution is set by policy rather than inherited from whatever a geocoder happened to emit. It sits under Precision Standards in Epi-Mapping, within the broader Spatial Epidemiology Fundamentals & Data Standards section.

Problem Context & Constraints

The naive approach — let coordinates flow through the pipeline at whatever precision the geocoder returns, then format them at export with an f-string — fails on three counts that are specific to surveillance data, not generic GIS.

First, over-precision leaks identity. A geocoder routinely returns eight decimal places of longitude. At a typical mid-latitude, the eighth decimal resolves to roughly a millimetre; the sixth to about 0.1 m. A point published at six or more decimals can be reverse-geocoded back to a single front door, which is exactly the singling-out that HIPAA Safe Harbor and GDPR pseudonymization are written to prevent. Formatting at export does not help: by then the over-precise float has already been written into spatial-join keys, intermediate Parquet files, and the audit trail, so the disclosed resolution and the analyzed resolution no longer match.

Second, Python’s native round() is the wrong tool. It applies banker’s rounding (round-half-to-even) and operates scalar-by-scalar in a Python loop. On a few hundred thousand cases that is slow, but the real defect is topological: rounding the x and y of a vertex independently, with a half-to-even mode that differs from the rest of the pipeline, can detach a point that was meant to coincide with a polygon boundary or another case’s location. Deterministic spatial joins then silently miss matches.

Third, a single global “round to N decimals” rule under-specifies the control. The ground distance subtended by a decimal place of longitude is latitude-dependent — about 1.1 m per five decimals at the equator, but only ~0.8 m at latitude 45°. Choosing N therefore requires reasoning about the study extent’s latitude and the analytical scale you actually need, not copying a number from a tutorial.

The fix is to treat decimal precision as a vectorized, CRS-asserted ingestion gate that runs once, immediately after coordinate extraction and before any Coordinate Reference Systems for Public Health reprojection, spatial index, or join — and to record the precision it applied for audit.

Choosing N: ground resolution per decimal place

The displacement a given decimal place can introduce is $111{,}320 \cdot \cos(\phi) \cdot 10^{-n}$ metres in longitude and $\approx 110{,}540 \cdot 10^{-n}$ metres in latitude, where $\phi$ is the latitude of the study area and $n$ the decimal place. The table below uses equator values as the worst case; multiply longitude figures by $\cos(\phi)$ for your extent.

The disclosed precision must always be coarser than or equal to the precision analyzed internally: you may compute hotspots at 5 decimals and publish at 3, never the reverse.

Prerequisites

# Python 3.11+
# geopandas==0.14.4
# shapely==2.0.4
# pyproj==3.6.1
# numpy==1.26.4
#
# Input state required before this step runs:
#   - A GeoDataFrame of case locations in EPSG:4326 (WGS84 geographic).
#     Decimal-degree precision control is only meaningful in a geographic CRS;
#     projected CRSs require metre-based snapping, handled separately.
#   - Geocoder positional accuracy known per record (or a conservative global
#     bound), so the chosen tolerance is never finer than the input's own error.
#   - Coordinate extraction complete; this gate runs before reprojection,
#     spatial indexing, joins, or any clustering method.

Two input facts must be settled before writing any code. The CRS must be geographic (EPSG:4326) so that “decimals” maps to ground distance through the formula above; if your cases arrived projected, snap in metres and convert the requirement accordingly. And the target tolerance must be coarser than the geocoder’s own error — rounding to 1.1 m resolution when the geocoder reports ±50 m manufactures false certainty that biases downstream confidence intervals.

Step-by-Step Solution

The function below replaces per-row formatting with a single deterministic, NumPy-backed rounding pass. It asserts the CRS, drops null/empty geometries that would propagate NaNs, restricts to Point geometries (polygon and line vertices need shapely.ops.transform), rounds x and y together, rebuilds geometry with the CRS re-asserted, and stamps the precision it applied onto the frame’s metadata for the audit trail.

import logging
from typing import Optional

import geopandas as gpd
import numpy as np

logger = logging.getLogger(__name__)


def enforce_coordinate_precision(
    gdf: gpd.GeoDataFrame,
    precision: int = 5,
    audit_tag: Optional[str] = None,
) -> gpd.GeoDataFrame:
    """Enforce deterministic decimal precision on WGS84 point coordinates.

    Validates CRS, filters invalid geometries, applies vectorized rounding,
    and attaches compliance metadata for the audit trail. Point geometries
    only; for polygons/lines use shapely.ops.transform with a rounding fn.
    """
    if gdf.crs is None or gdf.crs.to_epsg() != 4326:
        raise ValueError(
            "Input must be EPSG:4326 (WGS84) for decimal-degree precision control."
        )

    # Drop null/empty/invalid geometries before rounding so NaNs do not propagate.
    valid = gdf.geometry.notna() & gdf.geometry.is_valid & ~gdf.geometry.is_empty
    dropped = int((~valid).sum())
    if dropped:
        logger.warning("Dropping %d invalid/empty geometries before rounding.", dropped)
    gdf = gdf[valid].copy()

    # Restrict to points; polygon/line vertex rounding needs shapely.ops.transform.
    if not (gdf.geom_type == "Point").all():
        raise ValueError(
            "Point geometries only. Use shapely.ops.transform for polygons/lines."
        )

    # Round x and y together, deterministically, with NumPy's fixed half-even mode
    # applied uniformly (not Python's scalar round) so vertices stay coincident.
    coords = np.column_stack([gdf.geometry.x.to_numpy(), gdf.geometry.y.to_numpy()])
    rounded = np.round(coords, decimals=precision)

    gdf["geometry"] = gpd.points_from_xy(rounded[:, 0], rounded[:, 1])
    gdf = gdf.set_crs("EPSG:4326")

    # Stamp precision provenance for the audit trail.
    gdf.attrs["coordinate_precision"] = precision
    gdf.attrs["rounded_record_count"] = len(gdf)
    if audit_tag:
        gdf.attrs["compliance_tag"] = audit_tag

    logger.info(
        "Applied precision=%d to %d points (tag=%s).", precision, len(gdf), audit_tag
    )
    return gdf

The key choices are deliberate. np.round is vectorized and applies one rounding mode across the whole array, so coincident vertices stay coincident; rounding x and y from the same column_stack keeps the pair atomic; and set_crs (not to_crs) re-asserts the geographic CRS on the rebuilt geometry without reprojecting. The precision and record count land in gdf.attrs, which survives into the GeoParquet sidecar metadata when you serialize — see Converting Shapefiles to GeoJSON for Epi Pipelines for carrying these attributes through format conversion.

Validation & Edge Cases

Precision enforcement must be paired with a post-condition check, because silent topology degradation is the failure that does not raise.

1. Displacement exceeds tolerance. After rounding, confirm no point moved farther than the resolution the chosen precision permits. A point that jumps more than 0.5 * 10**(-precision) degrees indicates the input was already coarser, or a units mix-up.

import numpy as np

orig = np.column_stack([before.geometry.x, before.geometry.y])
new = np.column_stack([after.geometry.x, after.geometry.y])
max_shift_deg = np.abs(orig - new).max()
tol = 0.5 * 10 ** (-precision)
assert max_shift_deg <= tol + 1e-12, (
    f"Rounding shifted a point by {max_shift_deg:.3e}° > tolerance {tol:.3e}°"
)
# Log line on a clean run:
# INFO Applied precision=5 to 48213 points (tag=county-q2-release).

2. Polygon and line inputs silently centroided. If facility footprints or jurisdictional boundaries reach this function, the Point-only guard raises rather than collapsing them to centroids. Round their vertices instead:

from shapely.ops import transform

def round_vertices(geom, precision=5):
    return transform(lambda x, y, z=None: (round(x, precision), round(y, precision)), geom)

3. NaN propagation from empty geometries. A single empty geometry left in the frame turns .x/.y into NaN and poisons the whole rounded array. The valid mask removes them up front; the WARNING Dropping N invalid/empty geometries log line is the signal that your upstream geocoder is emitting failures you should reconcile against the case count, not ignore.

Compliance Notes

For regulatory defensibility, the precision gate must leave a record that ties the disclosed resolution to a policy decision. Log, per run: the input CRS (must be EPSG:4326), the applied precision, the record count before and after the invalid-geometry drop, the geocoder accuracy bound that justified the chosen N, and the audit_tag naming the release tier (e.g. public-3dp vs dua-5dp). The HIPAA Safe Harbor boundary — coarsen to the point where a coordinate can no longer single out an individual — should be asserted in code, not assumed: a public release path forces precision <= 3, while a 4–5 decimal path is permitted only behind a data use agreement. Persist these fields in the same metadata record that captures de-identification lineage; the schema for that lineage is defined in Building a HIPAA-Compliant Spatial Metadata Schema, and the surrounding obligations in Compliance Mapping Frameworks.

Frequently Asked Questions

Related Topics

• Precision Standards in Epi-Mapping — the parent guide on tolerance bands, CRS-aware rounding, and precision manifests.
• Coordinate Reference Systems for Public Health — enforce the canonical CRS this gate assumes before any distance calculation.
• Building a HIPAA-Compliant Spatial Metadata Schema — where the precision and de-identification provenance gets recorded.