Probability: The Language of Uncertainty

Bayes rule as pointwise multiplication, the Gaussian as a shape you can drag, and the two parameterizations — moments and canonical — that every filter in this book trades between.

The problem: two numbers in, one number out

In Chapter 1 we put Rusty in a ten-cell corridor with doors at cells 1, 4 and 5, gave it a door detector that fires correctly 60% of the time at a door and 20% of the time at a blank wall, and told it nothing about where it started. The detector fired.

The reader was handed two numbers — a prior of $0.1$ per cell, a likelihood of $0.6$ or $0.2$ — and one answer: $0.1875$ at each door, $0.0625$ everywhere else. What was never justified was the recipe. Why multiply rather than average? Why is the answer not $0.6$? And what exactly is the division at the end doing?

Here is the whole computation, in one line per step. Multiply the prior by the likelihood, cell by cell:

Those ten numbers sum to $3(0.06) + 7(0.02) = 0.32$, not to one, so they are not yet a distribution — they are the right shape at the wrong scale. Divide by the total:

That is Bayes rule, and the $0.32$ is the only part that needs a name. It is the probability of seeing what we saw, averaged over everywhere we might have been, and it is a constant with respect to the cell index — which is why the whole normalization can be deferred to the end and written as a single symbol $\eta = 1/0.32$. Thrun, Burgard and Fox use $\eta$ this way throughout, and so does this book: compute the shape, normalize later.

Now make the corridor continuous. Rusty's odometry says it has travelled about 5 m, give or take 2 m. Its wall-range sensor says 6.5 m, and it is a better sensor. Multiplying ten numbers by ten numbers becomes multiplying one curve by another, and the widget below does exactly that, live.

Three things in that figure are worth more than the algebra that justifies them later.

The posterior is narrower than both inputs, always. Not narrower than the average of the two — narrower than the better of the two. Fusing an excellent estimate with a mediocre one still improves the excellent one, because two independent opinions genuinely carry more information than one.

The posterior mean is not a midpoint. It sits between the two means, but pulled toward whichever curve is sharper. When they are equally sharp it lands exactly halfway; when the sensor is ten times sharper, it lands almost on the sensor. This is precision weighting, and the weight is one over the variance.

In canonical form the arithmetic is addition. The ledger under the plot tracks two numbers, $\Omega = 1/\sigma^2$ and $\xi = \mu/\sigma^2$, and the fusion is nothing but $\Omega_1 + \Omega_2$ and $\xi_1 + \xi_2$. That observation looks like a curiosity here. In Chapter 6 it becomes the information filter, and in Chapter 15 it becomes the reason modern SLAM back-ends never store a covariance matrix at all.

Building intuition: a Gaussian is a shape you can drag

One distribution dominates this book, and it earns that place for three reasons that are all theorems: the product of two Gaussians is Gaussian, the marginal of a Gaussian is Gaussian, and a linear map of a Gaussian is Gaussian. A belief that starts Gaussian and only ever meets those three operations stays Gaussian forever, described by a mean vector and a covariance matrix and nothing else. That is the entire premise of the Kalman filter.

So the covariance matrix deserves to be understood as geometry before it is understood as algebra. Play with it first.

The figure shows the same distribution three ways at once, and each way answers a different question. The cloud answers what would I see if I drew from it — and notice that the 500 dots are never re-sampled: they are one fixed set of standard-normal numbers being pushed through a changing matrix, which is all that $x = \mu + Lz$ means. The ellipses answer where does the density live, at one, two, and 2.4477 standard deviations. The gray axes answer which directions are special: they are the eigenvectors of $\Sigma$, and along them the two coordinates stop interfering with each other.

The misconception this kills is a stubborn one. It is tempting to read the ellipse's half-widths as the standard deviations of $x_a$ and $x_b$. They are not, and the blue marginal curves painted along the two edges prove it: as $\rho$ sweeps and the ellipse swings through 45 degrees, those curves do not move by a pixel. Correlation changes how the two coordinates covary; it does not change how uncertain either one is on its own. What it does change is how much one of them can teach you about the other — and that is a different question, taken up in Derivation 6.

The mathematics

Notation

Every symbol below is Thrun-compatible; the manifold operators wait for Chapter 3.

