Skip to content

Evaluation

Trajectory evaluation metrics.

Definitions follow the conventions established by Sturm et al., "A Benchmark for the Evaluation of RGB-D SLAM Systems" (IROS 2012) and by the TUM VI benchmark (Schubert et al., IROS 2018), because those are the numbers a reader will compare against. They are implemented here from the published definitions.

Absolute pose error (APE / "ATE") For every estimate pose the reference trajectory is interpolated at the estimate's timestamp and the position difference is taken. APE is the norm of that difference. RMSE over all poses is the headline number. The alignment mode is reported alongside every value, because "ATE = 0.09 m" is meaningless without it.

Relative pose error (RPE) Error in the relative transformation between two poses separated by a fixed time delta seconds or a fixed travelled distance delta metres. The translation error is ||(Q_i^-1 Q_j) - (P_i^-1 P_j)|| on the translational part; the rotation error is the geodesic angle angle(R_i^T R_j).

Drift Reported two ways, because they answer different questions:

``drift_pct_path_length``
    ``100 * APE(t) / path_length_gt(t)``. "How much error per metre
    travelled" -- comparable across sequences of different size.
``drift_pct_time``
    ``100 * APE(t) / t``. "How much error per second" -- what a
    navigation system actually cares about for a time-to-fix estimate.
Cross-check

tests/test_trajectory_io.py contains a TUM VI regression test that reads the published room1/512/16 estimate and ground truth and asserts ATE against the benchmark's own figure. That is a real validation of this implementation against a third party, not a self-consistency check.

Summary of a set of scalar errors.

Absolute trajectory error under a specific alignment.

Error growth relative to travelled distance and elapsed time.

Least-squares slope of error against time, in metres per second.

Reported over the second half of the sequence only. A single slope over the whole run is dominated by the initial convergence transient, which says nothing about steady-state drift.

RMSE / mean / median / std / min / max of a scalar error sequence.

Greedy nearest-timestamp association, the TUM benchmark convention.

For every estimate pose the closest reference pose within max_difference_s is selected. This is what evaluate_ate.py from the TUM RGB-D and TUM VI benchmark tools do, and reproducing it is what allows a direct comparison with published numbers.

Returns index arrays into estimate and reference. An estimate pose may be associated with the same reference pose as its neighbour; that is the behaviour of the original tool and is kept for comparability.

Why not interpolate? Interpolating the reference at the estimate's timestamp is the better estimator when the reference is smooth and sampled densely, and it is the default here. Association is offered because it is what published numbers are computed with, and because it never invents a reference pose. See docs/metrics.md.

Fit an alignment, returning it with the reference stream and the mask.

Only the positions are used to fit the transform. Fitting on 6-DoF poses would require a pose-graph solve for a result that, for a rigid sensor platform, is identical to the least-squares rigid fit on positions; that is the standard choice in the SLAM evaluation literature and it is stated here so the number can be reproduced.

The second element is a :class:Trajectory, not a position array: under association="nearest" the caller needs ref.metadata["index"] to map back to the estimate's own samples, which an array would not carry. The annotation previously claimed np.ndarray and was simply wrong.

ATE for one alignment choice.

alignment="none" reports the raw error with no transform applied, which is the honest number when the estimate's frame is already the reference frame (the case for a GNSS-aided filter). "rigid" and "similarity" fit and apply a similarity transform first.

The four ATE variants reported for every run.

none -- no alignment, raw error in the estimator's own frame. rigid -- rigid alignment over the whole sequence. rigid_start -- rigid alignment fitted on the first 20% only. The strictest of the four, and the one that cannot absorb accumulated drift into the fit. This mirrors the TUM VI ate_start metric. similarity -- alignment with a free scale, for scale-free estimates.

The rigid_start result is fitted with kind="rigid" but is labelled "rigid_start", because an ATE reported as "rigid" while having been fitted on 20% of the sequence is a different number under the same name.

RPE at a fixed time offset or a fixed travelled distance.

mode="time" pairs every pose with the one delta seconds later. mode="distance" pairs every pose with the first pose at least delta metres further along the reference path.

A single global alignment cannot be applied to a relative error: the definition is invariant to it, and the error that matters is the change in relative pose. Poses beyond the end of the sequence are excluded rather than clamped, because clamping would compare a pose with itself.

Drift statistics from an ATE result and the reference path.

Uses the rigid_start-style error curve: the alignment must have been fitted on the start of the sequence, otherwise "drift" measures the alignment, not the estimator. The caller passes the ATE variant it wants.

Seconds from the end of window until the error settles below a threshold.

window is (start_s, end_s) in trajectory time. hold_s requires the error to stay below the threshold for that long, so a single lucky sample does not count as recovery. Returns None if it never recovers.

Per-outage error growth, using baseline_error_m as the reference.

For each outage window the report contains the error at the moment the outage begins, the peak error inside it, the error when it ends, and the growth across the window. Growth measured against the pre-outage error rather than against zero is what makes a GNSS outage interpretable: an outage of a well-localised system should show a rise from a small baseline, not a large absolute number.

