Monomial-Gamma Hamiltonian MCMC

class pints.MonomialGammaHamiltonianMCMC(x0, sigma0=None)

Implements Monomial Gamma HMC as described in [1] - a generalisation of HMC as described in [2] - involving a non-physical kinetic energy term.

Uses a physical analogy of a particle moving across a landscape under Hamiltonian dynamics to aid efficient exploration of parameter space. Introduces an auxilary variable – the momentum (p_i) of a particle moving in dimension i of negative log posterior space – which supplements the position (q_i) of the particle in parameter space. The particle’s motion is dictated by solutions to Hamilton’s equations,

\[dq_i/dt = \partial H/\partial p_i, dp_i/dt = - \partial H/\partial q_i.\]

The Hamiltonian is given by,

\[H(q,p) = U(q) + K(p) = -log(p(q|X)p(q)) + \Sigma_{i=1}^{d} ( -g(p_i) + (2/c) \text{log}(1 + \text{exp}(cg(p_i))))\]

where d is the dimensionality of model, U is the potential energy and K is the kinetic energy term. Note the kinetic energy is the ‘soft’ version described in [1], where,

\[g(p_i) = (1 / m_i) \text{sign}|p_i|^{1 / a}\]

To numerically integrate Hamilton’s equations, it is essential to use a sympletic discretisation routine, of which the most typical approach is the leapfrog method,

\[\begin{split}p_i(t + \epsilon/2) &= p_i(t) - (\epsilon/2) dU(q_i)/ dq_i\\ q_i(t + \epsilon) &= q_i(t) + \epsilon d K(p_i(t + \epsilon/2))/dp_i\\ p_i(t + \epsilon) &= p_i(t + \epsilon/2) - (\epsilon/2) dU(q_i + \epsilon)/ dq_i\end{split}\]

The derivative of the soft kinetic energy term is given by,

\[d K(p_i)/dp_i = |p_i|^{-1 + 1 / a}\text{sign}(p_i) \times \text{tanh}(c|p_i|^{1/a}\text{sign}(p_i) / {2 m_i}) / {a m_i}\]

In particular, the algorithm we implement follows eqs. (4.14)-(4.16) in [2], since we allow different epsilon according to dimension.

Extends SingleChainMCMC.

References

[1]

Towards Unifying Hamiltonian Monte Carlo and Slice Sampling Yizhe Zhang, Xiangyu Wang, Changyou Chen, Ricardo Henao, Kai Fan, Lawrence Cari. Advances in Neural Information Processing Systems (NIPS)

[2] (1,2)

MCMC using Hamiltonian dynamics Radford M. Neal, Chapter 5 of the Handbook of Markov Chain Monte Carlo by Steve Brooks, Andrew Gelman, Galin Jones, and Xiao-Li Meng.

a()

Returns a in kinetic energy function.

ask()

See SingleChainMCMC.ask().

c()

Returns c in kinetic energy function.

divergent_iterations()

Returns the iteration number of any divergent iterations.

epsilon()

Returns epsilon used in leapfrog algorithm.

hamiltonian_threshold()

Returns threshold difference in Hamiltonian value from one iteration to next which determines whether an iteration is divergent.

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.

leapfrog_step_size()

Returns the step size for the leapfrog algorithm.

leapfrog_steps()

Returns the number of leapfrog steps to carry out for each iteration.

mass()

Returns mass m in kinetic energy function.

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.

scaled_epsilon()

Returns scaled epsilon used in leapfrog algorithm.

set_a(a)

Sets a in kinetic energy function.

set_c(c)

Sets c in kinetic energy function.

set_epsilon(epsilon)

Sets epsilon for the leapfrog algorithm.

set_hamiltonian_threshold(hamiltonian_threshold)

Sets threshold difference in Hamiltonian value from one iteration to next which determines whether an iteration is divergent.

set_hyper_parameters(x)

The hyper-parameter vector is [leapfrog_steps, leapfrog_step_size, a, c, mass].

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.

set_leapfrog_step_size(step_size)

Sets the step size for the leapfrog algorithm.

set_leapfrog_steps(steps)

Sets the number of leapfrog steps to carry out for each iteration.

set_mass(m)

Sets mass m in kinetic energy function.

tell(reply)

See pints.SingleChainMCMC.tell().