Metropolis Adjusted Microcanonical Hamiltonian Monte Carlo
Authors: Jakob Robnik
Reuben Cohn-Gordon
Uros Seljak
Venue: NeurIPS 2025 poster ()
Keywords: MCMC
Primary area: probabilistic_methods
TL;DR: We develop Metropolis adjusted scheme for microcanonical dynamics and an automatic tuning scheme for its hyperparameters.

Abstract: Sampling from high dimensional distributions is a computational bottleneck in many scientific applications. Hamiltonian Monte Carlo (HMC), and in particular the No-U-Turn Sampler (NUTS), are widely used, yet they struggle on problems with a very large number of parameters or a complicated geometry. Microcanonical Langevin Monte Carlo (MCLMC) has been recently proposed as an alternative which shows striking gains in efficiency over NUTS, especially for high-dimensional problems. However, it produces biased samples, with a bias that is hard to control in general. We introduce the Metropolis-Adjusted Microcanonical sampler (MAMS), which relies on the same dynamics as MCLMC, but introduces a Metropolis-Hastings step and thus produces asymptotically unbiased samples. We develop an automated tuning scheme for the hyperparameters of the algorithm, making it applicable out of the box. We demonstrate that MAMS outperforms NUTS across the board on benchmark problems of varying complexity and dimensionality, achieving up to a factor of seven speedup.

Metropolis Adjusted Microcanonical Hamiltonian
Monte Carlo
Jakob Robnik
Reuben Cohn-Gordon
Physics Department, University of California at Berkeley, Berkeley, CA 94720, USA bettobrrobnikbbekkeee edu
Physics Department, University California Berkeley, Berkeley, CA 94720, USA
reubenharry&gaill.co
Uroš Seljak
Physics Department,
University of California at Berkeley
and Lawrence Berkeley National Laboratory, Berkeley, Berkeley, CA 94720, USA usel jeeljak@berkeeey edu
Abstract
Sampling from high dimensional distributions is computational bottleneck in many scientific applications. Hamiltonian Monte Carlo (HMC), and in particular the No-U-Turn Sampler (NUTS), widely used, yet they struggle prob- lems with a very large number of parameters or a complicated geometry. Micro- canonical Langevin Monte Carlo (MCLMC) has been recently proposed alternative which shows striking gains in efficiency over NUTS, especially for high-dimensional problems. However, produces biased samples, with bias that is hard to control in general. We introduce the Metropolis-Adjustte Micro- canonical sampler (MAMS), which relies on the same dynamics as MCLMC, but introduces a Mietropolis-Hastinne samples. We develop automated tuning scheme for the hyperparameters of the algorithm, making it applicable out of the box. We demonstrate that MAMS outperforms NUTS across the board on benchmark problems of varying complexity and dimensionality, achieving up to factor of seven speedup.
step and thus produces asymptotically unbiased
Introduction
Drawing samples from given probability density p(x)), for E Rd has applications in a wide range of scientific disciplines, from Bayesian inference for statistics (Strumbelj et al. 2024 Carpenter et al. 2017) to biology (Gelman and Rubin, 1996). statistical physics (Janke 2008). quantum mechanies (Gattringer and Lang. (2100) and cosmology (Campagne et al.) 2023). From the perspective of practitioner these fields, what often desirable is black-box algorithm, the sense of taking as input an unnormalized density and returning samples from the corresponding distribution. An important special case is where the density is differentiable (either analytically or by automatic differentiation (Griewank and Walther 2008))
Markov Chain Monte Carlo Metrroplli al. 1953 Hastings 1970 broad class of methods suited to this task, which construct a Markov Chain az 1 such that æ; are samples from the target distribution p. When the density's gradient is available, Hamiltonian Monte Carlo (HMC) (Duane et al.] 1987, Neal 2011) is leading method. In particular, the No-U-Turn sampler is black-box version of HMC, where users do not need to manually select the hyperparameters. has
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

