Deterministic Langevin Monte Carlo with Normalizing Flows for Bayesian Inference
Authors: Richard D.P. Grumitt
Biwei Dai
Uros Seljak
Venue: NeurIPS 2022 Accept (2022)
Keywords: probabilistic methods
Bayesian Inference
Normalizing Flows

Abstract: We propose a general purpose Bayesian inference algorithm for expensive likelihoods, replacing the stochastic term in the Langevin equation with a deterministic density gradient term. The particle density is evaluated from the current particle positions using a Normalizing Flow (NF), which is differentiable and has good generalization properties in high dimensions. We take advantage of NF preconditioning and NF based Metropolis-Hastings updates for a faster convergence. We show on various examples that the method is competitive against state of the art sampling methods. 

Deterministic Langevin Monte Carlo with Normalizing Flows for Bayesian Inference
Richard D.P. Grumitt Department of Astronomy, Tsinghua University, Beijing 100084, China
Biwei Dai Physics Department, University of California Berkeley, CA 94720, USA
Uroš Seljak
Physics Department, University of California and Lawrence Berkeley National Laboratory Berkeley, CA 94720, USA
Abstract
We propose a general purpose Bayesian inference algorithm for expensive likeli- hoods, replacing the stochastic term in the Langevin equation with deterministic density gradient term. The particle density is evaluated from the current particle positions using Normalizing Flow (NF), which is differentiable and has good generalization properties in high dimensions. We take advantage of NF precondi- tioning and NF based Metropolis-Hasting updates for faster convergence. We show on various examples that the method is competitive against state of the art sampling methods.
Introduction
The task of Bayesian inference is to determine the posterior p(æly) of d parameters æ given the likelihood p(y[a=) of some data y and given the prior p(z), such that the posterior is given by Bayes theorem as p(zly) = p (())((p(((p() While we have access to the joint p(æ, y) (xzx), the normalization p(y) is generally unknown. Bayesian posterior inference can be related to the general task of finding a stationary distribution in an external potential by defining Uxx)= - In (y))z) x In p(x) L(x) + P(x), where L(x) is the negative log likelihood and P(x) the negative log prior. common approach to this task is to use a particle based Langevin equation, which solves stochastic differential equation (SDE) for the particle evolution (Langevin Monte Carlo, LMC). related method Hamiltonian Monte Carlo (HMC) [Duane et al.) 1987), which uses particle positions and velocities to evolve the particles. In these cases the long time equilibrium solution is that of the target distribution, the posterior p(zly). LMC and HMC are two of the most popular gradient based implementations of Markov Chain Monte Carlo (MCMC), but for finite step size they require Hetropolis-Hastingg (MH) adjustment to become unbiased. In this case they have theoretical convergence guarantees the stationary distribution (under mild ergodicity assumptions), but suffer from chain element correlations, which may render the convergence to be very slow.
In this paper propose Deterministic Langevin Monte Carlo (DLMC) with Normalizing Flows (NF). uses a random initialization from the prior, but the subsequent evolution is deterministic, replacing the stochastic velocity term in the Langevin equation with deterministic density gradient term. To avoid the need to solve the Fokker-Planck equation directly we evaluate the density defined by the current particle positions via NF, and then use its gradient update the first order particle Langevin dynamics. NFs are differentiable and have good generalization properties high dimensions.
36th Conference on Neural Information Processing Systems (NcurIPS 2022).

---

Deterministic sampling methods have until now been limited to particle based methods. Two examples are Stein Variational Gradient Descent (SVGD) [Luu and Wang 2016], and interacting particle solutions of Maoutsa et al. 2020), which proposes particle based scores for the gradient of the log-density. However, these particle based methods scale poorly to high dimensions. The main novelty of our contribution is the introduction of NFs to the class of deterministic methods for evaluating stationary distributions: we use NFs for density estimation and gradient evaluation, as well as for MH adjustment to correct for imperfect NF density estimation and accelerate convergence.
Our primary target application to expensive likelihoods (wall clock time of seconds or more). An example are inverse problems in scientific applications, where often one must solve an expensive Ordinary or Partial Differential Equation (ODE or PDE) forward model, where the computational cost can be minutes or hours. In such applications the computational cost of the NF (which is of order seconds in our experiments) is not dominant. For fast model evaluations standard MCMC methods suffice and our target are problems where the cost of standard methods is prohibitive. The computational cost of forward model is unrelated the complexity of the posterior distribution: often the posteriors are simple Gaussian distributions even when the computational cost is high. In such settings the main goal of posterior analysis sis is to minimize the number of likelihood calls to reach some prescribed precision on the posterior.
DLMC has the following advantages in comparison to several other MCMC samplers:
DLMC takes advantage of gradient information of the target distribution, but without stochastic noise, which slows down the particle evolution towards the target LMC and HMC. As result, DLMC gradient based optimization can rapidly move the particles into the region of high posterior mass for fast convergence.
DLMC is not a sequential Markov Chain, and DLMC particles can be updated in parallel at each time step, which can take advantage of trivial machine parallelization of likelihood evaluations, reducing the wall-clock time compared standard MCMC.
DLMC particles initialized random draw from the prior, which covers the posterior. Each particle produces an independent sample in the limit of a large number of particles (N).
DLMC uses NF for density estimation, and takes advantage of NF preconditioning and NF based MH updates, for faster convergence and to help correct for imperfect density estimation.
DLMC can handle multimodal posteriors, in contrast to standard MCMC samplers.
2 Langevin and Fokker-Planck equations
The overdamped Langevin equation is a stochastic differential equation describing particle motion in an external potential and subject to a random force with zero mean,
á = =-== -U<))U(<>++ +
(1)
where (n)))) = 0 and (n = 2ōt5((t (').
1'). For simplicity we set the diffusion coefficient, temperature and mobility to unity. The Langevin equation is first order velocity equation which has a deterministic velocity vÜ(x) and a stochastic velocity n. Here the gradient operator is with respect the parameters of U (x) VU(x)= dU (a) /dx.
In practice we need to discretize the Langevin equation, and for any finite step size At the result biased distribution. An example is the Ornstein-Uhlet Uhlenbeck process, where U(x) is harmonic potential, exp(-U'(x)) (x)) = N(x; μ, E). Standard Langevin algorithm updates lead to the solution ((x) = (H, Σ - ArE 1/2) Wibisono 2018]. We can see that for a finite fixed stepsize At the solution converges to biased answer regardless of the number of time integration steps taken, or the number of samples used. A solution to this is to supplement the Langevin updates with the Haastios (MH) Adjustment Langevin Algorithm (MALA), where MH acceptance guarantees detailed balance [Roberts and Tweedie 1996]. A similar MH adjustment is also required for HMC (Duane al. 1987].
The I Langevin equation can be viewed as a particle implementation of the evolution of the particle probability density q(x(t)), which is governed by the deterministic Fokker-Planck equation,

---

continuity equation for the density,
((x(t)) + V. 0, J= )=- VV = ((
(2)
We defined ((x(()) = - -Ing(=(() and expressed current as density times velocity, where the two terms in the probability current J correspond the two velocity terms in the Langevin equation. One can see from equation] 2 that the stationary distribution (( L) = 0 is given by p(zly) = exp(-U(())/pp()) The Fokker-Planck equation PDE, while the Langevin equation is stochastic (SDE). Langevin diffusion provides particle description of the density g t) and in the large time limit both equations lead stationary distribution equal to the posterior p( (g)).
If we replace the stochastic velocity in the Langevin equation I with the deterministic velocity in equation 2, we obtain the deterministic Langevin equation (Maoutsa al. 2020 Song et al. 2020]]
(t) v - (z)((t)) - (*x(t)))
Equation 3 gives the dynamics of a particle representation of the density (x)()), which in the large time limit + oo converges to the same distribution as the solution equation [ This follows from the fact that the stationary solution of equation 3 (+) = 0 is given by V V(x(()) = (xx(()) + with с being constant independent of æ, which must thus be equal to log ply), since q(xx()) is normalized. There is stochastic noise in Equation 3| so it can be interpreted as deterministic analog of the Langevin equation, which is a first order equation for velocity.
The main difficulty in solving equation 3 is the instantancous density term ((x(t)) needed to obtain (at)). Here we solve this by evolving N particles together, using NFs to evaluate their particle density. While this makes the particles interacting, in the limit of large N the particle evolution is independent. Furthermore, we can use the NF to do MH adjustment
correct for imperfect NF
density estimation (Section 3.2) Upon time discretization of Equation 3 we obtain Algorithm [ A detailed discussion of the algorithm implementation is given in Appendix B. We initialize /= by providing initial samples drawn at random from the known prior distribution p(ær), similar to Sequential Monte Carlo [el Moral et al. [006]. This has several positive features: since the posterior proportional to the prior times the likelihood, the prior always covers the posterior. Moreover, when the likelihood is weakly informative the convergence of DLMC fast. We then apply DLMC gradient updates to obtain the particle position updates. When initializing from the prior the initial density distribution cancels out the prior in the target, and the initial gradient update is the log likelihood term VL(x)
Equation 3 can be interpreted as gradient based minimization of the time dependent objective (()(()) V(a(()). In this view we can optimize the objective using any optimization method. Initially, if we start the particles from the prior p(x)), DLMC is simply performing optimization of the target likelihood p)y) Later, as ((x( ~ p(zly), the gradient of the optimization objective becomes zero, and the particles stop moving. One can use any optimization method to move the particles, so if the gradient of U(x) is not available we can use gradient free optimization methods to do so. For gradient based optimization we can use many of the standard optimization methods [Nocedal and Wright 2006). such first order (stochastic or deterministic) gradient descent, momentum based methods, or (quasi) second order methods (e.g.. L-BFGS, Newton, conjugate gradient etc.). The algorithm can be stopped when the particle distribution becomes stationary. For example, we can stop when the estimated first and second moments of the posterior stop changing.
Given the formulation above, the performance of DLMC depends on the ability of the NF to approximate the particle density and the corresponding gradient at each time step. The addition of MH updates (Section 3.2p help to correct for imperfect density estimation. However, in situations where the NF fit is a poor match for the target the MH acceptance can be very low. Furthermore, in high dimensions (O(1))) we are limited by the training cost to use simple NF approximations. However, in this scenario DLMC can still be used to move particles quickly towards the typical set, acting an initialization for MCMC in the NF latent space. This is discussed further in Section 5.5
3 Normalizing Flow density estimation
Normalizing Flows provide powerful framework for density estimation and sampling (Dinh et al. 2017, Papamakarios et al. 2017 Kingma and Dhariwal] 2018 Dai and Seljak] [02T]. These models map the d-dimensional data to d- d-dimensional latent variables through sequence of invertible
3

---

transformations f f f1 o f2 o... o fL. such that x == f() and z is mapped to a base distribution π (z), which is chosen to standard Normal distribution NN,, (0, I). The probability density of data can be evaluated using the change of variables formula:
det (f((a) dæ
) (f(x)) II det
[0f((z) dz
)
((x) = v(1) =
((f(x)) =((===))
(4)
The Jacobian determinant of each transform Jr det( af(()))
ax: ) must be easy to compute for evaluating the density, and the transformation ft should be easy invert for efficient sampling. In this paper we use Sliced Iterative Normalizing Flow (SINF) [Dai and Seljak] [222]] which has been
shown to achieve better performance for small training data (below a few thousand particles), and considerably faster, than alternatives.
3.1 Normalizing Flow as preconditioner
While the core of DLMC is the NF based density gradient, the same NF can be used for additional convergence acceleration towards the target. The proposed algorithm can be slow when the condition number of the problem high. A typical solution to precondition with a covariance matrix, which mimics the global curvature of the problem. This becomes more difficult if the curvature is position dependent. NFs can address this problem by warping the space in which we perform the optimization. This warping maps spatially varying curvature to the latent space where the curvature is constant and isotropic [Parno and Marzouk] 2018 Hoffman al., 2019]. We can implement DLMC in the NF latent space where the target distribution pl y is transformed into
p( y) p(f [()))))((z).
(5)
where ((z) = det (2) 0z
is the absolute value of the Jacobian determinant. We can define
the potential latent space as
(f) J(z).
(6)
we can evaluate (U)x) using methods such as auto-differentiation and if the NF Jacobian is also differentiable we can also evaluate the gradient in latent space x as
dæ J(z)
VU(z) (U(x)
dz
In latent space the latent variables are distributed as T(2) exp(-V(z)) = N(O, V(z) = 2/2 + (d/2) (n(T) and (V 2. This gives the deterministic Langevin gradient descent updates as
xt + Ar) =t) VU(=(())A + (()Att.
is interesting to compare this deterministic update to the standard stochastic Langevin equation update, which in latent space is
(t At) =x({) - (=(()))++++
where 77'~~ ~ NO, 21A)). Both updates have the potential gradient term, which the deterministic drift velocity. The difference between the two is that DLMC also uses the current position of the sample z for a deterministic update of the density gradient (() which moves the particles radially outwards from the center of the latent space, while the stochastic Langevin algorithm adds Gaussian random variable ry' to the update. We find that performing DLMC updates in the NF latent space accelerates convergence.
3.2 Normalizing Flow for Metropolis- adjustment
The main requirement of DLMC that g(x(()) properly describes the density of samples a=1) at time /. This can only be ensured in the large N limit. Moreover, for non-convex optimization problems possible that the samples get stuck at local minima far from the global minimum. One can add Metropolis Metrsiolls-Hasttnns (MH) adjustment correct for imperfect density estimation,

---

Algorithm Deterministic Langevin Monte Carlo with Normalizing Flows
1: Input: initial samples at 0) drawn at random from prior p(ar), size N, access to L(xr)
P(x:) = (x) and its gradient.
2: Initial update:
3: for = 1 to N do
4: <Att) - x' 0) VC))
5: end for
6: while not converged do
7: t < + At
8: Run Normalizing Flow on current particle positions to obtain g(z(t)) = exp(-V(2(tt) and
map =(1) =000=[(am (x)() for
9:
N do
10:
if latent space update then
11:
Map to latent space: z'()) (t) = f(x
12:
DL update: (t + At) = =t) (VU(2'()) ~
13:
Inverse map: (t + At) = (z
14:
else
15:
DL update: (t + At) = (t) - VIU(>' (t)) ((=
16:
end if
17:
MH update: draw a sample ,' q()(()). replace x' with x' with probability r(x', x²)
18:
end for
19:
end while
20:
Output: set of partieles (t) distributed as target distribution p(æly)-
eliminating samples that get stuck far from the global minimum. This step is particularly necessary for disconnected multi-modal peaks of the target to be properly equilibrated. We note that whilst LMC (i.e., MALA) and HMC utilize local MH adjustment, they typically cannot equilibriate between separate posterior peaks without the use of methods such as annealing.
Specifically, at given / we can draw independent particles a' g(x(()), one for each existing particle æ². We compare the new particle against an existing particle decide whether we accept or reject it. Since these new samples are independent the MH acceptance rate is (e.g.. Albergo et al. 2019))
r(æ, .) min {
p( p(a p(æ; ())()))))))))))))))))))))))
(10)
Note that the normalization constant cancels out this ratio, so we can evaluate this with knowledge of (x) only.
Initially, when g(x(()) is very broad compared to p(æ1z) the acceptance rate of MH will be low, since most of the samples drawn will have low p(xly), so this process may only eliminate the worst performers, and replace them with samples the higher density region of p(æ y). As DLMC progresses q( (()())) becomes closer to p(aly) y)) and the acceptance increases and reaches r r(x, a)) = the NF learns the target perfectly. In high dimensions the MH acceptance can remain low. However, the MH update is still necessary for proper equilibration of multi-modal peaks.
4 Related work
NFs and transport maps for Bayesian posterior analysis have been explored in Parno and Marzouk [2018), Hoffman et al. 2 Arbel et al. [202]] These approaches use NFs as latent space preconditioners, using standard MCMC methods such as HMC, LMC or MH deployed in latent space, where the geometry may be more favorable for fast mixing. In Albergo et al. 2019 NFs are used as an MH transition proposal, giving independent samples. DLMC goes beyond these NF applications in that it also uses NFs for the density gradient in the deterministic Langevin equation.
DLMC is an ensemble method, that all of the particles inform the NF density evaluation, which is in turn used to update each particle position. There have been many other ensemble methods proposed in the literature, such as the Affine Invariant Ensemble Sampler [Goodman and Weare
2010 Differential Evolution MCMC (eer Braak 2006). DLMC differs from these methods in
5

---

W
kecb MI15
kdbb Et DL
CLANC
Tuch
-12
-15
(a) Double Gaussian
(b) Rosenbrock
(c) III conditioned Gaussian
Figure I: Panel (a): results on d = 100 Gaussian mixture after 90 DLMC iterations (blue), compared the true distribution (orange). We show 1d density projections and 2d partiele distributions along the first two parameters. Panel (b): same for Rosenbrock function in d = 32 after 50 DLMC iterations. Panel (c): particle distributions for the d = 1000 ill conditioned Gaussian after 28 DLMC iterations, along with the distribution following 5 subsequent latent space HMC iterations (green). The Gaussian mixture uses 500 particles, the Rosenbrock uses 1000, and the ill conditioned Gaussian uses 2000.
02
-12
-10
20
10
Mila
10
-10
(a) Neal's funnel (a = 00)
(b) σ
(c) 0.1
Figure 2: DLMC results on d = 101 hierarchical variance problem (funnel problem), compared to very long NUTS run. We show posteriors of / and 21. Left panel shows Neal's funnel: corresponds needed since the likelihood update is zero. Middle panel shows σ = 5: despite the large noise, the d = 100 latent variables are informative of a and the posterior is narrower than the prior. Right panel shows 0.1: here the posterior is narrowed further, and is close to a Gaussian.
to noise - co,
updates
that it evaluates the particle density explicitly, and uses the density and target gradient to update the particle evolution.
Two deterministic algorithms related DLMC that have recently been proposed are Stein Variational Gradient Descent (SVGD) [Liu and Wang 2016). which uses direct pair summation over all particles, and score based interacting particles Maoutsa al. 2020). inspired by generative score based models [Song et al. 2020). which use the particles to directly estimate the density gradient term. DLMC differs from these in using NFs instead of direct pair interactions, which we find requires fewer particles to accurately represent the particle density.
2017 that both approximate the DLMC is related Variational Inference (e.g., Blei However, DLMC does not use standard posterior with parametrized normalized function q (x).
variational optimization such as KL minimization to optimize q(x), and does not use (xx) for the posterior, which is reported by the particles instead. In our experiments we observe that the posterior distribution of DLMC particles agrees better with the exact posteriors than the samples from NF density. This is expected, since we use the MH step to correct the NE. DLMC thus contains components of both sampling and optimization: sampling as optimization in the space of distributions has been pursued in many directions since Jordan et al. [998]

---

169
10"
MUTB-luning
NUTS-huning
16a
NUTS
8-OLMC
S-DLMC
10
P-OLMC
P-LLIC
10"
16"
10
10"
10"
10
10
10"
10
10
10
10
10
10"
Wall Time (mins)
Wall Time (mins)
(a) Rosenbrock
(b) funnel
10a
NUTS-Juning
NUTS St--tunnng
NUTS
10"
NUTS
S-OLMC
DLMC
P-OLMC
10
10
10"
10
10
10"
10
10
10
Wall Time (mins)
Gradient Evaluations
(c) funnel 0.1
(d) sparse logistic regression
Figure 3: Bias-squared on the second moment for the d = 32 Rosenbrock function (panel (a)), 101 hierarchical funnel wall clock time (panels (b) and (c)) and for 51 sparse logistic regression versus number of likelihood gradient evaluations (panel (d)). In panels (a)-(c) we show serial (S-DLMC) and parallel (P-DLMC) versions. Solid circles denote the end of any burn-in phases.
MLTE
Figure 4: Results on d = 51 sparse logistic regression after 100 iterations (blue), compared to the distribution obtained from a very long NUTS run (orange). We show three representative variables, showing good agreement. The corresponding bias-squared is 0.01 (figure 4).
Examples
The main DLMC application is for expensive likelihoods. We can simulate this by synthetically creating expensive model evaluations, but will also compare the performance in terms of the number of likelihood evaluations. The implementation of DLMC requires NF and particle evolution. A discussion of implementation details can be found in Appendix B. For NF we SINF, which has [Dai and Seljak] 2021]. is fast, and iterative. The number of layers
very few hyper-parameter can be chosen based on cross-valle validvvai,,,onn where set aside 20% of the samples, and iterate until validation data start to diverge from the training data. However, for the d = 1000 Gaussian (Section 5.5) we fix L 5. Typical SINF training time of order seconds on CPU.
The implementation of particle evolution requires us to specify the learning rate for the particle updates. At each iteration, take Adagrad updates in the (U((tt) V(æ(t))) direction (Duchi et al. 2011]. We use learning rates between 0.001 and 0.1, with smaller learning rates being more robust for targets with complicated geometries such as funnel distributions.
In our plots we focus on the quality of marginal posteriors, because that typically the primary goal of Bayesian inference. To quantify the quality of the posteriors numerically we follow Hoffman et al.
7

---

2019 and show the squared bias of the second moment versus wall clock time or number of gradient evaluations. The precision of this number limited by the effective number of samples (ESS). Here we choose bias-squared of 0.01 as sufficient in terms of posteriors, and corresponds for a Gaussian distribution to an ESS of 200. The number of particles N is a hyperparameter that must exceed this ESS. Overall we find there is tradeoff between number of particles and number of iterations, such that starting with a larger number of particles does not always imply a higher overall computational cost, due to the improved NF density estimation resulting in fewer DLMC iterations.
The main baseline we compare against the No-U-Turn Sampler (NUTS) (Hoffman et al. 2014], an adaptive HMC variant implemented in the NumPyro library (Phan et al. 2019]. It requires tuning and burn-in, which need to be included in the computational budget. NUTS a standard baseline since typically outperforms other samplers such as LMC, SVGD and SMC. Where we include NUTS as a baseline, we use 500 tuning steps. NUTS has the number of leapfrog steps as a tunable parameter, and when it is one HMC becomes LMC/MALA [Girolami and Calderhead 2011].
In some of our experiments we compare default DLMC (denoted serial, S-DLMC) to parallel DLMC (P-DLMC), where for the latter we show the wall clock time based on the assumption that we can evaluate the likelihood gradient in parallel on N CPU cores. We assume likelihood gradient cost of 1 minute, and the cost of the NF itself (seconds) is negligible.
5.1 Gaussian mixture
Multimodal distributions are a failure mode for many standard MCMC samplers as they cannot move from one peak to another when the potential barrier between them is too high. They need to be augmented with an annealing procedure, where one slowly morphs the target distribution from the prior to the posterior. Annealed Importance Sampling [Neal [001] and Sequential MC (SMC) (Del Moral et al. 2006) are examples of such procedures. They require equilibrating the distribution at every temperature level using standard MCMC, and a result significantly slower than standard MCMC. DLMC uses a different strategy from annealing, but can also handle multi-modal distributions, as long as an NF can approximate it. As an example we take Gaussian mixture with 1/3 of the posterior mass in the first peak and 2/3 in the second peak, in d = 100.
We choose N- 500, drawn from a broad initial prior 2) The log of the prior posterior volume ratio is 415, which is a challenge for non- non-gradient based MCMCs, as they become very slow in such settings. We apply logit transforms unconstrained space. DLMC is gradient based and makes rapid progress towards the region of high posterior mass. The resulting DLMC posterior after 90 iterations for the first two variables is shown and is in a near perfect agreement with the true distribution.
all the variables such that we sample in an
in figure [
5.2 Rosenbrock function
The 32-d example comes from a popular variant of the multidimensional Rosenbrock function [Rosenbrock] 1960), which is composed of 16 uncorrelated 2-d bananas. The model log- likelihhood is In L - n/2 [9(x211122 22.)² /Q + (xz--1 1) .Q 0.1, n = 32, with Gaussian prior (0, on all the parameters. We assume g(x)) = x and use an initial burn-in phase with N = 10 for 50 iterations, before upsampling to N 1000. The posterior for the first two variables is shown in figure [ We see that DLMC perfectly matches NUTS. The quality of the posterior versus wall clock time is shown in figure 3 where we are idealizing a data assimilation process by assigning high computational cost of minute to the evaluation of g(x). This allows us to compare serial and parallel DLMC. S-DLMC achieves results that are competitive with NUTS. P-DLMC can improve the wall clock time by N = 1000.
5.3 Hierarchical Bayesian analysis
Hierarchical Bayesian analyses are very common in many applications. They often lead to funnel type posteriors, which are challenging for standard samplers. Moreover, forward models relating the variance of the latent space to some parameters of interest can be expensive. As an example, in cosmology data analysis applications we measure Fourier modes and want to determine their power spectrum, which in turn constrains cosmological parameters. The evaluation of the power spectrum as a function of cosmological parameters requires expensive PDE evolution of structure in the
8

---

universe as function of time [Lewis et al., 2000]. Here we choose simplified version where the power variance P g()) is determined by cosmological parameter and controls the variance of parameters 2,. Parameters 2; are not directly observable, and instead we observe their noisy version with likelihood of the data y. so a combined prior and likelihood are
p(0) = (0, (() = N(0, P = g(0)), p(y 24, 0) =
i i=
(11)
where o² is the measurement noise and g(f) is an expensive function that evaluates the power spectrum P as a function of 8. To connect to Neal's 's funnel we assume the result of the PDE is that the power spectrum depends exponentially on 8, with an approximate relation P g(f) cxp(#). with a minute computational cost of performing g The goal is the posterior of parameter marginalized over Zi. We use N 500.
In the limit o² + oo the setup corresponds to Neal's funnel, which is deemed problematic for many samplers such HMC (Neal et al. 2011]. Note however that the DLMC solution in the o² - oo is trivial: the data likelihood is non-informative, and after drawing samples from the prior of equation [] the likelihood updates of equation I1 are zero, so the DLMC stopping criterion (stationary posterior) is satisfied immediately.
To make the problem closer to actual applications in cosmology we analyze it with a finite noise o². such that the data are informative of the parameter A. We create a simulated data set with the true value Orrue = 0 and use d = 101. The resulting posterior for a 5 and o 0.1 is compared to NUTS in figure 2 We also show the prior (i.e., Neal's funnel) to highlight that the posterior is narrow compared to the prior. For = 5 we used 15,000 likelihood calls, which is orders of magnitude fewer likelihood evaluations than the NUTS baseline, which struggles to converge at all (figure without the help of reparametrization preconditioner (Hoffman et al. 2019].
To highlight the advantage of DLMC with parallelization we further show the performance of (embarrassingly) parellel DLMC, where we evaluate the likelihood update for all 500 particles in parallel, so each iteration is minute. This converges to the required precision b² 10 in 30 minutes. Sequential MCMC methods such as NUTS cannot be so efficiently parallelized because they require a hyperparameter tuning and burn-in phase (each of order 10 likelihood calls). This leads to 3-4 orders of magnitude larger wall clock time. SMC can take advantage of parallelization, but needs lot of temperature levels and a lot of steps to equilibriate at each temperature level (Del Moral et al.) 2006 Wu et al.) 2017). such that the overall wall clock time is still orders of magnitude larger than parallel DLMC.
5.4 Sparse logistic regression
Hierarchical logistic regression with a sparse prior applied to the German credit dataset is a popular benchmark for sampling methods (Dua and Graff 2017]. We use the example from Hoffman et al. [2099] with d = 51, applying a log-transform to bounded variables such that we sample in an unconstrained space. This a challenging example for NUTS, with samples requiring on average 125 leapfrog steps.
We run extended burn-in phase with 10 particles to shrink from the prior to posterior volume, using 500 iterations. We then upsample to 5000 samples, with DLMC converging to the correct solution within 10 iterations. Upsampling to lower sample sizes, such as 1000 particles, also converges to the correct solution, but requires more DLMC iterations such that the computational cost is not reduced. The resulting marginals are shown in figure 4 for three representative variables showing good agreement.
5.5 III conditioned Gaussian
To asses the performance of DLMC in high dimensions we consider an ill conditioned Gaussian target in d 1000. We construct the target covariance Σ by sampling its eigenvalues from Gamma 1), giving a condition number of ~ 4000. The prior and likelihood are then defined as
p(x) (00, Σx 100ºI), p(y (x) = N(x, Σg -E-- (Σー1ーΣ;
(12)
where the prior and likelihood covariance are chosen such that the posterior covariance is Σ.
In this high dimensional setting we are limited by the computational cost of the NE, and the ability of the NF to approximate the current particle density. However, use simple flow with 5 layers,

---

to move particles quickly towards the typical set with DLMC. In this experiment we evolved 2000 particles until satisfied the condition V Σ ii=1 (x (u(x)) S 1.1dd, where the sum is over the target dimensions. This virial condition is chosen because initially we have V > d, with the particle ensemble satisfying V = d at equilibrium [Leimkuhler and Matthews 2015).
For this problem we are unable to achieve low bias with DLMC alone. However, we can treat DLMC as an initialization for HMC i.e., we run DLMC until we reach the virial threshold and then perform parallel HMC with each particle in the NF latent space. DL.MC reaches the virial threshold in 28 steps. Running parallel HMC chains in latent space with 10 leapfrog steps and the step size targeting an acceptance rate of 0.65, we reach a bias-squared of 0.01 after HMC iterations. This corresponds to ~ 10 serial likelihood evaluations, and ~ 100 parallll evaluations. The resulting particle distributions from DLMC and the subsequent latent space HMC are shown in figure | By comparison, NUTS requires 10 serial likelihood evaluations to reach bias- bias-squared of 0.0I. For this example the serial cost of DLMC+HMC and NUTS is comparable, with DLMC+HMC providing performance improvement of ~ 102 when implemented in parallel. In this experiment have used a very simple implementation of HMC, and the convergence rate could likely be improved by using some ensemble adaptation scheme e.g., Hoffman and Sountsov 2022
5.6 Ablation studies and comparison to particle deterministic methods
In Appendix C we compare to kermel density estimation based DLMC and SVGD, both of which give inferior results, confirming that the NF the key ingredient for the success of deterministic methods.
In Appendix D we present results of ablation studies. We show DLMC without MH, DLMC without NF preconditioning, and pure NF MH without DLMC. In all cases the results are inferior to DLMC with preconditioning and with MH adjustment. We also show NF with MH without DL is inferior to DLMC.
Conclusions
We present new general purpose Bayesian inference algorithm we call Deterministic Langevin Monte Carlo with Normalizing Flows. The method uses the deterministic Langevin equation, which requires the density gradient, which we propose be evaluated using a Normalizing Flow determined by the positions of all the particles at given time. Particle positions are then updated in the direction of log target density minus log current density gradient. To correct for imperfect NF density estimation, we add an MH step. The process is repeated until the particle distribution becomes stationary, which happens when the particle density equals the target density.
In our experiments we find that DLMC significantly outperforms baselines such as LMC, HMC and SVGD, and even more in terms of wall-clock time when likelihoods can be evaluated in parallel on a multi-core computer. A possible limitation is that NF training may not learn the true target distribution for small N. While DLMC with MH does offer some robustness against imperfect density estimation, it can fail in situations where the MH acceptance rate is very low, such that it cannot reach the stationary distribution in finite number of steps. This is particularly apparent in O(10²) dimensions, where we are limited by the NF training cost and the fidelity of its gradient approximation. Here, convergence on the correct target can be achieved by running DLMC until the particle ensemble reaches the typical set, at which point we perform short MCMC run in latent space. Indeed, treating DLMC as an initialization for a latent space MCMC in this way can provide some additional robustness through the associated asymptotic guarantees, whilst avoiding the expensive tuning and burn-in process typically required for MCMC. It would nonetheless be interesting to explore what NF architectures, regularization, and training procedures safeguard against low MH acceptance rates.
Acknowledgements
This material is based upon work supported by the U.S. Department of Energy, Office of Science, Office of Advanced Scientific Computing Research under Contract No. DE-AC02-05CH1 1231 at Lawrence Berkeley National Laboratory to enable research for Data-intensive Machine Learning and Analysis. We thank Minas Karamanis, David Nabergoj and James Sullivan for useful discussions.
10

---

References
MS Albergo, G Kanwar, and PE Shanahan. Flow-based generative models for markov chain monte carlo in lattice field theory. Physical Review D, 100())003515,, 2019.
Michael Arbel, Alexander G. de G. Matthews, and Arnaud Doucet. Annealed flow transport monte carlo. In Marina Meila and Tong Zhang. editors, Proceedings of the 38th International Conference on Machine Learning. ICML 2021, 18-24 July 2021, Virtual Event, volume 139 of Proceedings of Machine Learning Research, pages 318-330. PMLR, 2021. URI. wtt.httt::/prrceeeeee .mlr
gtess/v139// arbellia...a.a. html
David M Blei, Alp Kucukelbir, and Jon D McAuliffe. Variational inference: A review for statisticians. Journal of the American statistical Association. 112(518):859 877, 2017.
Biwei Dai and Uros Seljak. Sliced iterative normalizing flows. In Marina Meila and Tong Zhang, editors, Proceedings of the 38th International Conference on Machine Learning, ICML 2021, 18-24 July 2021, Virtual Event, volume 139 of Proceedings of Machine Learning Research, pages 2352-2364. PMLR, 2021. URL htRphht/://proceeeiie
mlr www..ppees//ee Nossss/vv339 (1/33///da211 html
Pierre Del Moral, Arnaud Doucet, and Ajay Jasra. Sequential monte carlo samplers. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 68(3):41 14446,146,,,, 2006.
Laurent Dinh, Jascha Sohl-Dickstein, and Samy Bengio. Density estimation using real NVP. In 5th International Conference on Learning Representutions, ICLR 2017, Toulon, France, April 24-26, 2017, Conference Truck Proceedings. OpenReview.net, 2017. URL
net/
forum?id-HHppbnhi11.
Dheeru Dua and Casey Graff. UCI machine learning repository, 2017. URL http:://rrciiee..ieee uci edu/ml
Simon Duane, Anthony D Kennedy, Brian J Pendleton, and Duncan Roweth. Hybrid monte carlo. Physics letters B. 195(2):216-222, 1987.
John Duchi, Elad Hazan, and Yoram Singer. Adaptive subgradient methods for online learning and stochastic optimization. Journal of Machine Learning Research, 12(61) 2121122121-2159, 2011. URL
http://jmlr org/ wttpg//appsss/v11/ paapers/v12/uuchi1111 htcp//uc.hhhhnnnnnnana
Mark Girolami and Ben Calderhead. Riemann manifold langevin and hamiltonian monte carlo methods. Journal of the Royal Statistical Society: Series B (Statistical Methodology). 73(2):
2011.
Jonathan Goodman and Jonathan Weare. Ensemble samplers with affine invariance. Communications in applied mathematics and computational science, 5(1):65-880 2010.
Matthew Hoffmann, Pavel Sountsov, Joshua Dillon, Ian Langmore, Dustin Tran, and Srinivas Vasudevan. Neutra-liznn lizing bad geometry in hamiltonian monte carlo using neural transport. arXiv preprint w:Xiv:1903.03700 2019.
Matthew D. Hoffman and Pavel Sountsov. Tuning-free generalized hamiltonian monte carlo. In Gustau Camps-Valls Francisco J. R. Ruiz, and Isabel Valera, editors, Proceedings of The 25th International Conference on Artificial Intelligence and Statistics, volume 151 of Proceedings of Machine Learning Research, pages 7799-7813. PMLR, 28-30 Mar 2022. URL https:// proceedings mlr, ppwess/v11/hooffmmeemm.
html
Matthew D Hoffman, Andrew Gelman, al. The no-u-turn sampler: adaptively setting path lengths in hamiltonian monte carlo. J. Mach. Learn. Res., 15(1) (5(1):1593-16))) 2014.
Richard Jordan, David Kinderlehrer, and Felix Otto. The variational formulation of the fokker-planck equation. SIAM journal on mathematical analysis, 29(1):1-17, 1998.
Diederik P. Kingma and Prafulla Dhariwal. Glow: Generative flow with invertible 1x1 convolutions. In Samy Bengio, Hanna M. Wallach, Hugo Larochelle, Kristen Grauman, Nicolò Cesa-Bianchn and Roman Garnett, editors, Advances in Neurul Information Processing Systems 31: Annual Conference on Neural Information Processing Systems 2018, NewrIPS 2018, December 3-8, 2018, Montréal, Canada, pages 10236-10245, 2018.

---

Benedict Leimkuhler and Charles Matthews. Molecular Dynamics: With Deterministic and Stochastic Numerical Methods. Interdisciplinary Applied Mathematics, Springer, May 2015. ISBN 978-3- 319-16374-1. doi: 10.
A. Lewis, A. Challinor, and A. Lasenby. Efficient Computation of Cosmic Microwave Background Anisotropies in Closed doi: 10.1086/309179.
Models. ApJ, 538:473-476, August 2000,
Qiang Liu and Dilin Wang- Stein variational gradient descent: A general purpose bayesian in- ference algorithm. In Daniel D. Lee, Masashi Sugiyama, Ulrike Luxburg, Isabelle Guyon, and Roman Garnett, editors, Advances in Neural Information Processing Systems 29: Annual Conference on Neural Information Processing Systems 2016, December 5-10, 2016. Barcelona, Spain, pages 2370-2378, 2016. URL htattps://ppocceeiinngg
1 2 Absto0aa Abstract...rtt.... html
Dimitra Maoutsa, Sebastian Reich, and Manfred Opper. Interacting Particle Solutions Fokker- Planck Equations Through Gradient-Log--Dens. Estimation. Entropy, 22(8):802, July 2020. doi: 10.3390/e22080800.
Radford M Neal. Annealed importance sampling. Statistics and computing, 11(2):125-139, 2001.
Radford M Neal et al. Meme using hamiltonian dynamics. Handbook of markov chain monte carlo, 2(11):2, 2011.
J. Nocedal and S. J. Wright. Numerical Optimization. Springer, New York, 2nd edition, 2006.
George Papamakarios, lain Murray, and Theo Pavlakou. Masked autoregressive flow for density estimation. In Isabelle Guyon, Ulrkee von Luxburg, Samy Bengio, Hanna M. Wallach, Rob Fergus, S. V. N. Vishwanathan, and Roman Garnett, editors, Advances in Neural Information Processing Systems 30: Annual Conference on Neural Information Processing Systems 2017, December 4-9, 2017, Long Beach, CA, USA, pages
2017.
Matthew D Parno and Youssef Marzouk. Transport map accelerated markov chain monte carlo. SIAM/ASA Journal Uncertainty Quantification, 6(2):645-68) 2018.
Du Phan, Neeraj Pradhan, and Martin Jankowiak. Composable effects for flexible and accelerated probabilistic programming in numpyro. arXiv preprint arXiv:19/2.11554 2019.
Gareth O Roberts and Richard I. Tweedie. Exponential convergence of langevin distributions and their discrete approximations. Bernoulli, pages 341-3414 341-363, 1996.
HoHo Rosenbrock. An automatic method for finding the greatest or least value of a function. The Computer Jourmal, 175) 1960.
Yang Song, Jascha Sohl-Dickstein, Diederik Kingma, Abhishek Kumar, Stefano Ermon, and Ben Poole. Score-based generative modeling through stochastic differential equations. arXiv preprint 1rX1v:2011.13456, 2020.
Cajo JF Ter Braak. A markov chain monte carlo version of the genetic algorithm differential evolution: easy bayesian computing for real parameter spaces. Statistics and Computing, 16(3):239-244 2006.
Andre Wibisono. Sampling optimization in the space of measures: The langevin dynamics composite optimization problem. In Conference on Learning Theory, pages 2093-3027. PMLR, 2018.
Yuhuai Wu, Yuri Burda, Ruslan Salakhutdinov, and Roger Grosse. the quantitative analysis of decoder-based generative models, 2017. httpp::://aa...oo//aaaaaaae///ee
001.////////11/1022003103 04273
12

---

A Checklist
1. For all authors...
(a) Do the main claims made in the abstract and introduction accurately reflect the paper's
contributions and scope? [Yes]
(b) Did you describe the limitations of your work? [Yes]
(c) Did you discuss any potential negative societal impacts of your work? [N/A]
(d) Have you read the ethics review guidelines and ensured that your paper conforms to them? [Yes]
2. If you are including theoretical results...
(a) Did you state the full set of assumptions of all theoretical results? [N/A]
(b) Did you include complete proofs of all theoretical results? [N/A]
3. If you experiments..
(a) Did you include the code, data, and instructions needed to reproduce the main experi- mental results (either the supplemental material or as URL)? [Yes]
(b) Did you specify all the training details (e.g., data splits. hyperparameters, how they were chosen)? [Yes]
(c) Did you report error bars (e.g.. with respect the random seed after running experi- ments multiple times)? [Yes]
(d) Did you include the total amount of compute and the type of resources used (c.g., type of GPUs, internal cluster, or cloud provider)? [Yes]
4. If you are using existing assets (e.g., code, data, models) or curating/releasing new assets...
(a) If your work existing assets, did you cite the creators? [Yes]
(b) Did you mention the license of the assets? [No] German credit data taken from public
UCI repository.
(c) Did you include any assets either in the supplemental material URL? [Yes] Code included as supplemental material.
(d) Did you discuss whether and how consent was obtained from people whose data you're 're using/curating? [N/A]
(e) Did you discuss whether the data you are using/curating contains personally identifiable information or offensive content? [N/A]
5. If you used crowdsourcing or conducted research with human subjects...
(a) Did you include the full text of instructions given to participants and screenshots, applicable? [N/A]
(b) Did you describe any potential participant risks, with links Institutional Review Board (IRB) approvals, if applicable? [N/A]
(c) Did you include the estimated hourly wage paid to participants and the total amount spent on participant compensation? [N/A]
13


Reviews:
[{"review_id": "vPIoQV978EN", "reviewer": "Reviewer_QrcG", "summary": "The work proposes a Bayesian inference procedure constructed using a deterministic Langevin equation. The resulting algorithm uses a number of particles initialized at random, their dynamics then following the deterministic equation toward its equilibrium distribution, while at regular time intervals a normalizing flow is fit to the collection of particles and used as an independence proposal within Metropolis-Hastings to propose replacing each particle with a new draw. The method is demonstrated on a number of toy examples and standard data sets.\n", "questions": "* Before the Metropolis-Hastings step, aren't the particles only approximately distributed according to the target distribution, given time discretization to get them there, and applying the Metropolis-Hastings kernel will not leave them properly distributed, either? This is not the case with HMC or MALA, say, where the time-discretized move is itself the proposal (over the current state of the Markov chain), and the Metropolis-Hastings rule accepts or rejects it. Here it seems there is a particle drawn from the time discretization, and a proposal drawn from the normalizing flow fit to the collection of those particles, and we accept or reject the swap. Perhaps the authors could better explain how this yields samples from the target distribution, or clarify the claims (e.g. first paragraph of conclusion claims \"To make it asymptotically unbiased we add MH step...\").\n\n* Independence proposals rely on the proposal being a good approximation of the target distribution in order to have a reasonable acceptance rate, and so for the Markov chain to mix well. I can accept that a normalizing flow can achieve a good approximation with sufficient size, but I do not see any empirical results on this part of the algorithm, e.g. what is the average acceptance rate of this step in the examples?\n\n* Bottom p.7 there is a comment that, \"...we find there is a tradeoff between number of particles and number of iterations, such that starting with a larger number of particle does not always imply a higher overall computational cost.\" I can see how this might be the case, but can also seem the limits of it, as in e.g. the extreme case of reducing to one particle. A plot showing this tradeoff would be a good addition.\n\n* As above, specific numerical results on the swap moves (e.g. acceptance rate) would be a good addition.\n\n* As above, some clarity on the presentation of the parallel version.\n", "limitations": "No concerns.", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "There are some interesting aspects to this work, such as the use of the deterministic Langevin equation, and the empirical results are quite good, especially the ablation studies in the supplementary material that separate out the various components that come together in the final algorithm. There are some gaps, however, that are stopping me from recommending acceptance.\n\nWhile I think I understand the construction of the algorithm, it is not clear to me that this results in samples distributed according to the target distribution, even up to the approximations introduced (e.g. fitting the normalizing flow to the current stock of particles, time discretization). With that said, the empirical results suggest that the algorithm works as intended numerically; if there is a bias it may be small. I've added some questions for the authors on this topic below and hope they can help me to better understand their work.\n\nIt is unclear to me how the results for P-DLMC have been obtained (e.g. in Figure 3). The authors note that these are based \"...on the assumption that we can evaluate the likelihood gradient in parallel on N CPU cores\" (bottom p.7). Is there actually a parallel implementation here? In my opinion, unless there is a parallel implementation to provide empirical results, this ought to be an item for the discussion, not the experimental section. Otherwise, one is not accounting for the real-world challenges of parallel computing (fewer than N cores, possibly variable execution time between particles leading to load balancing issues, communication overhead in gathering all particles to update the normalizing flow, etc). \n\nIf it is a model rather than actual implementation, have the authors simply divided through by N? Or have they used a model where only some p portion of the sequential program can be parallelized, thus using a scaling factor such as (1 - p + p/N)? There is a comment about the assumption that evaluating the likelihood gradient takes about a minute, while the cost of the normalizing flow is in the seconds, but even that limits the maximum speedup to much less than N (say 60 seconds likelihood gradient parallelized down to essentially zero, 1 second for normalizing flow, that still limits speedup to factor of 61, i.e. 61 seconds down to 1 second, even when N >> 61). These are some of the challenges of scaling these sort of population/ensemble algorithms in practice (SMC has similar challenges).\n", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "l1sD92UFFL", "reviewer": "Reviewer_qZmF", "summary": "The authors propose a mechanism to produce unbiased samples from expensive posteriors given that it is easy to sample from the prior and possible (although expensive) to evaluate the likelihood. The approach proposes concrete particle dynamics based on NFs that approach the target distribution and corrects them using Metropolis-Hastings. ", "questions": "1. It would seem that a lot of the power of this approach is coming from the strong modeling ability of the NF. I am missing comparisons with simpler approaches that can leverage this power. For instance:\n\n(a) What if we simply run variational inference to optimize the normalizing flow? One can easily compute the gradients with respect to the parameters. (For the SINF case, it might be possible to modify its main algorithm to take an unnormalized log-pdf instead of a set of samples - but the gradients for variational inference can always be approximated stochastically in a NF). Then one can find the closest q(z) within the NF domain and then sample from it and correct using MH. This decouples the number of particles used for learning and the ones used for sampling: More independent samples can be generated on demand after training. And does not require any particle dynamics.\n\n(b) What if we only use SINF and MH? MH produces a new set of samples closer to the posterior, then we use those samples and SINF to update q(z) in closed form. Then we repeat. Again, we have eliminated the separate particle dynamics.\n\nWhile (a) and (b) might be significantly worse than the proposal of the paper, it'd be reassuring to know that this is the case, and that all the heavy lifting is not being done by the NF being able to learn the energy function. Otherwise, the proposed dynamics don't seem justified enough.\n\n2. The experiments in which the posteriors are calculated are pretty simple. No strong baselines for those models are used, just NUTS which of course is slow.\n\n3. The tested models are _very_ low-dimensional. The high dimensional setting might have very poor performance, this is not known.\n\n4. I expect multimodal likelihoods with broad priors to result in a lack of mixing - I doubt all modes will be discovered even in low dimensional cases, when they are spaced and sharp. Again this is not tested even in a 2-dimensional experiment with many gaussians that are well separated. \n\n\nMinor comments:\n- Parentheses are used for bibliography citations, should be square brackets, not to be confused with cited equations.\n- In the MH update in Algorithm 1, time indices are missing.\n\n[1] Black Box Variational Inference. Rajesh Ranganath et al.\n\n\n", "limitations": "Addressed in the previous section.\n\n\n", "rating": 7, "confidence": 4, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths:\n- Technically correct\n- Clear presentation\n- Original dynamics\n- Uses the gradient information of the target distribution \n- Produces multiple uncorrelated samples\n\nWeaknesses:\nExperimental validation, in terms of:\n- Complexity of the models sampled from\n- Strength of the baselines\n- Ablations", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "KnlnpQk4OQA", "reviewer": "Reviewer_yr5Z", "summary": "The authors of this work propose using the deterministic “dual” of the langevin equation for simulating from a target posterior. This dual presents difficulties since it requires estimating the score of the FPK equations solutions, to do this the authors propose using normalizing flows, furthermore , the authors “precondition” the flows in latent space via simulating what seems to be latent Gaussian dynamics at equilibrium, this “preconditioning” improves/stabilizes the overall procedure.", "questions": "Comments/Suggestions:\n\n1. The deterministic Langevin Equation is well known and used in the diffusion score matching literature [1] (see Eq X) and should be acknowledged appropriately and maybe briefly discussed. Notice the velocity term coincides with the score. \n2. Equation 2 is broken U(x) should be U(x(t)) in order to be consistent, however this is also wrong, the FPK equation is defined pointwise for all x \\in R^d  that is its defined as q(x,t) or q_t(x) , it's not only defined pathwise on x(t), this notation is inaccurate. The f(x(t)) notation should be reserved for the deterministic langevin equation, but the FPK equation is defined pointwise. Across the whole paper the usage of q(x), q(x(t)) and x vs x(t) in general is quite inconsistent.\n3. From flows setup you would expect p(z|y) to be Gaussian but this is never made clear, instead you refer to \\pi(z)  being Gaussian which was never explicitly defined in the context of p(x|y) ? let's assume that in your notation \\pi(z) = p(z|y) = N(0,1) then indeed the score of q_t(z) is infact -z (i.e. V(z)=z) however this implies that z(t) evolves entirely according to linear Gaussian dynamics , not only that but if the score is indeed -z this implies that q_t(z) = N(z; 0,1) which is at equilibrium at does not depend on time whatsoever. What is the point of simulating gaussian noise at equilibrium in latent space ?  Formally I cannot see how this can obtain finite time estimates of q_t(x)  (t< \\infty) \n4. To clarify the previous point a bit further if indeed V(z)=z then this implies directly that -ln q_t = || z||^2/2 which implies directly that q_t is an independent Gaussian where q_t(z) solves the FPK equation for the density implied by the ODE dynamics in Eq 8. Remember U(z) and V(z) cannot be chosen independently, given U(z)  then grad(-U(z)) implies an FPK equation  where V(z) depends on the solution to that equation. One can't just go  oh let V(x) = ln Gaussian without first proving that q_t(z) = Gaussian.\n5. Continuing the line of reasoning dz_t =grad(-U(z))dt +dW_t should have as its steady state N(0,1)  if z_0 where samples from \\pi(z)=p(z|y) then in fact V(z) = z can easily be derived, but again you'd be just simulating an OU process at equilibrium … What is the advantage of this? Can you formally prove an advantage of doing this? Or at least motivate it conceptually ?\n\n\nQuestions on methodology:\n\n1. In your algorithm you don't have any pseudocode for the z(t) dynamics ? Why is this ? I found this quite confusing.\n2. For estimating q_t(x) you just say run normalizing flows  on x(t) ? so here you just do density estimation with normalizing flows correct ? (you learn flows in the standard way) . It really feels this is the essence to the method, all that is being done is running normalizing flows to estimate q_t(x) and then V(x) and then some adhoc stabilizing noise procedure is done with a lack of motivation. \n3. How is the preconditioning used then ? you do some extra gaussian dynamics in latent space with the learned flows to stabilize things a bit ?  (Please clarify and give more motivation)\n4. I can't seem to find an ablation for the effect of preconditioning/simulating the Gaussian, what effect does the adhoc preconditioning routine have? You can still use NF to fit the score without this step.\n", "limitations": "Overall there are clear counter points in this work that lack detailed discussion:\n\n1. When going from Langevin to the deterministic version, more approximation is required, in  sense it is a harder problem, whilst this allows for a hybrid version between ULA / normalizing flows, the authors do not motivate well why one would want to pursue such a problem (which requires more difficult approximation schemes).\n2. The whole derivation of the latent dynamics in the flow/gaussian space seems lacking and also is motivated. At the end of the day the equations you have seem to simulate Gaussian dynamics at equilibrium, this is either not particularly useful or wrong (if somehow they are out of equilibrium in practice). \n\n\n\n[1] Song, Y., Sohl-Dickstein, J., Kingma, D.P., Kumar, A., Ermon, S. and Poole, B., 2020. Score-based generative modeling through stochastic differential equations. arXiv preprint arXiv:2011.13456.", "rating": 4, "confidence": 4, "soundness": 3, "presentation": 2, "contribution": 2, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Pros:\n1. A sound methodology for simulating the deterministic “dual” of the Langevin equation using normalizing flows.\n2. Sound and successful experimental procedure.\n\nCons:\n\n1. The overall guarantees of the method (theoretically) are not clear and there's no acknowledgement if the method is mostly heuristic. \n2. Furthermore, intuition behind the method's preconditioning step is also not presented clearly. There are some potential conceptual/theoretical issues illustrated by the reviewer in point 3 (with a potential contradiction).\n3. Lack of motivation when moving from langevin to its stochastic counterpart.\n", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}]

Author rebuttal:

Decision: Accept