Random variables, joints, and conditionals

A random variable $X$ takes values according to a distribution. If the value space is discrete we write $p(X = x)$, abbreviated $p(x)$, with $\sum_x p(x) = 1$ and $p(x) \ge 0$. If it is continuous we assume a probability density $p(x)$ with $\int p(x)\,dx = 1$. A density is not a probability: it can exceed one. Only its integral over a region is a probability, and this book will use "probability", "density", and "distribution" interchangeably wherever no confusion can result.

The joint $p(x, y)$ is the density of both events together. The conditional is defined, whenever $p(y) > 0$, as

which is worth reading as a statement about renormalization: slice the joint at $Y = y$, and scale the slice so it integrates to one again. That is the definition the entire chapter hangs on.

$X$ and $Y$ are independent when $p(x, y) = p(x)\,p(y)$, equivalently $p(x \mid y) = p(x)$: knowing $Y$ tells you nothing about $X$. They are conditionally independent given $Z$ when

Finally, the theorem of total probability reconstructs a marginal from a conditional by integrating out what you introduced:

Bayes rule

The quantity $\eta^{-1} = p(z)$ has a name worth remembering: the evidence. The filter treats it as housekeeping, but it is a live diagnostic: a surprisingly small value means the measurement was unlikely under everything you currently believe. Chapter 11 rejects outliers with it and Chapter 12 uses it to notice that a robot has been kidnapped.

Expectation and covariance

The expectation is the probability-weighted average, and it is linear:

Linearity holds whether or not the components of $X$ are independent, which is why it will survive every approximation in this book. The covariance is the expected outer product of the deviation from the mean:

Three properties follow immediately and get used constantly. $\Sigma$ is symmetric, since the outer product is. It is positive semi-definite, since for any vector $a$ we have $a\T \Sigma a = \E[(a\T(X - \mu))^2] \ge 0$ — a variance, and variances cannot be negative. And its diagonal entries are the variances of the individual coordinates while its off-diagonal entries measure how they move together. A covariance with a negative eigenvalue is not an unlucky covariance; it is a bug, and Chapter 6 spends real effort making sure filters never produce one.

The Gaussian

In one dimension,

and in $n$ dimensions, with $\mu \in \R^n$ and $\Sigma$ symmetric positive-definite,

The second is a strict generalization of the first. Everything interesting lives in the exponent's quadratic form: the prefactor exists only to make the integral one, which is precisely why $\eta$ can absorb it and why so much of what follows is bookkeeping about quadratics.

It is worth being explicit about the working definition this book uses: a Gaussian is a quadratic in the exponent. Multiplying densities adds quadratics; conditioning fixes some variables in a quadratic; a linear map substitutes into a quadratic. All three stay quadratic, so all three stay Gaussian, and the next five derivations are just careful applications of that one fact.

The product of two Gaussians

Canonical form, where multiplication is addition

Linear transformations

Iso-density contours are ellipses

The Mahalanobis distance introduced there is worth naming on its own, because it is the only sensible way to answer "is this measurement consistent with what I expected?". Euclidean distance cannot answer it: 30 cm of error is nothing along an axis where the filter is uncertain and catastrophic along one where it is not. Mahalanobis distance measures in standard deviations, and its square is the $\chi^2$ statistic you compare against a table.

Marginals and conditionals

Partition a jointly Gaussian vector as $x = (x_a, x_b)$, with

There are two entirely different ways to get down to $x_a$ alone, and confusing them is the single most common error in a first estimation course.

Entropy

Uncertainty deserves a scalar, if only so that a robot can be asked which action would teach me the most. Shannon entropy is that scalar:

in bits, or with $\ln$ for nats. A uniform belief over $n$ cells has $\log_2 n$ bits; a certain belief has zero. For continuous variables the same expression with an integral is the differential entropy, and it can be negative — a density concentrated in a region narrower than one unit has $H < 0$. That is not a contradiction, only a reminder that differential entropy measures spread relative to the units you chose, and that only differences of differential entropy are physically meaningful. Information gain, which is a difference, is perfectly well behaved.

Samples: the third representation

So far a distribution has been a formula or a pair of moments. There is a third option, and half this book runs on it: a bag of numbers drawn from the distribution. The law of large numbers says that for any well-behaved $f$,

