diff --git a/src/pyrecest/filters/abstract_particle_filter.py b/src/pyrecest/filters/abstract_particle_filter.py index 8f3b8d113..add006788 100644 --- a/src/pyrecest/filters/abstract_particle_filter.py +++ b/src/pyrecest/filters/abstract_particle_filter.py @@ -8,6 +8,7 @@ array, hstack, isfinite, + max as backend_max, ndim, ones_like, random, @@ -230,10 +231,18 @@ def predict_nonlinear_nonadditive(self, f, samples, weights): raise ValueError("Noise weights must be finite.") if not bool(all(weights >= 0.0)): raise ValueError("Noise weights must be nonnegative.") - weight_sum = sum(weights) - if not bool(isfinite(weight_sum)) or not bool(weight_sum > 0.0): + if weights.shape[0] == 0: raise ValueError("Noise weights must have positive finite total mass.") - weights = weights / weight_sum + weight_scale = backend_max(weights) + if not bool(isfinite(weight_scale)) or not bool(weight_scale > 0.0): + raise ValueError("Noise weights must have positive finite total mass.") + scaled_weights = weights / weight_scale + scaled_weight_sum = sum(scaled_weights) + if not bool(isfinite(scaled_weight_sum)) or not bool( + scaled_weight_sum > 0.0 + ): + raise ValueError("Noise weights must have positive finite total mass.") + weights = scaled_weights / scaled_weight_sum n_particles = self.filter_state.w.shape[0] noise_samples = random.choice(samples, n_particles, p=weights) @@ -275,7 +284,7 @@ def update_model( """Update using a reusable particle measurement model.""" if not hasattr(measurement_model, "likelihood"): raise TypeError( - "Particle-filter measurement models must expose a likelihood callable." + "Particle measurement models must expose a likelihood callable." ) return self.update_nonlinear_using_likelihood( diff --git a/tests/filters/test_particle_filter_extreme_nonadditive_weights.py b/tests/filters/test_particle_filter_extreme_nonadditive_weights.py new file mode 100644 index 000000000..b17bc8a72 --- /dev/null +++ b/tests/filters/test_particle_filter_extreme_nonadditive_weights.py @@ -0,0 +1,32 @@ +from unittest import mock + +import numpy as np +import numpy.testing as npt + +from pyrecest.backend import array, to_numpy +from pyrecest.distributions import LinearDiracDistribution +from pyrecest.filters.euclidean_particle_filter import EuclideanParticleFilter + + +def test_nonadditive_prediction_normalizes_extreme_finite_noise_weights(): + particles = array([[10.0], [20.0]]) + samples = array([[1.0], [2.0]]) + backend_dtype = to_numpy(array([1.0])).dtype + max_weight = np.finfo(backend_dtype).max + weights = array([max_weight, max_weight / 2.0]) + particle_filter = EuclideanParticleFilter(n_particles=2, dim=1) + particle_filter.filter_state = LinearDiracDistribution(particles) + + with mock.patch( + "pyrecest.filters.abstract_particle_filter.random.choice", + return_value=samples, + ) as choice_mock: + particle_filter.predict_nonlinear_nonadditive( + lambda particle, noise: particle + noise, + samples, + weights, + ) + + normalized_weights = choice_mock.call_args.kwargs["p"] + npt.assert_allclose(to_numpy(normalized_weights), [2.0 / 3.0, 1.0 / 3.0]) + npt.assert_allclose(to_numpy(particle_filter.filter_state.d), [[11.0], [22.0]])