Denoising Diffusion Samplers
Authors: Francisco Vargas
Will Sussman Grathwohl
Arnaud Doucet
Venue: ICLR 2023 poster (2023)
Keywords: diffusion models
importance sampling
monte carlo
variational inference

Abstract: Denoising diffusion models are a popular class of generative models providing state-of-the-art results in many domains. One adds gradually noise to data using a diffusion to transform the data distribution into a Gaussian distribution. Samples from the generative model are then obtained by simulating an approximation of the time-reversal of this diffusion initialized by Gaussian samples. Practically, the intractable score terms appearing in the time-reversed process are approximated using score matching techniques. We explore here a similar idea to sample approximately from unnormalized probability density functions and estimate their normalizing constants. We consider a process where the target density diffuses towards a Gaussian. Denoising Diffusion Samplers (DDS) are obtained by approximating the corresponding time-reversal.  While score matching is not applicable in this context, we can leverage many of the ideas introduced in generative modeling for Monte Carlo sampling. Existing theoretical results from denoising diffusion models also provide theoretical guarantees for DDS. We discuss the connections between DDS, optimal control and Schr\"odinger bridges and finally demonstrate DDS experimentally on a variety of challenging sampling tasks.

Published as a conference paper at ICLR 2023
DENOISING DIFFUSION SAMPLERS
Francisco Vargas Will Grathwohl & Arnaud 2 University of Cambridge, DDeepMind
ABSTRACT
Denoising diffusion models are popular class of generative models providing state-of-the-an results in many domains. One adds gradually noise to data using a diffusion to transform the data distribution into Gaussian distribution. Samples from the generative model are then obtained by simulating an approximation the time-reversal of this diffusion initialized by Gaussian samples. Practically, the intractable score terms appearing in the time-reversed process are approximated using score matching techniques. We explore here similar idea to sample ap- proximately from unnormalized probability density functions and estimate their normalizing constants. We consider process where the target density diffuses to- wards a Gaussian. Denoising Diffusion Samplers (DDS) are obtained by approxi- mating the corresponding time-reversal. While score matching is not applicable this context, we can leverage many of the ideas introduced in generative modeling for Monte Carlo sampling. Existing theoretical results from denoising diffusion models also provide theoretical guarantees for DDS. We discuss the connections between DDS, optimal control and Schrödinger bridges and finally demonstrate DDS experimentally on variety of challenging sampling tasks.
INTRODUCTION
Let be a probability density on Rd of the form
7(x) Z
x ((x)
Z
(x)dxx Bdd
(1)
where y IR - R+ can be evaluated pointwise but the normalizing constant Z is intractable. We are here interested in both estimating Z and obtaining approximate samples from π.
A large variety of Monte Carlo techniques has been developed address this problem. In par- ticular Annealed Importance Sampling (AIS) (Neal 200T) and its Sequential Monte Carlo (SMC) extensions (Del Moral al. 2006) are often regarded as the gold standard to compute normalizing constants. Variational techniques are popular alternative to Markov Chain Monte Carlo (MCMC)
and SMC where one considers flexible family of easy-to-sample distributions q whose parameters are optimized by minimizing a suitable metric, typically the reverse Leiblec discrepancy
KL(4"I|>> (m). Typical choices for q" include mean-field approximation (Wainwright Jordan 2008) or normalizing flows (Papamakarios et al., 2021). To be able to model complex variational distri- butions, often useful to model (xr) as the marginal of an auxiliary extended distribution; i.e.
(x) f u)du. As this marginal is typically intractable, 0 is then learned by minimizing discrepancy measure between g (xx u) and extended target (x. u) (x))p((ula) where p" (u(z)) is an auxiliary conditional distribution (Agakov Barber 2004).
Over recent years, Monte Carlo techniques have also been fruitfully combined to variational tech- niques. For example, AIS can be thought of a procedure where qa(x,, u) is the joint distribution of Markov chain defined by sequence of MCMC kernels whose final state is while p (x. u) is the corresponding AIS extended target (Neal) (001). The parameters 0 of these kernels can then be optimized by minimizing (L p²) using stochastic gradient descent (Wu al. 2020 Geffner &
Domke 2021
Thin et al.) 2021 Zhang et al. 2021 Doucet et al. 2022 Gelfner & Domke 2022).
Instead of following AIS-type approach to define a flexible variational family, follow here an approach inspired by Denoising Diffusion Probabilistic Models (DDPM), a powerful class of
Work done while DeepMind

---