Convenience: timestamps, aligned positions and per-pose errors.

Is the estimator's uncertainty honest?

A pose estimate without a trustworthy uncertainty is not much use for flight planning. An aircraft is allowed to be lost; it is not allowed to be lost quietly. This module measures whether the covariance an estimator reports actually describes the error it makes, and provides the machinery to correct it when it does not.

The claim this project makes is narrow and falsifiable: a filter tuned for accuracy is not thereby calibrated, and the gap is measurable, correctable, and distribution-dependent. Each function below exists so that the claim can fail loudly rather than quietly.

What "calibrated" means here

If the reported covariance is correct, the squared Mahalanobis distance of the error, e^T P^-1 e, follows a chi-square distribution with as many degrees of freedom as the state dimension. That gives two independent ways to check:

  • the mean normalised error should equal the dimension;
  • the coverage of the claimed ellipsoid should equal P(chi2_dof <= r**2).

Coverage is the more useful of the two for a planner, because it is directly interpretable: it answers "how often is the aircraft inside the box I drew?" The per-axis sigma and the coverage are not interchangeable, and :mod:navkit.eval.statistics quantifies by how much. Every coverage number in this module is stated against the chi-square reference, never against the nominal percentage it was derived from.

Three levels of rigour

CoverageReport Descriptive. What fraction of epochs fell inside the claimed ellipsoid, with a Wilson interval so that a coverage measured on 400 epochs is not reported to more precision than it has.

inflation_factor Corrective. The scalar by which the covariance must be multiplied for observed coverage to reach the target. One global scalar, so it cannot fix a filter that is miscalibrated only in some directions or only in some regimes -- which is the usual case. Report it; do not apply it silently.

conformal_radius Distributional. A finite-sample bound on the normalised error derived from an order statistic of a calibration set, valid at a stated confidence under exchangeability with that set. This is the only one of the three carrying a guarantee a reader can check, and the guarantee is conditional. The behaviour under distribution shift is measured by :func:shift_sensitivity rather than assumed away.

Per-epoch normalised error between an estimate and a reference.

Attributes:

Name Type Description
squared ndarray

e^T P^-1 e per epoch. Under a correct covariance this is chi-square with dof degrees of freedom, so its mean should be dof.
radius ndarray

sqrt(squared), the quantity whose distribution is chi with dof degrees of freedom. This is the object the conformal bound is built on.
dof int

Dimension of the state being scored.
epochs int

How many epochs contributed. Zero-length input is allowed and yields an empty series, so callers can pass a filtered sequence without guarding.

Empirical q-quantile of the normalised radius.

Observed versus expected coverage at one per-axis sigma.

Observed minus expected. Negative means the filter is overconfident.

Calibration verdict from a two-sided test on the coverage.

The question is whether the expected coverage lies inside the Wilson interval of what was observed. Comparing the observed and expected values directly, as a sign test would, flags almost every correct filter as miscalibrated, because an observed coverage essentially never equals its expectation exactly.

Coverage of an estimator at a range of per-axis sigmas.

Everything measured about one estimator's uncertainty on one dataset.

Verdict from the bulk of the error distribution, via NEES.

The coverage verdict above is a tail test: it asks whether errors ever escape the k-sigma ellipsoid, and it uses a Wilson interval because an observed coverage essentially never equals its expectation exactly. That alone is not enough to describe a filter, because the two diagnostics detect different failure modes and can disagree.

A filter can sit slightly above its covariance almost everywhere (bulk NEES above dof, meaning marginally overconfident) while never once producing an error outside the 2-sigma ellipsoid (coverage at nominal). Calling that "underconfident" because the tail looks clean, or "overconfident" because NEES is high, would each discard half the evidence. Both are reported.

The test is the standard chi-square one: NEES * n is distributed as chi2(n * dof) when the filter is correctly calibrated. Polarity follows the mean: a NEES above dof means the real errors are larger than the covariance admits, so the filter is overconfident; below it, the filter is underconfident. That is the same sense as :attr:CoveragePoint.gap, where too few covered epochs also means overconfidence.

Verdict from the tail, i.e. from coverage at 2 sigma per axis.

Combined verdict, naming a disagreement rather than hiding it.

bulk and tail disagreeing is itself the finding: it means the error distribution is misshapen in a way neither number describes on its own, and that is a different defect from a simple miscalibration.

True only when the bulk and the tail both agree with expectation.

This is stricter than either test alone, which is the intent: a filter counts as calibrated only if it is neither systematically overconfident in the bulk nor producing tail escapes.

Wilson score interval for a binomial proportion.

Preferred over the normal approximation, which is badly behaved at the coverage values that matter here: a filter that claims 95% coverage and achieves 0% is exactly the case where a Wald interval produces nonsense.

Score an estimate against a reference given its position covariance.

