Measuring Geocoding Match-Rate and Precision in Python

TL;DR: Match-rate is the fraction of ground-truth cases a geocoder resolves within a distance threshold; pair it with median and 95th-percentile haversine error and a coverage figure, compute all four with a numpy metric function, and roll them up per provider with a pandas groupby. These are the numbers the regression corpus workflow gates on.

Match-rate on its own is a misleading headline. A geocoder that resolves 99% of inputs but places each result 300 metres off scores a high match-rate at a loose threshold and a terrible one at a tight threshold — the number is meaningless without the distance it was measured at, and without the error distribution behind it. This page defines the four metrics that together describe geocoding precision against a ground-truth set, and implements them in pure Python with numpy and in a vectorized pandas rollup.

The Four Metrics

Metric Definition Reads as Guards against
Match-rate (hit rate) Fraction of cases whose returned point is within hit_threshold_m of truth Overall accuracy at a stated tolerance Silent quality loss
Median error (m) 50th percentile of haversine error over resolved cases Typical positional accuracy Central drift
p95 error (m) 95th percentile of haversine error over resolved cases Tail behaviour on hard cases A degrading minority the median hides
Coverage Fraction of cases the provider resolved to any coordinate at all Recall, independent of accuracy A provider that returns nothing

The distinction between match-rate and coverage matters: coverage counts any answer, match-rate counts a correct answer. A provider can have high coverage and low match-rate (it always answers, often wrongly) or the reverse (it answers rarely but accurately). Report both, and always alongside median and p95 error so the shape of the distribution is visible, not just its centre.

The Metric Function

The function below takes arrays of returned coordinates and truth coordinates plus per-case thresholds, computes haversine error, and reduces to the four metrics. It is pure numpy with explicit type hints and treats an unresolved case (a NaN returned coordinate) as covered=False rather than a zero-distance hit.

from __future__ import annotations

from typing import TypedDict

import numpy as np

_EARTH_RADIUS_M = 6_371_000.0


class Metrics(TypedDict):
    n: int
    match_rate: float
    median_error_m: float
    p95_error_m: float
    coverage: float


def haversine_m(
    lat1: np.ndarray, lon1: np.ndarray, lat2: np.ndarray, lon2: np.ndarray
) -> np.ndarray:
    """Vectorized great-circle distance in metres between aligned coordinate arrays."""
    p1, p2 = np.radians(lat1), np.radians(lat2)
    dphi = np.radians(lat2 - lat1)
    dlmb = np.radians(lon2 - lon1)
    a = np.sin(dphi / 2) ** 2 + np.cos(p1) * np.cos(p2) * np.sin(dlmb / 2) ** 2
    return 2 * _EARTH_RADIUS_M * np.arcsin(np.sqrt(a))


def compute_metrics(
    got_lat: np.ndarray,
    got_lon: np.ndarray,
    truth_lat: np.ndarray,
    truth_lon: np.ndarray,
    threshold_m: np.ndarray,
) -> Metrics:
    """Compute match-rate and precision metrics for one provider against truth.

    Unresolved cases must arrive as NaN in got_lat/got_lon; they count against
    coverage and match-rate but are excluded from the error percentiles.

    Args:
        got_lat, got_lon:     Returned coordinates (NaN where unresolved).
        truth_lat, truth_lon: Ground-truth coordinates.
        threshold_m:          Per-case hit threshold in metres.

    Returns:
        A Metrics dict for the whole array.
    """
    got_lat = np.asarray(got_lat, dtype=float)
    resolved = ~np.isnan(got_lat) & ~np.isnan(np.asarray(got_lon, dtype=float))
    n = int(got_lat.size)
    if n == 0:
        return Metrics(n=0, match_rate=0.0, median_error_m=float("nan"),
                       p95_error_m=float("nan"), coverage=0.0)

    error = np.full(n, np.nan)
    error[resolved] = haversine_m(
        got_lat[resolved], np.asarray(got_lon)[resolved],
        np.asarray(truth_lat)[resolved], np.asarray(truth_lon)[resolved],
    )
    hit = resolved & (error <= np.asarray(threshold_m, dtype=float))
    resolved_err = error[resolved]

    return Metrics(
        n=n,
        match_rate=float(hit.sum() / n),
        median_error_m=float(np.median(resolved_err)) if resolved_err.size else float("nan"),
        p95_error_m=float(np.percentile(resolved_err, 95)) if resolved_err.size else float("nan"),
        coverage=float(resolved.sum() / n),
    )

