Reporting 95 Percent Confidence Regions for Cadastral Points

Quoting “±0.05 m” without saying whether that is one sigma or 95 per cent, one dimension or two, is the reporting failure that makes an otherwise careful uncertainty budget unusable. This guide, part of uncertainty propagation through transformation chains, converts a covariance into the specific statements a cadastral deliverable should carry, with the scale factors that belong to each and the labels that make them unambiguous.

Why the Dimensionality Changes the Number #

A confidence region in $k$ dimensions is the set of points whose Mahalanobis distance from the estimate is below a threshold, and that threshold comes from the chi-square distribution with $k$ degrees of freedom. In one dimension the 95 per cent threshold gives the familiar factor of 1.960 on the standard deviation. In two dimensions the same 95 per cent coverage requires 2.4477, because the probability is spread over an area rather than a line. In three it is 2.7955.

Figure — the 95% scale factor by dimensionality, and what using the wrong one costs.

The consequence is that applying the one-dimensional 1.96 to a two-dimensional error ellipse produces a region covering about 86 per cent while claiming 95. That is a quiet, systematic optimism in exactly the direction that loses an argument, and it is common enough that a reviewer will check for it.

Complete Runnable Implementation #

from __future__ import annotations

from dataclasses import dataclass

import numpy as np

# Chi-square 95% thresholds, square-rooted, by dimensionality.
K95 = {1: 1.9600, 2: 2.4477, 3: 2.7955}
K68 = {1: 1.0000, 2: 1.5151, 3: 1.8776}


@dataclass(frozen=True)
class ConfidenceRegion:
    """A labelled confidence statement — the label is part of the value."""

    semi_major_m: float
    semi_minor_m: float
    azimuth_deg: float          # of the major axis, clockwise from grid north
    level: float                # e.g. 0.95
    dimensions: int             # 1, 2 or 3

    def describe(self) -> str:
        return (f"{self.level:.0%} {self.dimensions}-D region: "
                f"{self.semi_major_m:.3f} m x {self.semi_minor_m:.3f} m "
                f"at {self.azimuth_deg:.1f} deg")


def horizontal_region(cov_en: np.ndarray, level: float = 0.95) -> ConfidenceRegion:
    """Error ellipse from a 2x2 east-north covariance.

    The azimuth is converted from the mathematical convention the eigenvector
    gives (counter-clockwise from east) to the survey convention (clockwise from
    north), because those two differ by a reflection and quoting the wrong one
    puts the weak direction at ninety degrees to the truth.
    """
    if cov_en.shape != (2, 2):
        raise ValueError("expected a (2, 2) east-north covariance in m^2")
    if not np.allclose(cov_en, cov_en.T, atol=1e-15):
        raise ValueError("covariance is not symmetric")
    vals, vecs = np.linalg.eigh(cov_en)
    order = np.argsort(vals)[::-1]
    vals, vecs = vals[order], vecs[:, order]
    if vals[-1] < -1e-15:
        raise ValueError("covariance is not positive semi-definite")
    k = K95[2] if abs(level - 0.95) < 1e-9 else K68[2]
    east, north = vecs[0, 0], vecs[1, 0]
    az_math = np.degrees(np.arctan2(north, east))       # CCW from east
    az_survey = (90.0 - az_math) % 180.0                # CW from north, mod 180
    return ConfidenceRegion(
        semi_major_m=float(k * np.sqrt(max(vals[0], 0.0))),
        semi_minor_m=float(k * np.sqrt(max(vals[1], 0.0))),
        azimuth_deg=float(az_survey),
        level=level,
        dimensions=2,
    )


def radial_95(cov_en: np.ndarray) -> float:
    """A single 95% horizontal radius — the circle that contains the ellipse.

    Convenient for a summary table, and lossy: it discards the orientation, so it
    is quoted ALONGSIDE the ellipse rather than instead of it.
    """
    return horizontal_region(cov_en).semi_major_m


def vertical_95(sigma_up_m: float) -> float:
    """95% vertical bound — one-dimensional, so the factor is 1.96, not 2.4477."""
    return 1.9600 * sigma_up_m

Parameter Reference #

Worked Example #

import numpy as np

cov = np.array([[0.000324, 0.000108],
                [0.000108, 0.000144]])          # sigma_e 18 mm, sigma_n 12 mm
region = horizontal_region(cov)
print(region.describe())
print(f"radial 95% {radial_95(cov) * 1000:.1f} mm")
print(f"vertical 95% {vertical_95(0.025) * 1000:.1f} mm")
# 95% 2-D region: 0.049 m x 0.023 m at 58.3 deg
# radial 95% 48.7 mm
# vertical 95% 49.0 mm