with an error that shrinks like $\sigma_f/\sqrt{N}$ by the central limit theorem. Samples have one enormous advantage over moments — they can represent any shape, including the three-peaked belief of Chapter 5 that no Gaussian can express — and one enormous disadvantage, which is that $1/\sqrt{N}$ is a miserable rate.

Both halves of that widget are worth internalizing before Chapter 8 turns samples into a filter. Ten samples look nothing like a bell curve and yet nothing is wrong; the estimate is unbiased at every $N$ and merely noisy. And getting one more decimal digit of accuracy costs a hundred times the samples, forever, in every dimension. That is why particle counts are argued about, why importance sampling exists, and why nobody runs a naive particle filter over a six-dimensional pose.

Implementation in Rust

The type carries an invariant — the covariance is symmetric positive-definite — so the constructor is the only place that can create one, and it proves the invariant by factoring. The factor is then cached, and every query below is a triangular solve rather than an inverse.

use nalgebra::{Cholesky, Const, SMatrix, SVector};
use rand::Rng;
use rand_distr::StandardNormal;

/// A Gaussian in moments form: p(x) = N(x; mean, cov).
///
/// The Cholesky factor is cached because it is the only cubic-cost object here.
/// Given L with cov = L Lᵀ, the density, the Mahalanobis distance, the entropy
/// and the sampler are all O(n²) — and none of them ever forms cov⁻¹, which is
/// slower, less accurate, and silent when it goes wrong.
#[derive(Clone, Debug)]
pub struct Gaussian<const N: usize> {
    mean: SVector<f64, N>,
    cov: SMatrix<f64, N, N>,
    chol: Cholesky<f64, Const<N>>,
}

/// The one way to fail: a covariance that is not positive-definite. Usually a
/// filter bug three steps upstream, which is exactly why we refuse it here.
#[derive(Clone, Copy, Debug, PartialEq, Eq)]
pub struct NotPositiveDefinite;

impl<const N: usize> Gaussian<N> {
    pub fn new(
        mean: SVector<f64, N>,
        cov: SMatrix<f64, N, N>,
    ) -> Result<Self, NotPositiveDefinite> {
        // Symmetrize first. A covariance update that is symmetric on paper
        // drifts by a few ULPs in floating point, and Cholesky is unforgiving.
        let cov = (cov + cov.transpose()) * 0.5;
        let chol = Cholesky::new(cov).ok_or(NotPositiveDefinite)?;
        Ok(Self { mean, cov, chol })
    }

    pub fn mean(&self) -> &SVector<f64, N> { &self.mean }
    pub fn cov(&self) -> &SMatrix<f64, N, N> { &self.cov }

    /// ln|Σ| = 2 Σᵢ ln Lᵢᵢ. Free, given the factor — and it never overflows the
    /// way a product of eigenvalues does for a 30-dimensional SLAM state.
    pub fn ln_det(&self) -> f64 {
        2.0 * self.chol.l().diagonal().iter().map(|d| d.ln()).sum::<f64>()
    }

    /// d²(x) = (x − μ)ᵀ Σ⁻¹ (x − μ) by one triangular solve (Algorithm above).
    pub fn mahalanobis2(&self, x: &SVector<f64, N>) -> f64 {
        let d = x - self.mean;
        d.dot(&self.chol.solve(&d))
    }

    pub fn ln_pdf(&self, x: &SVector<f64, N>) -> f64 {
        let ln_2pi = std::f64::consts::TAU.ln();
        -0.5 * (self.mahalanobis2(x) + self.ln_det() + N as f64 * ln_2pi)
    }

    pub fn pdf(&self, x: &SVector<f64, N>) -> f64 {
        self.ln_pdf(x).exp()
    }

    /// H = ½ ln((2πe)ⁿ |Σ|), Derivation 7. Note μ is absent: entropy says how
    /// uncertain the robot is, never where it thinks it is.
    pub fn entropy(&self) -> f64 {
        let ln_2pi_e = (std::f64::consts::TAU * std::f64::consts::E).ln();
        0.5 * (N as f64 * ln_2pi_e + self.ln_det())
    }

