M9-T3 — Resampling, Degeneracy and Impoverishment

Learning objectives

After this lesson, you should be able to:

• Explain weight degeneracy and why it is unavoidable without resampling.
• Compare multinomial, stratified and systematic resampling by the variance they add.
• Explain sample impoverishment and why it makes the cloud overconfident.
• Apply the remedies: a resampling threshold, jitter (roughening), and a sufficient number of particles.

The intuition

The weights are products of likelihoods. After a few readings from an accurate sensor, a particle that is a little off the truth has a likelihood several times smaller at every step, so its weight shrinks geometrically. In the experiment of the left panel the effective sample size of 200 particles falls from 29 at step 1 to 1.8 at step 10: after ten readings, two particles are doing the work of 200. This is degeneracy, and it is unavoidable: it follows from the product form of the weights.

Resampling is the cure: replace the particles by copies of the heavy ones. But copies are identical. After each resampling the set contains fewer distinct ancestors, and if the motion model adds little noise the copies stay almost identical. The right panel tests this on a robot at rest (the filter's motion model has no noise): after 40 resamplings only 3.2 of the 100 original particles are ancestors of the cloud, the cloud is a few clumps, and its spread, 3.5 cm, is less than half the real error of its mean (7.5 cm; the exact posterior has a standard deviation of 4.7 cm). The filter is overconfident. This is sample impoverishment.

Adding a small random jitter to each copy (roughening) restores diversity. With 1 cm of jitter the 100 particles are all distinct, the cloud's spread, 5.7 cm, matches its real error, 5.7 cm, near the exact 4.7 cm. With 5 cm the cloud is wider than the truth (12.6 cm): the jitter inflates it. The jitter must be small compared with the posterior width.

Resampling trades one kind of waste for another: it spends particles where the probability is, and it spends diversity. Resample only when you need to, and keep a little noise in the motion model.

Definitions

Degeneracy

The weights become so uneven that $N_\mathrm{eff}$ falls to about 1. Without resampling the variance of the weights increases with time.

Impoverishment

After repeated resampling few distinct particles remain: the cloud has the number of particles of a good filter and the diversity of a bad one.

Number of distinct ancestors

The number of initial particles that have a descendant in the current set: a direct measure of diversity.

Resampling threshold

Resample only when $N_\mathrm{eff}<\theta N$, typically $\theta=\tfrac12$. Fewer resamplings mean less impoverishment.

Roughening (jitter)

Add small Gaussian noise to the copies after resampling, with a standard deviation scaled to the spread of the cloud and to $N^{-1/d}$ for a $d$-dimensional state.

Derivation: why the weights degenerate and why copies lose diversity

Degeneracy

Without resampling the weight of particle $i$ after $t$ steps is the product of its likelihoods, $w_t^{[i]}\propto\prod_{s\le t} p(z_s\mid x_s^{[i]})$.

1. $$\log w_t^{[i]}=\sum_{s\le t}\log p\bigl(z_s\mid x_s^{[i]}\bigr)+\text{const}$$

   The log-weight is a sum of one term per step.

2. $$\mathrm{Var}\bigl(\log w^{[i]}_t\bigr)\approx t\,\sigma^2_\ell$$

   The terms of different particles differ randomly (each particle explores a different path), with variance $\sigma_\ell^2$ per step, so the spread of the log-weights grows linearly with time.

3. $$\frac{\sum w}{\max w}\sim e^{\sigma_\ell\sqrt t\cdot c}\quad\Rightarrow\quad N_\mathrm{eff}\to1$$

   The ratio of the best weight to the mean grows exponentially in √t; the weights become a heavy-tailed (log-normal) set dominated by the maximum.

The sharper the sensor, the larger $\sigma_\ell$, and the faster the collapse: here 5 steps to lose 97% of the sample size.

Impoverishment as genetic drift

Resampling is a branching process: each particle has a random number of children. Without any added noise the number of distinct ancestors is non-increasing, and for a fixed number of particles it falls to 1 eventually, as in the Wright–Fisher model of population genetics: with $N$ individuals the expected time until all descend from one ancestor is of order $N$ generations. The only things that restore diversity are noise in the motion model and jitter.

1. $$\text{ancestors}(t+1)\le\text{ancestors}(t)$$

   Copying cannot create new ancestors.

2. $$\sigma_\mathrm{cloud}^2\ \text{shrinks by a factor}\ \bigl(1-\tfrac1N\bigr)\ \text{per resampling}$$

   The sample variance of N draws with replacement from a set is, on average, $(1-1/N)$ times that of the set (the resampling adds the noise of a bootstrap, which is a contraction without jitter).

3. $$\text{with jitter }\sigma_j:\quad \sigma^2_\mathrm{cloud}\to\bigl(1-\tfrac1N\bigr)\sigma^2_\mathrm{cloud}+\sigma_j^2$$

   Jitter balances the contraction at a spread $\sigma_j\sqrt N$ if nothing else acted; in a filter the readings shrink the cloud as well, so the final spread is smaller. Too much jitter still leaves the cloud wider than the posterior, as the 5 cm run shows.

Degeneracy is a property of the weights, impoverishment of the copies. Resampling converts the first into the second, so use it sparingly and give the particles a reason to differ.

A worked example

Three programs.

1. Degeneracy. 200 particles track a robot that moves 20 cm per step (5 cm of noise) with a position sensor of 10 cm noise, with no resampling. The effective sample size is 28.9, 20.5, 6.4, 1.8, 3.0 and 2.0 at steps 1, 2, 5, 10, 20 and 30. The error of the estimate stays small for a few steps and is 22 cm at step 30, since two particles cannot follow the robot.
2. Resampling schemes. One fixed set of 100 weighted particles is resampled 4 000 times. The standard deviation of the resampled mean around the weighted mean is 0.066 for multinomial, 0.038 for stratified and 0.030 for systematic resampling: the latter two add half as much noise, or less.
3. Impoverishment. A robot at rest is observed 40 times with a 30 cm sensor, resampling each time. With no jitter 3.2 ancestors remain of 100, the cloud's spread is 3.5 cm and the error of the mean 7.5 cm. With 1 cm jitter all 100 particles are distinct, the spread is 5.7 cm and the error 5.7 cm (the exact posterior standard deviation is 4.7 cm). With 5 cm the spread is 12.6 cm and the error 9.4 cm.

import math
from robotprob import (Rng, gaussian_pdf, normalize_weights, effective_sample_size, resample_multinomial,
                       resample_systematic, resample_stratified, mean)

# 1. degeneracy: track a robot with a sharp position sensor, never resampling
def track(seed, n=200, steps=30, scheme='none', thr=0.5, sig_m=0.05, sig_z=0.1, jitter=0.0):
    rng = Rng(seed)
    x = 0.0
    ps = [rng.gauss(0.0, 1.0) for _ in range(n)]
    ws = [1.0 / n] * n
    ids = list(range(n))                                  # the ancestor of each particle
    log = []
    for t in range(1, steps + 1):
        x += 0.2 + rng.gauss(0.0, 0.05)
        z = x + rng.gauss(0.0, sig_z)
        ps = [p + 0.2 + rng.gauss(0.0, sig_m) for p in ps]
        ws = normalize_weights([w * gaussian_pdf(z, p, sig_z) for w, p in zip(ws, ps)])
        neff = effective_sample_size(ws)
        if scheme != 'none' and (scheme == 'always' or neff < thr * n):
            idx = (resample_multinomial if scheme == 'always' else resample_systematic)(rng, ws)
            ps = [ps[i] + (rng.gauss(0.0, jitter) if jitter else 0.0) for i in idx]
            ids = [ids[i] for i in idx]
            ws = [1.0 / n] * n
        est = sum(w * p for w, p in zip(ws, ps))
        log.append((t, x, est, neff, len(set(ids))))
    return log

print('no resampling, 200 particles: step  N_eff  error')
lg = track(2)
for t, x, est, neff, u in lg:
    if t in (1, 2, 5, 10, 20, 30):
        print(' ', str(t).rjust(4), str(round(neff, 1)).rjust(7), str(round(est - x, 3)).rjust(8))

# 2. the three resampling schemes: the estimate of the mean of a fixed weighted set, repeated
rng = Rng(11)
base = [rng.gauss(0.0, 1.0) for _ in range(100)]
w0 = normalize_weights([gaussian_pdf(1.0, p, 0.7) for p in base])
exact = sum(w * p for w, p in zip(w0, base))
print('resampling a fixed set of 100 weighted particles, 4000 times: spread of the estimated mean')
for name, fn in (('multinomial', resample_multinomial), ('stratified', resample_stratified), ('systematic', resample_systematic)):
    r = Rng(5)
    est = [mean([base[i] for i in fn(r, w0)]) for _ in range(4000)]
    print(' ', name.ljust(12), 'std', round(math.sqrt(mean([(e - exact) ** 2 for e in est])), 4))

# 3. impoverishment: a robot standing still, a filter that trusts its motion model (no motion noise)
def still(seed, jitter, n=100, steps=40):
    rng = Rng(seed)
    x = 2.0
    ps = [rng.gauss(0.0, 2.0) for _ in range(n)]
    ws = [1.0 / n] * n
    ids = list(range(n))
    for t in range(steps):
        z = x + rng.gauss(0.0, 0.3)
        ws = normalize_weights([w * gaussian_pdf(z, p, 0.3) for w, p in zip(ws, ps)])
        idx = resample_systematic(rng, ws)
        ps = [ps[i] + (rng.gauss(0.0, jitter) if jitter else 0.0) for i in idx]
        ids = [ids[i] for i in idx]
        ws = [1.0 / n] * n
    mu = mean(ps)
    sd = math.sqrt(mean([(p - mu) ** 2 for p in ps]))
    return len(set(ids)), mu, sd, len(set([round(p, 9) for p in ps]))
print('a robot at rest, 100 particles resampled every step for 40 steps, 200 runs each:')
print('  jitter   ancestors left   distinct particles   rms error of the mean   mean spread of the cloud   exact posterior std')
true_sd = 0.3 / math.sqrt(40)
for jit in (0.0, 0.01, 0.05):
    R = [still(seed, jit) for seed in range(1, 201)]
    print(' ', str(jit).rjust(6), str(round(mean([r[0] for r in R]), 1)).rjust(14), str(round(mean([r[3] for r in R]), 1)).rjust(19),
          str(round(math.sqrt(mean([(r[1] - 2.0) ** 2 for r in R])), 3)).rjust(22), str(round(mean([r[2] for r in R]), 3)).rjust(24), str(round(true_sd, 3)).rjust(21))

no resampling, 200 particles: step  N_eff  error
     1    28.9   -0.043
     2    20.5    -0.01
     5     6.4    0.085
    10     1.8    0.017
    20     3.0   -0.063
    30     2.0   -0.219
resampling a fixed set of 100 weighted particles, 4000 times: spread of the estimated mean
  multinomial  std 0.0656
  stratified   std 0.0379
  systematic   std 0.0295
a robot at rest, 100 particles resampled every step for 40 steps, 200 runs each:
  jitter   ancestors left   distinct particles   rms error of the mean   mean spread of the cloud   exact posterior std
     0.0            3.2                 3.2                  0.075                    0.035                 0.047
    0.01            3.8               100.0                  0.057                    0.057                 0.047
    0.05            5.7               100.0                  0.094                    0.126                 0.047

In robotics

A practical resampling policy

• Resample when $N_\mathrm{eff}<N/2$ (not every step).
• Use systematic or stratified resampling.
• Keep a small noise in the motion model even when the robot is still (a floor on the process noise), or add roughening of a fraction of the cloud's spread.
• Monitor the number of distinct particles or the cloud's spread against the error in tests.

Sharp sensors are the hard case

An accurate sensor (a lidar, a tag reader) gives a likelihood much narrower than the prior cloud, so most particles are discarded at every step and N_eff collapses at once. Remedies: inflate the sensor noise in the likelihood (a modest tempering, treating the sensor as less sharp than it is, which also covers model error), use more particles, or use a better proposal that moves particles towards the reading before weighting (next lesson and the optimal proposal).

Why the number of particles must grow with the dimension

In a high-dimensional state the posterior occupies a tiny part of the prior's volume and the weights degenerate at once. This is why the plain particle filter is used for 3-D poses and not for full maps: SLAM filters factor the problem (the map is conditioned on the path, Module 10).

Reporting an estimate

Use the weighted mean and the covariance of the cloud only if it is unimodal; check the number of distinct particles and the fraction of resamplings. A cloud made of a handful of identical particles looks confident and is not.

Common misconceptions

Key takeaways

1. Without resampling the weights degenerate: log-weights spread linearly in time, and N_eff fell from 29 to 2 in ten steps with a sharp sensor.
2. Resampling cures degeneracy but removes diversity: copies are identical, and the number of distinct ancestors only decreases.
3. Systematic and stratified resampling add much less noise than multinomial (0.030 and 0.038 against 0.066).
4. An impoverished cloud is overconfident: 3.5 cm claimed for a 7.5 cm error. A small jitter (1 cm) restored agreement; a large one (5 cm) over-spreads the cloud.
5. Resample on N_eff < N/2, keep noise in the motion model, use enough particles, and watch the number of distinct particles.

References and further reading

Each source is marked free or paid. Free to read does not mean free to reuse: the course resources page lists every source with its licence.

• S. Thrun, W. Burgard and D. Fox, Probabilistic Robotics, MIT Press, 2005, sections 4.3.4 (resampling, variance and the low-variance sampler) and 4.3.5 (practical considerations). Paid book.
• R. Douc and O. Cappé, “Comparison of Resampling Schemes for Particle Filtering”, Proceedings of the 4th International Symposium on Image and Signal Processing and Analysis, 2005. Paid or institutional access.
• S. Särkkä and L. Svensson, Bayesian Filtering and Smoothing, 2nd edition, 2023 , chapter 11 (particle filtering, resampling). Free to read as the authors’ online PDF; the book is copyrighted.
• C. Stachniss, “Particle Filter and Monte Carlo Localization” (video lecture) , covers the low-variance sampler. Free to watch.

Knowledge check

Three questions. Write an answer, then reveal the explanation. Your answers stay in your browser.

What’s next

With the algorithm and its failure modes in hand, the last lesson of the module puts it to work for Monte Carlo localisation of a robot in a map: with a laser-like sensor, global localisation, and recovery from kidnapping by injecting random particles.

Related content

Nothing below is required to finish this lesson.

• Before
  M9-T2: The Particle Filter Algorithm

  Derive the sequential importance resampling filter from the Bayes filter and implement its three steps.

• Lab
  M9-L1: Monte Carlo Localization

  Rover is placed in a known room and must find itself from range readings to four beacons. Implement the weights and the resampler, then make the filter survive a kidnapping.

Other courses on this site treat neighbouring topics from their own angle:

• SLAM M3-T8: Particle Filters and FastSLAM
• SLAM M3-L4: Explore Particle-Filter Localization