Hybrid Monte Carlo methods for sampling probability measures on submanifolds

Abstract

Probability measures supported on submanifolds can be sampled by adding an extra momentum variable to the state of the system, and discretizing the associated Hamiltonian dynamics with some stochastic perturbation in the extra variable. In order to avoid biases in the invariant probability measures sampled by discretizations of these stochastically perturbed Hamiltonian dynamics, a Metropolis rejection procedure can be considered. The so-obtained scheme belongs to the class of generalized Hybrid Monte Carlo (GHMC) algorithms. We show here how to generalize to GHMC a procedure suggested by Goodman, Holmes-Cerfon and Zappa for Metropolis random walks on submanifolds, where a reverse projection check is performed to enforce the reversibility of the algorithm for any timesteps and hence avoid biases in the invariant measure. We also provide a full mathematical analysis of such procedures, as well as numerical experiments demonstrating the importance of the reverse projection check on simple toy examples. Mathematical setting and resultsWe first describe in Section 2.1 the target probability measures we are interested in sampling. The core of our analysis is presented in Section 2.2, where we show how to rigorously formalize the reversibility properties of discretizations of the Hamiltonian dynamics on a submanifold for possibly large timesteps. One crucial ingredient in this analysis is the Lagrange multiplier function, which associates to a state on the submanifold and another state close to the submanifold (obtained by an unconstrained step of the dynamics) a new state on the submanifold. Examples of such Lagrange multiplier functions are provided in Section 2.3. Once the discretization of the Hamiltonian part of the dynamics is clear, GHMC schemes follow in a straightforward way; see Section 2.4. Geometric settingLet us first introduce the measures we are interested in sampling, namely the target measure (2) below, as well as the probability measure (4) on an extended space, which admits (2) as a marginal. We only provide the most essential objects needed for our analysis. For more details and motivations for the definitions and results given here, we refer the reader to the standard reference textbooks [1,4], and also to [22, Chapter 3.3.2] for a self-contained presentation. Target measureWe denote by d the dimension of the ambient space, and consider a submanifold M ⊂ R d defined as the zero level set of a given smooth function ξ :For example, the function ξ encodes molecular constraints or reaction coordinates in molecular dynamics, or is a summary statistics in Approximate Bayesian Computations. Let M ∈ R d×d be a fixed symmetric positive definite matrix, which is interpreted as the mass tensor below. The ambient space R d is endowed with the scalar product q,q M = q T Mq. One could take for simplicity M = Id but it is sometimes useful to consider non identity mass matrices for numerical purposes [13]. The function ξ is assumed to be smooth on a neighborhood of M in R d and such thatis an invertible matrix for...

Coulomb Gases Under Constraint: Some Theoretical and Numerical Results

“…An issue when combining a Metropolis-Hastings rule with a projection on a submanifold is that reversibility may be lost, which introduces a bias. A recent strategy has been to introduce a reversibility check in addition to the standard acceptionrejection rule, which makes the HMC scheme under constraint reversible [52,35]. Note that [53] proposes an interesting alternative to the scheme used here, which is however not compatible with a Metropolis selection procedure in its current form.…”

“…Note that [53] proposes an interesting alternative to the scheme used here, which is however not compatible with a Metropolis selection procedure in its current form. We thus present the algorithm as written in [35], with some simplifications and adaptations to our context, for which we introduce next the constrained Langevin dynamics.…”

“…When considering linear statistics, we prove in Section 3.2 that conditioning P n boils down to modifying the confinement potential V . In Section 3.3 we consider the case of quadratic statistics, which amounts to modifying the interaction kernel g. In order to validate our theoretical results and explore cases in which explicit solutions are not available, we also propose a method for sampling the law of Y n for a fixed n. Based on the Hamiltonian Monte Carlo (HMC) method used in [11] for sampling Gibbs measures associated to Coulomb and Log-gases, we describe and implement the generalized Hamiltonian Monte Carlo algorithm proposed in [35] for sampling probability measures on submanifolds. The method, detailed in Section 4.1, shows remarkable performance, and the results presented in Section 4.2 are in agreement with the theory.…”

Ergodic SDEs on submanifolds and related numerical sampling schemes

“…Besides of sampling on the entire space, the Metropolis-adjusted samplers on submanifolds, using either MCMC or Hybrid Monte Carlo, have been considered in several recent work [9,33,45,34]. Reversible Metropolis random walk on submanifolds has been constructed in [45], which is then extended in [34] to develop Hybrid Monte Carlo algorithms by allowing non-zero gradient forces in the proposal move. The numerical schemes in these work are unbiased and therefore large time step-sizes can be used in numerical estimations.…”

Sampling from manifold-restricted distributions using tangent bundle projections

“…For example, Riemannianmanifold versions of HMC and MALA (Girolami and Calderhead, 2011;Xifara et al, 2014) can make use of the pullback metric induced on parameter space, which is specified completely by the first derivative of the map. Manifold-constrained HMC methods (Brubaker et al, 2012;Lelièvre et al, 2018) also work on the tangent bundle of the manifold, and hence require the first derivative to define a basis for each tangent space. More generally, any sampling of some probability distribution on a submanifold of Euclidean space can be informed by its density with respect to the Hausdorff measure on the submanifold, which may be computed using the first derivative as well (Diaconis et al, 2013).…”