    /// x = μ + L z, Derivation 4 read backwards.
    pub fn sample<R: Rng + ?Sized>(&self, rng: &mut R) -> SVector<f64, N> {
        let z = SVector::<f64, N>::from_fn(|_, _| rng.sample(StandardNormal));
        self.mean + self.chol.l() * z
    }

    /// Push forward through y = A x + b (Derivation 4). The const generics are
    /// doing real work: a 2×3 A applied to a 3-D state can only produce a 2-D
    /// Gaussian, and getting it wrong is a compile error rather than a panic.
    pub fn transform<const M: usize>(
        &self,
        a: &SMatrix<f64, M, N>,
        b: &SVector<f64, M>,
    ) -> Result<Gaussian<M>, NotPositiveDefinite> {
        Gaussian::new(a * self.mean + b, a * self.cov * a.transpose())
    }

    /// (μ, Σ) → (ξ, Ω). ξ comes from a solve; Ω has to be materialized, because
    /// that is what the canonical form *is*. Pay the O(n³) once, then add.
    pub fn to_canonical(&self) -> Canonical<N> {
        Canonical { xi: self.chol.solve(&self.mean), omega: self.chol.inverse() }
    }
}

/// p(x_a | x_b = β) for a 2-D joint (Derivation 6), written in scalars so the
/// Schur complement is visible. The gain `k` moves the mean *and* shrinks the
/// variance; in Chapter 6 it acquires a name and a subscript.
pub fn condition_second(g: &Gaussian<2>, beta: f64) -> Gaussian<1> {
    let (mu, s) = (g.mean(), g.cov());
    let k = s[(0, 1)] / s[(1, 1)];
    Gaussian::new(
        SVector::from([mu[0] + k * (beta - mu[1])]),
        SMatrix::from([[s[(0, 0)] - k * s[(0, 1)]]]),
    )
    .expect("the Schur complement of an SPD matrix is SPD")
}

The canonical form is the same distribution asking a different question, and its entire implementation is Derivation 3.

use nalgebra::{Cholesky, SMatrix, SVector};
use super::gaussian::{Gaussian, NotPositiveDefinite};

/// Canonical (information) form: p(x) ∝ exp(−½ xᵀΩx + xᵀξ).
///
/// Moments form answers "where is it and how wide". Canonical form answers
/// "how much do I know". They carry identical information; what differs is
/// which of the two Bayes-filter operations is cheap.
#[derive(Clone, Debug, PartialEq)]
pub struct Canonical<const N: usize> {
    /// ξ = Σ⁻¹ μ
    pub xi: SVector<f64, N>,
    /// Ω = Σ⁻¹ — sparse for a pose graph, dense for its covariance.
    pub omega: SMatrix<f64, N, N>,
}

impl<const N: usize> Canonical<N> {
    /// The Bayes product: exponents add, so the parameters add. O(n²), and
    /// there is genuinely nothing else to it.
    pub fn product(&self, other: &Self) -> Self {
        Self { xi: self.xi + other.xi, omega: self.omega + other.omega }
    }

    /// Back to moments: one factorization, one solve, one inverse. O(n³).
    pub fn to_moments(&self) -> Result<Gaussian<N>, NotPositiveDefinite> {
        let chol = Cholesky::new(self.omega).ok_or(NotPositiveDefinite)?;
        Gaussian::new(chol.solve(&self.xi), chol.inverse())
    }
}

Randomness gets its own module for one reason: so that there is exactly one place in the workspace where a generator is constructed, and it takes a seed.

use rand::{rngs::SmallRng, Rng, SeedableRng};

/// Every demo, figure and test in this book starts here. Widgets display it.
pub const BOOK_SEED: u64 = 0x5EED_2026;

/// The only RNG constructor in the workspace. `thread_rng()` appears nowhere:
/// a figure that cannot be reproduced is not evidence, and a flaky test is a
/// test that will eventually be deleted rather than fixed.
pub fn rng(seed: u64) -> SmallRng {
    SmallRng::seed_from_u64(seed)
}

