Algorithmic warm starts for Hamiltonian Monte Carlo

Matthew S. Zhang; Jason M. Altschuler; Sinho Chewi

cs.DSMar 24, 2026

Algorithmic warm starts for Hamiltonian Monte Carlo

Authors: Matthew S. Zhang, Jason M. Altschuler, Sinho Chewi

Organizations: UToronto · UPenn · Yale

Abstract

Generating samples from a continuous probability density is a central algorithmic problem across statistics, engineering, and the sciences. For high-dimensional settings, Hamiltonian Monte Carlo (HMC) is the default algorithm across mainstream software packages. However, despite the extensive line of work on HMC and its widespread empirical success, it remains unclear how many iterations of HMC are required as a function of the dimension $d$ . On one hand, a variety of results show that Metropolized HMC converges in $O(d^{1/4})$ iterations from a warm start close to stationarity. On the other hand, Metropolized HMC is significantly slower without a warm start, e.g., requiring $Ω(d^{1/2})$ iterations even for simple target distributions such as isotropic Gaussians. Finding a warm start is therefore the computational bottleneck for HMC. We resolve this issue for the well-studied setting of sampling from a probability distribution satisfying strong log-concavity (or isoperimetry) and third-order derivative bounds. We prove that non-Metropolized HMC generates a warm start in $\tilde{O}(d^{1/4})$ iterations, after which we can exploit the warm start using Metropolized HMC. Our final complexity of $\tilde{O}(d^{1/4})$ is the fastest algorithm for high-accuracy sampling under these assumptions, improving over the prior best of $\tilde{O}(d^{1/2})$ . This closes the long line of work on the dimensional complexity of MHMC for such settings, and also provides a simple warm-start prescription for practical implementations.

Figures & tables

Panel of Figure 1 — Figure 1: The convergence of MHMC is heavily dependent on the step size $h$ . Large step sizes $h\asymp d^{-1/4}$ classically lead to fast convergence from a warm start, but can get stuck in cold starts due to very low acceptance probability ( left ). Small step sizes $h\asymp d^{-1/2}$ fix that issue but lead to slow movement, requiring at least $1/h\asymp d^{1/2}$ steps to traverse the space, let alone mix ( right ). Illustrated for the simple target $\pi=\mathcal{N}(0,I)$ in dimension $d=10^{4}$ , with “cold start” initialization at the mode $x_{0}=0$ . Reproducibility details: MHMC is repeatedly integrated for $T=1$ unit of continuous time via $1/h$ leapfrog steps of size $h$ , as is standard. Similar qualitative phenomena are observed for other settings.

Panel of Figure 2 — Figure 2: A key algorithmic insight is that unadjusted HMC rapidly escapes cold starts using large step sizes $h\asymp d^{-1/4}$ . This algorithm quickly leads to iterates which, if used as an initialization for MHMC, would have high acceptance probability ( right ). This motivates our two-phase algorithmic proposal: escape the cold start via unadjusted HMC, then exploit the warm start using MHMC. Removing the Metropolis filter bypasses the issue of low acceptance probability at cold starts ( left ) without requiring the use of small step sizes which result in slow movement ( middle ).

Figure 3: Schematic diagram of the one-shot coupling in the proof of Theorem 5.9 . Top and bottom (in black): the two OHO processes $(\bar{x}_{0},\bar{p}_{0})\mapsto(\bar{X}_{h}^{\mathsf{OHO}},\bar{P}_{h}^{\mathsf{OHO}})$ and $(x_{0},p_{0})\mapsto(X_{h}^{\mathsf{OHO}},P_{h}^{\mathsf{OHO}})$ are coupled to use the same Gaussian noise increments $\bar{B}_{1},\bar{B}_{2}$ . Diagonal (in blue): different noise increments $B_{1},B_{2}$ enable an auxiliary OHO process $(\bar{x}_{0},\bar{p}_{0})\mapsto(X_{h}^{\mathsf{OHO}},P_{h}^{\mathsf{OHO}})$ that starts from one process and ends at the other. To achieve this interpolation, the auxiliary noise increments $B_{1},B_{2}$ are uniquely determined as a function of $x_{0},p_{0},\bar{x}_{0},\bar{p}_{0},\bar{B_{1}},\bar{B_{2}}$ . Orange: $B_{1}$ in the first “O” step is uniquely determined so that the resulting momentum $P_{h}^{{\operatorname{target}}}$ will enable the correct $x$ -coordinate $X_{h}^{\mathsf{OH}}$ after the “H” step. $B_{2}$ is then uniquely determined so that the $p$ -coordinate $P_{h}^{\mathsf{OHO}}$ will match after the final “O” step.

