From de9c023b04b2e2c94045f1698d64f1ce81294fed Mon Sep 17 00:00:00 2001 From: "Michael J. Williams" Date: Tue, 27 Oct 2026 11:43:47 +0100 Subject: [PATCH 1/7] ENH(mutation): add global independence proposals --- CHANGELOG.md | 7 + .../mutation_configuration.rst | 41 +++++ src/sequince/steps/__init__.py | 2 + src/sequince/steps/mutation/__init__.py | 2 + src/sequince/steps/mutation/kernels.py | 78 +++++++++ .../test_independent_proposal.py | 161 ++++++++++++++++++ 6 files changed, 291 insertions(+) create mode 100644 tests/test_mutation/test_independent_proposal.py diff --git a/CHANGELOG.md b/CHANGELOG.md index d9094f64..bccb9714 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -6,6 +6,13 @@ The format follows [Keep a Changelog](https://keepachangelog.com/en/1.1.0/). ## [Unreleased] +### Added + +- Independence Metropolis-Hastings proposals using `IndependentProposalKernel` + and `ProposalDistribution`, with cached proposed log densities. Proposals must + describe mutation coordinates; physical-space priors require physical-space + mutation. + ## [0.1.0b1] - 2026-10-23 First public beta release. diff --git a/docs/user_guide/sampling_options/mutation_configuration.rst b/docs/user_guide/sampling_options/mutation_configuration.rst index f03c5f7c..69f28755 100644 --- a/docs/user_guide/sampling_options/mutation_configuration.rst +++ b/docs/user_guide/sampling_options/mutation_configuration.rst @@ -140,6 +140,47 @@ Use these steps to compose mutation updates: The adaptation wrappers support eager execution only. Place conditional steps outside them, for example ``MutateIfResampled(mutation_step)``. +Global independence proposals +----------------------------- + +:class:`~sequince.steps.mutation.kernels.IndependentProposalKernel` draws +proposals independently of each particle's current position. It accepts a +:class:`~sequince.proposals.ProposalDistribution`: ``sample_and_log_prob`` draws +the requested batch using the supplied ``rng`` and ``xp``, and ``log_prob`` +evaluates one log density per particle. Both methods must describe the same +distribution. Given a proposal named ``prior_proposal`` and a target named +``target``: + +.. code-block:: python + + from sequince.steps import ( + IndependentProposalKernel, + MetropolisHastingsMutation, + ) + + mutation = MetropolisHastingsMutation( + target, + IndependentProposalKernel(prior_proposal), + n_steps=10, + transform=None, + ) + +Acceptance includes ``log q(current) - log q(proposed)``. The kernel caches the +sampled log densities, so it only evaluates current-particle densities when +computing that correction. Samples must match the input particle shape, and +log densities must have shape ``(n_particles,)``. The proposal should cover the +target's support, and its distribution must remain fixed during mutation. + +The proposal must describe the coordinates supplied to the kernel. Use +``transform=None`` for physical-space prior proposals; when configuring a +prebuilt sampler, also set ``particle_transform=None``. Transformed proposals +must sample and evaluate densities in the transformed coordinates, including +the appropriate change of variables. A physical-space prior cannot be passed +unchanged to a latent-space kernel. + +The mutation-kernel factory accepts ``"independent-proposal"`` with a +``proposal`` object in ``kernel_kwargs``. + Prebuilt samplers ----------------- diff --git a/src/sequince/steps/__init__.py b/src/sequince/steps/__init__.py index c45fcdca..14b7984d 100644 --- a/src/sequince/steps/__init__.py +++ b/src/sequince/steps/__init__.py @@ -38,6 +38,7 @@ ESJDMetric, ESJDScaleWeightPolicy, FixedLengthMutation, + IndependentProposalKernel, KernelSample, MaskedGaussianRandomWalk, MetropolisBlock, @@ -124,6 +125,7 @@ "ErynMutation", "FixedLengthMutation", "MaskedGaussianRandomWalk", + "IndependentProposalKernel", "KernelSample", "MetropolisHastingsMutation", "RedBlueMetropolisHastingsMutation", diff --git a/src/sequince/steps/mutation/__init__.py b/src/sequince/steps/mutation/__init__.py index cf708b05..63009baa 100644 --- a/src/sequince/steps/mutation/__init__.py +++ b/src/sequince/steps/mutation/__init__.py @@ -37,6 +37,7 @@ ComplementaryEnsembleKernel, DifferentialEvolutionKernel, GaussianRandomWalk, + IndependentProposalKernel, KernelSample, MaskedGaussianRandomWalk, MixtureMutationKernel, @@ -98,6 +99,7 @@ "RedBlueMetropolisHastingsMutation", "MiniPCNMutation", "MaskedGaussianRandomWalk", + "IndependentProposalKernel", "KernelSample", "MetropolisBlock", "MixtureMutationKernel", diff --git a/src/sequince/steps/mutation/kernels.py b/src/sequince/steps/mutation/kernels.py index d3c228c2..39dd50b5 100644 --- a/src/sequince/steps/mutation/kernels.py +++ b/src/sequince/steps/mutation/kernels.py @@ -17,6 +17,8 @@ if TYPE_CHECKING: from orng import RandomGenerator + + from ...proposals import ProposalDistribution else: # pragma: no cover - runtime fallback from orng import RandomGenerator # type: ignore @@ -321,6 +323,81 @@ def evaluate_log_transition_ratio( ) +@dataclass(frozen=True) +class IndependentProposalKernel: + """Independence MH proposals in the coordinates supplied to the kernel. + + Parameters + ---------- + proposal : ProposalDistribution + Distribution defining samples and densities in mutation coordinates. + A physical-space prior is suitable only without a coordinate transform. + Sample batches must have exactly the input particle shape; densities + must have shape (n_particles,). Sampled densities are retained locally. + """ + + proposal: ProposalDistribution + + def sample(self, particles: Array, rng: RandomGenerator) -> KernelSample: + """Draw a batch and retain its log densities for the MH correction.""" + xp = namespace(particles) + values, density = self.proposal.sample_and_log_prob( + particles.shape[0], rng=rng, xp=xp + ) + values = asarray(values, xp=xp, like=particles, dtype=particles.dtype) + density = asarray( + density, xp=xp, like=particles, dtype=particles.dtype + ) + if values.shape != particles.shape: + raise ValueError("Proposal samples must match the particle shape.") + if density.shape != (particles.shape[0],): + raise ValueError( + "Proposal densities must have shape (n_particles,)." + ) + return KernelSample(values, {"log_proposal": density}) + + def log_transition_ratio( + self, from_particles: Array, to_particles: Array + ) -> Array: + """Return ``log q(current) - log q(proposed)`` per particle.""" + if from_particles.shape != to_particles.shape: + raise ValueError("Proposal samples must match the particle shape.") + return self.log_transition_ratio_from_sample( + from_particles, + KernelSample( + to_particles, + {"log_proposal": self.proposal.log_prob(to_particles)}, + ), + ) + + def log_transition_ratio_from_sample( + self, from_particles: Array, sample: KernelSample + ) -> Array: + """Evaluate the correction using the cached proposed log densities.""" + if from_particles.shape != sample.particles.shape: + raise ValueError("Proposal samples must match the particle shape.") + xp = namespace(from_particles) + current = asarray( + self.proposal.log_prob(from_particles), + xp=xp, + like=from_particles, + dtype=from_particles.dtype, + ) + proposed = asarray( + sample.auxiliary["log_proposal"], + xp=xp, + like=from_particles, + dtype=from_particles.dtype, + ) + if current.shape != (from_particles.shape[0],) or ( + proposed.shape != current.shape + ): + raise ValueError( + "Proposal densities must have shape (n_particles,)." + ) + return current - proposed + + @dataclass(frozen=True) class GaussianRandomWalk(SymmetricMutationKernel): """Isotropic Gaussian random-walk mutation kernel.""" @@ -882,6 +959,7 @@ def log_transition_ratio_from_sample( MUTATION_KERNELS: dict[str, type] = { "gaussian-random-walk": GaussianRandomWalk, + "independent-proposal": IndependentProposalKernel, "mixture": MixtureMutationKernel, } diff --git a/tests/test_mutation/test_independent_proposal.py b/tests/test_mutation/test_independent_proposal.py new file mode 100644 index 00000000..84824209 --- /dev/null +++ b/tests/test_mutation/test_independent_proposal.py @@ -0,0 +1,161 @@ +from __future__ import annotations + +from typing import Any + +import numpy as np +import pytest +from array_api_compat import device as get_device +from orng import RandomGenerator + +from sequince import ParticleState +from sequince.array_api import namespace, to_numpy +from sequince.steps import ( + IndependentProposalKernel, + MetropolisHastingsMutation, + StepContext, +) +from sequince.steps.mutation.kernels import ( + evaluate_log_transition_ratio, + get_mutation_kernel, +) +from sequince.typing import Array +from sequince.utils.random import resolve_rng + + +class NormalProposal: + def sample( + self, + n_samples: int, + context: StepContext | None = None, + *, + rng: RandomGenerator | int | None = None, + xp: Any | None = None, + ) -> Array: + return resolve_rng(rng).normal(size=(n_samples, 2)) * 2 + + def log_prob( + self, particles: Array, context: StepContext | None = None + ) -> Array: + xp = namespace(particles) + return -xp.sum(particles**2, axis=-1) / 8 + + def sample_and_log_prob( + self, + n_samples: int, + context: StepContext | None = None, + *, + rng: RandomGenerator | int | None = None, + xp: Any | None = None, + ) -> tuple[Array, Array]: + values = self.sample(n_samples, rng=rng, xp=xp) + return values, self.log_prob(values) + + +def test_independence_backend(backend, backend_factory) -> None: + case = backend_factory(backend) + particles = case.particles + kernel = IndependentProposalKernel(NormalProposal()) + sample = kernel.sample(particles, case.rng(5)) + ratio = evaluate_log_transition_ratio(kernel, particles, sample) + assert sample.particles.shape == particles.shape + assert sample.particles.dtype == particles.dtype + assert get_device(sample.particles) == get_device(particles) + assert case.matches_namespace(sample.particles) + np.testing.assert_allclose( + to_numpy(ratio), + to_numpy( + kernel.proposal.log_prob(particles) + - kernel.proposal.log_prob(sample.particles) + ), + rtol=1e-6, + atol=1e-6, + ) + with pytest.raises(ValueError, match="shape"): + kernel.log_transition_ratio(particles, sample.particles[:1]) + + +def test_support_and_validation() -> None: + class Bounded(NormalProposal): + def log_prob(self, particles, context=None): + return np.where(particles[:, 0] >= 0, -particles[:, 0], -np.inf) + + kernel = IndependentProposalKernel(Bounded()) + ratio = kernel.log_transition_ratio( + np.array([[-1.0, 0], [1, 0]]), np.array([[1.0, 0], [-1, 0]]) + ) + np.testing.assert_array_equal(ratio, [-np.inf, np.inf]) + with pytest.raises(ValueError, match="shape"): + kernel.sample(np.zeros((3, 1)), RandomGenerator(seed=5)) + + class BadDensity(NormalProposal): + def log_prob(self, particles, context=None): + return np.zeros((len(particles), 1)) + + with pytest.raises(ValueError, match="densities"): + IndependentProposalKernel(BadDensity()).sample( + np.zeros((3, 2)), RandomGenerator(seed=5) + ) + + +@pytest.mark.parametrize("shape", [(5,), (1, 2), (5, 2, 3)]) +@pytest.mark.parametrize("dtype", [np.float32, np.float64]) +def test_independence_particle_shapes(shape, dtype) -> None: + class ShapedProposal(NormalProposal): + def sample(self, n_samples, context=None, *, rng=None, xp=None): + return rng.normal(size=(n_samples, *shape[1:])) + + def log_prob(self, particles, context=None): + return ( + -np.sum( + np.reshape(particles, (len(particles), -1)) ** 2, axis=1 + ) + / 2 + ) + + particles = np.ones(shape, dtype=dtype) + kernel = IndependentProposalKernel(ShapedProposal()) + sample = kernel.sample(particles, RandomGenerator(seed=3)) + assert sample.particles.shape == shape + assert sample.particles.dtype == dtype + ratio = evaluate_log_transition_ratio(kernel, particles, sample) + assert ratio.shape == (shape[0],) + assert ratio.dtype == dtype + + +def test_cached_density() -> None: + class CountedProposal(NormalProposal): + def __init__(self) -> None: + self.calls = [] + + def log_prob(self, particles, context=None): + self.calls.append(particles) + return super().log_prob(particles, context) + + proposal = CountedProposal() + independent = IndependentProposalKernel(proposal) + particles = np.ones((10, 2)) + rng = RandomGenerator(seed=22) + sample = independent.sample(particles, rng) + evaluate_log_transition_ratio(independent, particles, sample) + assert len(proposal.calls) == 2 # Sampled density plus current density. + assert proposal.calls[-1] is particles + + +def test_known_target(context_factory) -> None: + rng = RandomGenerator(seed=44) + state = ParticleState(rng.normal(size=(1500, 2)) * 4) + step = MetropolisHastingsMutation( + lambda x: -np.sum(x * x, axis=-1) / 2, + IndependentProposalKernel(NormalProposal()), + n_steps=100, + ) + state, _, _ = step(state, context_factory(), rng) + np.testing.assert_allclose(np.mean(state.particles, axis=0), 0, atol=0.10) + np.testing.assert_allclose(np.var(state.particles, axis=0), 1, atol=0.12) + + +def test_factory() -> None: + proposal = NormalProposal() + kernel = get_mutation_kernel("independent-proposal", proposal=proposal) + assert isinstance(kernel, IndependentProposalKernel) + assert kernel.proposal is proposal From c88868f38563b15c1f4d80db9b71465d3e8442c9 Mon Sep 17 00:00:00 2001 From: "Michael J. Williams" Date: Tue, 27 Oct 2026 11:50:25 +0100 Subject: [PATCH 2/7] BUG(typing): annotate backend namespace predicates --- src/sequince/array_api.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/sequince/array_api.py b/src/sequince/array_api.py index 80881b0e..34b268c2 100644 --- a/src/sequince/array_api.py +++ b/src/sequince/array_api.py @@ -1,5 +1,6 @@ from __future__ import annotations +from collections.abc import Callable from typing import Any import numpy as _np @@ -23,7 +24,7 @@ from ._scipy_compat import enable_scipy_array_api as _enable_scipy_array_api from .typing import Array -IS_NAMESPACE_FUNCTIONS = { +IS_NAMESPACE_FUNCTIONS: dict[str, Callable[[Any], bool]] = { "numpy": is_numpy_namespace, "torch": is_torch_namespace, "jax": is_jax_namespace, From 252d01ce525fdfcd3d4debad797b44472459ae58 Mon Sep 17 00:00:00 2001 From: "Michael J. Williams" Date: Tue, 27 Oct 2026 11:56:41 +0100 Subject: [PATCH 3/7] DOC: rename global proposals section --- docs/user_guide/sampling_options/mutation_configuration.rst | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/docs/user_guide/sampling_options/mutation_configuration.rst b/docs/user_guide/sampling_options/mutation_configuration.rst index 69f28755..d424c3bb 100644 --- a/docs/user_guide/sampling_options/mutation_configuration.rst +++ b/docs/user_guide/sampling_options/mutation_configuration.rst @@ -140,8 +140,8 @@ Use these steps to compose mutation updates: The adaptation wrappers support eager execution only. Place conditional steps outside them, for example ``MutateIfResampled(mutation_step)``. -Global independence proposals ------------------------------ +Global proposals +---------------- :class:`~sequince.steps.mutation.kernels.IndependentProposalKernel` draws proposals independently of each particle's current position. It accepts a From a248e6996a22aa12aef325fc8c9d85e9e0b44082 Mon Sep 17 00:00:00 2001 From: "Michael J. Williams" Date: Tue, 27 Oct 2026 12:00:31 +0100 Subject: [PATCH 4/7] DOC: condense global proposal documentation --- .../mutation_configuration.rst | 44 +------------------ 1 file changed, 2 insertions(+), 42 deletions(-) diff --git a/docs/user_guide/sampling_options/mutation_configuration.rst b/docs/user_guide/sampling_options/mutation_configuration.rst index d424c3bb..c559e461 100644 --- a/docs/user_guide/sampling_options/mutation_configuration.rst +++ b/docs/user_guide/sampling_options/mutation_configuration.rst @@ -45,7 +45,8 @@ Each step takes a target callable or :class:`~sequince.target.Target` object. MH steps take a proposal through the ``kernel`` argument. See :mod:`sequince.steps.mutation.kernels` for Gaussian random walks, masked -updates, mixtures, and differential-evolution kernels. Red-blue MH requires a +updates, global proposals, mixtures, and differential-evolution kernels. +Red-blue MH requires a :class:`~sequince.steps.mutation.kernels.ComplementaryEnsembleKernel`. Common configuration @@ -140,47 +141,6 @@ Use these steps to compose mutation updates: The adaptation wrappers support eager execution only. Place conditional steps outside them, for example ``MutateIfResampled(mutation_step)``. -Global proposals ----------------- - -:class:`~sequince.steps.mutation.kernels.IndependentProposalKernel` draws -proposals independently of each particle's current position. It accepts a -:class:`~sequince.proposals.ProposalDistribution`: ``sample_and_log_prob`` draws -the requested batch using the supplied ``rng`` and ``xp``, and ``log_prob`` -evaluates one log density per particle. Both methods must describe the same -distribution. Given a proposal named ``prior_proposal`` and a target named -``target``: - -.. code-block:: python - - from sequince.steps import ( - IndependentProposalKernel, - MetropolisHastingsMutation, - ) - - mutation = MetropolisHastingsMutation( - target, - IndependentProposalKernel(prior_proposal), - n_steps=10, - transform=None, - ) - -Acceptance includes ``log q(current) - log q(proposed)``. The kernel caches the -sampled log densities, so it only evaluates current-particle densities when -computing that correction. Samples must match the input particle shape, and -log densities must have shape ``(n_particles,)``. The proposal should cover the -target's support, and its distribution must remain fixed during mutation. - -The proposal must describe the coordinates supplied to the kernel. Use -``transform=None`` for physical-space prior proposals; when configuring a -prebuilt sampler, also set ``particle_transform=None``. Transformed proposals -must sample and evaluate densities in the transformed coordinates, including -the appropriate change of variables. A physical-space prior cannot be passed -unchanged to a latent-space kernel. - -The mutation-kernel factory accepts ``"independent-proposal"`` with a -``proposal`` object in ``kernel_kwargs``. - Prebuilt samplers ----------------- From aa1ed1572c64f2dc5d3e317413e7e59e098cb238 Mon Sep 17 00:00:00 2001 From: "Michael J. Williams" Date: Tue, 27 Oct 2026 12:08:21 +0100 Subject: [PATCH 5/7] DOC: shorten global proposal changelog entry --- CHANGELOG.md | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index bccb9714..ca420e06 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -8,10 +8,8 @@ The format follows [Keep a Changelog](https://keepachangelog.com/en/1.1.0/). ### Added -- Independence Metropolis-Hastings proposals using `IndependentProposalKernel` - and `ProposalDistribution`, with cached proposed log densities. Proposals must - describe mutation coordinates; physical-space priors require physical-space - mutation. +- Global proposals with `IndependentProposalKernel` + ([#139](https://github.com/sequince-dev/sequince/pull/139)). ## [0.1.0b1] - 2026-10-23 From 90138536e3d1a8fd28cbbf8d43662a75e40164f6 Mon Sep 17 00:00:00 2001 From: "Michael J. Williams" Date: Tue, 27 Oct 2026 12:12:47 +0100 Subject: [PATCH 6/7] DOC(mutation): cite independence Metropolis reference --- src/sequince/steps/mutation/kernels.py | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/src/sequince/steps/mutation/kernels.py b/src/sequince/steps/mutation/kernels.py index 39dd50b5..0d39c7e4 100644 --- a/src/sequince/steps/mutation/kernels.py +++ b/src/sequince/steps/mutation/kernels.py @@ -334,6 +334,11 @@ class IndependentProposalKernel: A physical-space prior is suitable only without a coordinate transform. Sample batches must have exactly the input particle shape; densities must have shape (n_particles,). Sampled densities are retained locally. + + Notes + ----- + This is the independence Metropolis sampler described in + `Tierney (1994), Section 2.3.3 `_. """ proposal: ProposalDistribution From 68b5172190d142880873f064ed9da6b7f74f52d2 Mon Sep 17 00:00:00 2001 From: "Michael J. Williams" Date: Tue, 27 Oct 2026 12:18:32 +0100 Subject: [PATCH 7/7] DOC(mutation): remove independence sampler reference --- src/sequince/steps/mutation/kernels.py | 5 ----- 1 file changed, 5 deletions(-) diff --git a/src/sequince/steps/mutation/kernels.py b/src/sequince/steps/mutation/kernels.py index 0d39c7e4..39dd50b5 100644 --- a/src/sequince/steps/mutation/kernels.py +++ b/src/sequince/steps/mutation/kernels.py @@ -334,11 +334,6 @@ class IndependentProposalKernel: A physical-space prior is suitable only without a coordinate transform. Sample batches must have exactly the input particle shape; densities must have shape (n_particles,). Sampled densities are retained locally. - - Notes - ----- - This is the independence Metropolis sampler described in - `Tierney (1994), Section 2.3.3 `_. """ proposal: ProposalDistribution