/// Thrun et al., Table 5.4 — how you sampled a normal in 2000. Kept because
/// the book discusses it; never called, because `rand_distr::StandardNormal`
/// is exact, faster, and has tails.
///
/// Twelve draws from U(−1,1) sum to variance 12·(1/3) = 4, so the multiplier
/// that yields standard deviation `b` is 1/2. Support is bounded at ±6b, which
/// is the defect that matters: this sampler cannot produce the outlier that
/// breaks your filter, so it will never warn you that one exists.
pub fn sample_normal_distribution(b: f64, rng: &mut impl Rng) -> f64 {
    let sum: f64 = (0..12).map(|_| rng.random_range(-1.0..1.0)).sum();
    0.5 * b * sum
}

A worked example you can check by hand

Rusty's odometry gives a prior $\Normal(5.0,\, 4.0)$ — 5 m, standard deviation 2 m. The wall-range sensor gives a likelihood $\Normal(6.5,\, 1.0)$. Convert both, add, convert back:

Sanity-check the two claims from the widget against these numbers. The posterior variance $0.8$ is below both $4.0$ and $1.0$. The posterior mean $6.2$ sits $80\%$ of the way from the prior to the sensor, because the sensor is four times as precise: $\tfrac{1}{1} / (\tfrac{1}{4} + \tfrac{1}{1}) = 0.8$.

Now the 2-D act. Take $\mu = 0$ and

Conditioning on $x_b = 1$ gives, from Derivation 6, a mean of $1.9$ and a variance of $4 - 1.9^2/1 = 0.39$ — the same $0.39$, because in 2-D $\sigma_{a|b}^2 = \det\Sigma / \sigma_b^2$ and $\sigma_b^2 = 1$ here. The marginal standard deviation was $2.0$; the conditional one is $\sqrt{0.39} = 0.6245$. One scalar observation of a correlated quantity removed 69% of the uncertainty without measuring $x_a$ at all.

In entropy terms, the correlation is worth

The tests that pin those numbers

#[cfg(test)]
mod tests {
    use super::*;
    use crate::prob::sample::{rng, BOOK_SEED};
    use approx::assert_relative_eq;
    use nalgebra::{Matrix2, SMatrix, SVector, Vector2};

    /// The chapter's 1-D worked example, to the printed digit. If this test
    /// ever disagrees with the prose, the test is what settles it.
    #[test]
    fn worked_example_ch02_fuse_two_sensors() {
        let prior = Gaussian::<1>::new(SVector::from([5.0]), SMatrix::from([[4.0]])).unwrap();
        let likelihood = Gaussian::<1>::new(SVector::from([6.5]), SMatrix::from([[1.0]])).unwrap();

        let fused = prior.to_canonical().product(&likelihood.to_canonical());
        assert_relative_eq!(fused.omega[(0, 0)], 1.25, epsilon = 1e-12); // 0.25 + 1.00
        assert_relative_eq!(fused.xi[0], 7.75, epsilon = 1e-12);         // 1.25 + 6.50

        let posterior = fused.to_moments().unwrap();
        assert_relative_eq!(posterior.mean()[0], 6.2, epsilon = 1e-12);
        assert_relative_eq!(posterior.cov()[(0, 0)], 0.8, epsilon = 1e-12);
    }

    /// The 2-D worked example: determinant, Schur complement, and the entropy
    /// that correlation is worth.
    #[test]
    fn worked_example_ch02_correlated_pair() {
        let sigma = Matrix2::new(4.0, 1.9, 1.9, 1.0);
        let joint = Gaussian::<2>::new(Vector2::zeros(), sigma).unwrap();

        assert_relative_eq!(sigma.determinant(), 0.39, epsilon = 1e-9);

        let conditional = condition_second(&joint, 1.0);
        assert_relative_eq!(conditional.mean()[0], 1.9, epsilon = 1e-9);
        assert_relative_eq!(conditional.cov()[(0, 0)], 0.39, epsilon = 1e-9);

        let uncorrelated =
            Gaussian::<2>::new(Vector2::zeros(), Matrix2::new(4.0, 0.0, 0.0, 1.0)).unwrap();
        // ½ ln(4.0 / 0.39) = −½ ln(1 − ρ²), the nats that correlation is worth.
        let gain = uncorrelated.entropy() - joint.entropy();
        assert_relative_eq!(gain, 1.163_951_45, epsilon = 1e-8);
    }