Figure 4: Schematic diagram for the proof of Theorem 5.5 . Top and bottom (in black): the two processes $(\bar{X}_{kh},\bar{P}_{kh})$ and $(X_{kh},P_{kh})$ that we seek to compare. Each horizontal arrow is an OHO step. Diagonal (in blue): the auxiliary process $(X_{kh}^{\mathsf{aux}},P_{kh}^{\mathsf{aux}})$ we construct to apply the shifted composition rule. The auxiliary process interpolates between the original two processes at initialization (top left) and termination (bottom right). The auxiliary process updates by first shifting the momentum (grey, downwards arrow) and then performing an OHO step (black, right arrow).

Explore similar work

Sep 7, 2026

Smoothed Picard Hamiltonian Monte Carlo

We develop a new low-accuracy sampler, called smoothed Picard Hamiltonian Monte Carlo, which combines Gaussian smoothing, Picard iteration, and higher-order discretization. For a log-concave target

π\propto \exp(-V)

in dimension

d

satisfying

0 \prec αI \preceq \nabla^2 V \preceq βI

, with condition number

κ:= β/α

, smoothed Picard HMC returns a sample with

\sqrt α\,W_2(\cdot,π) \le \varepsilon

using

\widetilde O(κ^2 + κ^{7/6} d^{1/6}/\varepsilon^{1/3})

gradient queries. We also prove stronger

W_q

bounds, and then develop an algorithmic framework, the recursive warm start generator, to upgrade these

W_q

bounds to stronger divergence guarantees. This produces a warm start for the proximal bouncy particle sampler, introduced in a companion work, leading to a high-accuracy log-concave sampler with complexity

\widetilde O((κ^{7/6} d^{1/6} + κ^{1/2} d^{1/4})\mathrm{polylog}(1/\varepsilon))

Fan Chen, Sinho Chewi, Jianfeng Lu +1

Department of Electrical Engineering and Computer Science, Massachusetts Institute of Technology. · Department of Statistics and Data Science, Yale University. · Department of Mathematics, Duke University. +1

Markov Chain Monte Carlo Log-Concave Sampling

Ask PDF

Jul 14, 2026

Accelerated Mixing Time of Randomized Hamiltonian Monte Carlo

We show the Randomized Hamiltonian Monte Carlo (RHMC) algorithm has accelerated mixing time guarantees for sampling from log-concave probability distributions. RHMC proceeds by repeatedly simulating the continuous-time Hamiltonian dynamics for some random integration times, and resetting the velocity to be an independent Gaussian random variable between each simulation. We show that when the target distribution is log-concave and satisfies an

α

-Talagrand inequality (for example, if the target distribution is

α

-strongly log-concave), if we use a random integration time from either the triangular or the exponential distribution with mean

Θ(α^{-1/2})

, then RHMC converges exponentially fast in KL divergence, and the total integration time to reach error

\varepsilon

in KL divergence scales as

O(α^{-1/2} \log(\varepsilon^{-1}))

. We also show that when the target distribution is log-concave, if we use a sequence of random integration times from the triangular distribution with exponentially increasing means, then the total integration time to reach error

\varepsilon

in KL divergence scales as

O(\varepsilon^{-1/2})

. Our analysis relies on a bound on the average KL divergence along Hamiltonian dynamics, which is inspired by an analogous result on accelerated optimization methods based on Hamiltonian dynamics.

Siddharth Mitra, Vishwak Srinivasan, Xiuyuan Wang +1

Department of Computer Science, Yale University. · Department of Electrical Engineering and Computer Science, MIT.

Markov Chain Monte Carlo Log-Concave Sampling

Ask PDF

Feb 24, 2026

Counterdiabatic Hamiltonian Monte Carlo

Hamiltonian Monte Carlo (HMC) is a state of the art method for sampling from distributions with differentiable densities, but can converge slowly when applied to challenging multimodal problems. Running HMC with a time varying Hamiltonian, in order to interpolate from an initial tractable distribution to the target of interest, can address this problem. In conjunction with a weighting scheme to eliminate bias, this can be viewed as a special case of Sequential Monte Carlo (SMC) sampling \cite{doucet2001introduction}. However, this approach can be inefficient, since it requires slow change between the initial and final distribution. Inspired by \cite{sels2017minimizing}, where a learned counterdiabatic term added to the Hamiltonian allows for efficient quantum state preparation, we propose Counterdiabatic Hamiltonian Monte Carlo (CHMC), which can be viewed as an SMC sampler with a more efficient kernel. We establish its relationship to recent proposals for accelerating gradient-based sampling with learned drift terms, and demonstrate on simple benchmark problems.

Reuben Cohn-Gordon, Uroš Seljak, Dries Sels

University of California, Berkeley · New York University

Markov Chain Monte Carlo

Ask PDF