Published as a conference paper at ICLR 2023
generative models (Sohl-Dickstein et al. 2015 Ho et al. 2020 Song al., 2021c). In this context, one adds noise progressively to data using a diffusion to transform the complex data distribution into Gaussian distribution. The time-reversal of this diffusion can then be used to transform Gaussian sample into a sample from the target. While superficially similar to Langevin dynamics, this process mixes fast even in high dimensions as it inherits the mixing properties of the forward diffusion (De Bortoli et al., 2021 Theorem 1). However, as the time-reversal involves the derivatives of the logarithms of the intractable marginal densities of the forward diffusion, these so-called scores are practically approximated using score matching techniques. If the score estimation error is small, the approximate time-reversal still enjoys remarkable theoretical properties (De Bortoli 2022] Chen
et al. 2022
Lee et al. 2022).
These results motivate us to introduce Denoising Diffusion Samplers (DDS). Like DDPM, we con- sider forward diffusion which progressively transforms the target into Gaussian distribution. This defines an extended target distribution p(z, u) = (=)p((ulx) DDS are obtained by approx- imating the time-reversal of this diffusion using a process of distribution q' u). What distin- guishes DDS from DDPM is that we cannot simulate sample paths from the diffusion we want to time-reverse, used to approximate the score terms.
cannot sample its initial state from æ. Hence score matching ideas cannot be
We focus on minimizing KL(g" llp). equivalently maximizing an Evidence Lower Bound (ELBO), as in variational inference. We leverage a representation of this KI. discrepancy based on the intro- duction of a suitable auxiliary reference process that provides low variance estimate of this objective and its gradient. We can exploit the many similarities between DDS and DDPM to leverage some of the ideas developed in generative modeling for Monte Carlo sampling. This includes using the probability flow ordinary differential equation (ODE) (Song et al. 2021c) to derive novel normaliz- ing flows and the use of underdamped Langevin diffusions as a forward noising diffusion (Dockhorn) et al. 2022). The implementation of these samplers requires designing numerical integrators for the resulting stochastic differential equations (SDE) and ODE, However, simple integrators such as the standard
scheme do not yield valid BLBO in discrete-time. So as to guarantee one obtains a valid ELBO, DDS relies instead on an integrator for an auxiliary stationary refer- ence process which preserves its invariant distribution as well as an integrator for the approximate rive-seversal inducing a distribution absolutely continuous w.r.t. the distribution of the discretized reference process. Finally compare experimentally DDS to AIS, SMC and other state- of- the-art Monte Carlo methods on a variety of sampling tasks.
2 DENOISING DIFFUSION SAMPLERS: CONTINUOUS TIME
to gain insight on the structure of the time- reversal we want to approximate. We introduce C = C( (0. TJR²)) the space of continuous functions from [0, 7) Rd and ((C) the Borel sets on C. We will consider in this section path measures which are probability measures on (C, B(C),, see Léonard (2014a) for formal definition. Numerical integrators are discussed in the following section.
We start here by formulating DDS in continuous-
2.1 FORWARD DIFFUSION AND ITS TIME- REVERSAL
Consider the forward noising diffusion given by an Amnstiin-Uhlenbee (OU) process
dæ; d + a 22,,dB.
Fg n.
(2)
where (Br) (Be) is d-dimension Brownian motion and t 3 is non-decreasing positive function. This diffusion induces the path measure P on the time interval [0, TJ and the marginal density of xr is denoted pt. The transition density of this diffusion given by (1)()a(20) =
(01 v1 AeTo, o2 X,T) where X 1 exp(-2 f B,ds).
From now on, we will always consider a scenario where 1-" 3,ds > so that (r(x)) x NeT: 0, IT). We can thus think of 2 as transporting approximately the target density to this Gaussian density.
This referred to Variance Preserving diffusion by Song et al. 2021c
2

---

Published as a conference paper at ICLR 2023
From (Haussmann & Pardoux
1986,, its time-reversal (v.).e(a.7] (IT-1)eE(0.)) where equality
here in distribution, satishes
dy 13 20 (o
yo PT,
where (Wi)rejo.T] is another d-dimensional Brownian motion. By definition this time-reversal starts from yo ~ (y())) s (yo; 0, 1) and is such that yr ~ л. This suggests that if we could approximately simulate the diffusion (3) then we would obtain approximate samples from и.
However, putting this idea in practice requires being able to approximate the intractable scores (V (z))(()) To achieve this, DDPM rely
matching techniques (Hyvärinen 2005 Vincent 2011) as it easy to sample from 2 This impossible in our scenario as sampling from (2 requires sampling 2n which is impossible by assumption.
2.2 REFERENCE DIFFUSION AND VALUE FUNCTION
In our context, it is useful to introduce a reference process defined by the diffusion following (2) but initialized at v (zo) N (x0; 0, o² rather than (20)) thus ensuring that the marginals of the resulting path measure prel all satisfy 2t (x+) = N(x;; 0, o²f). From the chain rule for KL for path measures (see e.g. (Léonard 2014b] Theorem 2.4) and (De Bortoli et al. 2021 Proposition 24)), Thus P can be identified as the path measure minimizing the KL discrepancy w.r.t. peo over the set of path measures Q with marginal (=())) x<zg) at time
one has KL( Q| preff KL Ib + -ug [KL ( ((- [-
i.c. arg min ( go = T).
time-reversal representation of pref is given by (yo)eEJo.T (Dr-thee satisfying
dyt = -- --3 + a 23r
yo ~ po
(4)
As Br 20 V In pf ( (v:) ) 3r-1
rewrite the time-reversal
3 of
dy: = 3 - 20 (3:))d + +o 23 yo ~ PT:
(5)
where (x(x) (z))//p (z) is so-called value function which can be shown to satisfy Kol- mogorov equation such that or (z,) = (00 the expectation being w.r.t. peet.
2.3 LEARNING THE TIME-REVERSAL
To approximate the time-reversal induced by
3 of P, consider path measure Q whose time-reversal is
dy = Br-clu + 20 - (. + 23 yo ~ N o2
(6)
that yt -r- To obtain 2t) (z)), we parameterize so (1, by neural network whose parameters are obtained by minimizing
kL()")1P) oZ
(0 (KL(( (-Bb))
(7)
(T, o²T)(lrr) + o T [[ Br- 1, v (
where we have used the chain rule for KI. then Girsanov's theorem (see e.g. Klebaner 2012). This expression of the KI. is reminiscent of the expression obtained in (Song et al.) 2021b Theorem 1) in the context DDPM. However, the expectation appearing in 7 is here w.r.t. Q" and not w.nt. P and we cannot get rid of the intractable scores (V In pe using score matching ideas.
Instead, taking inspiration from
( we reparameterize the time-reversal of Q using
dy -BT- 20 [T dt + 23
yo ~o0,000²²
(8)
ie. (t. = so (t,) - In (x) = so x) + x/02. This reparameterization allows us to express KL (P) in compact form.
Proposition 1. The Radon-Nikodynn derivative d apme (MIo.T) satisfies under Q
dQd dppret
In
0² S Br Br-
1dt + - 2 -t TdW
(9)
3

---

Published as a conference paper at ICLR 2023
From the identity KL(@"I KL(L" + y
[In
Po (ur) polur)
)1
it follows that
KL(O2)|) = o (o² B Br- - Me 2 +
in
(10)
T(yr)
For minimizing 10 approximate samples from и can obtained by simulating 8 and returning yT ~g.
We can obtain an unbiased estimate of Z via the following importance sampling identity
((r)
dp
(Ur; 0, o2 d@ (JOoTT)),
(11)
where the second term can be computed directly from (9) and V/,, is obtained by solving 88.
2.4 CONTINUOUS-TIME NORMALIZING FLOW
To approximate the log likelihood of DDPM, it was proposed by [Snn et al. (2021c) to rely on the probability flow ODE. This is ODE admitting the same marginal distributions (Pr)teju.21 as the noising diffusion given by dæ. - Be { + o² V In pe(3)) Jdt. We here this ODE design continuous-time normalizing flow to sample from For O such that fo(t, x) V In (x (obtained
by minimization of the KL in (0)). we have Ia + o²V In Pa (x(x) as folt,, ((t,,)). æ). So is possible to sample approximately from by integrating the following ODE over []
dyt- = o² - L, )t,, If denote by o the distribution of then yr o is an approximate sample from We can use this normalizing flow to obtain an unbiased estimate of Z using importance sampling 2 = ) / (yr) for yr ~ % Indeed, contrary to the proposal of induced by 8 we can compute pointwise %6 using the instantaneous change of variables formula (Chen et al. 2018) such that In 7 (yr) In (wo: 0, o, o² - (T -1, )rl and where (ye arises from the ODE 12) The expensive divergence term can be estimated using the Hutchinson estimator Grathwohl et al. 2018: Song et al. 2021c).
~ o²),
(12)
2.5 EXTENSION TO UNDERDAMPED LANGEVIN DYNAMICS
For DDPM, it has been proposed extend the original state x Rd by momentum variable E Rd. One then diffuses the data distribution augmented by a Gaussian distribution on the mo- mentum using an underdamped Langevin dynamics targeting Wx; o2I) (m, O, M) where M is positive definite mass matrix (Dockhorn et al.) 2002). It was demonstrated empirically that the resulting scores are smoother, hence easier to estimate and this leads to improved performance. We adapt here this approach to Monte Carlo sampling: see Section B for more details.
We diffuse 7 (2x m) = π (=) 0, .)) using an underdamped Langevin dynamics, i.e.
Ft mtdt. -- dm, 02 df Sem,dt + /23.M/ (Zo, mo) ~ n.
dæ, = M
(13)
The resulting path measure on [0, T] denoted and the marginal of (I,, mi) is denoted ne-
Here, the reference process preff is defined by the diffusion 13 initialized using zo (0 o²I), mo ~ (O, M).. This ensures that ni (x,, mi) = A(r;; o² (" 0, M) for all t and the time-reversal process (μ,, Nz)te(o,1] (i c) (where equality is in distribu-
tion) of this stationary diffusion satisfies
dy: = - dn,
dt Br-mn, + , o²
(14)
Using manipulations identical Section 2.2 the time-reversal of P can be also be written for (1, m) = (I, m)m / ( (x, m) as dy: M -n and
dn: 3/ di ddd 23r v., In ør d ne)d M23 -
(15)
To approximate P. we will consider parameterized path measure @" whose time reversal is defined for (yo, no) ~ (yo; 0, 12 0, M) by dy, = -M -n and
y/1 dn. " dl dn 2B3 M (o 1, ne)dt +
(16)
We can then express the KL of interest a compact way similar to Proposition I see Appendix B.2 A normalising flow formulation is presented in Appendix B.3
4

---

Published as a conference paper at ICLR 2023
Algorithm DDS Training
1: Input: > y R - Rt, (3x) E (R+)*, E R2, X>
2: for 1...., training iterations do
3:
for = -...., N do # Iterate over training batch size.
4:
Sample yo NO 1) and set ro = 0
5: for k O,..., K 1 do # Iterate over integration steps.
6:
AK- = VI ak-k
7:
gi+1, = / K-KKAK + 20 A fo- k, 3k to taakkanooaakkkk... Ekm
iid. NO
rk+I, = Tk,n + 20 2R [| - ki Uk.,n) 112
8:
ax-x
end for 9:
Q+0 0 1 En=1 (rK,n + In (N(uk.;
10:
11: end for
12:
end for
13: Return:
3 DENOISING DIFFUSION SAMPLERS: DISCRETE-TIME
We introduce here discrete-time integrators for the SDEs and ODEs introduced in Section 2 Con- trary to DDPM, we not only need an integrator for Qs but also for prf to be able to compute an approximation of the Radon-Nikodym derivative cl / /dpref Additionally this integrator needs to be carefully designed as explained below to preserve an ELBO. For sake simplicity, we consider constant discretization step á > 0 such that K T/8 is an integer. In the indices, we write k for t kō to simplify notation; c.g. Zy for Økg-
3.1 INTEGRATOR FOR pref
Proposition T shows that learning the time reversal in continuous-time can be achieved by min-
imising the objective KL(Q") KL(Q"I[P) given in 10 This expression is obtained by using the identity
kL(Q" (P) = ( [preft + wo [In (Po (ur)/Pu(Ur)) This is also equivalent to maximiz-
ing the corresponding ELBO, Eo [In 2). for the unbiased estimate of Z given in 11
An Eurer-Maruyama (EM) discretisation of the corresponding SDEs is obviously applicable. How- ever, it is problematic as established below.
Proposition 2. Consider integrators of ppe and @ leading to approximations 1td (xo.k) and (IO.K). Then (L + E yx [In (p (UK)//(KK)) is guaranteed to be non-negative and Egg [In 7] an ELBO for 7 (π(yx) (Vo.k) ()):k))((A (yK;; 0, 0g (vo:k) if one has PK (μк) = o (vK). The EM discretisation of piel does not satisfy this condition.
A direct consequence of this result is that the estimator for In Z that uses the EM discretisation be such that Bo [In 2) > In Z observed empirically in Table 5 To resolve this issue, derive in the next section integrator for proe which ensures that pta y ppof (y) for all k.
3.2 ORNSTHIN-UHLANBEER
We integrate exactly
4 This is given by yo (N(0, and
uk+1 V1 ak- /OKKKEEK
Ek iid. N(OO I),
(17)
for = 1-exp (-2
p6 11k 3.ds) This defines the discrete- time reference process. We propose
the following exponential type integrator (ee Bortoli 2022 for initialized using yo (O(
1/++1 11-ak- -KKKK + 20 / ax - So(K k. yx) + o /OKK
(18)
Henceforth, it assumed that ag can he computed exactly.
5

---

Published as a conference paper at ICLR 2023
196 traiiing oi locs (wreophers stap4 owr ovaluss
005 mrw.ing loss benespherr, stape over
2700
0330
3200
1.00
0222
Figure 1: Training loss per hyperparameter: PIS (left) vs DDS (right).
where EN iid. (V, Equations
17 and 18 define the time reversals of the refer- ence process (x):k) and proposal (To.k) For its time reversal, we write q (w0..) = (yo: 0, o)
II k= of 11k (vK ++1 (y abusing slightly notation, and similarly for (yo:x).
We will be relying the discrete-time counterpart of (9)
Proposition 3. The log density ratio between qº (yo.k) and pr (yo.kc satisfies for yo.k qa (yo:k)
(((00.k) p (Jo:K)
K
K
A
In
20t
fo yK 111 20
folk,
-k) T (19)
Okk
Vark
k=1
k=1
where dE= : - VT Ork and 6k defined through (189 is such that &k N(O, 1). In particular, one obtains from KT((T"1pp) KKL(("" + Egg In ((
)] that
X2 (k. yk xxl + (
(L Enp
20
(20)
xyxx
We compute an unbiased gradient of this objective using the reparameteriization trick and the JAX software package (Bradbury et al.) 2018). The training procedure is summarized in Algorithm i Appendix D Unfortunately, contrary to DDPM, & K is not available in closed form for h: < K
so we can neither mini-batches over the time index & without having to simulate the process until the minimum sampled time nor reparameterize F& Ho et al.; (2020). Once we obtain the parameter 9 minimizing (0 DDS samples from q" using (8). The final sample μK has a distribution % approximating by design. By using importance sampling,
obtain unbiased estimate of
the normalizing constant Z (7(yK) (yo.k))
O, (yo.k)) for yo:k qa (Mo. k). Finally Appendices and B.6 extend this approach to provide a similar approach to discretize the underdamped dynamics proposed in Section 2.5 In this context, the proposed integrators rely on a
leapfrog scheme (Leimkuhler & Matthews
2016).
3.3 THEORETICAL GUARANTEES
Motivated by DDPM, bounds on the total variation between the target distribution and the distribu- tion of the samples generated by a time-discretization of approximate time reversal of forward noising diffusion have been first obtained in (De Bortoli et al., 202T) then refined in (De Bortoli 2022 Chen et al., 2022, Lee et al. 2022) and extended the Wasserstein metric. These results are directly applicable to DDS because their proofs rely on assumptions on the score approximation error but not on the way these approximations are learned, i.e. via score matching for DDPM and reverse KL for DDS. For example, the main result in Chen et al. (2022) shows that if the true scores are L-Lipschitz, the L2(p:) error on the scores is bounded and other mild integrability assumptions then DDS outputs samples e-close in total variation to in ((L2d/e²) time steps. As pointed out by the authors, this matches state-of-the-aar complexity bounds for Langevin Monte Carlo algorithm for sampling targets satisfying log-Sobole inequality without having to make any log- concavivyvty assumption on म. However, the assumption on the approximation error for score estimates is less realistic for DDS than DDPM as we do not observe realizations of the forward diffusion.

---

Published as a conference paper at ICLR 2023
funnel
losp
ion
Inz estimate
extimate
130
130
1350
330
230
209
Staps
100
200
Steps
30
30
Figure 2: In Z estimate (median plus upper/lower quartiles) as a function of number of steps K a) Funnel, b) LGCP, e) Logistic lonosphere dataset. Yellow dotted line MF-VI and dashed magenta is the gold standard.
3.4 INTEPRETATION, RELATED WORK AND EXTENSIONS
DDS as KL and path integral control. The reverse KL we minimize can be expressed as
KL p) = E.g [In
(N 0, 0, T(Fo)
(21)
Po- (22-) (t,)
This objective is a specific example of KL control problem (Kappen et al. 2012). In continuous- time, (10p corresponds to a path integral control problem; see e.g. (Kappen & Ruiz 2016)
Heng et al. (2020) use KL control ideas so as to sample from a target π and estimate Z. How-
ever, their algorithm relies on pn (x): being defined by a discretized non-homogeneeus Langevin dynamics such that px (xx) is typically not approximating known distribution. Additionally it
approximates the value functions (Φx) using regression using simple linear/quadratic functions. Finally it relics cability of this methodology to a restricted class of models.
good initial estimate of p(æ0:k) obtained through SMC. This limits the appli-
Connections to Schrödinger Bridges. The Schrödinger Bridge (SB) problem (Léonard 2014a De Bortoli et al., 202T) takes the following form in discrete time. Given a reference density (T"( we want to find the density p3h (xo:k) s.t. pth arg min, (KL(q]lp"))
qo Auu 9K
μx} where μo, P.K are prescribed distributions. This problem can be solved using iterative propor- tional fitting (IPF) which is defined by the following recursion with initialization p² = pref
pdn arg min (KL(g) p22=-1) go =oot, p2m+1 = arg min (KL())p²²") : = pK). (22)
Consider the SB problem where (z)) = π (2o). uu (#k) = (x 0, o² and the time-reversal
of p (œ0:k) is defined through 17) In this case, pd p corresponds to the discrete-time version of the noising process and pd to the time-reversal of but initialized at μи instead of pк. This
the process DDS is approximating. As Pк μK for K large enough, we have approximately pth Ps ~ P2. We can thus think of DDS as approximating the solution to this SB problem.
Consider now another SB problem where (o(2)) =(x0), px (xx) ôg (x and pre (zo:k)
5x(0x) II ((xx; 80 I), i.e. o is a pinned Brownian motion running backwards in time.
This SB problem was discussed in discrete-time in (Beghi 1996) and in continuous-tin time in (Föllmer)
1984 Dai Pra 1991 Tzen & Raginsky 2019). In this case, p² (x = T ((xo))²= (x (r)) is modified "noising" process that transports π to the degenerate measure do and it easily shown that ybb = p² Sampling from the time-reversal of this measure would generate samples from и starting from do- Algorithmically, Zhang al.) (2021) proposed approximating this time-reversal by some projection on some cigenfunctions. In parallel, Barr et al. (202)). Vargas et al. (202T) and Zhang & Chen 2022) approximated this SB by using a neural network parameterization of the gradient the logarithm of the corresponding value function trained by minimizing a reverse KI. We will adopt the Path Integral Sampler (PIS) terminology proposed by Zhang & Chen (2022) for this approach.
DDS can thus be seen an alternative to PIS which relies on reference dynamics corresponding to an overdamped or underdamped OU process instead of a pinned Brownian motion. Theoretically the drift of the resulting tive-seversal for DDS is not steep as for PIS (see Appendix A.2) and empirically this significantly improves numerical stability of the training procedure; see Figure I The use of pinned Brownian motion is also not amenable the construction of normalizing flows.

---

Published as a conference paper at ICLR 2023
brownian
In2 estimate
Inz astimate
Ind
224
Steps
Steps
Steps
Figure 3: In Z estimate as a function of number of steps K a) Logistic Sonar dataset, b) Brownian motion, e) NICE, Yellow dotted line is MF-VI and dashed magenta is the gold standard.
Forward KL minimization. Jing al. 20022, diffusion models ideas are also use to sample unnormalized probability densities. The criterion being minimized therein is the forward KI. as in DDPM, that is KL(l)lg"). As samples from p are not available, an importance sampling approxi- mation of based on samples from de is used to obtain estimate of this KL and its gradient w.r.. The method is shown to perform well in low dimensional examples but is expected to degrade significantly in high dimension the importance sampling approximation will be typically poor in the first training iterations.
4 EXPERIMENTS
We present here experiments for Algorithm In our implementation, So follows the PIS-GRAD network proposed in (Zhang & Chen 2022 (o(k, x)) = NN (K, x; 0) + (N;(k; o V In (x)) Across all experiments we use two layer architecture with 64 hidden units each (for both networks), as in Zhang & Chen (202), with the exception of the NICE (Dinh et al. (014) target where we use 3 layers with 512, 256, and 64 units respectively. The final layers are initialised to in order to
make the path regularisation term null. We 0t 1/2 1/2
Omax cos? (3 1-k/K+s with 0.008 as I+a
in (Nichol & Dhariwal, 2021). We found that detaching the target score stabilised optimization in both approaches without affecting the final result. We adopt this across experiments, an ablation of this feature can be seen in Appendix C.9.1
Across all tasks we compare DDS to SMC (Del Moral et al., 2006 Zhou al. 2016). PIS (Barr ett al., (2200) Vargas et al. 202T Zhang & Chen 2022). and Mean Field-Vi (MI- VI) with a Gaus- sian variational distribution. We also compare DDS to AIS Neal 2001) and optimized variants of
AIS using score matching (MCD) (Doucet et al. 2022 Geftner & Domke (022) for two standard Bayesian models. Finally we explore a task introduced in (Doucet et | (022) that uses a pre-trained normalising flow as a target. Within this setting we propose a benchmarking criterion that allows us to assess mode collapse in high dimensions and explore the benefits of incorporating inductive biases into fg. We carefully tuned the hyper- typer-parameters of all algorithms (e.g. step size, diffusion coefficient, and such), details can be found in Appendix C.2 Finally training time can be found in Appendix C.5 Additional experiments for the normalizing flows presented in Appendix C.11 and for the underdamped approach in Appendix CC.22 We note that these extensions did not bring any benefit compared to Algorithm
4.1 BENCHMARKING TARGETS
We first discuss two standard target distributions which are often used to benchmark methods; see e.g. (Neal 2003 Arbel et al. 2021 Heng et al. 2020 Zhang & Chen 2022). Results are presented Figure
Funnel Distribution:
This 10-dimensional challenging distribution is given by 7 = 0, Tx) where of = (Neal 2003).
Zhang & Chen 2022) inadvertently considered the more favourable oy scenario for PIS but used 0f for other methods. This explains the significant differences between their results and ours.
8

---

Published as a conference paper at ICLR 2023
Log Gaussian Cox process: This model arises in spatial statisties (Möller al. 1998). We use a = M x M = 1600 grid, resulting in the unnormalized target density 7(z) =
N( K. Iep exp(z.y a exp((z))).
4.2 BAYESIAN MODELS
We explore standard Bayesian models and compare them with standard AIS and MCD bench- marks (See Appendix (C8) in addition to the SMC and VI benchmarks presented so far. Results for lonosphere can be found in the 3rd pane of Figure whilst SONAR and Brownian are in Figure
Bernouli(sigmoiddl u.)). This Bayesian logistic model is evaluated on two datasets, Ionosphere (d 32) and Sonar (d = 61).
Logistic Regression: We set a ~ (O, o (),
Brownian Motion: We consider a discretised Brownian motion with a Gaussian observation model and latent volatility as time series model, d = 32. This model, proposed in the software package developed by Sountsov et al. (2020). is specified in more detail in Appendix C.3
4.3 MODE COLLAPSE IN HIGH DIMENSIONS
Normalizing Flow Evaluation: Following Doucet et (2022) we train NICE (Dinh et al. 2014)
on down-sampled d = 14 x 14 variant of MNIST (LeCun & Cortes 2010) and use the trained model as our target. As we can generate samples from our target, we evaluate the methods samples by measuring the Sinkhorn distance between true and estimated samples. This evaluation criteria allows to asses mode collapse for samplers in high dimensional settings. Results can be seen in Table T and the third pane of Figure 3
Method
DDS
PIS
SMC
MCD
In Z Estimate
3,2004 + 0.645
-3.933 + 0.754
-4.255 + 2.043
-6.25
W
0.00 (Terue: f)
658.079
658.778
750.245
NA
Table I: Results on NICE target, we performed 30 runs with different seeds for each approach. For
MCD we used the results from (Doucet et al. 2022).
5 DISCUSSION
We have explored the use DDPM ideas to sample unnormalized probability distributions and estimate their normalizing constants.
The DDS in Algorithm [ is empirically competitive with state-of-the-ar SMC and numerically more stable than PIS. This comes at the cost of non-negligible training time compared SMC. When accounting for it, SMC often provide better performance on simple targets. However, in the challenging multimodal NICE example, even a carefully tuned SMC sampler using Hamiltonian Monte Carlo transitions was not competitive to DDS. This is despite the fact that DDS (and PIS) are prone to mode dropping any method relying on the reverse KL.
We have also investigated normalizing flows based on the probability flow ODE as well as DDS based on an underdamped dynamics. Our experimental results were disappointing both cases in high dimensional scenarios. We conjecture that more sophisticated numerical integrators need to be developed for the normalizing flows to be competitive and that the neural network parameteriza- tion used in the underdamped scenario should be improved and better leverage the structure of the logarithmic derivative of the value functions.
Overall DDS are class of algorithms worth investigating and further developing. There has much work devoted to improving successfully DDPM over the past two years including, among many others, modified forward noising mechanisms (Hoogeboom & Salimans 2022), denoising diffusion implicit models (Song et al. 2021a) and sophisticated numerical integrators (Karras et al. 2022). highly likely that some of these techniques will lead to more powerful DDS. Advances KL and path integral control might also be adapted to DDS to provide better training procedures; see e.g.
Thalmeier al. 2020).

---

Published as a conference paper at ICLR 2023
REFERENCES
Felix V Agakov and David Barber. An auxiliary variational method. In Advances in Neural Infor- mation Processing Systems, 2004.
Michael Arbel, Alex Matthews, and Arnaud Doucet. Annealed flow transport Monte Carlo. In International Conference on Machine Learning, 2021.
Ariel Barr, Willem Gispen, and Austen Lamacraft. Quantum ground states from reinforcement learning. In Mathematical and Scientific Machine Learning, 2020.
Alessandro Beghi. On the relative entropy of discrete-time Markov processes with given end-point densities. IEEE Transactions on Information Theory, 42(5): 1533,1552-5 1996.
James Bradbury, Roy Frostig, Peter Hawkins, Matthew James Johnson, Chris Leary, Dougal Maclaurin, George Necula, Adam Paszke, Jake VanderPlas, Skye Wanderman- Milne, and Qiao Zhang. JAX: composable transformations of Python+NumPy programs, 2018. URL http:
//github. com/ cww.///o/////aaaaaooooo
Ricky TQ Chen, Yulia Rubanova, Jesse Bettencourt, and David Duvenaud. Neural ordinary differential equations. Advances in Neural Information Processing Systems, 2018.
Sitan Chen, Sinho Chewi, Jerry Li, Yuanzhi Li, Adil Salim, and Anru R. Zhang. Sampling is as easy as learning the score: theory for diffusion models with minimal data assumptions. arXiv preprint
arXiv:22209.11225. 2022.
Paolo Dai Pra. A stochastic control approach to reciprocal diffusion processes. Applied Mathematics and Optimization, 23(1):313 329, 1991.
Valentin De Bortoli. Convergence of denoising diffusion models under the manifold hypothesis. Transactions on Machine Learning Research, 2022.
Valentin De Bortoli, James Thomton, Jeremy Heng, and Arnaud Doucet. Diffusion Schrödinger bridge with applications to score-based generative modeling. In Advances in Neural Information Processing Systems, 2021.
Pierre Del Moral, Arnau Doucet, and Ajay Jasra. Sequential Monte Carlo samplers. Journal of the
Royal Statistical Sociery: Series B (Statistical Methodology), 68(3):411 114411446, 2006.
Laurent Dinh, David Krueger, and Yoshua Bengio. Nice: Non-linear independent components esti-
mation. arXiv preprint 110.55141..8516 2014.
Tim Dockhorn, Arash Vahdat, and Karsten Kreis. Score-based generative modeling with critically- damped Langevin diffusion. In International Conference on Learning Representations, 2022.
Arnaud Doucet, Will Grathwohl, Alexander G de G Matthews, and Heiko Strathmann. Score-based diffusion meets annealed importance sampling. In Advances in Neural Information Processing Systems, 2022.
Hans Föllmer. An entropy approach to the time reversal of diffusion processes. Lecture Notes in Control and Information Sciences, 69: 156-1,663, 1984.
Tomas Geffner and Justin Domke. MCMC variational inference via uncorrected Hamiltonian an- nealing. Advances in Neural Information Processing Systems, 2021.
Tomas Geffner and Justin Domke. Langevin diffusion variational inference. arXiv preprint 0rX10:2208.07743, 2022.
Will Grathwohl, Ricky TQ Chen, Jesse Bettencourt, Ilya Sutskever, and David Duvenaud. Ffjord: Free-form continuous dynamics for scalable reversible generative models. In International Con- ference on Learning Representations, 2018.
Klaus Greff, Raphaël Lopez Kaufman, Rishabh Kabra, Nick Watters, Christopher Burgess, Daniel Zoran, Loic Matthey, Matthew Botvinick, and Alexander Lerchner. Multi-object representation learning with iterative variational inference. International Conference on Machine Learning, 2019.
10

---

Published as a conference paper at ICLR 2023
Ulrich G Haussmann and Etienne Pardoux. Time reversal of diffusions. The Annals of Probability, 14(3):1188-122) 1986.
Jeremy Heng, Adrian N Bishop, George Deligiannidis, and Arnaud Doucet. Controlled sequential Monte Carlo. The Annals of Statistics, 48(5):2904 2929, 2020.
Matteo Hessel, David Budden, Fabio Viola, Mihaela Rosca, Eren Sezener, and Tom Hennigan. Optax: composable gradient transformation and optimisation, in jax., 2020.
Jonathan Ho, Ajay Jain, and Pieter Abbeel. Denoising diffusion probabilistic models. In Advances in Neural Information Processing Systems, 2020.
Emiel Hoogehoom and Tim Salimans. Blurring diffusion models. arXiv preprint 22055,20990555 2022.
Aapo Hyvärinen. Estimation of non-normalized statistical models by score matching. The Journal of Machine Learning Research, 6:695-709, 2005.
Bowen Jing, Gabriele Corso, Jeffrey Chang, Regina Barzilay, and Tommi Jaakkola. Torsional diffu- sion for molecular conformer generation. arXiv preprint arXiv:2206.01729, 2022.
Hilbert Kappen and Hans Christian Ruiz. Adaptive importance sampling for control and inference. Journal of Statisticul Physics, 162(5);1244-120) 2016.
Hilbert J Kappen, Vicenç Gómez, and Manfred Opper. Optimal control as a graphical model infer- ence problem. Machine learming, 17221119--182, 2012.
Tero Karras, Miik Aittala, Timo Aila, and Samuli Laine. Elucidating the design space of diffusion- based generative models. Advances in Neural Information Processing Systems, 2022.
Diederik P. Kingma and Jimmy Ba. Adam: A method for stochastic optimization. In International Conference on Learning Representations, 2015.
Diederik Kingma and Max Welling. Auto-encoding variational Bayes. arXiv preprint 13126613 14, 2013.
Fima C Klebaner. Introduction Stochastic Calculus with Applications. Imperial College Press, 2012.
Yann LeCun and Corinna Cortes. http://aam..ecun..co///edd///// 2010.
MNIST
handwritten digit database.
URL http:/ /:tp://yann..mm http:/aaannnnnnnn.nmmmm lecun. www.ccoccceooooooooo com/exdb/
mni st.
Yann LeCun, Léon Bottou, Yoshua Bengio, and Patrick Haffner. Gradient-based learning applied to
document recognition. Proceedings of the IEEE, 86(11):228--3324 1998.
Holden Lee, Jianfeng Lu, and Yixin Tan. Convergence of basre-based generative modeling for general data distributions. arXiv preprint arXiv:2209.122 12381, 2022.
Ben Leimkuhler and Charles Matthews. Molecular Dynamics with Deterministic and Stechastic Numerical Methods. Springer, 2016.
Christian Léonard. A survey of the Schrödinger problem and some of its connections with optimal transport. Discrete and Continuous Dynamical Systems-Series A, 34(4):1533-1574, 2014 4a.
Christian Léonard. Some properties of path measures. Séminaire de Probabilités XLVI, pp- 207-230. Springer, 2014b.
Luke Metz, Niru Maheswaranathan. Ruoxi Sun, C Daniel Freeman, Ben Poole, and Jascha Sohl- Dickstein. Using thousand optimization tasks to learn hyperparameter search strategies. arXiv preprint arXiv:2002. 1887, 2020.
Jesper Møller, Anne Randi Syversveen, and Rasmus Plenge Waagepetersen. Log Gaussian Cox processes. Scandinavian Journal of Statistics, (533):55-4822 1998.

---

Published as a conference paper at ICLR 2023
Radford M Neal. Annealed importance sampling, Statistics and Computing, 11(2):125-133 2001.
Radford M. Neal. Slice sampling. The Annals of Statisties, 31:705-767, 06 2003.
Alexander Quinn Nichol and Prafulla Dhariwal. Improved denoising diffusion probabilistic models. In International Conference on Machine Learning, 2021.
George Papamakarios, Eric Nalisnick, Danilo Jimenez. Rezende, Shakir Mohamed, and Balaji Lak- shminarayanan. Normalizing flows for probabilistic modeling and inference. The Journal of Machine Learning Research, 22(1):2617-2680 2021.
Paavo Parmas, Carl Edward Rasmussen, Jan Peters, and Kenji Doya. Pipps: Flexible model-based policy search robust to the curse of chaos. In International Conference on Machine Learning, 2018.
Jascha Sohl-Dickstein, Eric Weiss, Niru Maheswaranathan, and Surya Ganguli. Deep unsupervised learning using nonequilibrium thermodynamics. In International Conference Machine Learn- ing, 2015.
Jaaming Song, Chenlin Meng, and Stefano Ermon. Denoising diffusion implicit models. In Interna- tional Conference on Learning Representations, 2021a.
Yang Song, Conor Durkan, lain Murray, and Stefano Ermon. Maximum likelihood training of score- based diffusion models. In Advances in Neural Information Processing Systems, 2021b.
Yang Song, Jascha Sohl-Dickstein, Diederik P Kingma, Abhishek Kumar, Stefano Ermon, and Ben Poole. Score-based generative modeling through stochastic differential equations. In Interna- tional Conference on Learning Representations, 2021c.
Tommi Sottinen and Simo Särkkä. Application of Girsanov theorem to particle filtering of discretely observed continuous-time non-linear systems. Bayesian Analysis, 3(3):555-555 5844, 2008
Pavel Sountsov, Alexey Radul, and contributors. Inference gym, 2020. URL lnttps://yypi. http/project/nnfeeee
Dominik Thalmeier, Hilbert J Kappen, Simone Totaro, and Vicenc Gómez. Adaptive smoothing for path integral control. Journal of Machine Learning Research, 21:1-37, 2020.
Achille Thin, Nikita Kotelevskii, Alain Durmus, Eric Moulines, Maxim Panov, and Amaud Doucet. Monte Carlo variational auto-encoders. International Conference on Machine Learning, 2021.
Belinda Tzen and Maxim Raginsky. Theoretical guarantees for sampling and inference in generative models with latent diffusions. In Conference on Learning Theory, 2019.
Francisco Vargas, Andrius Ovsianas, David Fernandes, Mark Girolami, Neil Lawrence, and Nikolas Nüsken. Bayesian learning via neural Schrödinger-Fölller flows. arXiv preprint arXiv:2111.10510. 2021.
Paul Vicol, Luke Metz, and Jascha Sohl-Dickstein. Unbiased gradient estimation in unrolled com- putation graphs with persistent evolution strategies. In International Conference on Machine Learning, 2021.
Pascal Vincent. A connection between score matching and denoising autoencoders. Neural Compu- tation, 1367):166111674, 2011.
Martin J. Wainwright and Michael I. Jordan. Graphical Models, Exponential Families, and Varia- tional Inference. Foundations and Trends@ in Machine Learning, 1(1-2):1-305, November 2008.
Hao Wu, Jonas Köhler, and Frank Noé. Stochastic normalizing flows. Advances in Neural Information Processing Systems, 2020.
Benjamin Zhang, Tuhin Sahai, and Youssef Marzouk. Sampling via controlled stochastic dynamical systems. In / (Still Can't 1 Believe It's Not Better! NeurlPS 202/ Workshop, 2021.
12

---

Published as a conference paper at ICLR 2023
Guodong Zhang, Kyle Hsu, Jianing Li, Chelsea Finn, and Roger Grosse. Differentiable annealed im- portance sampling and the perils of gradient noise. In Advances in Neural Information Processing Systems, 2021.
Qinsheng Zhang and Yongxin Chen. Path integral sampler: stochastic control approach for sam- pling. In International Conference on Learning Representations, 2022.
Yan Zhou, Adam M Johansen, and John AD Aston. Toward automatic model comparison: an adaptive sequential monte carlo approach. Journal of Computational and Graphical Statistics, 25 (3):701-726, 2016.
A APPENDIX
A.1 PROOF OF PROPOSITION
The Radon-Nikodym derivative expression directly follows from application of Girsanov theorem (see e.g. Klebaner (2012). The path measures P and peef are induced by two diffusions following the same dynamics but with different initial conditions so that
P(w)= = ((W)
T(VT)
(wl(()))))
012r;
TUTT)
[w) 0y
Note we ap- ply Girsanov theorem to the time reversals of the path measures Q" and P using that the Radon- Nikodym between two path measures and their respective reversals are the same. As W2 is a Brow- nian motion under Q", the final expression (10) for KL (LP))) follows.
So it follows directly that KL(Q" PP) = KL(Q") [ + EE vn In
Pb
( (ur)
polor)
A.2 COMPARING DDS AND PIS DRIFTS FOR GAUSSIAN TARGETS
The optimal drift for both PIS and DDS can be expressed in terms of logarithmic derivative of the value function plus additional terms:
bDDs(ar (ar, 1) - (r 20 In ( (ris( () a o²V, In eT(((s)
(23)
where &e is the corresponding value function for DDS or PIS respectively.
Recall that In (x)) In (x) In (x) For a target (2x) NZ:;H, p, pf (Z2 pr brd , I) and Prjo(=a = -N(x,; a,zp, b,I), we obtain
(x(x)) / π ((20))11((4/00)dd = NNx; x b.l)
so
V -a) b,I) -x 0141 + (b,e) (6 )ol)
(24)
Corollary 1. For DDS, we have pafed O, lلet 021, a: -1- - 0² X with X
1-exp(-1e So 3,ds). So, for example, for o = Br = and D 1, we have
V (l))(( = eexp(-t), bDDS(Z, t) = x + 24 exp(p(((11 t)).
(25)
Corollary 2. For PIS, we have a reference process which is is pinned Brownian motion running
buckwards in time with μt = 0, best = (T I)I, aj = T by = (T=)) For = 1, obtain
T
(p1s(E, t) = -aa + b)) (x aT- +
(26)
In particular bris (x., 1) (x- 11) as 1. + 0.
13

---

Published as a conference paper at ICLR 2023
Hence for PIS the drift function explodes close the origin compared to DDS. This explosion holds for any target и as it is only related to the fact that pt concentrates to óg as - T. This makes harder to approximate for neural network. Additionally when discretizing the resulting diffusion, this means that smaller discretization steps should be used close to
A.3 PROOF OF PROPOSITION 2
Proof. Let us express the discrete time version of the reference process as
K-1
ndt (yo.k) = pb =p)=(wo) (wo)
II
ref Pk+
(27)
k=
and denote by pk the corresponding marginal density of yk which satisfies
Pk+1 (yk+1) =
PR+11R (16+1198)
(28)
The backward decomposition of this joint distribution is given by
(yo.k) PK (yK)
II R=
Pklk+1 (uk(yA++1)
(29)
where
pef Pk+1|8 (M+++ 1y (P)))
Ph 1++ (3x1BM++11)
(30)
Pk+1 (3k++1)
If our chosen integrator induces a transition kernel Ph+ 1/A (Me (18)) which is such that pk (uk)
=
pb (μK), then
T(KK)
K-1
pre (yo:k) Po (yK)
ग (BK)
II
Pkjk+1 (yn (8++1)
(31)
k=
is a valid (normalised) probability density. Hence it follows that
In (o((v))
KL(g" (Jo:k (vo.k)
(32)
x((x)
T TYKK) Pi (UK)
=KL
(v) ( (Ho.k)
2 0.
If the integrator does not preserve the marginals we have that
=(KK)
(UK) T TyK)
K-1
ref (Jo:K) PU (UK) (uk)
II
ref Pkik+1 (144)))11)
(33)
This not a probability density and thus objective is no longer guaranteed to be positive and consequently the expectation of our estimator of In Z will not be necessarily a lower bound for Z. Finally simple calculations show that the Euler discretisation does not preserve the invariant distribution of pref in DDS.
A.4
PROOF OF PROPOSITION 3
We have
In (vo:k) In dd (vo.Rc) p(yo:x) (xo:k) by construction,
+ N10000,0011
(34)
T(Y))
Now
we have from /I aukkK-k 02a and from
17 that pff (VK- -+1 UK (K
we obtain dh-1 (WK-6+1VKKk)
N(yK vI Ondi- 20 ork oxk) So (k, UK 1), o0
14

---

Published as a conference paper at ICLR 2023
It follows that
q (Ho:k)
In
(35)
(vo-
9x (vo)
= In
In
P (
k+1 WK pei 15-11x -
k=1
yg 2a 2
v1 CkyKK-K -xl?
1yKK++1 VI - - 20 (1- vI f yK yk)
where we have exploited the fact that gK (yo) 1o (vo) = (Noo0 0, 0².²
Now using ovao
- --------------a -oaKK- - 20 - vI 11--0)) cxa) (k; UK UK- we can rewrite
35 as
d (go.k) pro (yo:k)
In
(36)
K
[llek + 20 (1-v1 (1-vI-ar)
(ki yK -k]11 -||x111
2 x=11
VOk
vI aa..
K
VI atk)
20
|lfo(k, yK +1 +2
fo VK DK- ck.
O'k
Vak
1=1
k=
follows directly from and (66). Note we have for 8 << that 20 21--v--au
Now 19)
0² 3 and 2o (I=VT-a)
25 as expected from
Vaa
Finally the final expression 20 of the KL follows now from the fact that Eg (fo yK -x) Ten]
Eg Egg (ko =k) c =
Ege (fo (ki =x) Eg [cx bullk=lll 0.
A.5 ALTERNATIVE KULLBACK LEIBLER DECOMPOSITION
KL decomposition similar in spirit to the one developed for DDPM (Ho al. 2020) can be be derived. It leverages the fact that
Ph Pk- (1
N N(xai (o), 02
(37)
for ja 20)
1-6k
2o
oor (1-00-))
1-ût
lks
3x
1-āx-1
1-ax
BR.
Proposition 4. The reverse Kullback-Leiie discrepancy KL (pp) satisfies
KL(g"llp) E,g [[x((qg((=))llxxxxeekk)) KL KL 2n / 2k))
+ (20l
N( Bexx + o2 xx), 0² the terms (T,(f (#x- > (xx- xK)) KL between Gaussian distributions and can be calculated analytically.
So for (
We found this decomposition to numerically unstable and prone to diverging in our reverse KL setting.
15

---

Published as a conference paper at ICLR 2023
Proof.
The reverse KL can be decomposed follows
In o (T0:K)
p(ro.)) ((oox))) Eip In (xoo))
KL(G"llp)
Egg
Egp In 9 (2
p(011xx<00) )
Egg [[m])
L(0) 1e [xo)]
where
((0) = Egg In
o (((c)rk))) p( (p(01xx1000 1200)
Egg [[gx((k)) +
Σ In th (= (
K
PkoA- (aTA 1 ]
Now using the identity for k 2 2
Pk ( xx = Pa P 10 (x1)
= Ph (1= 1 2a Pa (Ix (xo);
we can rewrite L(D) as
l (0k-1
Pk- P In
2012 (Io (50021)))
L( = Eg0 In qu (xx) +
Σ
In
Pk- 1j0, (IX- 1 2z2)
PI (2) (zo)
2=
= Egg In qK (xx)
K
th-1jk (I 1
PK al (x (200) (20)
+
In
110, Th ( 100 xA
+ 60p1 (20
A=
K
= E,g [KL( (IX) (2)))))()))))))))))))))()))))))) +
KL 111 (xx 1 ( xo, xx))
Σ
k=
+ In 661 (Iol x1)
The result now follows directly.
B UNDERDAMPED LANGEVIN DYNAMICS
In the generative modeling context, it has been proposed to extend the original state E IRd by mo- mentum variable m E Rd One then diffuses in this extended space using an underdamped Langevin dynamics (Dockhorn et al. 2022) targeting Wax; 0, o² I [m:((;;)))))) 0, MM. It was demonstrated empir- ically that the resulting scores one needs to estimate are smoother and this led to improved per- formance. We adapt this approach to Monte Carlo sampling. This adaptation is trivial and in particular requires to design carefully numerical integrators.
B. CONTINUOUS TIME
We now consider an augmented target distribution π (x:) 0, .)) where M positive defi- nite mass matrix. We then diffuse this extended target using the following underdamped Langevin dynamics, i.e.
dr, = M 11 M.mm,,. m,d,
(38)
dm,
3t
dt 3.mdt V28,M22
where Ey ~ π. mg ~ (0, M). The resulting path measure on [0. 7] is denoted again P. From Haussmann Pardoux 1986), the time-reversal process is also diffusion satisfying
dy: - M In,t,
(39)
dn, 2 dl + mnd + -.8T,MV Iump- mild +
16

---

Published as a conference paper at ICLR 2023
for (Yo, nu) "y where Ni denotes the density of (Iz, m²) under 388
Now consider a reference process per on [0, I] defined by the forward process 38 initial- ized using zo (0, o²), mo (N0, In this case one can check that nt 21 me) N (z,;;0 0, o2 0, M) and the time-reversal process per satisfies
dy: = - n,dt,
(40)
Yt dt + Br-en 28r-inidt + 1
dng =
o²
Y1 dt mndt + - r²
as vn (y, v, (nVW 0, (m -M-In. Hence follows that the time-reversal 39 of can be also be written as
dy: - -
(41)
V1 dt nnt 2 In lnor-ely,, ne)dt /2 dn. r²
where de m²) == n,(1,, m)/mt (x, (x,m).
To approximate P, consider parameterized path measure Q" whose time reversal is defined for (yo, no) ~ (No: 0, ()No O, M) by dy = -M 'nddt and
1/1 dt dr + Mo 1,, mi)dt + 2Br dn. Φ²
(42)
B.2 LEARNING THE TIME-REVERSAL THROUGH KL MINIMIZATION
To approximate P, we will consider parameterized diffusion whose time reversal is defined by
dy: = --M- -M- -n,tt,
(43)
yt dt dnd + fo 1, m4)dt / M22 dm; σ²
for (Vo, na) (yo: O, (mo: 0, M) inducing a path measure @ on the time interval [0, 7).
We can now compute the Radon-Nikodys derivative between Q and pref using an extension of Girsanov theorem (Theorem A.3. in Sottinnn & Särkkä (2008))
dº dpiel
1/" V2 (T t, ne) dw, (" ||3T- - Ve, ny ne)l dt
In
(44)
To summarize, we have the following proposition.
Proposition 5. The Radon-Nikodym derivative
d apa (y no.r) satisfies under Q"
d de -1[-
[[Br Br- Br t. ve, dt) I' /2 fo(T m)) TdW. (45)
In
From KL(C"1 pIr + Eqo
In
(ir ((r)))))
] it follows that no holv,, .T)
t. KL(O")|P) Eg U. By mm d1 (
(Nlur:0,
(46)
The second term on the r.h.s. of 446 follows from the fact that In(p5 (yr, nT) po(yp, In(V( o² TTTT))((T)).
17

---

Tt dt B.MV m. In (1, me)dt, Published as a conference paper at ICLR 2023
B.3 NORMALIZING FLOW THROUGH ORDINARY DIFFERENTIAL EQUATION
The following ODE gives exactly the same marginals 1/r as the SDE 38 defining P
dæ, 1 m,dt,
(47)
dm, x²
dl Benydl B.MV my In (1 mi)dl,
Thus if we integrate its time reversal from 0 to T starting from (yo. 90) Tpr
dy -M
- n,dt,
(48)
dn; o²
dt BT-MMV vn In In ne)dt,
then we would obtain at time T sample (UT, (T) x (yT) ( 0, M).
In practice, this suggests that once we have learned an approximation so- (1, xr m) Vn In th m) by minimization of the ELBO then we can construct proposal using
dy: = --M-1nndh
(49)
Yt dt + - (T - 9, ny) dn,- o²
for (yo, no) (No: 0, 2T) 0, .). The resulting sample (yr, nr) 7j will have distribu- tion close to т x (0, M). Again it possible to compute pointwise the distribution of this sample to perform importance sampling correction.
B.4 FROM CONTINUOUS TIME TO DISCRETE TIME
We now need to come up with discrete-time integrator for the rimerreversal of pre given by 440 and the time-reversal of Q" given by (43 Let us start with (40). We split it into the two components
dy = - M nn,
dn, = o2
Yt
dt,
(50)
and
dy:= = 0. dn, - - + 2Br M
(51)
We will compose these transitions. To obtain (2h+1, 1, nR+1) from (y,, næ), we first integrate (50) To do this, consider the Hamiltonian equation dy, = M M-'nndt dn, = - preserves N(y; 0, o.o22))((1 0, .M) as invariant distribution. We can integrate this ODE exactly interval of length 5 and denote its solution ((y,n): see Section B.5 for details. We use its inverse --1(g ny. = Pmp n) o D n.) m (y. n) where Pmip (v, n) = (y,-n) so that (yk+1, n(,) ) = (yk, nx). Then we integrate exactly (51) using
Wt dt which
n /1 al-ak-urg + Vax-IM
(52)
We have thus design transition of the form td (y 11+1, n' |yke, nx)) = (g
(())(nx+1 /1 alkk--a-ktte -ak-kkkiii a ak- -.).
(53)
Now to integrate (43). split it in three parts. We first integrate (50) using exactly the same integrator (1k+1, Tt Φ "yk, n.g. Then we integrate
dn. 2 t., Vt, nt)dt
(54)
using
n = n/e + 2(1 /1 1 1a-aa--- M - k, yh+1 nk,,,
(55)
where we abuse notation and write fo(K - k, yk+1, n',) instead of fo ((K(- k)8, yk+ n4). Finally, we integrate the OU part using (52) but replacing n/ by with n². So the final transition we get
q" (yk M441, n/e nll Uk. n4) = S (1n+1, n +2(1 1ea1--1-an- -- (KK- (11)
x v -1aa----aae -n ax
However, we can integrate ni analytically to obtain
(y N441, nilVe, v2)) = 6 (x) (u-
ok xN(ne415 V1 ox ax 2(1 ox- So(K -k. -k (56)
18

---

Published as a conference paper at ICLR 2023
B.5 EXACT SOLUTION OF HAMILTONIAN EQUATION
Consider the Hamiltonian equation and M = I, the solution of dæ. = m,dt, dm,
2t dt
(57)
r²
can be exactly written as (m), (2 defined through It zng cos(l/o) mno sin(l/o), Mn sin(l/a) mg (1/0)
zo
(58)
This is the so-called harmonic oscillator.
Now the inverse of the Hamiltonian flow satisfies 4 (z, m) = mmp (I, m) o (xx,, m)) Pmp (z, m)
(see e.g. Leimkuhler & Matthews (2016) so we can simply integrate
y/ dt,
dy: n,tt,
dn.
(59)
using (y++1, m() $7 "uk, m4).
This gives
yk+1 = yk cos((To) - ngo sin(T/o)
(60)
n4g = (- yk sin(t/oO (T/o) ne (ts(t/o))
y) sin(T/o) + n& cos(t/o). σ
(61)
Obviously, if we have M oz then
me dl,
3t dl, dm,
dæ.
(62)
a
0²
then we use reparameterizatiio m m4/a and
2t dt.
dæ; = m,
dm
(63)
o²
and the solution is as above with being replaced by 0..a. So for ci = ng a. we have
yk+1 = yk cascoos////
- n - (iaesas()a-
sin(t/ + nu cos(t/ cos(t/(Vao)) (Voo)), 744 Vao
yk
so as R2 = ang and similarly níg = añk we finally obtain
Uk+1 yA cos(T cos((T((aaa)) n Va sin(T/(Vaa))
(64)
n1 Va yn sin(t/ +
(65)
B.6 DERIVATION OF UNDERDAMPENED KL IN DISCRETE TIME
In this section we derive the discrete-time KL for the underdamped noising dynamics.
Proposition 6. The log density ratio Ir ( (yo:k, no:k, 16.K) (/o no.k mk,k)) for Yo:K, Po:k 1o.K (-) equals
K
Ek fo yK -k-1, 1 TeA [ hb, yK -k11 R W + 112
2
(66)
Vou
1=1
VI an VI ak) and ex is obtained as a function of n441 and ni when describing the integrator for Qº. In particular, have Ek idd. N(00 I) and one obtains
where KR
K
KL(g"!lp) = =2E,0
[1 yK 90 k+1 n - + (N(ok:0, x(Ux)
(67)
ag
19

---

Published as a conference paper at ICLR 2023
The integrator has been designed that the ratio between the transitions of the proposal (56) and the reference process so cancel and
53 is well-defined as the deterministic parts are identical in the two transitions
o" (yk+1 1k+1 m/g lkk 12R)
(68)
(y6+1, Nk+1, y/s 11ps))
(m VI- a (n/ + + vI- MKok) k. 70)111)))):a---40
(n u MM)
The calculations to compute the Radon-Nikodyn derivative are now very similar to what we did in the proof of Proposition
We have
(yo.K, no:k, mo:k) In
(69)
prd (1o:K TAD:K, mlu.))
4K (x0. no) In
K
In
PK (yo, no)
Pk
[ k+1
-
Inkk++1 an(n' + 2(1 x
where we have exploited the fact that 9x (Yo, no) po no) NNo: 0, 0² 11:1NN(())) 0, M)
Now let us introduce
M-
Cl:
(nk -+1 VI
-21 hat
(70)
Va
hence
M 1/2
(nk T4K- VI -vocnnkkk)))
Vak
21-1/2
=6 +
VI ax
(K
Vak
2//1/2
=c +
VI ax vI fo yk -
Vak
We can rewrite (69)
as
qº (yo.k nork v/xx)
In
(71)
(yo:K, no:k. , mo:KK)
K
2/2
|lea +
- -1- o
||eu|]?
Vou
[l,
UK- -kk (1-0
(1 (1-0
4M
2
4M1
VI- a1
vI ax) UK ktl, n'x k)) TK
Vak
Σ fo (k, y
213e1 /I ap)
2/ ar vI Oik
fo (k, yk -+t n x))
a4k
where we note that TT-mnl
v TT
B 0ra/2 / ~ BasÓ.
20

---

Published as a conference paper at ICLR 2023
C EXPERIMENTS
In this section incorporate additional experiments and ablations as well as further experi- mental detail. For both PIS and DDS we created grid with d K T 0.05 and values of T {3.4, 6.4, 12.8, 25.6) and corresponding number of steps K E {64, 128, 128, 256). When tuning PIS we explore different values of which in practice, when discretised, amounts to the same effect changing For DD$ we apply the cosine schedule directly on the number of steps however wC ensure that Ee Cr Qmax such that we have similar scaling in PIS. Finally both PIS and DDS where trained with Adam (Kingma Ba 2015) at most 11000 iterations, although in most cases most both converged in less than 6000 iterations. Across all experiments used sam- ple size 300 to estimate the ELBO at train time and of 2000 to for the reported IS estimator of Z for all methods.
Finally for the ground truth where available we use the true normalizing constant as is the case with the Funnel distribution and the Normalising Flows used in the NICE Dinh et al. (2014) task. For the rest of the tasks we follow prior work (Zhang & Chen 2022; Arbel et al. 2021) use a long run SMC chain (1000 temperatures, 30 seeds with 2000 samples each).
C.1 NETWORK PARAMETRISATION
In order to improve numerical stability we re-scale the neural network drift parametrisation as fol- lows:
folk,y) = 2- fo(k. = NN, y: 0) (N: O Vln((y).
This done such that the reciprocal term oK in the DDS objective is cancelled as the term reaches very small values in the boundaries causing the overall objective to be large. Finally as
as it follows that our proposed re-parametrisation converges the same A 11------ ag 2
SDE as before, but now with stabler and simpler updates:
1k+1,n = 11-11aaaea a + +o'ax-nFolK' k, yn,n) + aoaa EK, iid N(O.1),
(72)
Tk+1, Tk, +2 2 o'ak-=lllo(KK hi 11?
(73)
C.2 TUNING HYPER-PARAMETES
For both DDS and PIS we explore grid of 25 hyper-parameter values. For PIS search for the best performing value of the volatility coefficient over 25 values, depending on the task we vary the end points of the drift however we noticed that PIS leads numerical instabilities for y >> 4 thus we never search for values larger than this. For DDS we searched across 5 values for and values for amaxx this led total of 25 combinations we explored for each experiment. Finally for SMC we searched over 3 of its different step sizes leading to total of 5 125.
C.3 BROWNIAN MOTION MODEL WITH GAUSSIAN OBSERVATION NOISE
The model is given by
Qinn ~ LogNormal(0., 2),
Globa LogNormal(0, 2),
31 ~ NO0, ainn):
2; N(I;-11 Qiun),
i= ...2... i= 1,....333 30.
~ N(x;, Grbb);
The goal is to perform inference over the variables ainns aabs and {ar given the observations (v:) 10 (22))020
21

---

Published as a conference paper at ICLR 2023
PIS drift magnitude
DDS drift magnitude
30
25
20
0
10 20
30
40
50
60
Step Number
Figure 4: Magnitude of the learnt neural net approximation of the drift V, In (x)) (see (23
function of
Training Times per Iteration seconds (LGGC)
64
128
256
512
DDS
0.104 + 0.0004
0.205
0.0004
0.406
0.0004
0.809 0.0004
PIS
0.104 = 0.0004
0.205 + 0.0003
0.406 + 0.0005
0.808 0.0006
Table 2: Training time per ELBO gradient update on 300 samples.
C.4 DRIFT MAGNITUDE
In Figure 4 we compare the magnitudes of the NN drifts learned by PIS and DDS scalar target A (6, o²) target. The drifts where sampled on randomly evaluated trajectories from each of the learned samplers. We observe the PIS drift to have large magnitudes with higher variance more prone to numerical instabilities, whilst the OU drift has a notably lower magnitude with significantly less variance. We conjecture this is the reason that PIS becomes unstable training time for large number of steps across several of our PIS experiments, many of which did not converge due numerical instabilities. For simplicity in this experiment we set ag be uniform.
C.5 TRAINING TIME
We evaluate the training times of DDS and PIS on an 8-chip TPU circuit and average over runs. Table 2 shows running times training iteration. In total we trained for 11000 iterations thus the total training times are 19 minutes; 37 minutes; hour and 14 minutes; hours and 26 minutes for I 64, 128, 256, 512 respectively.
C.6 VARIATIONAL. ENCODER IN LATENT SPACE
Following Arbel et al. (2021) we explore estimating the posterior of pre-trained Variational Auto Encoder (VAE) (Kingma Welling 2013). We found all approaches to reach a very small error for In Z very quickly however as seen with PIS across experiments results became unstable for large T and thus many iypaeepasameters when tuning the volatility for PIS led to the loss exploding.
C.7 GENERATED IMAGES
In this section we provide some of the generated samples for the normalising flow evaluation task. In Figure 0 we observe how both PIS and DDS able to generate more image classes than SMC due to mode collapse. Out the 3 approaches we can see that DDS mode collapses the least whilst SMC generates the higher quality images. Using neural network with inductive biases such
22

---

Published as a conference paper at ICLR 2023
-100
InZ Eztemtatee
-1105
-111.5
D05
R5
-12200
SMC
54
125
256
512
Number of Stapp
Figure 5: Results on pretrained VAE from Arbel et al. (2021)
63
X-
9
SV
D
SS
63
AQ
Drg
Figure 6: a) DDS, b) PIS c) SMC. We can see how SMC mode collapses significantly as it only seems have 5 image classes. For PIS identified classes whilst for DDS we found out of 10.
as a convolutional neural network (LeCun et al. 1998) should improve the image quality over the simple feed-forward networks we have used, believe this is the reason behind the low quality in the images as both minimising the path KI. and sampling high quality images a difficult task for small feed-forward network.
C.8 COMPARISON TO AIS AND MCD
We compare our results to the AIS and MCD baselines presented in Geffner & Domke (022)) Results can be seen in Tables and 3] we can see that both DDS and PIS outperform LDV across all number steps and tasks with DDS outperforming PIS on higher values of due to inconssitent and unstable behaviour presented by PIS.
C.9 EULER-MAYURAND VS EXPONENTIAL INTEGRATOR
In this section we demonstrate empirically how using the Euler Mayurama integrator directly on the DDS objective leads to overestimating In Z. We compare on the Funnel dataset for which we know In Z be and the lonosphere dataset which is small and for which we have a reliable estimate of In Z via long run SMC. Results in Table show case how for lower values of & 64, 128
Tonosphere
Sonar
ULA
MCD
UIIA
LDVI
DDS
PIS
ULA
MCD
UIIA
LDVI
DDS
PIS
64
-113.8
-112.5
-112.8
-112.1
-11.7
-111.6
-115.3
-I11I1
-111.9
-109.7
-109.4
-18.99
:: = 128
-113.1
-12.2
-11.3
-111.9
-1116
-112.0
-113.5
-1102
-110.6
-109.1
-108.9
-109.3
b 256
-112.7
-112.1
-12.11
-1117
-1116
-113.0
-112.1
-109.7
-109.7
-108.
-108.9
-108.9
Table 3: Comparison of DDS and PIS to AIS based approaches on logistic regression models.
23

---

Published as a conference paper at ICLR 2023
Brownian Motion
ULA
MCD
UHA
LDVI
DDS
PIS
k: 64
-0.7
0.1
0.1
-0.5
1.0
1.0
4: 128
-0.3
0.2
0.4
0.7
1.1
1.0
k 256
-0.1
0.5
0.6
0.9
1.0
0.7
Table 4: Comparison of DDS and PIS AIS based approaches on Brownian motion time series target.
Funnel
Ionosphere
Exponential Integrator
Mayurama
Ground Truth Exponential Integrator
Euler Mayurama
Ground Truth
0.2006
0.059
2.425 0.266
-111.093 + 0.169
-106.364
0.371
k 128
-0.176
0.099
2.011 =: 0.532
-111.587 1: 0.136
-111.158 1: 0.399
00.221
0
-111.575
-112:271
-111.560
k = 512
-0.176 + 0.068
-0.154 0.066
- -111.582 = 0.136
-112.912 + 0.353
Table 5: Comparing EM Exponential integrators on Funnel and Ionosphere datasets. We can how EM significantly overestimates In Z.
the -uler-Maruyama (EM) based approach significantly overestimates In Z and overall does not perform well.
C.9.1 DETACHED GRADIENT ABLATIONS
In this section we perform an ablation over our proposed modification of the PIS-Grad network architecture. We compare the same architecture with and without the gradient attached. We found that detaching gradient led to a favourable performance as well as more stable training. Results are presented in table and were carried out for K = 128 and using the best found diffusion coefficient from each task.
We chose to explore this feature it is known that optimization through unrolled computation graphs can chaotic and introduce exploding/vanishing gradients (Parmas et al.) 2018 Metz al.) 2020) which work (Greff et al.) al., 2019) where detached scores are provided as features to the recogmition networks. Alternative approaches exist based smoothing the chaotic loss-landscape (Vicol et al., 2021 but we leave this for future work.
numerically introduce bias. This has been successfully applied in related prior
C.10 COSINE BASED DECAY FOR PIS
Unlike DDS, with PIS it is not immediately clear how to schedule the step-sizes in similar manner to Nichol & Dhariwal 2022). The simplest approach would be to schedule 8 (the discretisation step). We perform ablation with the cosine schedule
the discretisation step with PIS. Results are displayed Table C.T0, we can see that overall the cosine discretisation decreases or preserves the performance of PIS with the exception of the Sonar data-set in which it increases. Whilst there may be a way to improve PIS with a bespoke discretisation we were unable to find such this work. Additionally
from Table C0 how DDS improves significantly across every task when using the cosine squared schedule compared to uniform.
C.11 PROBABILITY FLOW ODE
Figures 7 ladd X show samples obtained from the probability flow ODE of a trained DDS model. We can see that the probability flow ODE is able to perfectly sample from uni-modal Gaussian, matters became more challenging with multi-modal distributions. Initially we discretised the ODE (2) using the same type of integrators as the SDEs to obtain Wk+1 = y/k + - /I ax-x)fo( fo (k) for Yo (0; I21). Unfortunately under this discretisation ODE to become stiff with more complex distributions such as the mixture of Gaussians, resulting in strange effects in the samples (matching the modes but completely wrong shapes). By using Heun
found the probability flow
24

---

Published as a conference paper at ICLR 2023
Methed
Funnel
LGCP
Somar
lonosplacee
Brownian
NICE
PIS Cos decay
-0.256
0.084
501 0000
0.881
-110.025 d: 0.0720
-108.843 1: 00.2:2
-111.618
0.127
1,035 0.097
-4.212 22: 0.620
PIS miiti
-0.268
0.07
502
002
1: 0.957
-110.025 1. 0.070
-109.319 1: 0.356
-111.002
0.145
0,580 1: 145
-3.965
1: 0.3310
DDS
-0.176 0.008
497 944 0.994
-110.012
8.071
-108 9033 226
-111.587
0.136
1,477 0. 156
-3.400
-1 0.566
-0.275
0.127
1.338
-10.246 0.188
-100.356
0..66
0.161
0.475 0.166
-7.047
0.681
Table 6: Ablation for cos scheduling.
Method
Pannel
LGCP
VAE
Sonar
lonoosphere
Brownian
NICE
PIS Grad Detuch
-0 268 1 0.07
502 002
0.957
-10,025
0.0701
-109.349
0.356
-11 602 1.454
0.900 0.145
-3.965 1. 0.340
Grad
273 0.798
-100033-
-109.577
-5.703 0.490
Table 7: Results for PIS grad with and without detaching the gradient.
integrator as proposed in Karras al. 2022) we were able to obtain better results, however also found we needed to increase the network size to improve the results as well as train with a learning rate decay of 0.99. With these extra features we were able simulate a probability flow ODE that roughly matched the marginal densities of the SDE. However even with these fixes we still found the probability flow based estimator of In Z to drastically overestimate, this is not entirely surprising as in discrete time the expectation of this estimator is not guaranteed to be an ELBO.
C.12 UNDERDAMPED OU RESULTS
We perform some additional experiments with the underdamped OU reference process. Similar to the damped setting parametrise the network as:
(k x, p) NN1 (k. x. p; + NN2( In T(x).
(74)
Unfortunately, unlike in DDS and PIS we cannot directly aid the update/proposal for with V In (z), thus this naive parametrisation is not the ideal inductive bias.
Results in Figures 9 and 10 show some experiments for this approach. The results perform worse than DDS for a standard OU reference process and PIS. However they are still better than VI. We believe that future work exploring better inductive biases for the network fo fk....,)) could narrow the performance gap in this approach. Additionally we highlight that in Figure | (Funnel Distribu- tion) we can see the underdamped approach overestimates In Z for k 64. We believe this may be due to numerical error when generating a sample.
C. BOx PLOTS FOR MAIN RESULTS
For completeness in this section we also present our results via box plots be seen in Figures and [2]
Figure 7: a) Target distribution (MoG), b) DDS SDE samples, e) trained probability flow ODE.
25

---

Published as a conference paper at ICLR 2023
Figure 8: a) Target distribution (66,, 1). b) DDS SDE c) Trained DDS probability flow ODE.
D DENOISING DIFFUSION SAMPLERS IN DISCRETE-TIME: FURTHER DETAILS
D.1 FORWARD PROCESS AND ITS TIME REVERSAL
In discrete-time, we consider the following discrete-time "forward" Markov process
K
p(zo.k) = (
PkJ ( k=1
II
(75)
we use Paja (I = Wa v1 axtx- a The parameters (o are such that (*(*x) = W( 022
Following
17
K
We propose here to sample approximately from by approximating the ancestral sampling scheme for 75) corresponding to the backward decomposition
K
Pk- 11-11 (1x- with Pk 1lk (#x- (xx) Pk- Pk-1( 1)P
p( pK (xx) II
Pk (xrk
1=1
(76)
where po (xo) (xro). we could thus sample FK ~ pK(.) then Tk-1 ~ Pk- (1 (-12a) for
K K,..., then 2i would indeed be sample from æ. However as the marginal densities (p&) available in closed-form, we cannot implement exactly this ancestral sampling procedure.
As we have by design (x( ss N(=k;;, O²!²), we can simply initialize the ancestral sam- pling procedure by Gaussian sample. The approximation of the backward Markov kemnlss Pa ila (#x-11 (ma) is however more involved.
D.2 REFERENCE PROCESS AND BELLMAN RECURSION
We introduce the "simple" reference process o (xo:k) defined by
K
pd (xo.k) = Nwoo 0, o² (T((EE))))) Noxo 0, 21) II
Pkjk P
(77)
which is designed by construction to admits the marginal distributions p (x+) = (axx 0, 02 for all can be easily verified using the chain rule for KI. that the extended target process p is the distribution minimizing the reverse (or forward) KL discrepancy w.r.t. of over the set of path measures q(ZO:K) with marginal qo (xo) = (xx)) at the initial time, i.e.
arg min KL 40
The intractable backward Markov densities Pы-1jA one wants to approximate can also be rewritten pe as twisted versions of the tractable backward Markov densities p-- by some value functions (oh )
26

---

Published as a conference paper at ICLR 2023
+233
SBC
Figure 9: I estimate function of number of steps K Funnel, b) LGCP, c) Logistic lonosphere dataset. Yellow Dotted line MF-VI and dashed magenta is the gold standard.
-3915
10:
-304
128
355
212
512
Figure 10: In Z estimate as function of number of steps a) Logistic Sonar dataset, b) Brownian motion, NICE. Yellow dotted line is MF-VI and dashed magenta is the gold standard.
Proposition 7. We have for k = K
De-1]6(2a-11000 Φx-1 1 D (#R=1
for Φb(xK) = Pk ()(()))
(78)
Φ (x*)
pk p (xx)
The value functions (dxx k= satisfy a forward Bellman type equation
T(ZO)
δx(IK)) =
@o(xo))
(79)
Po (cr)
follows that
=f f (o((z) Pojk (2ol
0x (xx)
(80)
Proof. To establish
778.. we use Bayes' rule
Pk (Ir PPNj-- 1
Pk 1|k (I
(x(x
Φ (x 1 )pk_ (I 1)ppkk 1 xx-1)
@u (I )p (IR-1)
2 (r )pr 1/1 (0x (xnk
(xR)
where we have used the fact that Pa (1xx) the (xxx )p (xx)) poot (x ( = P
and Bayes' rule again. Now we have f P (Jo 1| [ )dx) = 1 for any wa. so follows directly from the expression of this transition kernel that the value function satisfies
Tx0))
(x( 0x-1(xx-1)PF
@o(xo)
(81)
p (xo))
By iterating this recursion, follows that
Φn (x(x)) = ff dol(a)) pola (2ol
(82)
27

---

Published as a conference paper at ICLR 2023
-233
-1.0
MCC
204
PMC
24
Figure 11: In Z estimate function of number of steps K a) Funnel, b) LGCP, Logistic lonosphere dataset. Yellow Dotted line MF-VI and dashed magenta is the gold standard.
-1225
-1875
4322
4324
D33
BMC
Figure 12: In Z estimate as function of number of steps K a) Logistic Sonar dataset, b) Brownian motion, c) NICE. Yellow dotted line is MF-VI and dashed magenta is the gold standard.
D.3 APPROXIMATING THE BACKWARD KERNELS
We want approximate the backward Markov transitions Pk-1|K- From 78 this would require
only to approximate the value functions Monte Carlo estimates based on (80) are high variance but also sample from the resulting twisted kernel which is difficult. However, if we select 83 x O, then we have on (x) se (x) and by Taylor expansion of V In db. around vx we obtain
Pk- Pk-1/A (3k-1 lan)
P (Ix- (xx)) exp(In 1 (xx In (xx))
~ N vI axox + o² - vI as ) In
(83)
By differentiating the identity
80 w.r.t. Ik. we obtain the following expression for de
Vln /
V Finp (xol(a)) Paja (xol=x)d
(84)
Alternatively, we can also use Vn = VInpA(<x) + = where (In) can be easily shown satisfy
V (Inp( v Pe Pojk (xp Iy dx).
(85)
For DDPM, the conditional expectation in (85) is reformulated as the solution to a regression prob- lem, leveraging the fact that we can easily obtain samples from (z
æ is the data distribution in this context. In the Monte Carlo sampling context considered here, we could sample approximately from Pou (org(Tro [zx)) oc T(xu)Ps o (2x//20)
xx)) using MCMC so to approximate the condi- tional expectations in (84) and (85) but this would defeat the purpose of the proposed methodology.
We propose instead to approximate V In dhs by minimizing a suitable criterion.
In practice, we will consider distribution gdt (xo:k) approximating p(Fo: of the form
(x(u)) N(xx: o² II gi 1|2 (
(86)
Le. select (xx) (x(x)) 0, 0 as (k( M, O2I. We want ( approximate Pk-1je Pe- (II-1 (x)), i.c. fo(k, x VIn Inspired by (8).
will consider an approximation of the form
(1x- NN auea + t
ofas))
(87)
where fo(k; 2R) is differentiable w.r.t. 0.
28

---

Published as a conference paper at ICLR 2023
D.4 HYPERPARAMETERS
D.4.1 FITTED HYPERPARAMETERS
To aid reproducibility we report all fitted hyperparameters for each of our methods and PIS across all experiments in Tables 8 and
D.4.2 OPTIMISATION HYPERPARAMETERS
As mentioned in the experimental section across all experiments modulo the Funnel we use the Adam optimiser with a learning rate of 0.0001 with no learning decay and 11000 training iterations, for the rest of the optimisation parameters use the default settings provided by the Optax library (Hessel et al. 2020) which are b = 0.9, b2 0.999, = 10-8 naming as per Kingma & Ba 2015
From the github repository Zhang Chen (2022) we were only able find hyperparameters reported for the Funnel distribution. In order to first reproduce their results we used the a learning rate 0.005 and learning rate decay of 0.95 as per their implementation, their results were initially not reproducible due to bug in setting of 1 despite comparing methods at the less favourable values of o 3. For 0 = 1 we were able to reproduce their results. However we report results at Øf = 3 as this is the traditional value used for this loss. As other optimisation configuration files were reported we used the more conservative learning rate 0.0001 since PIS was very unstable for 0.005 with decay 0.95 across many of our tasks. Finally would like clarify that the exact same optimiser settings where used for both PIS and DDS in order to ensure a fair comparison.
D.4.3 DRIFT AND GRADIENT CLIPPING
We follow the same gradient clipping as in Zhang Chen (2022) that is
NN2(k; @ ceppppppplaaiaeeeoo 102,, -00)
-10*, 10a)
(88)
29

---

Published as a conference paper at ICLR 2023 3
0. 25.
1.65
0 0.666 = 3.85, m 0.0000
3.85, 0,,000
= 1.0
1.0.
255
0 == 0.6 = 3-8.88850
0 0.6 0.6.a 3.85,
1.0
UDMP
0.6,0
UDMP
05.75
11.00
PS
0 01253
0.00
0.:16
2.125
d= 13. 23
PIS
d
1.5
8 = 01688, a 1.463
0 0.688, a 0.688
nice
DDS
-03.a 1.075
a = 0..3 = 01688
163 48,m 22
Table 9: Fitted hyperparameters for sonar, brownian, vac and nice experiments.
0.118 48,m
DDS
Table 8: Fitted hyperparameters for funnel, lgep and ion experiments.
= 1.11.. a 25, m 0.4
14.0 2.5, 0.4
===11 1.4. = 4.5, m = 0.4
UDMP
1,7,0 4.5, m
0.000
0.0408
PIS
UDMP
2.55
1.8
0= 1.068
0.712
== 0.579
d 0.416
bnowniee
11.00
0.1.
DDS
PIS
UDMP NA
0 = 21, = 1.50
0.75
= 010000
1,500
NA
PIS 01253
1122521
a 0.466
1-221,0
= 2.1a
0 21,0
lgcp
DDS
0.80.. = 1870
0.81 22
0 = 1.8,0 = 1,67, m 0.9
188,a 3,7, 0.0
= 1.075 a 2. = 0.9
0,088.. a 3.7, m
DDS
Vie
1.7 22
9.1
UDMP
1.7,0
2.0 m
15
0 = 1.068
0.742
0 0,253
UDMP
PIS
251
a = 3,025, a = 1,075
0,0875
= 1.85 = 0.3
5 463, a 0.3
PIS
660
0.75
1,075,
Tunnel
DDS
013
0.3
DDS
= 64
128
K = 256
64
5/2
K
X
30


Reviews:
[{"review_id": "ehZtBRQknUr", "reviewer": "Reviewer_cpvn", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 1, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "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": "The paper describes sampling (MCMC parlance) as a diffusion problem. In the forward direction, the target distribution diffuses to Gaussian noise, while the reverse direction of interest allows diffusing from noise to the target distribution. As can be seen and is highlighted by the authors, this problem is classical in MCMC.\nThey develop a sampler which is inspired by recent advances in diffusion based models. Building a reference diffusion process that diffuses from noise to noise, they describe some particular criterion they optimize to train their sampler (eq 10), basically based on minimizing the KL divergence between the distribution it yields and the reference distribution.\nThe same story is then told again for different kinds of flavours, involving  augmenting the original state by a momentum variables and also studying the discrete case, which is of interest for implementations.\nExperiments are numerous and mostly done on small data, but apparently matching what is usually done in the field.", "summary_of_the_review": "Probably correct theoretical paper about the connections between diffusion models and MCMC. Hard to follow. It doesn't look like great care was taken to help the practitionner use the approach.", "strength_and_weaknesses": "The paper is definitely timely, considering the hype around diffusion models, and it reads pretty rigorous. I guess it may be of interest to the very few readers able to understand it.\nOn the downside, I must say that I had great difficulty following it, because my math is not strong enough. As a practitionner, I am left with some frustrating feeling that I probably don't really even know how to use that sampler in practice, the whole thing being quite burried in the maths.\nMy guess is: this paper is very nice for a strong theoretical connection between diffusion models and samplers, but I can't really assess it.\n\nSome attempts of comments on the go:\n* below (2), I guess you mean \\lambda_T, not \\alpha_T\n* maybe it is trivial, but I probably missed a proof for what is claimed in the first paragraph of 2.2\n* between (5) and (6), isn't there a minus sign \"-\" mismatch ?\n* \"where we have use\" -> \"used\"\n* It looks like only the equation below (8) $f_\\theta(t,x) = ....$ is necessary, not really (8), which looks like just a copy of (5)\n* the discussion in 2.4 and what it means as a difference with the section before is difficult for me to understand\n* above (14), I don't know what \\pi(x; 0, \\sigma^2 I) means. This section looks like some extension of section 2, and I must confess I just gave up following the overwhelming amount of maths at this point.\n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "2: The contributions are only marginally significant or novel."}, "extra_text": {}}, {"review_id": "UDR7XbIAKrm", "reviewer": "Reviewer_bg9B", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "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": "The paper introduces Denoising Diffusion Sampler. The task at hand is to sample from an un-normalized target density with an unknown normalization constant; and also to estimate said constant. Inspired by Denoising Diffusion generative models, the paper proposes a reverse-time SDE of a process that diffuses the target distribution to the Gaussian distribution. The authors propose an approximate learning objective so that to get rid of intractable scores of the marginals of the process. The authors evaluate their method on benchmarks, comparing it to other sampling methods.", "summary_of_the_review": "My issue with the paper lies within its empirical evaluation and comparative advantage of the proposed method as it relates to other sampling techniques; and also within the clarity of the presentation.\n\nUPDATE: after the discussion period, I am willing to increase the score.", "strength_and_weaknesses": "Strength: the method is connected to Denoising Diffusion generative models, which gained much popularity.\n\nWeakness: I have concerns regarding the empirical evaluation and the clarity of the presentation. The figures which are supposed to represent comparison against other sampling methods (mainly SMC and path integral sampler) are barely discernible. But it can be seen that the DDS seems to work on par with PIS. It is also not clear the training times that are mentioned to be non-negligible compared to SMC. Besides, the experiments seem to be relatively small or medium scale. \n\nRegarding the style of the presentation it is sometimes obscured by heavy-weight formulae, e.g. in section 3.2. I think proper formatting should remedy that. Tables do not seem to be referenced anywhere, and Figures 2, 3 neither. Sections' 3.1, 3.2 titles contain the same typo: \"Ulhenbeck\".", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "Sqq1xizeRXU", "reviewer": "Reviewer_Vb4W", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "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": "The present paper proposes a novel method to approximately sample from unnormalized density functions and estimate their normalizing factors. Motivated by recent successes of denoising diffusion models, the main idea is to reverse an Ornstein-Uhlenberg (OU) process which diffuses the target density into an approximate normal distribution. This requires an approximation to the score of this process which is learned using a neural network by minimizing the reverse Kullback-Leibler distance. Moreover, extensions to underdamped Langevin dynamics and the connection to Schrödinger brides and stochastic optimal control are presented. Finally, the performance of the method is demonstrated by several numerical experiments.", "summary_of_the_review": "The problem tackled in the paper is of high interest, e.g., in Bayesian statistics and computational sciences. Furthermore, it seems to be a promising direction to transfer methods from denoising diffusion models to improve methods in stochastic optimal control used for sampling from unnormalized densities. The paper provides a good starting point in terms of theoretical contributions, which, however, to a large extend directly follows from the theory developed for denoising diffusion models. The method used for the numerical experiments only differs from existing methods by using a different reference SDE (OU process vs. pinned Brownian motion). This indeed provides improved numerical stability and better estimates of normalizing constants for the considered examples. However, the presentation of the paper could be improved and a full picture of the benefits and drawbacks of the proposed method seems to be lacking. ", "strength_and_weaknesses": "**Strenghts:**\n\nThe paper shows how to transfer several ideas developed for denoising diffusion models to the task of sampling from unnormalized densities, where the denoising score matching objective cannot be applied. This represents an interesting research direction which is of high interest to practitioners. To this end, the paper derives a suitable objective for learning the drift of the reverse-time process using the Kullback-Leibler divergence on the path-space. Further, the forward process is extended to underdamped Langevin dynamics and the corresponding discrete-time analogous are presented.\nDifferent from previous methods in stochastic optimal control, this leads to a reference process given by an overdamped or underdamped OU process instead of a pinned Brownian motion, which improves the numerical stability of the training procedure (especially for larger number of steps), yields better estimates for the normalizing constants, and allows to construct reverse-time ODEs.\n\n**Weaknesses:**\n\nComparing the objectives of the proposed method and already existing methods (e.g., the so-called path integral sampler), there is basically only one main difference, i.e., the use of an OU-process instead of a pinned Brownian motion. Judging from the numerical experiments, this can improve training stability and estimation of normalizing constants. It remains unclear whether these benefits stem from the SDE type (variance preserving OU vs. variance exploding Brownian motion), which could also easily be adapted for the path integral sampler, or the initial Gaussian distribution instead of a fixed starting point. \n\nThe paper also seems to lack evidence whether these changes systematically improve sample quality (very minor relative improvement in Table 1 and hardly visible differences in Figure 6) and further experiments in this direction would be valuable. Moreover, there are other methods to sample from unnormalized densities using diffusion models, see, e.g., https://arxiv.org/abs/2206.01729 (Section 3.6. and 4.5), which could be compared to the proposed method and which propose to use effective sample size as additional evaluation metric.\n\nWhile the theory developed in the first sections of the paper is interesting to read, a couple of ideas actually do not benefit the numerical performance in practice, e.g., underdamped Langevin dynamics and the reverse-time ODE (continuous-time normalizing flow), see Appendix C.8 and C.9. I think this might distract the reader from the practically relevant parts and one could, for instance, move corresponding sections to the appendix.\n\nFinally, it would be great if the paper included a discussion why the exponential type integrator is considered in the discrete-time setting. There are better theoretical guarantees as compared to the Euler-Maruyama scheme, however, the latter seems to be used for the standard implementation of the path integral sampler. Also, it would be interesting if the adjoint method could still be employed to compute memory-efficient gradients for the objectives.\n\n**Minor issues:**\n1) The proofs of the statements and additional material in the appendix should be referenced in the main paper.\n2) 'Logarithmic derivatives of the intractable marginal densities': Should this rather be 'derivatives of the logarithms of the intractable marginal densities'?\n3) 'As its initial state $x$ is distributed according to $\\pi$': The problem is rather that we cannot sample from the density $\\pi$.\n4) The time-reversal is written as $(y_t)\\_{t \\in [0,T]}=(x_{T-t})_{t \\in [0,T]}$, which could be interpreted as an equality in an almost sure sense. However, I think that this only holds in distribution.\n5) The considered set of path measures should be properly defined.\n6) While the chain rule for KL divergences can easily be verified for time-discrete settings, it would be beneficial to include a reference for the time-continuous case in order to see what regularity assumptions are needed.\n7) Figure 1 could be described in more detail.\n8) I guess the symbol $\\perp$ in the beginning of Section 4 should denote the inverse. It feels counter-intuitive to detach the only part of the loss which contains information on the target density and it would be interesting to have more insights regarding this choice.\n9) It would be beneficial for the reader to provide more information on the problem settings and metrics.\n10) For some settings the presented box plots are difficult to compare and it would help the reader to provide, for instance, relative errors.\n11) Typos:\n    - 'descent' is missing after 'gradient' on page 2.\n    - 'the use *of* underdamped diffusions' on page 2.\n    - $\\alpha_T$ should probably be $\\lambda_T$ on page 2.\n    - $\\beta_{T-t}$ is missing in the inline equation before (5).\n    - The last sentence in Section 2, the sentence starting with 'For example' in Section 3.3, and the second sentence in Appendix C.4 seem to be grammatically incorrect.\n    - There is a missing full stop in the last paragraph on page 5.\n    - 'At' should be 'As' in the beginning of Section 3.2.\n    - Wrong references: 'Figures C.7 and C.7' in Appendix C.7. and 'Table C.7.1' in Section C.7.1. \n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "2: The contributions are only marginally significant or novel.", "empirical_novelty_and_significance": "2: The contributions are only marginally significant or novel."}, "extra_text": {}}, {"review_id": "K5quhU_nqzw", "reviewer": "Reviewer_SEx8", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "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": "Authors propose novel theoretical framework based on Girsanov formula to learn sampler based on reversing diffusion process from non-normalised density function.", "summary_of_the_review": "The proposed approach definitely interesting and I am eager to see whether it can be successfully applied for non-convex optimization. I believe that this result is worth publishing due to novel approach for sampling. ", "strength_and_weaknesses": "Strength: novel theoretical approach to sampling from non-normalised density\n\nWeakness: absent experimental comparison with Langevin Monte-Carlo methods as they aim to solve the same problem.", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}]

Author rebuttal:

Decision: Accept: poster