    /// Round-trip: moments → canonical → moments is the identity. Property-test
    /// this over random SPD matrices; the fixed case is the regression guard.
    #[test]
    fn canonical_round_trip_is_identity() {
        let sigma = Matrix2::new(4.0, 1.9, 1.9, 1.0);
        let g = Gaussian::<2>::new(Vector2::new(1.0, -2.0), sigma).unwrap();
        let back = g.to_canonical().to_moments().unwrap();
        assert_relative_eq!(*back.mean(), *g.mean(), epsilon = 1e-10);
        assert_relative_eq!(*back.cov(), *g.cov(), epsilon = 1e-10);
    }

    /// The 95% ellipse had better contain 95% of the samples. If it does not,
    /// either the sampler or the χ² quantile is wrong — and every covariance
    /// ellipse printed in this book is lying by the same amount.
    #[test]
    fn ellipse_coverage_is_what_it_claims() {
        let g = Gaussian::<2>::new(Vector2::zeros(), Matrix2::new(4.0, 1.9, 1.9, 1.0)).unwrap();
        let mut r = rng(BOOK_SEED);
        const M: usize = 100_000;
        const CHI2_95: f64 = 5.991_464_547_107_98; // −2 ln 0.05

        let inside = (0..M)
            .filter(|_| g.mahalanobis2(&g.sample(&mut r)) <= CHI2_95)
            .count();
        assert!((inside as f64 / M as f64 - 0.95).abs() < 5e-3);
    }
}

Putting it together

The runnable artifact for this chapter is one example binary that does both acts and prints the numbers the prose claims.

use nalgebra::{Matrix2, SMatrix, SVector, Vector2};
use pr_core::prob::{rng, Gaussian, BOOK_SEED};

fn main() {
    // Act I — one dimension, by hand-checkable arithmetic.
    let prior = Gaussian::<1>::new(SVector::from([5.0]), SMatrix::from([[4.0]])).unwrap();
    let likelihood = Gaussian::<1>::new(SVector::from([6.5]), SMatrix::from([[1.0]])).unwrap();
    let posterior = prior
        .to_canonical()
        .product(&likelihood.to_canonical())
        .to_moments()
        .unwrap();

    println!("prior       N(mu = {:.3}, var = {:.3})", prior.mean()[0], prior.cov()[(0, 0)]);
    println!("likelihood  N(mu = {:.3}, var = {:.3})", likelihood.mean()[0], likelihood.cov()[(0, 0)]);
    println!("posterior   N(mu = {:.3}, var = {:.3})", posterior.mean()[0], posterior.cov()[(0, 0)]);
    println!("entropy     {:.4} nats  (prior {:.4})", posterior.entropy(), prior.entropy());

    // Act II — two dimensions, checked against 10 000 seeded draws.
    let joint = Gaussian::<2>::new(Vector2::zeros(), Matrix2::new(4.0, 1.9, 1.9, 1.0)).unwrap();
    let mut r = rng(BOOK_SEED);
    let draws: Vec<_> = (0..10_000).map(|_| joint.sample(&mut r)).collect();

    let mean = draws.iter().sum::<Vector2<f64>>() / draws.len() as f64;
    let cov = draws
        .iter()
        .map(|x| (x - mean) * (x - mean).transpose())
        .sum::<Matrix2<f64>>()
        / (draws.len() - 1) as f64;

    println!("sample cov  [{:.3} {:.3}; {:.3} {:.3}]", cov[(0, 0)], cov[(0, 1)], cov[(1, 0)], cov[(1, 1)]);
    println!("entropy     {:.4} nats", joint.entropy());
}

prior       N(mu = 5.000, var = 4.000)
likelihood  N(mu = 6.500, var = 1.000)
posterior   N(mu = 6.200, var = 0.800)
entropy     1.3074 nats  (prior 2.1121)
sample cov  [3.984 1.891; 1.891 0.995]
entropy     2.3671 nats

