Metropolis-Adjusted Langevin Algorithm (MALA) MCMC

class pints.MALAMCMC(x0, sigma0=None)

Metropolis-Adjusted Langevin Algorithm (MALA), an MCMC sampler as described in [1].

This method involves simulating Langevin diffusion such that the solution to the time evolution equation (the Fokker-Planck PDE) is a stationary distribution that equals the target density (in Bayesian problems, the posterior distribution). The stochastic differential equation (SDE) given below ensures that if \(u(\theta, 0) = \pi(\theta)\), then \(\partial u / \partial t = 0\),

\[\mathrm{d}\Theta_t = 1/2 \nabla \; \text{log} \pi(\Theta_t) \mathrm{d}t + \mathrm{d}W_t\]

where \(\pi(\theta)\) is the target density and \(W\) is a standard multivariate Wiener process.

In general, the above SDE cannot be solved exactly and the below first-order Euler discretisation is used instead,

\[\theta^* = \theta_t + \epsilon^2 1/2 \nabla \; \text{log} \pi(\theta_t) + \epsilon z\]

where \(z \sim \mathcal{N}(0, I)\) resulting in a mean \(\mu(\theta^*) = \theta_t + \epsilon^2 1/2 \nabla \; \text{log} \pi(\theta_t)\).

To correct for first-order integration error that is introduced from discretisation, a Metropolis-Hastings acceptance probability is calculated after a step,

\[\alpha = \frac{\pi(\theta^*)q(\theta_t|\theta^*)}{\pi(\theta_t) q(\theta^*|\theta_t)}\]

where \(q(\theta_2|\theta_1) = \mathcal{N}(\theta_2|\mu(\theta_1), \epsilon I)\) and \(\theta^*\) is accepted with probability \(\text{min}(1, \alpha)\).

Here we consider a slight variant of the above method discussed in [1], which is to use a preconditioning matrix \(M\) to allow differing degrees of freedom in each dimension.

\[\theta^* = \theta_t + \epsilon'^2 1/2 \nabla \; \text{log} \pi(\theta_t) + \epsilon' z\]

leading to \(q(\theta_2|\theta_1) = \mathcal{N}(\theta_2|\mu(\theta_1), \epsilon')\).

where \(\epsilon' = \epsilon \sqrt{M}\) is given by the initial value of sigma0.

Extends SingleChainMCMC.

References

[1] (1,2)

Girolami, M. and Calderhead, B., 2011. Riemann manifold langevin and hamiltonian monte carlo methods. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 73(2), pp.123-214. https://doi.org/10.1111/j.1467-9868.2010.00765.x

acceptance_rate()

Returns the current (measured) acceptance rate.

ask()

See SingleChainMCMC.ask().

epsilon()

Returns epsilon which is the effective step size used in proposals.

in_initial_phase()

For methods that need an initial phase (see needs_initial_phase()), this method returns True if the method is currently configured to be in its initial phase. For other methods a NotImplementedError is returned.

n_hyper_parameters()

See TunableMethod.n_hyper_parameters().

name()

See pints.MCMCSampler.name().

needs_initial_phase()

Returns True if this method needs an initial phase, for example an adaptation-free period for adaptive covariance methods, or a warm-up phase for DREAM.

needs_sensitivities()

See pints.MCMCSampler.needs_sensitivities().

replace(current, current_log_pdf, proposed=None)

Replaces the internal current position, current LogPDF, and proposed point (if any) by the user-specified values.

This method can only be used once the initial position and LogPDF have been set (so after at least 1 round of ask-and-tell).

This is an optional method, and some samplers may not support it.

set_epsilon(epsilon=None)

Sets epsilon, which is the effective step size used in proposals. If epsilon not specified, then epsilon = 0.2 * diag(sigma0) will be used.

set_hyper_parameters(x)

The hyper-parameter vector is [epsilon].

The effective step size (epsilon) is step_size * scale_vector.

See TunableMethod.set_hyper_parameters().

set_initial_phase(in_initial_phase)

For methods that need an initial phase (see needs_initial_phase()), this method toggles the initial phase algorithm. For other methods a NotImplementedError is returned.

tell(reply)

See pints.SingleChainMCMC.tell().