been implemented in libraries like Stan (Carpenter et al. 2017) and Pyro (Bingham et al. 2019).
serving wide community of scientists.
NUTS is a powerful method, but there have been many attempts in the past decade to replace with algorithms that achive the same accuracy at a lower computational cost. One such method is Microcanonical Langevin Monte Carlo (MCLMC) (Robnik et al. (024);; it replaces the HMC dynamics with velocity-norm preserving dynamics, resulting in a method that is more stable to large gradients. Benchmarking in cosmology (Simon-Ontroy et al. 2025). Bayesian inference (Robnik
ct al. 2024 Sommer et al. 2024), and field theories (Robnik and Seljak 2024) suggests MCLMC is
a promising candidate to replace NUTS as go-to gradient-based based sampler.
The drawback of MCLMC is that it is biased: that is, the samples it produces do not correspond to samples from the true target distribution p. but rather to samples from a nearby distribution p. Although this bias is controllable, in many fields, asymptotically unbiased samplers
wanted, limiting the widespread utility of MCLMC. For HMC, this problem is resolved by the use of a Metropolis-Hastins (MH) step (Metropolis et al.) 1953 Hastings 1970), which accepts or rejects proposed moves x + x' of the Markov chain according to an acceptance min (1, e W (e' 2)). where
W((a p(x') ((>aa²) p(x)) q( (x)
But for MCLMC, the corresponding quantity W has not been previously derived. Moreover, an adaptation scheme for choosing the hyperparameters of the algorithm in the MH adjusted case is needed make the algorithm usable out of the box, so that can serve the same use cases as NUTS.
Contributions In this paper, we derive the acceptance probabilities for microcanonical dynamics with and without Langevin noise in Section and Section 4 respectively. Notably, W turns out be the energy eπror induced by discretization of the dynamics, as in HMC. We term the resulting sampler the Nietropolis-Adjuster Microcanonical Sampler (MAMS), and develop an automatic adaptation scheme (Section ()) to make MAMS applicable without having to specify hyperparameters manually. We test MAMS on standard benchmarks in Section and find that outperforms the state-of-the fff-hheari HMC with NUTS tuning by factor of two at worst, and seven at best. The algorithm is implemented in blackjax (Cabezas et al., 2024),, applicable
and is publicly available, together with documentation and tutorials] The code for reproducing numerical experiments is also available"
2 Related work
A wide variety of gradient- samplers have been proposed, including Metropolis Adjusted
Langevin trajectories (Riou-Durand and Vogrine 2022), generalized HMC (Horowitz) 1991a la Neal
2020), the Metropolis Adjusted Langevin Algorithm (MALA; Grenander and Miller ( (1994)) ). De-
terministic Langevin Monte Carlo (Grumitt et al.) (022), Nose-Hoover (Evans and Holian 1985
Leimkuhler and R
Riemannian HMC (Girolami and Calderhead] 2011). Magnetic HMC
(Tripuraneni et al. (017)) and the Barker proposal (Livingstone and Zanella [022]. Some of these
methods come with automatic tuning schemes that make them black-box, for example MALT (Riou-
Durand and Vogrine 2022 Kiou-Durande et al. (023). generalized HMC (Hoffman et al., 2021) and
HMC (Sountsov and Hottman [022] et al. (202T), but these schemes are designed for the
many-short-chaan MCMC regime (Sountsov et al.) [222]] Margossian al. 2024), which we do not consider here [ To our knowledge, NUTS remains the of art black-box method for selecting the trajectory length in HMC-like algorithms for the single chain regime.
The dynamics described by Equation 3 have been independently proposed several times. In computational chemistry, they were derived by constraining Hamiltonian dynamics to have fixed velocity norm (Tuckerman et al., 2001] MMnary et al., (003) and termed isokinetic dynamics. More recently, Steeg and Galstyan| (202)) proposed them as time-rescaling of Hamiltonian dynamics
htttps//bbaackjaax devs hitrrssggttuub.. http://gggttuubb.oooppiittii.. www.aaapppinnooooooooooooooooole bookaagooooooo/oooiiioaoo aagori tetstthmmmmmnnt html https://github
httpps///giccbbb.oom/euuubarrrrt www.///wwwwaarrr/eeeeoeee saaprr/ssamlrrrarrrrr hewpr//emaammllebbbbbaaarrrs benchmarks
$ Maan--sort-chhins approach is to run multiple short chains in parallel instead of single long chain. This
regime is interesting when parallel resources are available. However, it is often not applicable, either because one only has a single CPU, or because the parallel resources are needed elsewhere, for example, for parallelizing the model (Gattringer and Lang 2010 for performing multiple sampling tasks (Robuik al. 2024).

---

with non-standard kinetic energy and no momentum resampling. Robnik et al. (2024) observed that while HMC aims to reach a stationary distribution known statistical mechanies as the canonical distribution, it is also possible to target what is known as the microcanonical distribution, i.e. the delta function at some level set of the energy. The Hamiltonian H must then be chosen carefully to ensure that the position marginal of the microcanonical distribution is the desired target and one such choice is the Hamiltonian from Steeg and Galstyan (202I). They propose adding velocity resampling every steps or Langevin noise every step method to obtain ergodicity.
In all of these instances, MCLMC has been proposed without Metropolis- Hastings, and such, has been proposed as a biased sampler. This is the shortcoming that the present work resolves. We refer to our sampler as using microcanonical dynamics reference to previous work on MCLMC, In contrast, we will refer to the dynamics in standard HMC as canonical dynamics.
3 Technical Preliminaries
Hamiltonian Monte Carlo Let L be the negative log likelihood of p up a constant, i.e. p(æ)
e- fe If gradients (L are available and are sufficiently smooth, HMC is the gold standard proposal distribution. In HMC, each parameter I; has an associated velocity u, Parameters and their velocities evolve by set of differential equations
e-[(x) /Z, where 7 =
i V p(r)),
which are designed to have Pнмс(a, u) pl x ((() their stationary distribution. Here /N is the
standard normal distribution. Note that the marginal distribution f PPмM (x, u)du is equal p(x) the distribution we want to sample from. Thus, sampling from p(x) reduces to solving Equation ( In practice, the dynamics has be simulated numerically, by iteratively solving for x at fixed u and vice versa, and updating the variables by time step at each iteration. The discrepancy between this approximation and the true dynamics causes the stationary distribution to differ from the target distribution, but this can be corrected by MH; that is, we can use discretized Hamiltonian dynamics as a proposal q. Furthermore, to attain ergodicity, the velocities u must be resampled after every n steps.
The resulting algorithm has two hyperparameters: the discretization step size of the dynamics e and the trajectory length between each resampling L = ne. Choosing good values for these two hyperparameters is crucial (Beskos al., 2013 Neal 201T. Betancourt) 2018). and practical sampler must also provide a robust adaptation scheme for choosing them.
An alternative to HMC is microcanonical dynamics (Robnik et al. 2024 Tuckerman et al., 2001 Minary et al. 2003 Steeg and Galstyan 2021) defined by:
Microcanonical dynamics
ē= ==
u = uu T V log p(x)/(d 1)
(3)
where u has unit norm which is preserved by the dynamics. The proposed benefit is that the normalization of the velocity makes the dynamics more stable large gradients. When integrated exactly, these dynamics have PMCLMC (æt u) p(x) Mga-a (u) as a stationary distribution (see Appendix B. so that the marginal is still p(z)). Here Uga-1 is the uniform distribution the d - sphere. Robnik et al. (2024) propose using these dynamics without MH in order to approximately sample from p(æ)). In this case, the velocity partially resampled after every step, and the step size of the discretized dynamics is chosen small enough to limit deviation from the target distribution to acceptable levels.
While this algorithm works well in practice when the step size is properly tuned, the numerical integration error is not corrected, resulting in asymptotic bias which is hard to control. In
HMC this is solved by the MH step, which requires calculating W, as defined in Equation (I). In this case, W can be easily derived since the integrator is symplectic (volume preserving) and q(>]<<>) (are==')//r //((=///() |w) = 1. The integrator used for microcanonical dynamics is not symplectic, so it not immediately clear how to calculate W.
4 Metropolis adjustment for canonical and microcanonical dynamics
Both canonical and microcanonical dynamics can be numerically solved by separately solving the differential equation for the parameters æ, at fixed velocities u and vice versa. For a time interval e,
3

---

we refer to the position update A, u) and the velocity update as (x, u). The solution of the combined dynamics at time ne is then constructed by a composition of these updates:
4 T 1//n Pi/n
4 8/n-
(4)
+. = Be/2 A. Be/2-
(5)
This is known as the leapfrog (or velocity Verlet) scheme. A final time reversal map T(x, u) (x. -u), is inserted to ensure the map an involution, i.c. yp o id, where id the identity map. This is useful in the Metropolis step but does not affect the dynamics in any way, because a full velocity refreshment is performed after the Metropolis step, erasing the effect of time reversal.
=
Both HMC and MCHMC possess a quantity, which we refer to as energy, which conserved for exact dynamics, but only approximately conserved by the discrete update from Equation (5). The energy H is composed of two parts, a potential energy V and kinetic energy K. The position updates (ar u) Aa(ar, u) change the potential energy by
V(Ac(z)) - (z)) = - -log p(x') p(x))
(6)
while the velocity updates (x, u') = (xx u) change the kinetic energy by
1
((B.(z)) - ((z)
||e0112 ||uu|
(7)
2
for HMC and by
((B.(=)) - K(x) (d-)))=dd_)lo)))))) (d-1)loggcohhe 5 + e uish 6}
(8)
for MAMS. Here ||-|| is the Euclidean norm, e = - VL(()/ V[(x)] and 6 = [VL()))] /(d
1).
To derive the MH ratio, the key is to realize that the A and B updates are deterministic <<>'l=) 0(((z)) z'),
(9)
where the transition map is generated by a dynamical system (Fang et al. 2014) for = (x, u).
(() F(=(()).
(10)
Here, z(0) = z and 2(7) = (0z) = z'. The drift vector field F in canonical and microcanonical
dynamics can be read from Equations (2) and 33. respectively. For the former, it equals FA FA(æ, u)
(u, 0) during the A updates, and I'8( (z, u. = (0, VL(x)) during the B updates. For the latter, it
equals FA u) = (u, 0) during the A updates, and Fg (2, u) (0, (1 uu VV (x)) (d) (d/d/dd11) 1))
during the B updates. These fields can be used to explicitly solve for the A and B updates; the solutions are given in Appendix(B.3
Lemma 4.1. For proposals which are deterministic involutions generated by dynamical system of the form 10p. W in the MH acceptance probability [T] equals
p(z') T p(z)
L F(z(s)))
(z'', = loog
Proof. The first term comes from the first factor in Equation the ratio of transition probabilities is
I For the second term, observe that
q( (=11²') 5(x(2z) - z) 8(sp(z') - x) 0g (z) (( 2 2) ((p(2)) - z²) 8(2 - ((²')) Oz
de Oz (z))
where in the second step we have used reversibility, as well as standard properties of the delta function This last expression is the Jacobian determinant of the transition map y.. Finally, the second term of W in Lemma 4. follows from Eq. (1T) by the Abel-Jacobi-Liouvile identity. O
"Note that the MAMS dynamics of Equation 3 are Hamiltonian, for MAMS, K is not kinetic energy the standard sense. In Appendix (B.5 relationship between MAMS dynamics and Hamiltonian dynamics for which K actually is kinetic energy is given, justifying the name.
'Recall that 8(x ax) 6(x (a). and ((f())) E. 8(r a.) aa where the roots of
In case, f(z) d(z)= z'. that 6(0(2) 2) = 6(2 d(2')) De (d(z') 5(2 ((2')) (2))) de (z)
de

---

Lemma 4.2. For proposal of the form T for all k, W in both HMC and MAMS equals the total energy change of the proposal, that is, the sum of all the energy changes:
To Bex Ann o... An, where eN R and E R
((p(z),z) = Σ ((Ar((Z)1)) - V(zk-1) + ((B,, Aca (zk-1)) K [(Ae (zk-1))
where 2 Be, Anp .... Be, o Am (z) and energy changes are taken from Equations 6 and for HMC, and (6) and (8) for MAMS.
Proof. The position update A in both HMC and MAMS has a vanishing divergence: V FΛ
0u. 0, so during the position update, only the first term in Lemma 4.1 survives and W (z'z)
8m:
- -log u //p((((( (2, u) - lolloop')
AV from Equation 66
The velocity update in HMC has vanishing divergence V FB
_a 0, so during
Ou:
the HMC velocity update, only the first term of W in Lemma = ÁK from Equation The velocity update in MAMS, on the other hand, has non-zero divergence:
4.1
survives and W z)
-
llogp(x' u')d //((x, u) = - log
sinh ő + cosh 8(e u) cosh 8 + sinh 8(e u))
FB = [V(a)) e [[((()]]
d = = - dt
log[cosh + sinh d(e - =)).
In the first equality, we have used the divergence from Robnik and Seljak (2024), which we also re-derive in Appendix B.7 In the second equality, used the explicit form of the velocity update from Appendix B.3 namely Equation (28). So we find that the velocity update for MAMS has
B.7
V FB(z(s))ds = (d 1) log[ 8 + t--ushn sinh 8) AK W W (x), = /.
from Equation interested reader.
A more direct derivation of the above is provided in Appendix B.6 for the
W of a composition of maps of the form in Lemma 4.1 is a sum of the individual W, i.e., x(=))
W()p((())) =) W(p(x(=)), x(=)) ((x()), 2). This yields the formula in the statement
of the theorem.
This result shows how to perform Metropolis adjustment in MAMS: analogously to HMC, for every or A update that the algorithm performs in the proposal, the energy difference must be calculated and added to the cumulative energy difference W (B update energy difference is calculated by Equation 7 A update by Equa- tion (6). Final MH acceptance probability is min(1 e This is favorable result, be- cause both the HMC and MAMS numerical in- tegrators keep the energy error small, even over long trajectories (Leimkuhler and Matthews. 2015), implying that a high acceptance rate can be maintained. As an empirical illustration that the MH acceptance probability from Lemma 4.2 correct, Section 4 shows a histogram of 2 mil- lion samples from -i0-dimensional Gaussian (Ist dimension shown) using the MH-adjusted kernel (orange), given step size of 20. The ker- without MH adjustment is also shown (blue) and exhibits asymptotic bias.
Kerne without Mil best Kermll with Mi tea: Thus dedritbation
0.-
0.4
0.2
0.1
08
Figure 1: Histogram of MAMS samples with (or- ange) and without MH adjustment (blue) from 100- dim standard normal (Ist dim shown). Step size is chosen very large (e = 20) to highlight the bias the MH step removes.
5

---

5 Sensitivity to hyperparameters
The performance of HMC is known to be very sensitive to the choice of the trajectory length, and the problem becomes even more pronounced for ill-conditioned targets, where different directions may require different trajectory lengths for optimal performance (Neal 2011). This is further illustrated in Appendix D Two solutions to this problem are randomizing the trajectory length (Bou-Rabee and Sanz-Serna 2017) and replacing the full velocity refreshment with partial refreshment after every step, also known as the underdamped Langevin Monte Carlo (Horowitz) 1991b). Here, we will pursue both approaches with respect to MAMS.
Random integration length We randomize the integration length by taking T&g: 22 L/e integration steps to construct the k-th MH proposal. Here h2 can either random draws from the uniform distribution ((0, 1) or the k-th element of Halton's sequence, as recommended in Owen y: Hoffman et al. (202T). Other distributions of the trajectory length were also explored in the literature (Sountsov and Hoffman 2022) but with no gain in performance. The factor of two is inserted to make sure that we do L/E steps on average 6
Partial refreshment Partially refreshing the velocity after every step also has the effect of randomiz- ing the time before the velocity coherence is lost, and therefore has similar benefits to randomizing the integration length (Jiang] 2023). However, while the flipping of velocity, needed for the deterministic part of the update to be an involution, is made redundant by a full resampling of velocity, this is not the case for partial refreshment. This results in rejected trajectories backtracking some of the progress that was made in the previously accepted proposals (Riou-Durand and Vogrinc 2022). Skipping the velocity flip is possible, but it results in small bias in the stationary distribution (Akhmatskaya et al. 2009), and it is not clear that it has any advantages over full refreshment. Two popular solutions for LMC are to either use non-reversible MH acceptance probability in Neal (2020); Hoffman and Sountsov (2022) or to add a full velocity refreshment before the MH step as in MALT (Riou-Durany and Vogrine 2022). We will prove that both can be straightforwardly used with microcanonical dynamics, and then concentrate on the MALT strategy in the remainder of the paper.
We will generate the Langevin dynamics by the "OBABO" scheme (Leimkuhler and Matthews 2015). where BAB is the deterministic yPe map from Equation (5) and O is the partial velocity refreshment. In LMC, (w) = ciu + oz where is the standard normal distributed variable, C1 = e -/ and /1- - c3. Lpartial is parameter that controls the partial refreshment's strength and is comparable with HMC's s trajectory length L. With microcanonical dynamics, a similar expression that additionally normalizes the velocity has been proposed (Robnik and Seljak 2024):
c + 02 2/vd
Oc(u)= ciu + o2 2/Va
(11)
Denote by A(z', z) the energy error accumulated in the deterministic (sFe) part of the update. Note that for microcanonical Langevin dynamics, only the deterministic part of the update changes the energy, while in canonical Langevin dynamics, the O update also does but is not included in A.
Theorem 5.1. The Metstpolis-Hasttngs acceptance probability of the MAMS proposal q(>" z)) corresponding to TOBABO is min(1, -(('')))).
The generalized HMC strategy (Hoffman and Sountsov 2022) only uses the one-step proposal, so Lemma 4.2 shows that it can be generalized to the microcanonical update, simply by using the microcanonical energy instead of the canonical energy. The MALT proposal, on the other hand, consists of n LMC (or in our case, microcanonical LMC) steps and a full refreshment of the velocity. as shown in Algorithm|] 0
Theorem 5.2. {a:)iio defined in Alg is Markov chain whose stationary distribution is p(z).
Proofs of both theorems are in Appendix A
'More precisely, So [2uL/e] du # L/c, because of the ceiling function. In the implementation, we use the correct expression, which is n = [2yh L/e]] where y = -71
+1++11 and Y [ 1] is the integer part
of y. This follows from solving L/c B[nx] 1424... Y (u- x1(+1) (Y+ Y/2) for y.

---

Input:
negative log-density function C: _2 - R. initial condition æ0 Ed, number of samples N > 0, step size e > O, steps per sample L/e E N, partial refreshment parameter pprrial). The last three parameters can be determined automatically as in Section 0
Returns: samples (x") Jn=1 from p(z)) e e
for to
u ~ Mga-1
= - (a', u)
8 5+0 0
for + n do (z,z,
x' + (2)2)
+ (z))) 8 + 8+ (22 end for
draw a random uniform variable U ~ -0, 1)
if U e- then
x x++1 z= ===1200 [0]
else
x [+1 20
end if
end for
Algorithm 1: MAMS
Langevin
Automatic hyperparameter tuning
MAMS has two hyperparameters, stepsize and the trajectory length I., where L/je is the (average) number of steps in a proposal's 's trajectory. The Langevin version of the algorithm has an additional hyperparameter Lprrtial that determines the partial refreshment strength during the proposal trajecto- ries, i.c., the amount of Langevin noise. In addition, it is common to use a preconditioning matrix M to linearly transform the configuration space, in order to reduce the condition number of the covariance matrix. The algorithm's 's performance crucially depends on these hyperparameters, so we here develop an automatic tuning scheme. First, the stepsize is tuned, then the preconditioning matrix, and finally, the trajectory length. Lpaitall is directly set by the trajectory length, see Appendix [ We adopt schedule similar the one in (Robnik et al., 2024). where each of these three stages takes 10% of the total sampling time, so that tuning does not significantly increase the total sampling cost.
Stepsize We extend the argument from (Beskos et al. 2013 Neal 2011 to MAMS Appendix C showing that the optimal acceptance rate in MAMS is the same as in HMC, namely 65%. For HMC. larger acceptance rates have been observed to perform better in practice (Phan et al., 2019a), with some theoretical justification (Betancourt et al. 2015). We set the acceptance rate to 90%. In the
first stage of tuning, we use a stochastie optimization scheme, dual averaging (Nesterov 2009) from (Hoffman and Gelman 2014). to adapt the step size until desired acceptance rate achieved.
Preconditioning matrix In the second stage, we determine the preconditioning matrix. A simple choice of diagonal preconditioning matrix is obtained by estimating variance along each parameter.
Trajectory length Microcanonical and canonical dynamics are extremely efficient in exploring the configuration space, while staying on the typical set. Therefore, do not wish to reduce them to a comparatively inefficient diffusion process by adding too much noise, i.c., having too low L. On the other hand, we want to prevent the dynamics from being caught in cycles or quasi-cycles to maintain efficient exploration. Heuristically, we should send the dynamics a new direction at the time scale that the dynamics needs move different part of the configuration space, producing new effective sample (Robnik et al., 2024). This suggests two approaches for tuning L (Robnik et al. 2024].
The simpler is to estimate the size of the typical set by computing the average of the eigenvalues of the covariance matrix, which is equal to the mean of the variances in each dimension (Robnik et al.)
7

---

ESS #Gradient evaluaaiions
10
(Langestt) Hyperparameters gnt scarch
10
10
2x
10a Condition Number
Condition Number
10-1
ESS #Gradfent evailationn
10a
10
10->
101
10a
10a
10P
10a
10a
10'
Condition Number
Condition Number
Figure 2: Tuning performance Gaussians function of condition number. Gaussians are 100d with eigenvalues log- eig-unffomm distributed (top row) and outlier distributed (bottom). Left panels: the value of L from the automatic tuning algorithm is shown with a solid line, and the optimal L obtained by grid search (dashed), As can be seen, automatic tuning achieves close to optimal values. Right panels: ESS (for the worst parameter) shown tuning scheme (solid lines) achieves ESS, which is very close to the grid search results (dashed lines).
function of condition number. MAMS
MAMS scaling is similar to NUTS and goes as ESS x condition number 11/2 which is shown as grey lines in the background. However, MAMS has around four times better proportionality constant.
2024). With linearly preconditioned target, these variances and the estimate for the optimal L -- vd. We will this initial value. A more refined approach is to set be the same scale as the time passed between effective samples:
LALBA x E[time between effective samples] E [time between samples] Timt
L Tint (12)
The proportionality constant is of order one and will be determined numerically, based on Gaussian targets. Integrated autocorrelation time Tine the ratio between the total number of (correlated) samples in the chain and the number of effectively uncorrelated samples. It depends on the observable f(x) that are interested in and can be calculated as
Tint = 1
pelfl.
(13)
where
p B((((x(s)) E(JD((( +1) B[JD]
(14)
Var[f]
and harmonically average Tit [r:] over i. We determine the proportionality constant of Equation (12) in way that L equals the optimal L, determined by a grid search, for the standard Gaussian. We find proportionality constant of 0.3 for MAMS without Langevin noise and 0.23 with Langevin noise.
the chain autocorrelation function in stationarity. We take f(x,)
Experiments
We aim to compare MAMS with NUTS, which is the main competitive black-box sampler in the single chain setting, see Section 2 We also compare against MALA (i.c., HMC with one step per trajectory), which is also popular in some settings.

---

gradient calls Di < 0.01
Gaussian
Rosenbrock
10
NUTS
10
6 x 10
x 104
3 10
10a
10a
10
102
10a dimension
10a
dimension
Figure 3: Sampling performance function of the dimensionality of the problem. On the left, the problem is standard Gaussian, on the right, multiple independent copies of the Rosenbrock function banana- shaped target). Number of gradient calls convergence scales as d (grey lines in the background) for both MAMS and NUTS, but MAMS has a better proportionality constant.
Evaluation metric We follow Hoffman and Sountsov 2022 and define the squared error of the expectation value [x)(x)] as
6 ( = (Beampler f E(A))
(15)
Var[s]
and consider the largest second-moment error across parameters, b max sisd (xf)). because in our problems of interest, there is typically parameter of particular interest that has a significantly higher error than the other parameters (for example, a hierarchical parameter in many Bayesian models). For distributions which are product independent low-dimensional distributions (standard Gaussian, Rosenbrock function, and Cauchy problems), because all parameters should have the same error. For the Cauchy distribution, the second moment E[r?] diverges, so instead consider the expected value of - log p(z), i.c., the entropy of the distribution. 12 can be interpreted as an accuracy equivalent of 100 effective samples (Hoffman and Sountsov 2022
take average instead of the maximum
In typical applications, computing the gradients V log p(x) dominates the total sampling cost, so we take the number of gradient evaluations proxy of wall-clock time. For very simple models, gradient calls might not dominate the cost, and the exact implementation becomes important. Our implementation is efficient; for example, 1000 samples with L = 2 and stepsize = 1 for stochastic volatility take 1.5 seconds with MAMS on a single CPU and 4.8 seconds with NUTS. We do not report these numbers since they irrelevant for the more expensive models that the method meant to be applied to in practice. As (Hoffman and Sountsov] (022), we measure sampler's 's performance the number of gradient calls needed to achieve low error, 12 byax 0.01.
the numerical integration
We note that it is common to report effective sample size (ESS) per gradient evaluation instead, but both carry similar information, since can be interpreted 100/(ESS/ # gradient evaluations). Furthermore, we argue that the former is of primary interest and the latter is only used as a proxy for the former.
Scaling with the condition number and the dimensionality Figure compares MAMS with NUTS on 100-dimensional Gaussians with varying condition number. Two distributions of covariance matrix cigenvalues are tested: uniform in log and outlier distributed. Outlier distributed means that two eigenvalues are while the other eigenvalues are 1. For both samplers, the number of gradients low етrоr scales with the condition number i 3, but MAMS is faster by a factor of around 4. Figure compares MAMS and NUTS scaling with the problem's 's dimensionality. Both have the known d4 scaling law Neal 2011
albeit MAMS has a better proportionality constant.
Benchmarks Table I
compares MAMS with NUTS and MALA on set of benchmark problems, mostly adapted from the Inference Gym (Sountsov et al. 2020). Problems vary in dimensionality (36-2429), are both synthetic and with real data, and include distribution with a very long tail (Cauchy), bimodal distribution (Bimodal), and many Bayesian inference problems. Problem details

---

NUTS
MALA
MAMS
MAMS (Langevin)
MAMS (Grid Search)
Gaussian
19,652
11,010
3,249
3,172
3,121
Banana
95,519
140,524
14,078
14,818
15,288
Bimodal
210,758
> 10
139,418
136,770
123,295
Rosenbrock
161,359
> 10
94,184
103,545
93,782
Cauchy
171,373
824,429
110,404
155,963
87900
Brownian
29,816
597,119
13,528
15,232
14,015
German Credit
88,975
> 106
55,748
49,979
52,265
ItemResp
76,043
249,470
45,371
56,902
45,640
Stoch Vol
843,768
> 100
430,088
510,190
431,957
Funnel
> 10
> 107
2,346,899
1,765,311
1,013,048
Table 1: Number gradients calls needed to get the squared error on the worst second moment below 0.01. Lower is better; number of gradients is roughly proportional to wall clock time.
are in Appendix E.1 We do not include Bayesian neural networks because they lack an established ground truth (Izmailov et al. 2021), would require the use of stochastic gradients, which MAMS is not directly amenable to and because there is evidence that suggests that higher accuracy samples are not necessary for the good performance (Wenzel al., 2020), so unadjusted methods achieve superior performance (Sommer et al. (2024).
For all algorithms, MALA, the only remaining parameter to tnne is step size, which tuned by dual averaging, targeting 80% acceptance rate for NUTS and 57 47. for MALA (Neal 2011). For MAMS, we further tune using the scheme of Section 6 We take the adaptation steps as our burn-in, initializing the chain with the final state returned by the adaptation procedure. NUTS is run using the BlackJax (Cabezas et al.) 2024) implementation, with the provided window adaptation scheme. Table T shows the number of gradient calls the chain (excluding tuning) used to reach squared error of 0.01. To reduce variance in these results, we run at least 128 chains for each problem and take the median of the eтror across chains each step.
an initial run find a diagonal preconditioning matrix. For NUTS and
Results In all cases, MAMS outperforms NUTS, typically by a factor of two at worst and seven at best. MALA drastically underperforms relative to MAMS and NUTS. Using Langevin noise instead of the trajectory length randomization has little effect most cases, analogous to the situation for HMC (Jiang 2023 Riou-Durand et al. 2023). For Neal's Funnel, we need an acceptance rate of 0.99 for MAMS to converge. We were unable to obtain convergence for NUTS. This problem a known NUTS failure mode, so it is of note that MAMS converges.
To assess how successful the tuning scheme from Section 6 is finding the optimal value of L. we perform a grid search over L., by first performing a long NUTS run to obtain a diagonal of the covariance matrix and an initial I., and then for each new candidate value of L,, tuning step size by dual averaging with target acceptance rate of 0.9. In Table we show the number of gradients to low error using this optimal Performance is very close to optimal on benchmark problems.
8 Conclusions
Our core contribution is MAMS, out-of-the-box gradient-based sampler applicable in the same settings as NUTS HMC and intended as a successor to it. Our experiments found substantial performance gains for MAMS over NUTS in terms of statistical efficiency. These experiments were on problems varying in dimension (up to 10t) and included real datasets, multimodality, and long tails. This said, to reach the maturity of NUTS, the method needs to be battle-tested over many years on an even broader variety of problems.
We note that MAMS is simple to implement, with little code change from standard HMC. A promising future direction is the many-short-chains MCMC regime (Hoffman and Sountsov 2022 Margossian et al. 2024), since MAMS is not as control-fior heavy as NUTS and since MAMS with Langevin noise can have fixed number of steps per trajectory, making parallelization more efficient (Sountsoy al. 2024).
10

---

Acknowledgments and Disclosure of Funding
This material based upon work supported in part by the Grising-Simonn Foundation grant 2021- 3282 and by NSF CSSI grant award number 231 1559.
References
Akhmatskaya, E., N. Bou-Rabee, and S. Reich (2009, April). A comparison of generalized hybrid Monte Carlo methods with and without momentum flip. Journal of Computational Physics 228(6), 2256-2265.
Beskos, A., N. Pillai, G. Roberts, J.-M. Sanz-Serna, and A. Stuart (2013, November), Optimal tuning of the hybrid Monte Carlo algorithm. Bernoulli 19(5A), 1501-1534. for Mathematical Statistics and Probability.
Publisher: Bernoulli Society
Betancourt, M. (2018, July). arXiv:1701.024:: 02434 [stat].
A Conceptual Introduction to Hamiltonian Monte
Carlo.
Betancourt, M. J., S. Bymne and M. Girolami (2015, February). Optimizing The Integrator Step Size for Hamiltonian Monte Carlo. arXiv:141 [stat].
Bingham, E., J. P. Chen, M. Jankowiak, F. Obermeyer, N. Pradhan, T. Karaletsos, R. Singh, P. Szerlip, P. Horsfall, and N. D. Goodman (2019). Pyro: Deep universal probabilistic programming. Journal of machine learning research 20(28), 1-6.
Bou-Rabee, N. and J. M. Sanz-Serna (2017, August). Randomized Hamiltonian Monte Carlo. The Annals of Applied Probability 27(4). arXiv:1511.099: [math].
Cabezas, A., A. Corenflos, J. Lao, and R. Louf (2024). Blackjax: Composable bayesian inference in jax.
Campagne, J.-E., F. Lanusse, J. Zuntz, A. Boucaud, S. Casas, M. Karamanis, D. Kirkby, D. Lanzieri, Y. Li, and A. Peel (2023, April). JAX-COSMO: An End-to-End Differentiable and GPU Ac- celerated Cosmology Library. The Open Journal of Astrophysics 6, 10.21 23302//sttoo.22.05511 arXiv:2302.05166 [astro-ph].
Carpenter, B., A. Gelman, M. D. Hoffmann, D. Lee, B. Goodrich, M. Betancourt, M. Brubaker, J. Guo, P. Li, and A. Riddell (2017, January). Stan: A Probabilistic Programming Language. Journal of Statistical Software 76, 1-32.
Creutz, M. (1988, August). Global Monte Carlo algorithms for many-fermion systems. Physical Review D 38(4), 1238.. Publisher: American Physical Society.
Crooks, G. E. (1999, September). Entropy production fluctuation theorem and the nonequilibrium work relation for free energy differences. Physical Review E 60(3), 2721-2726. Publisher: American Physical Society.
Duane, S., A. D. Kennedy, J. Pendleton, and D. Roweth (1987, September). Hybrid Monte Carlo. Physies Letters B 195(2), 216-222.
Evans, D. J. and L. Holian (1985, October). The Nose- thermostat. Journal of Chemical Physics 83, 4069-4074. Publisher: AIP ADS Bibcode: 1985JChPh.88.4066
Evans, D. J. and D. J. Searles (1994, August). Equilibrium microstates which generate second law
violating steady states. Physical Review E 50(2), 1645-1648. Publisher: American Physical Society.
Evans, D. J. and D. Searles (2002, November). The Fluctuation Theorem. Advances in Physics 51(7), 1585. Publisher: Taylor & Francis _eprint: https://doi.or//10 0000000001830011115550
Fang, Y., J.-M. Sanz-Serna, and D. Skeel (2014, May). Compressible Generalized Hybrid Monte Carlo. The Journal of Chemical Physies 140(17), 174108. arXiv:1402.7107 [physies, stat].

---

Gattringer, C. and C. B. Lang (2010). Quantum Chromodynamics on the Lattice: An Introductory Presentation, Volume 788 of Lecture Notes in Physics. Berlin, Heidelberg: Springer Berlin Heidelberg.
Gelman, A. and D. B. Rubin (1996). Markov chain monte carlo methods in biostatistics. Statistical methods in medical research S(4), 339-355.
Girolami, M. and B. Calderhead (2011, March). Riemann Manifold Langevin and Hamiltonian Monte Carlo Methods. Journal of the Royal Statistical Society Series B: Statistical Methodology 73(2), 123-214.
Grenander, U. and M. I. Miller (1994, November). Representations of Knowledge in Complex the Royal Statistical Society Series B: Statistical Methodology 56(4), 549-581.
Systems. Journal the
Griewank, A. and A. Walther (2008, January). Evaluating Derivatives: Principles and Techniques of Algorithmic Differentiation, Second Edition. Society for Industrial and Applied Mathematics. Edition: Second.
Grumitt, R. D. P., B. Dai, and U. Seljak (2022, October). Deterministic Langevin Monte Carlo with Normalizing Flows for Bayesian Inference. arXiv:2205. 14240 [cond-mat, physicx:physies, stat].
Hastings, W. K. (1970, April). Monte Carlo sampling methods using Markov chains and their applications. Biometrika 57(1), 97-109.
Hoffman, M., A. Radul, and P. Sountsov (2021, March). An ACCP Scheme for Setting Trajectory Lengths in Hamiltonian Monte Carlo. In Proceedings of The 24th International Conference on Artificial Intelligence and Statisties, pp. 3907-3915. PMLR. ISSN: 2640-3498.
Hoffman, M. D. and A. Gelman (2014). The No-U-Turn Sampler: Adaptively Setting Path Lengths in Hamiltonian Monte Carlo. Journal of Machine Learning Research 15(47), 1593-1623.
Hoffman, M. D. and P. Sountsov (2022, May). Tuning-Free Generalized Hamiltonian Monte Carlo. In Proceedings of The 25th International Conference on Artificial Intelligence and Statistics, pp- 7799-7813. PMLR. ISSN: 2640-3498.
Horowitz, A. M. (1991a, October). generalized guided Monte Carlo algorithm. Physics Letters B 268(2), 247-252.
Horowitz, A. M. (1991b). A generalized guided monte carlo algorithm. Physics Letters B 268(2). 247-252.
Izmailov, P., S. Vikram, M. D. Hoffman, and A. G. G. Wilson (2021, July). What Are Bayesian Neural Network Posteriors Really Like? In Proceedings of the 38th International Conference on Machine Learning, pp. 4629-4640. PMLR. ISSN: 2640-3498.
Janke, W. (2008). Monte carlo methods in classical statistical physics. In Computational many- particle physics, pp. 79-140. Springer.
Jiang. Q. (2023, January). On the dissipation of ideal Hamiltonian Monte Carlo sampler. Stat 12(1), e629.
Leimkuhler, and C. Matthews (2015). Molecular Dynamies: With Deterministic and Stochastic Numerical Methods, Volume 39 of Interdisciplinary Applied Mathematics. Cham: Springer International Publishing.
Leimkuhler, B. and S. Reich (2004). Simulating hamiltonian dynamics. Number 14. Cambridge university press.
Leimkuhler, B. and S. Reich (2009, July). Metropolis adjusted Nosé-Hoover thermostat. ESAIM: Mathematical Modelling and Numerical Analysis 43(4), 743-755.
Livingstone, S. and G. Zanella (2022, April). The Barker Proposal: Combining Robustness and Efficiency Baadient-Based MCMC. Journal of the Royal Statistical Society Series Statistical Methodology 84(2), 496-523.
12

---

Margossian, C. C., M. D. Hoffman, P. Sountsov, L. Riou-Durand, A. Vehtari, and A. Gelman (2024, January). Nested RA: Assessing the Convergence of Markov Chain Monte Carlo When Running Many Short Chains. Bayesian Analysis -I(-1), 1-28. Publisher: International Society for Bayesian Analysis.
Metropolis, N., A. W. Rosenbluth, M. N. Rosenbluth, A. H. Teller, and E. Teller (1953, June). Equa- tion of State Calculations by Fast Computing Machines. The Journal of Chemical Physics (1(6), 1087-1092.
Minary, P., G. Martyna, and M. E. Tuckerman (2003, February). Algorithms and novel applications based on the isokinetic ensemble. II. Ab initio molecular dynamics. The Journal of Chemical Physics 1 118(6), 2527-2538.
Neal, R. M. (2011, May). MCMC using Hamiltonian dynamics. arXiv:1206.190) [physics, stat].
Neal, R. M. (2020, January). Non-reversibly updating a uniform [0.1] value for Metropolis ac- cept/reject decisions. arXiv:2001.1195) [stat].
Nesterov, Y. (2009, August). Primal-dual subgradient methods for convex problems. Mathematical Programming 120(1), 221-259.
Owen, A. B. (2017). A randomized halton algorithm in r. arXiv preprint arXiv:
Phan, D., N. Pradhan, and M. Jankowiak (2019a, December). Composable Effects for Flexible and Accelerated Probabilistic Programming NumPyro. arXiv:1912.1155: [stat].
Phan. D., N. Pradhan, and M. Jankowiak (2019b). Composable effects for flexible and accelerated probabilistic programming in numpyro. In Proceedings of the 2nd Symposium on Advances in Approximate Bayesian Inference. arXiv preprint 1tXiv:1912.115:5
Boou-Dd L., Sountsov, Vogrine, C. Margossian, and S. Power (2023, April). Adaptive Tuning for Metropolis Adjusted Langevin Trajectories. In Proceedings of The 26th International Conference on Artificial Intelligence and Statistics, pp. 8102-8116. PMLR. ISSN: 2640-3498.
Riou-Durand, L. and Vogrine (2022). Metropolis adjusted langevin trajectories: a robust alternative to hamiltonian monte carlo. arXiv preprint 2202..22200.2233
Robnik. J.. A. E. Bayer, M. Charisi, Z. Haiman, A. Lin, and U. Seljak (2024, October). Periodicity significance testing with null-signal templates: reassessment PIF's SMBH binary candidates. Monthly Notices of the Royal Astronomical Society 534(2), 1609-1620.
Robnik, J., R. Cohn-Gordon, and U. Seljak (2024). Black-box unadjusted Hamiltonian Monte Carlo. arXiv preprint arXiv:2412.0887
Robnik, J., G. B. De Luca, E. Silverstein, and U. Seljak (2024, March). Microcanonical Hamiltonian Monte Carlo. The Journal of Machine Learning Research 24(1), 311:14666-3111:14222
Robnik, J. and U. Seljak (2024, July). Fluctuation without dissipation: Microcanonical Langevin Monte Carlo. In Proceedings of the 6th Symposium on Advances in Approximate Bayesian Inference, pp. 111-126. PMLR. ISSN: 2640-3498.
Sevick, E. M., R. Prabhakar, R. Williams, and D. J. Searles (2008, May). Fluctuation Theorems. Annual Review of Physical Chemistry 59(1), 603-633, arXiv:0709.3885
[cond-mat].
Simon- Onfroy, H., F. Lanusse, and A. d. Mattia (2025, April). Benchmarking field-level cosmological inference from galaxy redshift surveys. arXiv:2504.20130 [astro-ph].
Skeel, R. D. (2009). What makes molecular dynamics work? SIAM Journal on Scientific Comput- ing 32),), 1363-1378.
Sommer, E., J. Robnik, G. Nozadze, U. Seljak, and D. Rügamer (2024, October). Microcanonical Langevin Ensembles: Advancing the Sampling of Bayesian Neural Networks.
Sountsov, P., C. Carroll, and M. Hoffman (2024, November). Running Markov Chain Monte Carlo on Modern Hardware and Software. arXiv:2411.04220 [stat].
13

---

Sountsov, P. and M. D. Hoffman (2022, May). Focusing on Difficult Directions for Learning HMC Trajectory Lengths. 11.15:2100.11576 [stat].
Sountsov, P., A. Radul, and contributors (2020). Inference gym.
Steeg, G. and A. Galstyan (2021, December). Hamiltonian Dynamics with Non-Newtonian Momentum for Rapid Sampling. arXiv:2111.022 1.24444 [physics].
Tripuraneni, N., M. Rowland, Z. Ghahramani, and R. Turner (2017, July). Magnetic Hamiltonian Monte Carlo. Proceedings of the 34th International Conference on Machine Learning, pp- 3453-3461. PMLR. ISSN: 2640-3498.
Tuckerman, M. E. (2023). Statistical mechanics: theory and molecular simulation. Oxford university press.
Tuckerman, M. E., Y. Liu, G. Ciccotti, and G. J. Martyna (2001), Non-hamiltonian molecular dynamics: Generalizing hamiltonian phase space principles to non-hamiltonian systems. The Journal of Chemical Physics 115(4), 1678-1702.
Wenzel, F., K. Roth, B. S. Veeling, Swiętkowski, L. Tran, S. Mandt, Snoek, Salimans, R. Jenatton, and S. Nowozin (2020, July). How Good is the Bayes Posterior in Deep Neural Networks Really? arXiv:2002.02405 [stat].
Strumbelj, E., A. Bouchard-Côté, J. Corander, A. Gelman, H. Rue, L. Murray, H. Pesonen, M. Plum- mer, and A. Vehtari (2024, February). Past, Present and Future of Software for Bayesian Inference. Statistical Science 39(1), 46-61. Publisher: Institute of Mathematical Statistics.
14

---

A Metropolis adjusted Microcanonical Langevin dynamics proofs
Denote by o(=1|=) the density corresponding the update and by q( 4(=112),
|2), the density correspond- ing to the single step proposal TOBABO. We will use shorthand notation for the time reversal: = Tz)) and denote by A z) the energy enror accumulated in the deterministic part of the update.
A.1 Proof of Theorem 5.1
Proof. For the MH ratio we will need
[
(16)
<<'12) ol
where we have used the delta function to evaluate the integral over Z'. We further simplify the numerator
(z Z dZ = / o(21 Z Z'Z / 0(21{(Z)
= - Oz2 of (Z dd = / olz)))((((())) (2) (x(ZZ)
az dz / 00 (0 of Z1z²
app(2) 2
dZ.
az
In the first step we have used that o y = o( y) and that time reversal an involution. In the second step, we have performed a change of variables from Z to Z (for which, the Jacobian determinant of the transformation is 1). In the third step used that sp Z o (Z). In the fourth step we change variables to o - instead of Z. In the last step we use that o x)) = o(x) y).
Since o only connects states with the same æ there is only one Z which makes the integral nonvanish- ing and we get
q =
0G((Z)
(=
az
as if there were no updates. The O updates also preserve the target density, so we see that the acceptance probability is only concerned with the BAB part of the update. In this case, the desired acceptance probability was already derived in Lemma
A.2 Proof of Theorem
Proof. Following a similar structure of the proof as in (Riou-Durand and Vogrinc 2022 we will work on the space of trajectories 20:1. = (zo..... 21) E M We will define kernel Q on the space of trajectories, with q as marginal-zo kernel. We will prove that Q is reversible with respect the extended density
P(zu.)) = II (2
and use it to show that g is reversible with respect to the marginal p(xr))
We define the Gibbs update, corresponding to the conditional distribution P(-|xo):
(( 5(2b) u0)Uss- II [(<(1=/_1)
The Gibbs kernel G is reversible with respect to P by construction. Built upon deterministic proposal of the backward trajectory
Zo.Z (Z,,, En),
(17)
introduce Metropolis update:
M (= PM 8(=b.1 20.L) + (1 PMH)6(z(61 20.1),
15

---

laxd @ -la-ixd-11 Both of these matrices have determinant =1,, the determinant of their where РмH (Zo.)) min(1, (1, e (- For 17 > 0, the distribution P admits density with respect
Lebesgue's measure. Therefore
P(Zo.)) dEo.
A(zo..)
(18)
P(zo.L)
0zo.L.
ensures that the Metropolis kernel M is reversible with respect to P.
Before proceeding with the proof, we express Equation (18) in simple, easy-to-compute form.
The Jacobian is 0.a. a @ 0z* where is the matrix of the permutation o(i) = L --2 and
Ozoi
Kronecker product also +11 and its absolute value is 1.
We get
P(Zo.))
p(žL) (Z,-1 z1)
II p(z,) I 9(z,-1
L(zoL)
(19)
P
p(z0)
I x1-1
Il p(z,-1) II g(z (21-1)
II
L
( (z))(zz))
Elen ((
((=; ( [=1
where A(z;,2i-1) is the energy error in step i, by Theorem 5.1 We are now position define the trajectory-space kernel:
Q GMG.
(20)
The palindromic structure of Q ensures reversibility with respect to P. Since the transition (-122o.)) (-((aa) only depends on the starting position zp E R and p(x) is the marginal II dz; defines marginally a Markov kernel
of P, we obtain that q<<///))) = on Rd, reversible with respect to p. In particular, the distribution {2:320 in Algorithm coincides with the distribution of Markov chain generated by q.
B Microcanonical dynamics
In this appendix, we establish relationship between the microcanonical dynamics of Equation and Hamiltonian system with energy from which it can derived by a time-rescaling operation. As well motivating the dynamics of Equation (3). this allows to show that W in Lemma for for the dynamics Equation (3) corresponds to the change in energy of the Hamiltonian system. also provide complete derivation of the form of W for microcanonical dynamics. Familiarity with the basics of Hamiltonian mechanics is assumed throughout.
B.1 Sundman transformation
We begin by introducing transformation Hamiltonian system known Sundman transform (Leimkuhler and Reich 2004)
S(F) = ((z(tt)F(((())
where w is any function 22 + R. Intuitively, this a z-dependent time rescaling of the dynamics. Therefore it surprising that:
Lemma B.I. The integral curves of S(F) are the same as of F (Skeel 2009
Proof. To this, first use 2G to refer the dynamics from field G, and posit that =s(F)())
=f(((s)), where di(a) = w(((s)). Then we see that
(3) dz dzz(1(s)) dt ds dt ds
dz dt F(<=(()))((z(()))
16

---

which shows that, indeed, (s(F) 2 8, where is function IR + R, which amounts to what we out to show.
However, note that the stationary distribution is not necessarily preserved, on account of the phase space dependence of the -inerescaliig, which means that in volume of phase space, different particles will move at different velocities.
B.2 Obtaining the dynamics of Equation
Consider the Hamiltonian system given by H = T V, with ( = (d - 1) log [I| and V(x) C(x). Here, II is the canonical momentum associated to position Then the dynamics derived from Hamilton's 's equations motion are:
d dt =
0H 8TT
(d
= F(z).
(21)
OH
-
Any Hamiltonian dynamics has p(=) x 8(H - C) as a stationary distribution, which can be sampled from by integrating the equations if ergodicity holds. As observed in (Steeg and Galstyan (202I) and Robnik et al. (024), the closely related Hamiltonian dlogg ||I|| + +(x) has the property that the marginal of this stationary distribution is the desired target, namely p(æ) xx -((x) However, numerical integration of these equations is unstable due to the | factor, and moreover, MH
adjustment is not possible since numerical integration induces error H, which would result in proposals always being rejected, due to the delta function.
Both problems can be addressed with Sundman transform and subsequent change of variables. To that end, we choose w z) [ (d-1) 1) (which corresponds, up to factor, to the weight r in Steeg and Galstyan (2021 and to w Robnik et (2024)). we obtain:
d
II v II| (d m]]
(22)
dl
Changing variables to u = /I/ tion
[|| we obtain precisely the microcanonical dynamics of Equa-
x [
BB Bu
dt
-1
wu y ((x))/(d 1)
where we have used that the Jacobian final change of variable only targets II.
IIIT
Note that B2 = S(F) since this
all
[
B.3 Diserete updates
For completeness, we here state the position and velocity updates of the canonical and microcanonical dynamies, which are obtained by solving dynamics fixed velocity for the position update and fixed position for the velocity update. For canonical dynamics, this amounts to solving
d
d B, - VV(x(()),
A. =(())
(23)
de
de
with initial condition Ag = (t) and Bo = u(t). These solution is trivial:
A. = xt) + eu(t)
B. u - ((xx)t))
(24)
'Here we follow Robnik et al (2024) and Steeg and Galstyan avoid the need for weighting scheme used in those papers.
2021 but our Hamiltonian differs by factor.
17

---

For microcanonical dynamics, one needs to solve
de A. u(()
de Be = -11--(1uuuu (l-=atuutt"")V(a(tt))(l))
(25)
with initial condition Ap = (t) and Bn = u(t). The veloeity equation is vector version of the Riccati equation (Steeg and Galstyan 2021 Denote = - ([(((t))(( - and replace the variable Be by ye such that
B. g ye
(26)
This is convenient, because the equation for B, is a nonlinear first-order differential equation, but equation for ye is linear second-order differential equation
d²
de U. (9g (g.
(27)
which is easy to solve and yields the updates
Be u(t) + (sinh 8 + e (coss 6 1))e
A, = (1() + (t)
(28)
cosh 5 zu() sinh
where 8= d /(d and e VV(==/
B.4 Obtaining the stationary distribution of Equation
We derive the stationary distribution of Equation
3 following the approach of Tuckerman (2023).
There, it is shown that for a flow F, if there is a g such that & log g -V P. and A is the conserved quantity under the dynamics, then p(z) x (A(2), where f is any function.
We note that V F gogg- = - ((x). Further, |u] is preserved by the dynamics if we initualize with woo = 1, as can easily be seen: d
v(()) के ((x), using Appendix B.7 in the first step. Therefore
2t (u-u) 2u-u- = 2u-(1- uI- (-V(((())(dd)))))( = u- 1
= 0.
Thus stationary distribution is:
p(x, u) x e [(x) [[[u|] - 1).
(29)
Importantly, because even the discretized dynamics are norm preserving, the condition 5(/uz - 1)
always satisfied, so that p(z)
p(=') is always well defined. This makes it possible to perform MH
adjustment, contrast to the original Hamiltonian dynamics as discussed in Appendix B.2
B.5 W as energy change
In the non-equilibrium physics literature, W (termed the dissipation function) is interpreted as work done on the system and the second term in Lemma 4.1 is the dissipated heat (Evans and Searles 1994 2002 Sevick al. 2008). W plays central role fluctuation theorems, for example, Crook's relation (Crooks) 1999) states that the transitions z + z' are more probable than z' + by factor e W(x'' = In statistics, this fact is used by the MH algorithm to obtain reversibility, or detailed balance, a sufficient condition for convergence to the target distribution.
Here we will justify why it can also be interpreted as an energy change in microcanonical dynamics.
Lemma B.2. W, calculated for the microcanonical dynamics over time interval [0, T1 is equal ΛE of the Hamiltonian (d (d-1) 1) log [I| + L(x) for an interval [s(0), s(T'), where is the time rescaling arising from the Sundman transformation w(z) = [I/l/(d 1).
Proof. Recall that for a flow field F:
W(zr (7)) = log p(zr(T)) (T
((z)(0)) (" V F.F(z((s))ds
(30)
18

---

Given the form of the stationary distribution induced by B, derived in Appendix B.4 we see that the first term of the work, log P()=(0))
(T))(T)) = ((=x(T)) - (())(0)) (( (T) ((
((xx(s(TT)) - ((z)(s(0))) which is equal to AV for an interval of time [s(0), (T)). As for the second term, observe that dK (I{()(()))) ds
0H
dIl
dt
dx dt
dII
u (-V((z))=-VV =
On
dt
ds
dt ds
di
which is preciscly the integrand of the second term.
This shows that W = AK+ = AE, where AE is the energy change of the original Hamiltonian, over the rescaled time interval s(T)) As we know, AE for the exact Hamiltonian flow, and indeed W O for the exact dynamics of Equation (3), which is to say that for the exact dynamics, no MH correction would be needed for an asymptotically unbiased sampler.
However, our practical interest is in the discretized dynamics arising from a Velocity Verlet numerical integrator. In this case, we wish to calculate W for Bu and Bx separately, and consider the sum, noting that W is an additive quantity with respect to the concatenation of two dynamics. Considering W with respect to only Ba. we see that the first term of W remains AV, since the stationary distribution gives uniform weight to all values of u of unit norm, and the dynamics are norm preserving. The second term vanishes, because Ba
As for Bu, since the norm preserving change in u leaves the density unchanged. the first term of W vanishes. Meanwhile, the second term A K. from the above derivation, since V B V₂ B₂ + Vn Bu = V, Bu- Thus, the full W is equal to AV + AK = A as desired. For HMC, it is easily seen that W for Fu is AV, and for Fu AT. Putting this together, we maintain the result of Lemma B.2 but now in setting where W is not 0 so that MH adjustment is of use.
V
B.6 Direct calculation of velocity update W
We here provide a self-contained derivation of the MH ratio for the velocity update from Equation (2) The MH ratio is a scalar with respect to state space transformations, i.e. it is the same in all coordinate systems. We can therefore select convenient coordinates for its computation. We will choose spherical coordinates in which e is the north pole and
u cos de + sin Of,
(31)
for some unit vector f. orthogonal to e. 1 is then a coordinate on the 5dl11 manifold. The velocity updating map from Equation (28).
sinh 5 + cos ((coh 8 - 1) cosh ö + cos & sinh 8 e sin 0
u +
(32)
cosh 8 + cos 0 sinh 5
sinh 8 + cos 0 cosh
cosh - cos 1 sinh +
cosh 5 + cos 1 sinh 0
(33)
can be expressed in terms of the variable:
sinh 8 + cos 0 cosh 5
sin 0 cosh 8 + cos s
cos H'
sin 0'
(34)
cosh + cos & sinh
The Jacobian of the @ + 1' transformation is di'
d0'
cos 8'
(35)
dø
d cos 0'
d0
cosh 8 + cos 0 sinh 8
and the density ratio is
p(0')
g(0')
sin 8'
d-2
(36)
p(0)
V
9(0) sin Ø
(cosh 5 + cos sinh 8)d-2+
where is the metric determinant on a 1 sphere. Combining the two together yields
W = (d - 1) log cosh + cos 0 sinh 8):
(37)
which is the kinetic energy from Equation
8
19

---

B.7 Direct calculation of the velocity update divergence
For completeness, we here derive the divergence of the microcanonical velocity update flow field F. We will use the divergence theorem, which states that the integral of the divergence of vector field over some volume S2 equals the flux this vector field over the boundary Q2 Here, flux is F n where is the unit vector, normal to the boundary.
We will use the coordinate system defined in Equation shell, centered around the north pole and spanning the 1 range [1), 0 + Au). The boundary of Ω are two spheres in d 2 dimensions with radia sin & and sin(0 Ad). Note that F is normal this boundary and flux is a constant on each shell. It is outflowing on the boundary which is closer to the north pole and inflowing on the other boundary.
31 and pick as the volume S thin spherical
Note that for A - O, we have that V F is a constant on Ω. The divergence theorem in this limit therefore implies
d d8 (ssin 2)
V V(S
()in 3)sin
(38)
where V(Sd-2) is the volume of the unit sphere in d-2 dimensions and we have used that the volume of re-dimensional sphere with radius r is V(S")r", By rearranging we get:
(Lin 0)d-2)
V F
(39)
(sin 0)
We have
F V[((x)
d-1 (1 uu")e,
(40)
[+((a]]
[V[((r)]
VVL(x)
FF
e(l wu
u)
sin D.
(41)
d d-1
d d-1
d
Inserting
FF in Equation
(39 yields
[V(x (d (sin 0) cos 0
V.
VLlle
(42)
d
(sin 0)d-2
C Optimal acceptance rate
The optimal acceptance rate argument for MAMS is analogous to the one in Neal (201)) We will use two general properties of the deterministic MH proposal:
1. The expected value of the MH ratio under the stationary distribution, Ez-ple W(=',=)), is
/ ((x)e--Wer ((''=))()
/ p(z) ( dz ||(z'| ((=1=))((=)
de (2)
d(p(())) 1.
0z
This is the Jarzynski equality. In statistical literature it was used by Neal (001 Creutz (1988) in the special case when is symplectic.
2. In equilibrium,
P(W > Olaccepted) = P(W Olaccepted) = 1
2
by the design of the MH algorithm Neal 2011). Since P(accepted[W 0) = 1, we have that
2 = P(w < 0laccepted) = P(W 0 (accepted) <
P(W < P(accopted)
P(accepted)
P(accepted) = 2P(W < 0).
20

---

Let us approximate the stationary distribution over W as (H, o²), as in Neal 2011 by the Jarzynski equality:
We then have
1 = - p(z)e 2(
(W (w-n>²/2 "20/200" W dw = = o² TM
(43)
V2T02
implying that 0² = 2p. By property (2) we then have
P(accept) = (- /1//H//2))) 24),
(44)
where & is the Gaussian cumulative density function. Denote by Kaccepted the number of accepted proposals needed for new effective sample. This corresponds to moving distance the order of the size of the typical set
K acceppeed Ne x vd.
(45)
since vd is the size of the standard Gaussian's 's typical set. The number of effective samples per gradient call then
Kacoepted oc <P(accept)
P(accept)
ESS / # gradient evaluations
(46)
Kiotal
Vd
The error of the MCHMC Velocity Verlet integrator for an interval of fixed length is (Robnik et al. 2024
0²/d x Cd
(47)
implying that 0² 2p x e e//d. Therefore
ESS / # gradient evaluations x n 1/ 4 p(2))//2)) d 1/4
(48)
so that the efficiency drops d-1/44 ESS is maximal at μ 0.41, corresponding to P(accept) 65%. From Equation (47) then see that the optimal stepsize grows d!/4 instead of the d! that would correspond to the unimpaired efficiency.
Note that this result is different if higher-order integrator is used. For example, when using a fourth order integrator o² /d x /d)d",, the optimal setting is 0.13 and P(accept) = 80%.
Empirically, find that even for second-order integrator, targeting a higher acceptance rate, of 90%, works well in practice; we use this our experiments.
D Optimal rate of partial refreshments
We will here demonstrate the sensitivity of HMC to trajectory length for ill-conditioned Gaussian distributions and show how Langevin dynamics alleviates this issue. We will then obtain the optimal rate of partial refreshments in the limit of large condition number. For MALT, Riou-Durand and Vogrine 2022) derive the ESS of the second moments in the continuous-time limit for the Gaussian
targets N (0,0) o):
1- p² ESS(B,T)
(49)
1+ p2
where
B sin
p -- (coswT ====ST
(50)
T is the trajectory length, B is the LMC damping parameter (y 233 in RRou--Duaam and Vogrinc (2022).
Figure shows ESS as a function of the trajectory length and the rate of partial refreshments, both for the standard Gaussian and for the Gaussian in the limit of very large condition number. For the standard Gaussian, HMC achieves better performance than MALT, but only if trajectory length is chosen well. For ill-conditioned Gaussian however, HMC drastically underperforms compared to MALT. Optimal MALT hyperpatameter settings for the ill conditioned Gaussians are 3 = 0.567 and T' 1.413. Therefore the optimal ratio of the decoherence time scales of the partial and full refreshment are
ppartial 1/3
= 1.25. T
(51)
(fall
We will use the same setting for Langevin MAMS, so Lpartill / L 1.25.
21

---

0.45
mifrorr: (ESS(o)) Stua.
ESS/Ona SeLL
Onak
T
2.3
0.15
0.1
04
B
Figure 4: Effective sample size in continuous time for MALT LMC on Gaussian targets. x-axis the LMC damping parameter, y-axis the trajectory length. no damping and no Langevin noise), I = is LMC with critical damping 3 = 1/ Left panel: isotropic Gaussian (0, max). Note that HMC achieves the optimal performance properly tuned, the only reason to introduce Langevin noise would be to potentially make the tuning easier. Right panel: extremely ill-conditional Gaussian with all scales (0, Tmoa] ESS along the worst direction is shown. HMC performs poorly as it cannot be tuned to all scales, damping of 30 = 0.57 performs best. Note that these results do not imply MALT having non-zero ESS in the infinite condition number limit: we only study continuous time MALT here.
is the HMC line (3 means there is
E Further experimental details
E.1
Benchmark Inference Models
We detail the inference models used in Section For models adapted from the Inference gym (Sountsov et al. 2020 we give model's inference gym name the parenthesis.
Gaussian is 100-dimensional with condition number 100 and eigenvalues uniformly spaced in log.
Banana (Banana) is two-dimensional, sanana-shaped target.
Bimodal:
A mixture of two Gaussians in 50 dimensions, such that
p(Φ) = (1 - (z0 /) + aN(x)u, o21,,
where 0.25, / (4, 0, ... 0) and = 0.6.
Rosenbrock a banana-shaped target in 36 dimensions. It is 18 copies of the Rosenbrock functions with Q 0Q=00 0. see O.
Cauchy product of 100 ID standard Cauchy distributions.
Brownian Motion (Brown anMot i ionnnknwwS lesMi s ngMi dd eObserva ticns) is a 32-dimensional hierarchical problem, where Brownian motion with unknown innovation noise is fitted the noisy and partially missing data.
Sparse logistic regression (GermanCredii tNumer www.esepppssoooeeeo stt.iiigsccceesssre s sicn) is 51-dimensional Bayesian hierarchical model, where logistic regression used
model the approval of the credit based on the information about the applicant.
Item Response theory (Synt het ic weittemReeppe tetteesssssseehrrrrr is 501-dimensional hier-
archical problem where students' ability is inferred, given the test results.
Stochastic Volatility is 2429-dimensional hierarchical non-Gaussian random walk fit the S&P500 returns data, adapted from numpyro (Phan et al. 2019b).
Neal's funnel (Neal 2011) is funnel shaped target with a hierarchical parameter 21 ~ N/(03)) that controls the variance of the other parameters 2; ~ No0,, e²/ for 2, 3, d. We take d = 20.
22

---

MAMS
MAMS 10)
MAMS (/ 10)
Gaussian
3249
3162
3169
Banana 14,078
13,098
13,469
Rosenbrock 94,184
91,620
94,533
Brownian 13,528
13,759
13,973
German Credit
55,748
59,777
57,118
ItemResp 45,371
43,840
450,40
Table 2: Number of gradients calls needed to get the squared enror on the worst second moment
below 0.01. MAMS is compared with MAMS where the initial step size has been multiplied (third column) or divided (fourth column) by factor 10. This demonstrates that the performance of the automatic hyperparameter adaptation is not sensitive to the initialization.
Brownian Motion
German Credit
101
10a
MAMS
109
10
10-1
10
10-2
10-3
10
0.0
0.2
0.4
0.6
0.8
1.0
0.0
0.5
1.0
1.5
2.0
2.5
Item Response Theory
Stochastic Volatility
104
109
10a
10-4
10-1
10-2
10a2
0.5
gradient calls 10s
# gradient calls 100
Figure 5: Sampling error, measured by b2 function of computational cost, measured by number of gradient calls. MAMS compared with NUTS on several benchmark problems, demonstrating it consistently outperforms NUTS. Gradient calls to b 0.01 dashed line is the number reported in Table
Ground truth expectation values E[x²] and Var[=*] = - E[=31>²] are computed analytically for the III Conditioned Gaussian, by generating exact samples for Banana, Rosenbrock and Neal's funnel and by very long NUTS runs for the other targets.
E.2 Sensitivity to the hyperparameter initialization
MAMS hyperparameters are automatically determined by the algorithm from Section 6 algorithm initialized with some stepsize Cina = 0.2 d and trajectory length Linit Vd. MAMS performance could in principle be sensitive to this initialization. In Table 2 we show that this is not the case.
However, the
E.3 Convergence curves
Figure 5 shows the convergence of the second moments function of gradient calls used for NUTS and MAMS. 12 b computed as Section Sampling error for both samplers steadily decreases, but MAMS is consistently faster.
23

---

Model
Metric
MAMS
NUTS
Banana
bmar
0.42%
1.52%
Bimodal Gaussian
Omar
6.34%
2.07%
Brownian Motion
0mar
0.36%
0.22%
Cauchy
Davy
5.37%
0.02%
German Credit
mmar
0.12%
0.09%
Item Response
b max
1.05%
0.17%
Rosenbrock
ban 0.06%
0.02%
Standard Gaussian
b aan 1.04%
0.13%
Stochastic Volatility
Uman 0.10%
0.09%
Table 3: Relative uncertainty associated with Table
E.4 Uncertainty in Table 1
Table shows the relative uncertainty in the results from Table T Uncertainty is calculated by bootstrap: for a given model, we produce a set chains (usually 128), and calculate the bias bmaz at each step of the chain. We then resample (with replacement) 100 times from this set, and compute our final metric (number of gradients to low bias) 100 times. We take the standard deviation this list of length 100 to obtain estimate of the error. This error, relative to the values in Table I is then reported in percent.
E.5 Computational architecture
The experiments were run on 128 CPU cores, where each core is 2x AMD EPYC 7763 (Milan)
CPU.
24

---

NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
The main claims made in the abstract and introduction accurately reflect the paper's contributions and scope. The paper clearly introduces the Metropolis-Aduused Microcanonical Sampler (MAMS), describing it as novel approach that outperforms the No-U-Turn Sampler (NUTS) across multiple benchmarks. It further claims theoretical derivations for the Metropolis adjustment and automated tuning, both of which are thor- oughly supported by the mathematical derivations and experimental results in the paper. The scope of the claims is aligned with the presented results, and the paper does not overstate its contributions.
Justification:
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No or
NA answer to this question will not perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results
be expected to generalize other settings.
is fine to include aspirational goals as motivation as long is clear that these goals are not attained by the paper.
2. Limitations
Question: Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification: Limitations of MAMS are very similar to that of NUTS as demonstrated in Experiments and discussed in the Conclusions.
Guidelines:
The answer NA means that the paper has limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors encouraged to create separate "Limitations" section in their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g., if the approach was only tested few datasets with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution is low images are taken in low lighting. Or speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness.
While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, worse outcome might be that reviewers discover limitations that aren't acknowledged in the paper. The authors should use their best
25

---

judgment and recognize that individual actions in favor of transparency play an impor- tant role in developing norms that preserve the integrity of the community, Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory assumptions and proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and a complete (and correct) proof?
Answer: [Yes]
The paper provides full set of assumptions and complete, rigorously derived proof's for all theoretical results, including the derivation of the Metropolis- Hastinggs ac- ceptance probability for the microcanonical dynamics (Section 4) and the calculation of the Jacobian and divergence of the velocity update (Appendix). Each theorem and lemma is clearly numbered and referenced in the main text, and the proofs are presented in the appendix, ensuring clarity and completeness. Assumptions for each result are explicitly stated.
Justification:
Guidelines:
The answer NA means that the paper does include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but they appear in the supplemental material, the authors are encouraged to provide short proof sketch to provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should properly referenced.
4. Experimental result reproducibility
Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that it affects the main claims and/or conclusions of the paper (regardless whether the code and data are provided or not)?
Question:
Answer: [Yes]
Justification: The paper provides comprehensive details necessary to reproduce the main experimental results. It clearly outlines the algorithms (including pseudocode in Algorithm []. the adaptation scheme, and the benchmarking methodology. The paper also specifies the hyperparameter settings and evaluation metrics. Additionally, the paper provides details about the benchmark problems, including dimensionality and problem structure (Section and Appendix D The implementation details, such the of BlackJax for NUTS and the code structure for MAMS, are also discussed, ensuring reproducibility.
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, a No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution is specific model and empirical evaluation, it may be necessary to either make possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g., in the case of large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
26

---

While NeurIPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) If the contribution is primarily new algorithm, the paper should make it clear how to reproduce that algorithm. (b) the contribution is primarily a new model architecture, the paper should describe the architecture clearly and fully.
(c) If the contribution is new model (e.g.. a large language model), then there should either way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome to describe the particular way they provide for reproducibility. the case of closed-source models, may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers have some path reproducing or verifying the results.
5. Open access to data and code
Question: Does the paper provide open access to the data and code, with sufficient instruc- tions to faithfully reproduce the main experimental results, as described in supplemental material?
Answer: [Yes]
Justification: MAMS is implemented in open source package BlackJax and comes with tutorial on how use it. The code for reproducing the experiments and implementations of NUTS and MALA are also provided as a GitHub repository.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurlPS code and data submission guidelines (https://nips.cm cc/
for more details.
While we encourage the release of code and data, we understand that this might not be possible, so "No" is an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to to reproduce the results. See the NeurIPS code and data submission guidelines (https
httpps.cc/pbblcc/uudddee
for more details.
The authors should provide instructions on data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if' applicable).
Providing as much information possible in supplemental material (appended to the paper) recommended, but including URLs to data and code permitted.
6. Experimental setting/details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: The paper specifies all the experimental details necessary to understand the results. It clearly describes the datasets used, including both synthetic and real-world bench- marks. Hyperparameters adaptation scheme is described in detail (Section 7 Experiments and Section 6 Adaptation). The choice of evaluation metrics explained, and the setup for the baseline methods, including NUTS with BlackJax, is described. The experimental
27

---

protocol, including the number of chains, initialization, and the method for calculating performance metic, is also detailed. These specifications ensure that the experimental setting is fully transparent.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented the core the paper to level of detail that is necessary to appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, as supplemental material.
7. Experiment statistical significance
Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Question:
Answer: [Yes]
Justification: we use bootstrap method to calculate standard deviation of results, which are fully reported in Appendix
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by error bars, confi- dence intervals, or statistical significance tests, least for the experiments that support the main claims of the paper.
The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.)
The assumptions made should be given (e.g.. Normally distributed errors).
It should be clear whether the error bar is the standard deviation or the standard error of the mean.
It is OK to report 1-sigma error bars, but should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative eтror rates).
If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text.
8. Experiments compute resources
Question: For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?
Answer: [Yes]
Justification: The paper provides information on the runtime of the algorithm on Stochastic volatility of the model comparison to NUTS and states the computational architecture used. That being said, the computational architecture is not really relevant, as the algorithm here is tested on simple models with cheap to compute gradient evaluations, but it is intended to be used on more involved models. Number of gradient evaluations is therefore the relevant metric, not the runtime our experiments, as we explain the Experiments section 7
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs well estimate the total compute.
28

---

The paper should disclose whether the full research project required more compute than the experiments reported the paper (e.g., preliminary or failed experiments that didn't make it into the paper).
9. Code of ethics
Question: Does the research conducted the paper conform, in every respect, with the NeurlPS Code of Ethics https::/neurrpp: http/publcc/Bhiiscuuieeee
Answer: [Yes]
The research conducted in this paper fully conforms with the NeurlPS Code of Ethics. It adheres the ethical principles, including transparency, integrity, and responsible use of computational resources. The paper does not involve any unethical data collection, misuse of computational power, other ethical concerns.
Justification:
Guidelines:
The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.
If the authors answer No, they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g.. there is a special consid- eration due to laws or regulations their jurisdiction).
10. Broader impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [NA]
Justification: The paper is on foundational research and not tied to particular applications.
Guidelines:
The answer NA means that there no societal impact of the work performed.
If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional unintentional) misuse of the technology.
If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification: The paper is foundational research and not tied to particular applications.
Guidelines:
29

---

The answer NA means that the paper poses no such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Are the creators or original owners of assets (c.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Question:
Answer: [Yes]
The paper properly credits the creators and original owners of all assets used, including code, data, and models. The relevant assets are cited with appropriate references, and their licenses and terms of use are respected.
Justification:
Guidelines:
The answer NA means that the paper does not use existing assets.
The authors should cite the original paper that produced the code package dataset.
The authors should state which version of the asset used and, possible, include URL.
The name of the license (e.g., CC-BY 4.0) should be included for each asset.
For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should provided.
assets are released, the license, copyright information, and terms use in the package should be provided. For popular datasets, paperswi phcooe thcode. com com/datasete
has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that re-packaged, both the original license and the license of the derived asset (if has changed) should be provided.
If this information is not available online, the authors are encouraged to reach out to the asset's creators.
13. New assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [Yes]
Justification: MAMS implementation BlackJax comes with a well-written documentation and tutorial on how to use it.
Guidelines:
The answer NA means that the paper does release new assets.
Researchers should communicate the details the dataset/code/moodel as part their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset is used.
At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsourcing and research with human subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well details about compensation (if any)?
30

---

Answer:
[NA]
Justification:
Paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does involve crowdsourcing nor research with human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included in the main paper.
According the NeurIPS Code of Ethies, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional review board (IRB) approvals equivalent for research with human subjects
Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained?
Question:
Answer: [NA]
Justification: Paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does involve crowdsourcing nor research with human subjects.
Depending on the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurlPS Code of Ethics and the guidelines for their institution.
For initial submissions, not include any information that would break anonymity (if applicable), such as the institution conducting the review.
16. Declaration of LLM usage
Question: Does the paper describe the usage of LLMs if is an important, original, or non-standard component of the core methods in this research? Note that if the LLM is used only for writing. editing, or formatting purposes and does not impact the core methodology, scientific rigorousness, originality of the research, declaration is not required.
Answer: [NA]
Justification: The core method development in this research does not involve LLMs as any important, original, or non-standard components.
Guidelines:
The answer NA means that the core method development in this research does not involve LLMs as any important, original, non stanndaddrrd components.
Please refer to our LLM policy for what should or should not be described.
wwtparrrisccc///ooooaaa cc/Conffreenite cc//oo/ferencce//222//11111
31


Reviews:
[{"review_id": "qwUg6YbSDr", "reviewer": "Reviewer_DLZV", "summary": "Microcanonical Langevin Monte Carlo (MCLMC) offers an efficient alternative to Hamiltonian Monte Carlo (HMC) on multiple target distributions. However, like for HMC, the produced samples converge to the correct distribution only for the continuous-time limit. To address this, HMC employs the Metropolis-Hastings test, which removes the bias by guaranteeing that the corresponding Markov chain preserves the target distribution. MCLMC has been missing such a test, and the current paper closes this gap. In particular, it designs the Metropolish-Hastings-Green test [1], taking the Microcanonical Langevin dynamics integrated in time as a proposal. Furthermore, the authors propose practical ways for choosing and tuning hyperparameters. Finally, they empirically study the algorithm's performance for various target distributions.\n\n[1] Green, Peter J. \"Reversible jump Markov chain Monte Carlo computation and Bayesian model determination.\" Biometrika 82, no. 4 (1995): 711-732.", "questions": "I have no questions or suggestions that would affect my evaluation. Personally, I'm wondering if the result of Theorem 5.1 is similar to the framework proposed in [4]. Not the evaluation of the MCLMC Jacobian, but the test of a non-measure-preserving involutive function.\n\n[4] Neklyudov, Kirill, Max Welling, Evgenii Egorov, and Dmitry Vetrov. \"Involutive MCMC: a unifying framework.\" In International Conference on Machine Learning, pp. 7273-7282. PMLR, 2020.", "limitations": "The authors study and discuss limitations of their work throughout the paper.", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 4, "clarity": 4, "significance": 3, "originality": 3, "strengths_and_weaknesses": "The paper presents a complete, novel, and relevant study. Especially, I would like to highlight the following points.\n1. All the results are stated very clearly, without overloaded notation, and with clear contributions.\n2. The motivation behind the proposed test is very clear and is additionally highlighted in Fig. 1.\n3. Besides their theoretical contribution, the authors provide practical guidelines, whose importance is hard to overestimate.\n4. The empirical study is well-motivated, and the choice of metrics, baselines, and benchmarks is very well explained.\nIn general, the paper meets the standards of a NeurIPS publication.\n\nIn my opinion, the paper does not have major weaknesses. However, I would suggest adding the applications and benchmarks that are more interesting for the NeurIPS audience. Indeed, in recent years, the community has transitioned from classical Monte Carlo benchmarks [2] (as considered in the current paper) to applied problems in the natural sciences [3]. I believe that this would significantly strengthen the paper.\n\n[2] Vargas, Francisco, Will Grathwohl, and Arnaud Doucet. \"Denoising diffusion samplers.\" arXiv preprint arXiv:2302.13834 (2023).  \n[3] Akhound-Sadegh, Tara, Jarrid Rector-Brooks, Avishek Joey Bose, Sarthak Mittal, Pablo Lemos, Cheng-Hao Liu, Marcin Sendera et al. \"Iterated denoising energy matching for sampling from Boltzmann densities.\" arXiv preprint arXiv:2402.06121 (2024).", "was_revised": true, "final_justification": "The paper has a solid theoretical contribution and a good empirical study meeting all the criteria for a NeurIPS publication. The rebuttal hasn't changed my evaluation.", "rebuttal": "We thank the reviewer for their thoughtful and encouraging comments. We especially appreciate the recognition of the clarity, practical relevance, and completeness of our contributions.\n\nWe agree that including more applied or domain-specific benchmarks, particularly those relevant to the natural sciences would further strengthen the impact of our work. We see this as a valuable direction for future research and plan to explore it in subsequent work.\n\nRegarding the connection to [4], indeed, **MAMS without Langevin noise** fits within the involutive MCMC framework presented in [4], as it involves a deterministic, involutive map. However, **MAMS with Langevin noise** falls outside this framework, as the Langevin proposal introduces stochasticity, which violates the determinism requirement in [4]'s construction. We have added an appropriate reference to [4]."}, {"review_id": "BytEBx6jYi", "reviewer": "Reviewer_K5dp", "summary": "The paper introduces a new version of mcmc method designed for efficient sampling from high-dimensional probability distributions. It is built on MCLMC by adding a MH correction step, removing the asymptotic bias. In addition , the paper proposes a fully automated hyperparameter tuning scheme for general use. The method is implemented in JAX and will be open-sourced.", "questions": "1. Can you provide results or qualitative comparisons to other gradient-based samplers?\n2. I would like the authors to highlight the novelty of the paper from either empirical or theoretical perspective.\n3. I wonder how sensitive is MAMS to poor step size or trajectory length initialization during tuning?", "limitations": "no. I don't think this work has potential negative societal impact.", "rating": 4, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "This paper presents MAMS, which applies a MH correction to MLMC to produce unbiased samples with efficiency advantages. Rigorous derivation of the acceptance rule is provided and implement a practical automatic tuning scheme. Empirical results show significant improvements over NUTS across a range of benchmarks.\n\nMy main concern is the novelty of the method. The core contribution seems to be an integration of existing ideas. In addition, the experimental comparison is restricted to NUTS.", "was_revised": true, "final_justification": "My primary questions and follow-up questions are addressed by the reviewers. The authors admit that the method requires much more evidence to be successor to NUTs. Given the overall novelty and theoretical contribution, I’m willing to raise my score.", "rebuttal": "We thank the reviewer for their thoughtful comments and for recognizing the technical soundness, clarity, and empirical strength of our work. We’d like to respond to the concern regarding novelty and the scope of comparisons.\n\n### Novelty and Significance\nOur main contribution is the development of a new black-box gradient-based sampler that is **asymptotically unbiased** and consistently outperforms NUTS — currently the most widely used black-box MCMC method (e.g., in Stan [1]). We believe this alone is a very significant result, with a direct impact for practitioners. From a theoretical standpoint, the derivation of the MH acceptance rule for our microcanonical sampler (with and without Langevin noise) is non-trivial and novel. These derivations are necessary to ensure unbiasedness. Additionally, we propose a new automatic hyperparameter tuning scheme — including for trajectory length — which is key to making the sampler practical. One reviewer remarked that \"best practices from the literature were incorporated at every level, including little tricks I wasn't familiar with...\" and highlighted the tuning scheme as a non-obvious, effective contribution. We believe this underscores the care and originality in the algorithm’s design.\n\n### Comparison to Other Methods\nRegarding the scope of baselines, we emphasize that our goal is to present a **general-purpose black-box sampling scheme**. While many gradient-based samplers exist in the literature, **very few include both a sampler and a tuning mechanism** that can be used out of the box — NUTS is the main such baseline. We therefore chose to focus on this meaningful and widely adopted comparator.\nAdding comparisons to methods that require significant hand-tuning would not be fair or informative in this context. That said, if the reviewer has specific methods in mind that meet the same black-box usability standard, we would be happy to discuss them.\n\n### Sensitivity to Initialization\nWe appreciate the question regarding sensitivity to step size and trajectory length initialization. We have now added a robustness experiment to the appendix: we repeated the tuning procedure using initial step sizes and $L$ values that were 10× larger and 10× smaller. The resulting tuned values and performance metrics (Table 1) were nearly unchanged, demonstrating that our scheme is **not sensitive to initialization**.\n\n### Conclusion\nWe hope this clarifies the novelty and significance of our contributions. In summary, the proposed method offers a new black-box unbiased algorithm that outperforms NUTS. Developing it required non-trivial derivation of the Metropolis-Hastings acceptance probability and a combination of tricks to make it black-box. The resulting algorithm is efficient, robust and applicable out of the box, and is intended to serve as a successor to NUTS. We believe these advances justify the paper's relevance and significance, and we respectfully ask the reviewer to consider revising their score in light of this clarification.\n\n[1] Stan: A Probabilistic Programming Language, Bob Carpenter, Andrew Gelman, Matthew D. Hoffman, Daniel Lee, Ben Goodrich, Michael Betancourt, Marcus Brubaker, Jiqiang Guo, Peter Li, Allen Riddell"}, {"review_id": "QNc2SmmxDE", "reviewer": "Reviewer_2xBS", "summary": "This paper designs a new black-box sampler leveraging recent developments with a new type of \"microcanonical\" dynamics for efficient exploration. \nBias is controlled by deriving the Metropolis Hastings acceptance rate for these dynamics.\nAn automatic hyper-parameter selection is designed, the analogue of NUTS for HMC, so that practitioners can use the sampler as a black-box without tuning. \nThe authors make the code available with a full featured sampler comparison suite.", "questions": "## Micro-canonical versus canonical \n\n> We refer to our sampler as using microcanonical dynamics in reference to previous work on MCLMC, In contrast, we will refer to the dynamics in standard HMC as canonical dynamics.\n\nIt's not critical but maybe you could give a little intuition to readers about these names. (E.g. most samplers are motivated by designing physical dynamics which match the target in thermal equilibrium, called the canonical ensemble in physics, while the proposed samplers were motivated by designing an energy preserving dynamics, or microcanonical ensemble in physics, that directly matches the target.) \nThe intro for Eq. 2 and 3 makes them look so similar it's not obvious where the difference comes from. (I like that the notation shows how effectively similar HMC and microcanonical are.) \n\n- The final presentation of the MH weight was not very readable. I guess one proceeds from Lemma 4.2, replacing expressions with Eq. 6 & 8. It goes through some intermediate notation about the steps. I get that the substep notation is important, to know which energy differences are incorporated, but it makes it hard to read. If it's possible, it would be nice to have a summary in the main text which summarizes the step and final weight, without intermediate notation. \n\n- Using NUTS for preconditioning is kind of disappointing, if the goal is to replace NUTS \n\n- Trajectory length. I didn't follow this part very well, though the intuition makes sense (to make trajectories long enough to reduce autocorrelation). Was this idea explored before, or is it similar logic to NUTS? My memory of NUTS is dim - obviously \"No U-Turn\" is the criteria, but I thought they determined that in a more dynamics-focused way, not based on autocorrelation. It just wasn't clear to me where your approach came from. (Obviously it works well, though!) \n\n\n- Footnote 2 \"Note that the dynamics... are not Hamiltonian\" was confusing, I think you just wanted to point out that the kinetic energy for MCHMC is not standard. It's still \"Hamiltonian dynamics\", with a non-standard kinetic energy \n\n- You mentioned involution a few times. Like other potential readers, I was not familiar with the term before I had read this nice paper on the useful properties of involutions in sampling.\nhttp://proceedings.mlr.press/v119/neklyudov20a/neklyudov20a.pdf", "limitations": "yes", "rating": 5, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 4, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strengths:\n\n- NUTS has been very influential because of its solid, hands-free sampling performance. Getting significant improvements with new dynamics while maintaining the automated, black-box usage is a significant advance. \n- Nice derivation of the MH correction for the microcanonical sampling dynamics, which to my knowledge has not appeared before in this literature. \n- I appreciate that the authors did not take for granted some of the developments that applied to HMC (like whether to do partial or full moment refresh, what acceptance rate is ideal), and investigated these aspects. \n- It seems like best practices from the literature were incorporated at every level, including little tricks I wasn't familiar with, like the step number randomization (Neal's original HMC paper mentioned the potential issue, but gave no solution). \n- The automatic hyper-parameter selection scheme, especially the trajectory length part, does not look trivial to come up with. The results comparing to grid search show that it is very effective. \n- It looks like it's being released with a pretty full-featured JAX sampling library, which is a great resource for the field. \n- Solid, well-thought out experiments\n\n\nWeaknesses: \n\nThe paper is technically heavy and presumes some deep background about samplers. I don't see any way around this, but I use the questions section to suggest areas where clarity may be improved.", "was_revised": true, "final_justification": "I maintained my original score to accept. \n- Authors did a good job with follow-ups. The extensive supplementary actually addressed some issues, and they have promised to clarify in text where possible.  \n- I didn't see any convincing arguments in the low score reviews. One reviewer found the paper hard to follow. I didn't find it hard to follow, but that may be because I was more familiar with papers in this specific area.", "rebuttal": "We sincerely thank the reviewer for the thoughtful and encouraging feedback. We address the points below.\n### Microcanonical vs. Canonical — Physical Background\nWe agree that it would be helpful to provide a clearer explanation of the physical intuition behind the equations 2 and 3 in the main text. While we discuss this in Appendix B, we now expand on it in the main text, as suggested.\n### MH weight\nThank you for pointing this out, we understand a clear statement of the final MH weight is valuable, especially if one wants to implement the algorithm. In simple terms, for every B or A update that the algorithm performs in the proposal, the energy difference must be calculated and added to the cumulative energy difference (B update energy difference is calculated by Equation 8, A update by Equation 6). Final MH acceptance probability is min(1, e^-total cumulative energy difference). We have clarified this in the paper and also moved the algorithm’s pseudocode to the main text.\n### Use of NUTS for Preconditioning\nWe agree that relying on NUTS for preconditioning may seem contradictory to our goal. In our experiments, we opted for NUTS-based preconditioning solely to isolate the sampler's ability to sample from its ability to precondition (which, granted, are strongly related). Preconditioning can be viewed as reparametrizing the parameter space, so we wanted to compare NUTS with MAMS on exactly the same problem, i.e., same target, same parametrization. \nIn practice, however, we recommend using MAMS for preconditioning as well which is in fact the default behavior in our released implementation. We have clarified this in the updated manuscript.\n### Trajectory Length Adaptation\nThe trajectory length adaptation used in MAMS is not related to NUTS, but instead follows the approach proposed in [1] for unadjusted MCLMC. We have adapted it here in a natural way to work with the Metropolized variant. We now state this more clearly in the text.\n### Equation 3 and Hamiltonian Dynamics\nTechnically, Equation 3 does not describe Hamiltonian dynamics in the usual sense, as the transformation is not volume-preserving (the Jacobian is not unity). However, it does represent reparametrized Hamiltonian dynamics, where the time has been rescaled such that the velocity has unit norm. This connection allows us to define the change in energy along these trajectories, as discussed in Appendix B. We have revised the footnote to make this distinction more precise.\n### Involutive MCMC Reference\nWe appreciate the pointer to Neklyudov et al. (2020). We now cite this reference in the relevant section and briefly discuss its connection. Specifically, MAMS without Langevin noise is compatible with the involutive MCMC framework, while the stochasticity introduced by Langevin noise moves it outside of the class of deterministic involutive proposals.\n\nOnce again, we thank the reviewer for the detailed and constructive suggestions. We believe these clarifications will improve the paper’s readability and accessibility.\n\n[1] Robnik, J., G. B. De Luca, E. Silverstein, and U. Seljak (2024). Microcanonical Hamiltonian Monte Carlo. The Journal of Machine Learning Research 24(1), 311:14696–311:14729"}, {"review_id": "j4I4TOC9ql", "reviewer": "Reviewer_rovM", "summary": "The paper proposes MAMS, a Metropolis-adjusted micro-canonical sampler that improves sampling efficiency in high dimensions. It corrects the bias in MCLMC while maintaining its efficiency and requires no manual tuning. MAMS also outperforms NUTS by achieving speedup empirically.", "questions": "- In eq (1), is the Markov chain $x\\rightarrow x'$ or $x'\\rightarrow x$?\n- The authors should define $d$ to be the dimension.\n- In eq (4), why do you need to flip the sign of the velocity after $n$-th leapfrog integration?\n- What is the explicit form of potential energy $V$ and kinetic energy $K$ in line 121? Are $V(x)=\\mathcal{L}(x)$ and $K(u)=\\frac{1}{2}||u||^2$? If so, these should be clearly defined.\n- In eq (9), does $\\delta(\\varphi(z)-z')$ mean $\\delta$ multiplied by $\\varphi(z)-z'$ or $\\delta$ evaluated at $\\varphi(z)-z'$? If the latter is true, $\\delta$ in line 125 should be defined as a function of $x$, i.e., $\\delta(x)$.\n- In line 122, why do you also update the velocity to $u'$ by the position map $A_{\\epsilon}$?\n- Why is $\\varphi$ defined differently in Lemma 4.2 compared to Equation (4)?\n- What is the explicit form of $\\delta(z',z)$ in line 207?\n- What is ESS in line 279? Is it defined or referred in this paper?\n- It would be more insightful if the authors provided convergence curves for all samplers involved in the experiment to enable a better comparison.", "limitations": "No. The authors claim that the limitation is mentioned in experiments and discussed in the conclusions, but it remains unclear what the specific limitation is. It would be helpful to present it more explicitly, preferably in a dedicated paragraph or section.", "rating": 4, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 2, "clarity": 1, "significance": 3, "originality": 3, "strengths_and_weaknesses": "**Strength:**\n- This paper studies a novel sampler-MAMS, which incorporates Metropolis-Hasting filter to MCLMC to correct bias and exhibits empirical acceleration over NUTS. It also introduces adaptive step size and trajectory length to enhance practicality and ease of use.\n\n- Extensive numerical experiments demonstrate the advantages of MAMS over NUTS.\n\n**Weakness:**\n- This paper is not well-written with several confusing and unclear notations (see Questions below). One of its main contributions—the proposed algorithm—is deferred to the appendix, which significantly hinders readability and understanding.\n- The theoretical contribution of this paper is marginal: while [1] provides a mixing analysis of NUTS for Gaussian target distributions, this paper lacks analysis of the mixing properties of MAMS.\n\n[1]Bou-Rabee, N., & Oberdörster, S. (2024). Mixing of the No-U-Turn Sampler and the Geometry of Gaussian Concentration. arXiv preprint arXiv:2410.06978.", "was_revised": true, "final_justification": "The authors convince me that the derivation of a valid Metropolis-Hastings rule for microcanonical dynamics is nontrivial. The experiment shows that it consistently outperforms NUTS in multiple settings. Thus I would like to raise my score to 4.", "rebuttal": "We thank the reviewer for their time and thoughtful feedback. We understand that the main concerns relate to the lack of theoretical mixing guarantees and aspects of clarity in notation and presentation. We address both below.\n\n### Theoretical Contribution and Significance\nWe respectfully disagree with the assessment that this work lacks substantial theoretical contribution. A key theoretical novelty of our paper is the **derivation of a valid Metropolis-Hastings acceptance rule for microcanonical dynamics**, enabling the first **asymptotically unbiased sampler** based on MCLMC. This derivation is nontrivial, especially for MAMS with Langevin noise, as a standard procedure needs to be generalized to non-volume preserving dynamics.\nWhile we do not include formal mixing-time bounds (as in [1]), our primary goal is **practical performance**. As NeurIPS emphasizes empirical relevance and impact across ML subfields, we focus on showing that **MAMS consistently outperforms NUTS** — the de facto standard black-box sampler — across benchmark problems. We also study its scaling with the dimensionality and condition number. This addresses a real need in practical sampling and represents an important advance for the applied ML and probabilistic programming communities.\nMixing analysis, especially in the microcanonical setting, is an interesting and challenging direction, and we agree that this is valuable future work. However, we believe that the current results — both theoretical and empirical — already justify publication. We ask the reviewer to consider this framing and the practical significance of our contribution in the NeurIPS context.\n\n### Presentation and notation\nWe appreciate the detailed list of clarity issues and have addressed them point-by-point below. We also note that other reviewers found the presentation strong (e.g., Reviewer DLZV: *\"All the results are stated very clearly, without overloaded notation, and with clear contributions.\"*). Nevertheless, we fully agree that papers must be accessible across subfields and have revised the manuscript accordingly:\n- x -> x’, we have clarified this in the text.\n- Although not explicitly stated, d can be seen to be the dimension in the first line of the introduction, where we state, $x \\in \\mathbb{R}^d$. We have explicitly stated that it is the dimension.\n- Velocity flip is necessary to make the proposal map $\\varphi$ an involution. Indeed its effect is erased by the full refreshment, but we still keep it for consistency. \n- MCLMC dynamics are a time-rescaled Hamiltonian dynamics, so the kinetic energy difference between two states z1 and z2 cannot be defined as K(z2) - K(z1), but is instead a nonlinear function K(z2, z1). This is why we have not directly defined K(z), but only its difference between the states. Note that it is only the difference that is relevant for the acceptance probability. We have now clarified this around Equations 6, 7 and 8. \n- $\\delta$ is a delta function so yes, it can be read as being evaluated at its argument (with all the usual subtleties of the delta “function”). We realize that the proof of Lemma 4.1 with delta functions might not appeal to the mathematically inclined audience, so we have now also included the reference to the other proofs from the literature, for example [1].\n- Thank you for catching this, it was a typo.\n- $\\varphi$ in Equation 4 is defined by a leapfrog update, in Lemma 4.2 it is slightly more general, allowing for other split-type integrators. Thank you for pointing this out, we have clarified it in the text.\n- Indeed this is a very sloppy notation, because $\\Delta$ does not actually depend on z and z’, but rather only on z’’ = O(z). We instead should have defined $\\Delta$(z’’) = K(BAB(z’’)) - K(AB(z’’)) + V(AB(z’’)) - V(B(z’’)) + K(B(z’’)) - K(z’’), we have now fixed this in the text.\n- Thank you for catching this, we have indeed not defined it. We define Effective sample size (ESS) as: number of samples / $\\tau_{int}$, from Equation 13. This is very standard, see for example [2]. \n- We agree this would be a nice illustration, we have added these curves to the appendix.\n\n[1] Neklyudov, Kirill, et al. \"Involutive MCMC: a unifying framework.\" *International Conference on Machine Learning*. PMLR, 2020.\n[2] Stan reference manual, Section 15.4: Effective sample size.\n\n### Limitations section\nFrom a practitioner’s standpoint, a limitation of MAMS compared to NUTS is that it has not been so widely tested. This is understandable, given that we have just proposed MAMS, while NUTS has been state-of-the-art for years. This limitation is highlighted in the discussion."}]

Author rebuttal:

Meta-review: The paper extends microcanonical Langevin Monte Carlo (MCLMC) that was proposed as an alternative to Hamiltonian Monte Carlo (HMC). The main and significant contribution is to add a Metropolis-Hastings correction step in order to ensure that the 
sampler converges to the target distribution. Adding such correction step requires to compute the ratio of the densities, which was not easy,  and the authors have shown how to do this. 

All reviewers agreed that this is a good contribution to the MCMC literature. The paper is well-written and the experiments show
the method can work better than NUTS.    

I do feel though that in the experimental comparison a couple of obvious (fully black box gradient-based samplers) should have been included and I will encourage the authors to add them for the final versions (reviewer K5dp also asked about this). 

The first obvious addition is standard preconditioned MALA: 
$$
q(x'|x) = N(x' | x - \frac{\sigma^2}{2} S \nabla L(x), \sigma^2 S)      
$$
where $S$ should be the same preconditioner used by NUTS and MAMS. Given $S$, $\sigma^2$ is further tuned 
to achieve acceptance rate $57.4$\%. Obviously this is completely black box.  This MALA is like HMC with L=1, and 
you could also add MAMS-MALA by just setting L=1 in MAMS and again tune the step size targeting acceptance rate  around 57.4. 
Both samplers require only one gradient evaluation per MCMC iteration.
Decision: Accept (poster)