particles.mcmc

MCMC (Markov chain Monte Carlo) and related algorithms.

Overview

This module contains various classes that implement MCMC samplers:

  • MCMC: the base class for all MCMC samplers;

  • GenericRWHM: base class for random-walk Hastings-Metropolis;

  • GenericGibbs: base class for Gibbs samplers;

  • PMMH, ParticleGibbs: base classes for the PMCMC (particle MCMC algorithms) with the same name.

For instance, here is how to run 200 iterations of an adaptive random-walk sampler:

# ...
# define some_static_model, some_prior
# ...
my_mcmc = BasicRWHM(model=some_static_model, prior=some_prior, niter=200,
                    adaptive=True)
my_mcmc.run()

Upon completion, object my_mcmc have an attribute called chain, which is a ThetaParticles object (see module smc_samplers). In particular, my_mcmc.chain has the following attributes:

  • theta: a structured array that contains the 200 simulated parameters;

  • lpost: an array that contains the log-posterior density at these 200 parameters.

See the dedicated notebook tutorial (on Bayesian inference for state-space models) for more examples and explanations.

Random walk Metropolis

Both GenericRWHM and PMMH rely on a random walk proposal; that is, given the current point theta, the proposed point is sampled from a Gaussian distribution, with mean theta, and some covariance matrix Sigma (or Sigma_t at iteration t the adaptive case, see below). Various parameters may be set to tune this type of proposal:

  • to run a standard random walk Metropolis algorithm, set adaptive=False, and specify the covariance matrix Sigma through parameter rw_cov. Otherwise, Sigma is set by default to the identity matrix (which could a terrible choice in certain problems).

  • to run an adaptive random walk Metropolis algorithm, where the matrix Sigma_t is progressively adapted to the past states of the chain, set adaptive=True. In that case, Sigma_t is set to a fraction of the running estimate of the covariance matrix of the target distribution, and rw_cov is used as a preliminary estimate for that target covariance. See parameter scale (in the documentation of GenericRWHM) for more details on the factor in front of the running estimate, and VanishCovTracker for how the running estimate is computed recursively.

By default, the adaptive version is used.

Beyond the bootstrap filter within PMMH

PMMH runs, at each iteration, a particle filter to approximate the likelihood of the considered state-space model. By default, a bootstrap filter is used (with default choices for the resampling scheme, and so on). It is possible to change the settings of this filtering algorithm, or run a different type of algorithm, as follows:

  1. You can set parameter fk_cls to a FeynmanKac subclass such as e.g. ssms.GuidedPF if you wish to use a guided filter (rather than a bootstrap filter). See module state_space_models for more details on these FeynmanKac subclasses derived from a given state-space model.

  2. You can use parameter smc_options (dict-like) to pass various parameters to class SMC when the algorithm is instantiated.

  3. You can even use parameter smc_cls to specify a different class for the

algorithm itself (a SMC subclass instead of SMC itself).

  1. Finally, if you need something even more general, you can also subclass

PMMH and redefine method alg_instance, which takes as argument theta (a dict-like object) and returns an algorithm (an instance of class SMC or one its subclasses).

Functions

msjd(theta)

Mean squared jumping distance.

Classes

BasicRWHM([niter, verbose, theta0, ...])

Basic random walk Hastings-Metropolis sampler.

CSMC([fk, N, ESSrmin, xstar])

Conditional SMC.

GenericGibbs([niter, verbose, theta0, ...])

Generic Gibbs sampler for a state-space model.

GenericRWHM([niter, verbose, theta0, ...])

Base class for random walk Hasting-Metropolis samplers.

MCMC([niter, verbose])

MCMC base class.

PMMH([niter, verbose, ssm_cls, smc_cls, ...])

Particle Marginal Metropolis Hastings.

ParticleGibbs([niter, verbose, ssm_cls, ...])

Particle Gibbs sampler (abstract class).

VanishCovTracker([alpha, dim, mu0, Sigma0])

Tracks the vanishing mean and covariance of a sequence of points.