The percentiles are computed only over resolved_err, so a provider is not rewarded for the low “error” of the cases it never answered. Match-rate and coverage, in contrast, use the full n as the denominator — an unanswered case is a failed case, not an excluded one.

Per-Provider pandas Rollup

When a single corpus run geocodes each case through several providers, store the results in long format — one row per (case, provider) — and let a groupby produce the comparison table in one pass. This is the shape that feeds a head-to-head vendor comparison as well as a single-provider gate.

import numpy as np
import pandas as pd

_EARTH_RADIUS_M = 6_371_000.0


def _haversine(df: pd.DataFrame) -> pd.Series:
    p1, p2 = np.radians(df["truth_lat"]), np.radians(df["got_lat"])
    dphi = np.radians(df["got_lat"] - df["truth_lat"])
    dlmb = np.radians(df["got_lon"] - df["truth_lon"])
    a = np.sin(dphi / 2) ** 2 + np.cos(p1) * np.cos(p2) * np.sin(dlmb / 2) ** 2
    return 2 * _EARTH_RADIUS_M * np.arcsin(np.sqrt(a))


def provider_metrics(df: pd.DataFrame) -> pd.DataFrame:
    """Build a per-provider metrics table from a long-format scored frame.

    Expects columns: provider, truth_lat, truth_lon, got_lat, got_lon,
    hit_threshold_m. Unresolved rows carry NaN in got_lat/got_lon.

    Returns:
        DataFrame indexed by provider with match_rate, median_error_m,
        p95_error_m, and coverage columns.
    """
    df = df.copy()
    resolved = df["got_lat"].notna() & df["got_lon"].notna()
    df["error_m"] = np.where(resolved, _haversine(df), np.nan)
    df["is_hit"] = resolved & (df["error_m"] <= df["hit_threshold_m"])

    return df.groupby("provider").agg(
        n=("is_hit", "size"),
        match_rate=("is_hit", "mean"),
        median_error_m=("error_m", "median"),
        p95_error_m=("error_m", lambda s: s.quantile(0.95)),
        coverage=("error_m", lambda s: s.notna().mean()),
    ).sort_values("match_rate", ascending=False)

Because is_hit is a boolean, mean() returns the match-rate directly, and error_m.notna().mean() gives coverage — the same denominator logic as the numpy function, expressed vectorized. The resulting table sorts providers by match-rate, ready to print in a CI log or hand to the comparing geocoding accuracy across providers methodology.

Edge Cases

  • All cases unresolved in a stratum. median and p95 return NaN when no case resolved. Handle NaN explicitly in the gate rather than letting it compare as False against a baseline — a NaN median means the provider answered nothing, which is itself a regression worth failing on.
  • Coordinate order swap. A transposed lat/lon produces uniform multi-thousand-kilometre errors across every case. Assert truth_lat and got_lat fall in [-90, 90] before scoring; a global error explosion is nearly always this, not a real regression.
  • Antimeridian and poles. The haversine formula is stable across the ±180° meridian and at the poles, but a provider that clips longitude to [0, 360) instead of [-180, 180) will show a false 40,000 km error near the antimeridian. Normalize returned longitude into [-180, 180) before computing distance.

Integration Note

These metrics are the measurement half of the regression corpus workflow; the sampling half — how the ground-truth cases and their per-case thresholds are assembled — is covered in building a stratified ground-truth address test set. The hit classification here assumes each returned coordinate is already vetted as trustworthy; that vetting, including confidence scoring and outlier flagging, is the job of validating geocoding accuracy and confidence scoring.