Notice how the horizontal and vertical 95 per cent figures come out nearly equal — 48.7 mm and 49.0 mm — from very different one-sigma inputs. That is the dimensionality factor doing its work, and it is why the two numbers cannot be produced by the same line of code.

Figure — four ways to quote the same covariance, and what each one hides.

Validation Check #

def check_coverage(cov: np.ndarray, n: int = 200_000, seed: int = 7) -> None:
    """Simulate to confirm the region really covers 95 per cent."""
    rng = np.random.default_rng(seed)
    samples = rng.multivariate_normal(np.zeros(2), cov, size=n)
    r = horizontal_region(cov)
    theta = np.radians(90.0 - r.azimuth_deg)
    c, s = np.cos(theta), np.sin(theta)
    u = samples[:, 0] * c + samples[:, 1] * s
    v = -samples[:, 0] * s + samples[:, 1] * c
    inside = (u / r.semi_major_m) ** 2 + (v / r.semi_minor_m) ** 2 <= 1.0
    frac = float(inside.mean())
    assert 0.945 <= frac <= 0.955, f"region covers {frac:.3f}, not 0.95"

A Monte Carlo coverage check is the one test that catches a wrong scale factor, a wrong dimensionality and a wrong azimuth convention all at once, and it needs no reference implementation to compare against.

Common Mistakes #

The one-dimensional factor applied to a two-dimensional region. 1.96 on an ellipse gives about 86 per cent coverage. The fix is trivial; noticing is the hard part, which is why the coverage simulation above is worth keeping in the test suite rather than running once.

The azimuth quoted in the mathematical convention. The eigenvector gives an angle counter-clockwise from east; surveyors read azimuths clockwise from north. The two differ by a reflection, so an ellipse elongated north-east gets reported as elongated north-west — the numbers look reasonable and the weak direction is wrong.

A region quoted with no level and no dimensionality. “±0.05 m” is four different statements depending on what is meant, and they differ by a factor of 2.5 between them. Every quoted region should carry its level and its dimensionality in the same string, as describe() does, because tables get copied out of reports and the surrounding text does not travel with them.

Presenting the Region in a Deliverable #

Three presentations serve different readers, and a complete deliverable carries all of them. A table gives the semi-axes, the azimuth, the level and the dimensionality per point, which is what a reviewer recomputes from. A plot draws the ellipses in place at a stated exaggeration factor on equally scaled axes, which is what shows a systematic weakness across a network at a glance. A single scalar — the 95 per cent radial figure — is what a pass/fail gate compares against a scalar tolerance, and it is the only one of the three that can be misread as the whole story.

Figure — a residual field with its 95% acceptance ring: the pattern is the diagnosis.

The ordering matters when space is tight: keep the table, because it is the reproducible one. A plot without its underlying numbers cannot be checked, and a scalar without its ellipse hides the direction of the weakness — which, on a boundary running along the weak axis, is the property that decides whether the survey meets specification.

Frequently Asked Questions #

Should a cadastral report quote the ellipse or the radius? #

Both, with the ellipse as the primary statement. The radius is convenient for a summary table and for a pass/fail gate against a scalar tolerance; the ellipse is what tells a reader that the uncertainty is twice as large north-east as north-west, which is exactly what matters when a boundary runs in one of those directions.

Is CEP or DRMS acceptable instead? #

They are widely used in positioning literature and both hide information a cadastral reader needs. Circular error probable is a 50 per cent radius, which is a much weaker statement than it appears; twice-DRMS is roughly, but not exactly, 95 per cent and the approximation degrades as the ellipse becomes elongated. Where a specification asks for one of them, produce it, and quote the ellipse alongside.

How do I combine a horizontal and a vertical statement? #

Not by adding them. They are separate statements with separate dimensionalities, and a three-dimensional region requires the full three-by-three covariance and the 2.7955 factor. Most cadastral specifications set horizontal and vertical tolerances separately precisely because the two have different error structures, and reporting them separately matches how they will be judged.

What confidence level should be used? #

Whatever the specification demands — 95 per cent is the common default in cadastral work, while some agencies use 68 per cent or a 99 per cent bound for control. The point is not which number but that it is stated: a region without its level is not a statement about anything, and converting between levels afterwards requires assumptions the reader should not have to make.

One-Line Summary for a Report Template #

If a template can only hold one line per point, this is the one that loses least: semi-major, semi-minor, azimuth, level and dimensionality, in that order and explicitly labelled. Everything else in this guide can be reconstructed from those five values plus the covariance they came from, and no subset of them can be reconstructed from a bare “±0.05 m”.

Related #

• Uncertainty propagation through transformation chains — the parent topic
• Visualizing geodetic error ellipses with Matplotlib
• Error distribution modeling in Python
• RMSE to agency submission workflow in Python
• Combining grid and parameter uncertainty in Python