Both trajectories are compared index by index, so the caller is responsible for resampling the reference onto the estimate's timestamps. Covariance blocks that are numerically singular are handled by symmetric regularisation rather than by dropping the epoch, because silently discarding the worst epochs is precisely how a miscalibrated filter gets to look well behaved.

Measure coverage at several per-axis sigmas.

Note that series.squared is the squared Mahalanobis distance, so a k-sigma-per-axis ellipsoid is the set squared <= k**2 * dof.

Scalar covariance inflation that reaches target_coverage.

Returns inf when the target is unreachable at any inflation, which happens when the filter is underconfident to a degree that no global scale can repair. Reporting infinity is the honest answer; clipping to a large finite number would disguise an unbounded calibration error as a merely large correction.

Finite-sample bound on the normalised error under exchangeability.

The bound is the order statistic at index ceil((n + 1) * confidence) / n of the calibration radii. The index, not the plain confidence quantile, is what makes the guarantee hold for finite n; using the plain quantile is a small and very common error that silently voids the guarantee. When the index exceeds n the target cannot be certified from this sample and the bound is reported as infinity rather than extrapolated.

How much of a calibrated guarantee survives a change of distribution.

The conformal bound is only valid for data exchangeable with the set it was built from. Applying it after a change in noise level, outage profile, or manoeuvre regime is exactly the situation a planner cares about, and the honest thing is to measure how far the realised coverage falls below the nominal one instead of carrying the nominal number forward.

Build a :class:CalibrationReport from a scored sequence.

Exact chi-square distribution functions for integer degrees of freedom.

Navigation metrics need the chi-square distribution in two places: to know what coverage a claimed ellipsoid should have, and to invert that relation when calibrating. Both are needed only at integer degrees of freedom (one per scalar axis, three for a position vector, four for a pose), and at integer degrees of freedom the distribution has a closed form.

For a chi-square variate with n degrees of freedom, X / 2 follows a gamma variate of shape n / 2, and the regularised lower incomplete gamma function P(a, t) satisfies the half-step recurrence

P(a + 1, t) = P(a, t) - t**a * exp(-t) / Gamma(a + 1)

So every integer n is reached from one of two exact seeds,

P(1/2, t) = erf(sqrt(t))            P(1, t) = 1 - exp(-t)

by applying the recurrence. That makes the coverage reference exact rather than tabulated, with no SciPy dependency and no interpolation error near the tails where the numbers actually matter.

The quantile function is obtained by bisection on the CDF. It is slower than a specialised implementation and completely unambiguous, which is the right trade for a function called a few thousand times per report.

P(X <= x) for a chi-square variate with dof degrees of freedom.

Inverse of :func:chi2_cdf by bisection.

q is a probability in [0, 1]. Bisection is used rather than a Newton iteration because the CDF is monotone and cheap to evaluate here, so there is no convergence or bracketing failure mode to guard against.

Coverage of an isotropic k-sigma-per-axis ellipsoid in dof dims.

This is where the usual shorthand breaks down. The familiar "3-sigma means 99.7%" is a statement about a one-dimensional Gaussian. Applied per axis to a three-dimensional position ellipsoid, the same 3-sigma threshold gives a squared Mahalanobis radius of 3 * 3**2 = 27, and a coverage of P(chi2_3 <= 27) -- a little over 0.999, not 0.997. Worse, the habit usually applied is "2 sigma means 95% per axis", which for three axes means P(chi2_3 <= 12) = 0.9926 of samples fall inside, while a reader who expects 95% is being told something twice as strict as they think.

Quoting the coverage instead of the per-axis sigma removes the ambiguity.

Per-axis sigma whose isotropic ellipsoid achieves coverage.

The inverse of :func:ellipsoid_coverage, and the form in which a calibrated bound should be stated.

Failure detection: turning an error curve into a list of events.

A metric is a number; an event is something an engineer has to act on. This module converts the per-pose error curve plus a set of thresholds into a list of discrete failure events with start, end, peak and duration.

Event types

threshold_exceeded The error stayed above an absolute threshold for longer than min_duration_s. The default threshold is scenario-independent and set in the config; a run that produces no event is a run that stayed inside the declared quality envelope, which is a meaningful statement.

divergence The error grew monotonically by more than divergence_growth_m over divergence_window_s without recovering. Monotonic growth is what separates divergence from a transient excursion, and it is the failure mode that matters most for a navigation system.

rejected_updates The estimator's innovation gate rejected a run of consecutive measurements. This is a filter-health event rather than a pose event: it says the filter stopped trusting its own measurements, which precedes divergence rather than following it.

Thresholds are explicit configuration, never hard-coded, and every event records the threshold that produced it so a report can be audited.

Quality thresholds for a run, all in metres and seconds.

One detected failure, with the threshold that produced it.

Detect failure events in an error curve.

t must be relative to the start of the sequence (the caller is responsible for that), so the event times are directly comparable to the scenario's outage windows.