Read the first entropy line against Derivation 7: $\tfrac{1}{2}\ln(2\pi e \cdot 4) = 2.1121$ before the measurement, $\tfrac{1}{2}\ln(2\pi e \cdot 0.8) = 1.3074$ after. The sensor was worth $0.80$ nats, or $1.16$ bits — and it would have been worth exactly the same number had it read $6.5$, $60$, or $-4$, because entropy does not care where the belief sits. The 2-D entropy, $2.3671$ nats, is $\ln(2\pi e) + \tfrac{1}{2}\ln 0.39$, and the $1.164$ nats separating it from the uncorrelated case is the correlation earning its keep.

The sample covariance agrees with $\Sigma$ to about half a percent at $N = 10^4$, which is the $1/\sqrt{N}$ rate of the Sampling Convergence widget doing exactly what it promised — and no better. Re-seed and the third digit moves; the second does not.

That is the whole toolkit. From here:

• Chapter 5 applies Bayes rule recursively, and the "compute the shape, normalize later" habit becomes line 3 of the Bayes filter.
• Chapter 6 industrializes Derivations 2, 4 and 6 into the Kalman filter, and the moments/canonical duality of Derivation 3 becomes the Kalman filter and the information filter side by side.
• Chapter 8 takes samples seriously as the belief itself.
• Chapter 11 gates data associations on the Mahalanobis distance of Derivation 5, using the very $\chi^2$ number that sizes the ellipses here.
• Chapter 24 turns the entropy of Derivation 7 into information gain, and lets the robot choose where to drive.

One last thing worth saying plainly: the purple curve in the Blob Multiplier is not an illustration of Canonical::product. It is the output of the TypeScript port of Canonical::product, fed the same numbers as the Rust test, agreeing to the last digit the screen can show. Where the prose, the figure and the code disagree in this book, the test is what settles it — and there is always a test.

Exercises

1. Foundation exerciseFDifficulty 2 of 3The product in n dimensions, twice

   Derive the product of two $n$-dimensional Gaussians in canonical form. It should take three lines. Then derive the same result in moments form, which will take a page and require the matrix inversion lemma. Conclude in one sentence why Chapter 6 offers two filters rather than one, and state which one you would choose for a robot with twenty sensors and one motion model.

2. Foundation exerciseFDifficulty 2 of 3Variance is positive, entropy is not

   Show that $\Cov[X] \succeq 0$ for every random vector with finite second moments, directly from the definition. Then specialize the Gaussian entropy formula to one dimension and find the values of $\sigma$ for which the differential entropy is negative. Explain in two sentences why a negative entropy is not a contradiction, and why the difference of two differential entropies is still meaningful when both are negative.

3. Foundation exerciseFDifficulty 3 of 3The evidence is a Gaussian

   Derivation 2 claims that the constant thrown away when two Gaussians are multiplied is exactly $\Normal(\mu_1; \mu_2, \sigma_1^2 + \sigma_2^2)$. Prove it by carrying the constants through Steps 1–4 instead of discarding them. Then explain why this makes an innovation-based outlier test and an evidence-based one the same test.

4. Conceptual exerciseCDifficulty 1 of 3Predict, then check: precision weighting

   In the Blob Multiplier, set the likelihood's variance to nine times the prior's. Before you look: what fraction of the distance from the prior mean to the likelihood mean will the posterior mean travel? Write the number down, then verify it. Now push the likelihood variance to its maximum and predict what the posterior converges to and why. Finally, find a setting where the posterior is narrower than the likelihood by less than one percent and say what that setting means physically.

   Hint

   The posterior mean is a weighted average with weights proportional to 1/σ².

5. Conceptual exerciseCDifficulty 2 of 3Predict, then check: the slice that does not widen

   Set $\rho = 0.95$ in Slice vs. Squash. Predict whether the conditional at $\beta = \mu_b$ is narrower, wider, or the same width as the conditional at $\beta = \mu_b + 2\sigma_b$. Most people get this wrong the first time. Verify with the widget, then explain the result using Derivation 6 by pointing at the one symbol that is missing from the conditional covariance. Finally, use the Gaussian Playground to explain why the ellipse's horizontal extent at $\rho = 0.95$ is not the width of the blue marginal curve underneath it.

6. Practical exercisePDifficulty 2 of 3Conditioning in general dimension

   Implement conditioning for a partitioned Gaussian<N>, returning the distribution of the unobserved block given values for the observed one. Rust's const generics cannot express Gaussian<{N - K}> on stable, so you will have to choose: fixed 2-block const parameters with an A + B == N assertion, a dynamically-sized DVector/DMatrix variant, or the nightly generic_const_exprs feature. Pick one, justify it in a comment, and property-test the result against brute-force numerical integration of the joint density on a 2-D grid.

7. Practical exercisePDifficulty 3 of 3Two samplers, one benchmark

   Implement sample_normal_distribution(b) from the Callout above and benchmark it against rand_distr::StandardNormal with Criterion at $10^6$ draws. Report throughput, and report the first four sample moments of each. Then do the part that actually matters: estimate $P(\lvert X \rvert > 5)$ with both samplers and compare against the exact $5.7 \times 10^{-7}$. Write two sentences on why a filter that gates outliers at $5\sigma$ must not be tested with the twelve-uniform sampler.

References

1. Thrun, S., Burgard, W., and Fox, D. (2005) Probabilistic Robotics. MIT Press.link to Probabilistic Robotics (opens in a new tab)

   Section 2.2 is the notation baseline for this chapter — η, the Gaussian, expectation, entropy. The canonical parameterization is deferred there to the information filter in Chapter 3; we pull it forward so the duality is a one-line observation later instead of a surprise.

2. Jaynes, E. T. (edited by G. L. Bretthorst) (2003) Probability Theory: The Logic of Science. Cambridge University Press.doi:10.1017/CBO9780511790423 (opens in a new tab)

   The argument that probability is extended logic rather than a theory of frequencies. Read Chapters 1–2 if the phrase 'the robot's belief' has ever felt like a category error.

3. Marsaglia, G. and Tsang, W. W. (2000) The Ziggurat Method for Generating Random Variables. Journal of Statistical Software 5(8), 1–7.doi:10.18637/jss.v005.i08 (opens in a new tab)

   The sampler behind rand_distr::StandardNormal, and the reason this book does not use Thrun's twelve-uniform approximation: exact tails, and faster.

4. Barfoot, T. D. (2024) State Estimation for Robotics, 2nd edition. Cambridge University Press.link to State Estimation for Robotics, 2nd edition (opens in a new tab)

   Chapter 2 covers the same Gaussian identities with more care about the linear-algebra details, including the block-inverse manipulations we sketch in Derivation 6 and prove in Appendix B. The author maintains a freely readable draft at the URL above.

5. Barfoot, T. D., Forbes, J. R., and Yoon, D. J. (2020) Exactly sparse Gaussian variational inference with application to derivative-free batch nonlinear state estimation. International Journal of Robotics Research 39(13), 1473–1502.doi:10.1177/0278364920937608 (opens in a new tab)

   A modern argument for living in canonical form: the inverse covariance of a batch estimation problem is block-tridiagonal, so storing Ω instead of Σ turns a dense O(n²) object into a sparse O(n) one. Derivation 3 is the toy version of the whole method.

6. Barfoot, T. D. (2020) Fundamental Linear Algebra Problem of Gaussian Inference. arXiv:2010.08022.link to Fundamental Linear Algebra Problem of Gaussian Inference (opens in a new tab)

   States precisely which entries of Σ you can recover cheaply once you have Ω, via the Takahashi recursions. The clearest available account of why the moments/canonical trade in Derivation 3 is not symmetric in practice.

7. Ortiz, J., Evans, T., and Davison, A. J. (2021) A visual introduction to Gaussian Belief Propagation. arXiv:2107.02308.link to A visual introduction to Gaussian Belief Propagation (opens in a new tab)

   Gaussians in canonical form, made interactive — and the closest thing in the literature to what this book is trying to be. Its factor-multiplication figures are the multi-variable sequel to the Blob Multiplier.

8. Placed, J. A., Strader, J., Carrillo, H., Atanasov, N., Indelman, V., Carlone, L., and Castellanos, J. A. (2023) A Survey on Active Simultaneous Localization and Mapping: State of the Art and New Frontiers. IEEE Transactions on Robotics 39(3), 1686–1705.doi:10.1109/TRO.2023.3248510 (opens in a new tab)

   Where the entropy of Derivation 7 goes: Section IV surveys the information-theoretic utilities — entropy reduction, D-optimality, log det Σ — that Chapter 24 uses to decide where a robot should drive next.