1 Introduction↩︎

Diffusion, or score-based generative models (DMs, SGMs) have shown remarkable performances in recent years [1]–[3], with applications ranging from computer vision [4], [5] and medicine [6] to natural language processing [7]. DMs aim to estimate a data distribution \(\mu\), with access only to a finite set of samples. They do this by constructing a stochastic process \((\overrightarrow{X}_t)_{t \in [0,T]}\) (the forward process), being the solution of a stochastic differential equation (SDE) over the time interval \([0, T]\) and initialized from the data distribution \(\mu\). Two classical instantiations of \((\overrightarrow{X}_t)_{t \in [0,T]}\) are either the \(d\)-dimensional Brownian motion or the Ornstein–Uhlenbeck process [8]. In our work, we adopt the latter, which admits the standard Gaussian, denoted by \(\gamma^d\), as stationary distribution. This construction defines a path measure connecting \(\mu\) to \(\gamma^d\), and thus an ideal generative model can be formed by taking a large value of \(T\) and considering the time-reversed (or backward) process associated with \((\overrightarrow{X}_t)_{t \in [0,T]}\), defined for any \(t \in [0,T]\) as \(\overleftarrow{X}_t := \overrightarrow{X}_{T - t}\). It can be shown that this too is a diffusion process, whose drift depends on the Stein scores \((t, x) \mapsto s(t, x)\) of the forward marginals [9], [10] which are unknown in practice. This score can be estimated using a family of neural networks \(s^{\theta} : (t, x) \mapsto s^{\theta}(t, x)\) parameterized by \(\theta \in \Theta\) [11] which can then be used to approximate the backward process with \((\overleftarrow{X}_t^{\theta})_{t\in [0,T]}\) using the network \(s^\theta\) in place of \(s\). In practice, the backward process is further approximated by initializing from \(\gamma^d\) and applying a numerical scheme to discretize the corresponding SDE. Formally, this is often framed as minimizing the Kullback-Leibler (KL) divergence between the data and model distributions.

Learning the score. Training and evaluating directly using the KL divergence is intractable as we do not have direct access to the density of the model distribution. To circumvent this, we typically rely on a tractable upper bound, given by, \[\label{eq:score-matching-bound} \text{\normalfont{{KL}}}\left(\mu | \overleftarrow{\mu}^{\theta}_T \right) \leqslant\frac{1}{4} \int_0^T \mathbb{E}_{} [ \varepsilon_t^\theta({\overrightarrow{X}}_t)^2 ] \text{\normalfont d}t + \text{\normalfont{{KL}}}\left(\overrightarrow{\mu}_T | \gamma^d \right),\tag{1}\] with \(\varepsilon_t^\theta(x)^2 := {\Vert s(t,x) - s^\theta(t,x) \Vert^2}\) and \(\overrightarrow{\mu}_T := \mathrm{Law}({\overrightarrow{X}}_T)\). The first term on the right-hand side is the score matching loss, defined as a weighted \(\mathrm{L}^2\) loss between the true and approximate score functions [11]. The second term is usually intractable, but does not depend on the learned score \(s^\theta\) and decays exponentially fast as \(T\) grows [12]. This loss function is further approximated using the denoising score matching loss, efficiently optimized using stochastic gradient methods [13], [14]. The theoretical justification for 1 is provided by [11], using a data processing inequality that bounds the KL between the marginals, \(\mu\) and \(\overleftarrow{\mu}^{\theta}_T\), by the KL between their path measures, which is further expressed via Girsanov’s theorem [15]. Complementing this perspective, [16]–[20] show that this bound can be understood as a variational evidence lower bound (ELBO) on the log-likelihood, grounding diffusion models within the variational inference framework.

The score-matching gap. Given that the score matching loss is an upper bound, with no guarantee of closely approximating the KL divergence, it might seem surprising that this has so effectively served as a theoretical foundation of the diffusion model framework. The use of a path measure data processing inequality produces an unpredictable looseness in the bound which we refer to as the score matching gap. This observation motivates our fundamental question:

The existence of this gap has non-trivial consequences. In particular, it renders the score matching objective an unreliable proxy for model evaluation; a lower score matching loss does not necessarily imply a better approximation. Consequently, evaluation typically utilizes sample-based metrics, which require many function evaluations and fail to disentangle the errors arising from score approximation with those arising from the discretization of the reverse process [1], [2].

Even without being formally defined, the score matching gap has been implicitly involved in the convergence bounds literature [21], motivated by the study of the strong generalization ability of DMs [22]. Such results consist in upper-bounding a discrepancy¹ \(\mathrm{D}\left(\mu, \overleftarrow{\mu}^{\theta}_T\right)\) by the score matching loss and other terms accounting for the initialization and discretization error. A common practice in this literature is to relate the error to a discretized version of \(\int_{[0,T]} \varepsilon_t^\theta({\overrightarrow{X}}_t)^2 \text{\normalfont d}t\) or assume an uniform bound of \(\varepsilon_t^\theta({\overrightarrow{X}}_t)^2\) over \([0,T]\), meaning implicitly that the score matching loss is the right objective to evaluate the method. Thus, the potential looseness of the score matching gap is discarded from the convergence analysis. However, it allows to retrieve discretized versions of 1 for various choices of \(\mathrm{D}\), including the KL divergence case [23]–[26] as an improvement of previous results for the Total Variation (TV) distance [27], [28]. To address certain limitations of KL and TV, Wasserstein convergence bounds have also been considered [29], often involving a uniform bound on the score matching loss [30]–[33].

Contributions. In this work, we provide a comprehensive analysis of the score matching gap with a focus on tightening score matching loss bounds like 1 . We begin 3 with a result suggesting that for \(T\) large, the score matching gap is tight in the worst case: it can only be tightened given additional information. In Section 3.1, we tighten 1 using a pseudo Lipschitz assumption on the score network \(s^\theta\), along with optional dissipativity and obtain, \[\begin{align} \label{eq:informal-forward-KL-bound} \text{\normalfont{{KL}}}\left(\mu | \overleftarrow{\mu}^{\theta}_T \right) \lesssim \int_0^T \lambda_T(t) \mathbb{E}_{} [ \varepsilon_t^\theta({\overrightarrow{X}}_t)^2 ] \text{\normalfont d}t + \mathrm{K}_T \;, \end{align}\tag{2}\] with \(\mathrm{K}_T := C_T \text{\normalfont{{KL}}}\left(\overrightarrow{\mu}_T | \gamma^d \right)\) for some \(C_T < 1\). The time weighting \(\lambda_T(t) < 1\) is decreasing in \(t\) (see Table ¿tbl:tab:decay-table? for explicit rates), prioritising the score matching loss near convergence. This also tightens, as a byproduct, the classical convergence bounds of [24], [26] in the time-continuous setting.

Inspired by the convergence analysis literature, we investigate the score matching gap for alternative topologies. In Section 4.1, we derive bounds of type 1 for the reversed KL divergence \(\text{\normalfont{{KL}}}\left(\overleftarrow{\mu}^{\theta}_T | \mu \right)\) with a time-decaying weighting that depends only on the concentration properties of the data distribution \(\mu\). We also extend our theory to Wasserstein distances in 4.2, where we show that under dissipativity and smoothness conditions, it holds that, \[\begin{align} \mathscr{W}_1(\mu, \overleftarrow{\mu}^{\theta}_T) \lesssim \int_0^T C_T(t) \mathbb{E}_{} [ \varepsilon_t^\theta({\overrightarrow{X}}_t) ] \text{\normalfont d}t + \mathscr{W}_1(\overrightarrow{\mu}_T, \gamma^d) \;, \end{align}\] where, again, \(C_T(t)\) is a decay function detailed in Table ¿tbl:tab:decay-table?.

Finally, Section 5 presents some experimental results on a low-dimensional setting showing that our time-decaying weightings provide better estimates of sample quality that are more useful for model evaluation.

Technical innovations. Our core technical innovation is to exploit the contraction properties of the backward processes, connecting the score-matching gap to ergodicity of the underlying diffusion process. For our KL bounds, we exploit an ‘entropy flow’ analysis (3.1) alongside logarithmic Sobolev inequalities [34]. These ideas appeared partially in the convergence bounds literature, focusing on time-discrete settings. In particular, [27], [30] used “\(\chi^2\)-flows” to reach TV and Wasserstein bounds, but only considered reversed \(\chi^2\) flows, without exhibiting an explicit decay like \(\lambda_T(t)\). For our Wasserstein bounds, we use a reflection coupling technique to tighten the time-dependence of the bounds. While using this technique is new in the context of DMs, some Wasserstein convergence bounds tightening the classical \(\mathcal{O}(T)\) factor also recently emerged [35], [36], as well as time-uniform results [29].

All omitted proofs can be found in the appendix.

Notation. We denote \(a \wedge b := \min(a,b)\). The set of Borel probability measures on \({\mathbb{R}^d}\) is written \({\mathcal{P}}({\mathbb{R}^d})\). Given \(\rho,\pi \in {\mathcal{P}}({\mathbb{R}^d})\), their Kullback-Leibler (KL) divergence is \(\text{\normalfont{{KL}}}\left(\rho | \pi \right) := \int \log (\text{\normalfont d}\rho / \text{\normalfont d}\pi) \text{\normalfont d}\rho\) and the relative Fisher information is \(\mathscr{I}\left(\rho | \pi \right) := \int \Vert \nabla \log (\text{\normalfont d}\rho / \text{\normalfont d}\pi) \Vert^2 \text{\normalfont d}\rho\). We use \(\lesssim\) for inequality up to an absolute constant. Given \(I \subset \mathbb{R}^d\), we note by \({\mathrm{C}}^{1,2}(I \times {\mathbb{R}^d})\) the functions \(f(t,x)\) that are \({\mathrm{C}}^1\) in \(t\) and \({\mathrm{C}}^2\) in \(x\). Given \(\mu,\nu \in {\mathcal{P}}({\mathbb{R}^d})\), \(\Pi(\mu,\nu) \subset {\mathcal{P}}({\mathbb{R}^d}\times {\mathbb{R}^d})\) is the set of couplings between \(\mu\) and \(\nu\). For a distance \(d\) on \({\mathbb{R}^d}\), the Wasserstein distance is \(\mathscr{W}_d (\mu,\nu) := \inf_{\pi \in \Pi(\mu,\nu)} \int d(x,y) \text{\normalfont d}\pi(x,y)\). For \(p\geqslant 1\), \(\mathscr{W}_p\) denotes \(\mathscr{W}_d\) with the \(p\)-norm \(d(x, y) = \left\Vert x-y \right\Vert_p\).

2 Setup and Technical Background↩︎

2.1 Background on Diffusion Models↩︎

In this paper, we consider the forward process to be an Ornstein-Uhlenbeck process on \({\mathbb{R}^d}\), i.e., \[\begin{align} \label{eq:forward-Ornstein-Uhlenbeck} \text{\normalfont d}{\overrightarrow{X}}_t = - {\overrightarrow{X}}_t \text{\normalfont d}t + \sqrt{2} \text{\normalfont d}\mathrm{B}_t \;, \quad {\overrightarrow{X}}_0 \sim \mu \;, \end{align}\tag{3}\] where \(\mu\) is the data distribution. This process admits the standard Gaussian \(\gamma^d := \mathrm{N}(0, \mathrm{I}_d)\) as invariant distribution. We denote \(\overrightarrow{\mu}_t := \mathrm{Law}({\overrightarrow{X}}_t)\), \(\overrightarrow{p}_t\) its Lebesgue density, for \(t>0\), and by \(\tilde{p}_t := \overrightarrow{p}_t/ \gamma^d\) the renormalized density.

Let \(T>0\) and define the time-reversal process \({\overleftarrow{X}}_t := {\overrightarrow{X}}_{T - t}\), which under mild conditions is shown to be a weak solution of the SDE [37], [38], \[\begin{align} \label{eq:true-backward-equation} \text{\normalfont d}{\overleftarrow{X}}_t = \left( -{\overleftarrow{X}}_t + 2 \nabla \log \tilde{p}_{T - t} ({\overleftarrow{X}}_t) \right) \text{\normalfont d}t + \sqrt{2} \text{\normalfont d}\overline{\mathrm{B}}_t\;, \end{align}\tag{4}\] with \({\overleftarrow{X}}_0 \sim \overrightarrow{\mu}_T\) and \((\overline{\mathrm{B}}_t)_{t\geqslant 0}\) a standard Brownian motion. Similarly, we denote \(\overleftarrow{\mu}_t := \mathrm{Law}({\overleftarrow{X}}_t)\) and by \(\overleftarrow{p}_t\) its Lebesgue density for \(t>0\). The score \(s(t,x) := 2 \nabla \log \tilde{p}_t(x)\) is usually estimated through a parametric family of score networks \(\{s^\theta : \mathbb{R}_+ \times {\mathbb{R}^d}\to {\mathbb{R}^d},~\theta \in \Theta \}\), where \(\Theta\) is a given hypothesis class. Given such a parameter \(\theta\in \Theta\), the backward process 4 is approximated by the following SDE, starting from \({\overleftarrow{X}}_0 \sim \gamma^d\) (instead of the unknown \(\overrightarrow{\mu}_T\)) and defined by \[\begin{align} \label{eq:inference-backward} \text{\normalfont d}{\overleftarrow{X}}_t^\theta = \{-{\overleftarrow{X}}_t^\theta + s^{\theta}(T-t, {\overleftarrow{X}}_t^\theta ) \} \text{\normalfont d}t + \sqrt{2} \text{\normalfont d}\mathrm{B}_t\;. \end{align}\tag{5}\] Popular discretization schemes to simulate 5 include the Euler-Maruyama and exponential Euler integrator schemes [39]. We denote \(\overleftarrow{\mu}^{\theta}_t := \mathrm{Law}({\overleftarrow{X}}_t^\theta)\) and \(\overleftarrow{p}^{\theta}_t\) its Lebesgue density. Throughout, we assume that \(s^\theta\) is locally bounded measurable.

In practice, \(\mu\) is unknown and the score has to be estimated from samples \((Z_1,\dots,Z_n) \sim \mu^{\otimes n}\). It has been shown by [13], [40] that \(\theta\) can be estimated by minimizing the denoising score matching (DSM) loss, \[\begin{align} \frac{1}{n} \sum_{i=1}^n \int \mathbb{E}_{} [ \Vert s^\theta (t, {\overrightarrow{X}}_t^{Z_i}) - 2 \nabla \log \tilde{p}_{t|0} ({\overrightarrow{X}}_t^{Z_i} | Z_i) \Vert^2 ] \text{\normalfont d}\varpi(t)\;, \end{align}\] where \(\varpi\) is a probability distribution on \([0,T]\), \({\overrightarrow{X}}_t^x\) corresponds to 3 initialized at \({\overrightarrow{X}}_0^x := x\) , and \(\tilde{p}_{t|0} (\cdot | x)\) is the conditional density of \({\overrightarrow{X}}_t\) given \({\overrightarrow{X}}_0=x\). This is motivated by the identity \(\nabla \log \tilde{p}_t({\overrightarrow{X}}_t) = \mathbb{E}_{} [ \nabla \log \tilde{p}_{t|0}({\overrightarrow{X}}_t | {\overrightarrow{X}}_0) ]\), coming from the Fisher’s identity [41]. Modern practices often involve a time-weighted version of this loss [4], [5].

2.2 Logarithmic Sobolev Inequalities↩︎

Let \(\Phi(x) := x\log(x)\), with \(\Phi(0) := 0\). Let \(\pi \in {\mathcal{P}}({\mathbb{R}^d})\), the entropy functional associated to \(\pi\) is defined as \(\mathrm{Ent}_\pi(f) := \mathbb{E}_{} \left[ \Phi(f(X)) \right] - \Phi(\mathbb{E}_{} \left[ f(X) \right])\), for \(X \sim \pi\) and \(f \in \mathrm{L}^1(\pi)\). Some of our contributions are based on the logarithmic Sobolev inequalities (LSI) associated with [eq:true-backward-equation,eq:inference-backward] [12], [42].

Let \(\mu,\nu \in {\mathcal{P}}({\mathbb{R}^d})\) such that \(\mu\ll\nu\), the LSI for \(\mu\) can be rewritten as \(2\text{\normalfont{{KL}}}\left(\nu | \mu \right) \leqslant\rho \mathscr{I}\left(\nu | \mu \right)\). The typical example is the standard Gaussian \(\gamma^d\), which satisfies the LSI with constant \(1\) [34]. LSIs classically appear in the ergodic theory of diffusion processes [12]. In particular, the fact that the invariant distribution of a given diffusion process (such as \(({\overrightarrow{X}}_t)_{t\geqslant 0}\)) satisfies an LSI is equivalent to the exponential convergence of the entropy along the associated semigroup. In our case, by reversibility of the invariant measure \(\gamma^d\), this translates to \[\begin{align} \label{eq:exponential-convergence-oub} \text{\normalfont{{KL}}}\left(\overrightarrow{\mu}_t | \gamma^d \right) \leqslant e^{-2t} \text{\normalfont{{KL}}}\left(\mu | \gamma^d \right) \;. \end{align}\tag{6}\] In our paper, we exploit these convergence properties of the forward and backward processes to obtain sharper bounds.

3 Theoretical Study of Score Matching Gaps↩︎

To motivate our work, we first investigate the question: can the score matching bound be tightened in general? In the proposition below, we provide a negative result.

Here, \(\nu\) represents any given distribution produced by a diffusion model and \(\mu\) is the target. The proposition shows that by increasing \(T\), the score matching gap can be made arbitrarily tight with the right choice of score. The construction is based on the solution to a dynamic Schrödinger bridge problem [43] and can be seen as a non-typical, irregular score function. This suggests that, in general, one cannot tighten the score matching gap without specifying particular score estimators. This motivates our direction of improving the score matching bound through the use of additional regularity assumptions. To this end, we first describe the entropy flow method for DMs.

3.1 Adapting the Entropy Flow Technique to DMs↩︎

6 is based on the celebrated De-Bruijn identity [44], \(\frac{\text{\normalfont d}}{\text{\normalfont d}t} \text{\normalfont{{KL}}}\left(\overrightarrow{\mu}_t | \gamma^d \right) = - \mathscr{I}\left(\overrightarrow{\mu}_t | \gamma^d \right)\) and follows directly from the LSI for \(\gamma^d\) and Grönwall’s lemma. Such “entropy flow” computations have also been used to bound the distance between the marginal distributions of Langevin processes [45], and in various subfields of machine learning [46]–[50]. The interest of this technique in our context is to provide natural upper bounds on the forward and reverse KL divergences between \(\mu\) and \(\overleftarrow{\mu}^{\theta}_T\).

When \(\mu\) has a Lebesgue density, the map \((t,x) \mapsto \overleftarrow{p}_t(x)\) is in \({\mathrm{C}}^{1,2}([0,T) \times {\mathbb{R}^d})\) and the spatial derivatives up to order \(2\) are bounded (see 8.2).

Background on Fokker-Planck equations is given in 7.1. Let us recall the notation, \[\begin{align} \label{eq:epsilon-definition} \varepsilon_t^\theta (x)^2 := \Vert 2 \nabla \log \tilde{p}_{t} (x) - s^{\theta} (t, x) \Vert^2\;, \quad x \in {\mathbb{R}^d}\;. \end{align}\tag{7}\] When stated, we make the following assumption, ensuring that several terms appearing in our theory are finite.

Note that the first part of the assumption only holds for \(t \in (0,T)\), in particular, it does not require that \(\mathscr{I}\left(\mu | \gamma^d \right)\) or that \(\mu\) is absolutely continuous w.r.t.\(\gamma^d\). The second part of the assumption is mild. In particular, square-integrability of \(\varepsilon_t^\theta ({\overleftarrow{X}}_t)\) is necessary for the score matching framework to be well-posed and the last part is typically implied by Hölder continuity of \(s^\theta\) when \(\mu\) has compact support.

The following proposition is an upper bound on the relative entropy between the true and estimated backward processes.

This result is an adaptation of the entropy flow method for DMs and follows from the Fokker-Planck equations satisfied by the backward processes (7.1). It is a natural generalization of the De-Bruijn identity for DMs. Whenever \(\overleftarrow{\mu}^{\theta}_t\) satisfies a LSI, i.e., the Fisher information can be replaced by a KL divergence, and 1 leads to a contraction, improving the bound.

In [sec:right-kl-bounds,sec:reversed-kl-bounds] we will apply respectively 1 to tighten the score matching gap via time-decaying weightings (see Table ¿tbl:tab:decay-table?, [fig:exponential-decay-non-reverse-KL,fig:exponential-decay-reversed]).

3.2 Tightening the bound under LSI↩︎

The next theorem is a generic bound on the forward KL \(\text{\normalfont{{KL}}}\left(\mu | \overleftarrow{\mu}^{\theta}_T \right)\) expressed in terms of the (possibly infinite) time-dependent LSI constant of \(\overleftarrow{\mu}^{\theta}_t\).

This theorem improves over 1 obtained by [11] using Girsanov’s theorem, which it recovers in the case of infinite LSI constants by choosing \(\alpha \equiv 0\) (roughly matching the analysis of [51]). However, as we always have \(\rho(0) = 1\) (as \(\overleftarrow{\mu}^{\theta}_0 = \gamma^d\)), we can systematically improve the time-dependence of 1 as soon as \(\rho\) is continuous.

We note that to benefit from an exponential decay 3, we need to choose \(\alpha\) such that \(\alpha(t) \in (0,1)\) at least on a subset of \([0,T]\) with non-zero Lebesgue measure. While this leads to improvements in the time-dependence of the bound, it might also lead to worsening the absolute constant \(1/4\) in front of the integral. In the following results, we make the choice \(\alpha \equiv 1/2\) to simplify the exposition. As we will see, this leads to non-trivial exponential decay in the bound, improving the time dependence at the cost of potentially worsening the absolute constants in the bound.

Now, our main challenge is to estimate the log-Sobolev constant of \(\overleftarrow{\mu}^{\theta}_t\), for all \(t\in [0,T]\). We note that 5 can be seen as a perturbation of the Ornstein-Uhlenbeck process 3 , with the perturbation drift \(s^\theta(t,\cdot)\). This observation allows us to exploit the rich literature on perturbative analysis of LSIs [52], [53] for Langevin processes. In the next theorem, we obtain score matching bounds based on a one-sided Lipschitz condition.

We refer to 4 for further upper bounds on the constant \(\mathrm{K}_T\) in 5. Interestingly, our theory links time-weightings improvements in score-matching bounds with regularity conditions on \(s^\theta\). We consider two particular cases: a one-sided Lipschitz constant uniformly bounded by \(M\) over \([0,T]\) and a weaker assumption of a one-sided Lipschitz constant \(M_t\) that can diverge in \(\mathcal{O}_{} \left( 1/t \right)\) as \(t \to 0^+\). In particular, if \(M < 1\), we obtain a time-uniform bound of order \(\mathcal{O}(\sup_t \mathbb{E}_{} [ \varepsilon_t^\theta({\overrightarrow{X}}_t)^2 ] )\). Moreover, note that the initialization error term, \(e^{-2T} \lambda_T(T) \text{\normalfont{{KL}}}\left(\mu | \gamma^d \right)\), is also improved over the classical analysis. In the idealized continuous-time setting, this result improves over the classical KL convergence analysis [26], [54]. Additionally, our proof technique can accommodate any formula for \(t \mapsto M_t\) (e.g. \(M_t= \mathcal{O}(1 / t^p), p > 0\)), even though it may lead to intricate expressions (see 8.4).

Time decay analysis. The map \(\lambda_T(t)\) is represented on 1 for various parameter values and are always decreasing and \(\lambda_T(0) = 1\). Qualitatively, this shows that the score approximation has more impact on the score matching gap for small forward times, i.e., low noise levels. This is to be compared with the standard scheduling practices when training DMs [5], [18]. Indeed, as noted by [4], the standard \(\epsilon\)-prediction objective implies to downweight the small noise scales. With our notation, a weighted DSM loss, \[\begin{align} \frac{1}{n} \sum_{i=1}^n \int_0^T \omega_t\mathbb{E}_{} [ \Vert s^\theta (t, {\overrightarrow{X}}_t^{Z_i}) - 2 \nabla \log \tilde{p}_{t|0} ({\overrightarrow{X}}_t^{Z_i} | Z_i) \Vert^2 ] \text{\normalfont d}t \end{align}\] is minimized with a finite dataset \((Z_1,\dots,Z_n) \sim \mu^{\otimes n}\), and \(t\mapsto \omega_t\) is increasing, i.e., it is smaller for small times scales. A byproduct of this procedure might be to reduce the overfitting at small noise scales and lead to a better score approximation at these small scales, as overfitting leads to generating only images close to the training dataset. While this procedure is initially motivated by reducing the variance of the objective [4], we argue that our theory may serve as a partial explanation of the good generalization error observed in practice, as it suggests that the small noise scales matter more for generalization than the large ones. We empirically investigate this in 5. Most importantly, our bounds decouple the time-weighting that should be used for training [4] and evaluation.

In 5, we require a small pseudo-Lipschitz constant (\(M < 1\)) to obtain an exponential decay \(\eta_T\) and, hence, a time-uniform bound. This is due to exponentially growing estimates of the LSI constant of \(\overleftarrow{\mu}^{\theta}_t\) when \(M>1\). As shown by [53], it is possible to improve these estimates by adding a dissipativity assumption, leading to the following condition.

which we exploit in the following theorem.

7 improves the classical ELBO of 1 and tightens the score matching gap. It requires the drifts of the backward equation to be one-sided Lipschitz (as in 13) and 2 adds a dissipativity-like condition, which we call "one-point dissipativity". Similar dissipativity conditions classically appear in the literature [53], [57].

4 Extending the Score Matching Gap to Alternative Topologies↩︎

We now investigate extensions of the score matching gap for the reverse KL and Wasserstein topologies. While losing the link with variational inference, this approach circumvents limitations of 3.2, in particular, the reverse KL results of 4.1 avoid regularity assumptions on \(s^\theta\). Alternatively, and inspired by the convergence analysis literature, Wasserstein distances studied in 4.2 avoid the absolute continuity constraints of KL-based metrics.

4.1 Reverse KL Bounds↩︎

Despite the fact that the KL divergence is not symmetric [58], the entropy flow toolbox is flexible enough to obtain a flow for the reverse divergence \(\text{\normalfont{{KL}}}\left(\overleftarrow{\mu}^{\theta}_t | \overleftarrow{\mu}_t \right)\). To state this result, we need an assumption, which is the reverse counterpart of 1.

We then have the following reverse KL bound.

Similarly, [27], [30] analyzed the time-derivative of the reverse \(\chi^2\)-divergence \(\chi^2\left(\overleftarrow{\mu}^{\theta}_t | \overleftarrow{\mu}_t \right)\) to get TV and Wasserstein convergence bounds. 2 may be seen as a generalization of their technique to forward and reverse KL divergences. Beyond yielding a stronger result, our approach provides explicit decay rates (see 8), which are hard to obtain in the \(\chi^2\) case.

The results of 3.2 involve estimating log-Sobolev properties of \(\overleftarrow{\mu}^{\theta}_T\), which requires regularity assumptions on \(s^{\theta}\). Interestingly, we can avoid such an assumption by considering the reverse KL divergence \(\text{\normalfont{{KL}}}\left(\overleftarrow{\mu}^{\theta}_{T} | \mu \right)\) which already appeared in the message passing literature [59] or in relation to Rényi variational inference [60]. To do so, we involve instead a log-Sobolev property for \(\mu\), as detailed below.

8 contains an explicit decay \(\bar{\lambda}\) (see 2) independent of \(s^{\theta}\). 8 shows that this decay is at least exponential, making the bound time-uniform, i.e., of order \(\sup_t \mathbb{E}_{} [ \varepsilon_t^\theta({\overleftarrow{X}}_{T - t}^\theta) ]\). Such LSI assumptions have already been considered in the literature [27], [30]. It allows us to estimate the LSI constant of \(\overleftarrow{\mu}_t\) (see 8.6).

As mentioned earlier, such a modified score matching gap cannot be linked to a classical ELBO [11], as minimizing the reverse KL divergence is not equivalent to maximizing the likelihood of the model. Despite this, we argue that the reverse KL divergence \(\text{\normalfont{{KL}}}\left(\overleftarrow{\mu}^{\theta}_T | \mu \right)\) also translates a facet of the problem by providing new insights on the evaluation of DMs: concentration properties of the data distribution can be exploited in the absence of exploitable regularity of the score network and we observe again that small noise scales are key to evaluate overfitting in DMs. 8 can also be seen as a reverse KL convergence bound for continuous-time DMs, which is new to the best of our knowledge. Interestingly, the term \(\mathbb{E}_{} [ \varepsilon_{ t}^\theta({\overleftarrow{X}}_{T - t}^\theta)^2 ]\) already appears in the literature on Wasserstein convergence bounds for DMs, but with a worse dependence on \(T\) [36], [65].

Note that 8 involves the score approximation w.r.t.\(\overleftarrow{X}^\theta_t\), which cannot be estimated from samples. We take a first step towards a more classical score approximation term involving \(X_t\) in 1 via Donsker-Varadhan’s change of measure [66].

4.2 Wasserstein Bounds↩︎

Wasserstein bound via concentration of measure. A corollary of 8 is a new bound in the \(\mathscr{W}_2\) distance. Indeed, it has been shown by [67], [68] that the LSI with constant \(\rho_0\) implies the Talagrand inequality \(\rho_0 \mathscr{W}_2 (\mu, \nu) \leqslant\text{\normalfont{{KL}}}\left(\nu | \mu \right)\), for all \(\nu \in {\mathcal{P}}({\mathbb{R}^d})\). Combined with 8, this proves that \[\begin{align} \label{eq:wasserstein-bound-via-talagrand} \frac{\mathscr{W}_2(\mu,\overleftarrow{\mu}^{\theta}_T)^2}{\rho_0} \leqslant&\bar{\mathrm{K}}_T + \int_0^{T} \frac{\bar{\lambda}(t)}{2} \mathbb{E}_{} [ \varepsilon_{ t}^\theta({\overleftarrow{X}}_{T - t}^\theta)^2 ] \text{\normalfont d}t \;. \end{align}\tag{9}\] This is of order \(\mathcal{O}((\rho_0^2 \sup_t \mathbb{E}_{} [ \varepsilon_{ t}^\theta({\overleftarrow{X}}_{T - t}^\theta)^2 ])^{1/2})\), thus time-uniform. However, this term still cannot be estimated from samples. In the rest of this section, we alleviate this issue by exploiting the reflection coupling technique.

Wasserstein Bounds via reflection couplings. Analogously to the KL case, we exploit Wasserstein contraction properties. Historically, exponential contraction in \(\mathscr{W}_p\) (\(p\in[1,+\infty)\)) were obtained for some Langevin SDEs [69], [70] via synchronous couplings. Major progress was achieved by [71], using reflection couplings [72] to obtain Wasserstein contraction rates under more generic assumptions (see 7.3). Inspired by these works, we rely on a time-dependent dissipativity at distance condition as defined in 2. Let us introduce the notation \[\begin{align} \kappa_b(r) := \inf_{\left\Vert x - y \right\Vert = r} \frac{\langle x - y, b(y) - b(x) \rangle}{\left\Vert x - y \right\Vert^2} \;. \end{align}\]

We observe that the \((L,K,R/2)\)-dissipativity at distance condition of 2 implies that \[\begin{align} \kappa_b (r) \geqslant\begin{cases} -L, & r\leqslant R, \\ K, & r \geqslant R \;. \end{cases} \end{align}\] This is instrumental in our proofs, where it is used to control the distance between the two coupled SDEs (see 9).

To obtain Wasserstein contraction, a key insight from [71] is to analyze the quantity \(f(r_t)\), where \(r_t\) is the distance between the coupled processes and \(f:\mathbb{R}_+ \to \mathbb{R}_+\) is a function constructed in the following lemma.

Additional details are gathered in 7.3. [71] showed that contraction in the distance \(\mathscr{W}_f\) induced by the distance \((x,y) \mapsto f(\left\Vert x - y \right\Vert)\) requires the existence of \(f\) as given by 3. For this reason, we focus in our main result on \(\mathscr{W}_f(\mu, \overleftarrow{\mu}^{\theta}_T)\). Note however that \(\mathscr{W}_f\) distance is equivalent to \(\mathscr{W}_1\). Indeed, we prove in 9 that \(\mathscr{W}_1 \leqslant\mathscr{W}_f \leqslant 2 e^{\frac{LR^2}{8}} \mathscr{W}_1\), with \(f\) given by 3,

Let us recall the notations \(b_t(x) := s(t,x) - x\) and \(b^\theta_t(x) := s^\theta(t,x) - x\), with \(s(t,x) := 2\nabla \log \tilde{p}_t(x)\). Since we wish to derive contractions for two SDEs with different drifts, we must make a technical but important modification to the method of [71]. In the present case, borrowing the methodology of [73], we couple the processes using a reflection coupling at distance and with a synchronous coupling when close. All details are rigorously stated in 9. This leads to the following theorem.

The main advantage of this result compared to 9 is that it directly involves the forward process \({\overrightarrow{X}}_t\), which can be estimated from samples as in the forward case. Similar to [sec:right-kl-bounds,sec:reversed-kl-bounds], we observe that the \(C_T\) factor comes from the regularity assumption on \(s^\theta\). As above, we observe that \(C_T\) is a decreasing function of \(t\), showing that the score approximation is more critical to the performance for small noise scales.

Using the 10 in appendix, we deduce the following corollary, which is a bound in the distance \(\mathscr{W}_1\).

This bound might also be interpreted as a convergence bound, in the idealized setting of a continuous-time backward process. Then, 2 has improved dependence on \(T\) compared to the recent results of [31], [35], [36]. From this convergence analysis perspective, the closest work from ours is [29], which, despite using a slightly different setting, can obtain time-uniform bounds when transferred to our setup. Their setup, proof technique, and assumptions differ from our work and may be seen as a complementary approach.

5 Empirical Study↩︎

We investigate the validity of our theoretical results and their potential practical utility on toy experiments and with a DDPM model on CIFAR-10. See 11 for further details regarding experiments.

5.0.0.1 Toy experiment.

We consider a toy setting where the data distribution is uniform on a unit circle embedded in \(\mathbb{R}^2\) with a training set of equally spaced points around the circle. We train a small DM and since the setting is low dimensional, we can directly estimate the KL divergence of model samples via a histogram approximation. In 3, we plot the test curves at different timesteps (i.e., the test denoising score matching loss for some times \(t\in[0,T]\), which is, up to an additive constant, an estimation of \(\mathbb{E}_{} [ \varepsilon_t^\theta ]\)), compared against the KL divergence over the training trajectory. Our results are gathered in 3. We find that the test curves at different timesteps are quite different: at larger timesteps, the test score matching error reaches its minimum further into training, or even never reaches it. Furthermore, we find that the test curves at smaller timesteps more closely track the KL divergence over the training trajectory, identifying a clear point of inflection. This agrees with our theory that smaller noise scales are more critical for model evaluation.

5.0.0.2 Evaluating regularity properties of the U-Net.

Next, we assess the validity of the Lipschitz and dissipativity assumptions used throughout our work in a high dimensional realistic setting. We use a DDPM U-Net model trained on CIFAR-10. In Figure 4, we compute the quantity, \[\label{eq:exp95kappa} \kappa(x, y) = \tfrac{\langle s(y, t) - s(x, t), x - y \rangle}{\|x-y\|^2},\tag{10}\] for different values of \(x\) and \(y\) and plot it against \(\|x - y\|\). This allows us to evaluate the value of the smoothness constant, \(M_t = \sup_{x, y} \kappa(x, y)\). We find that \(\kappa(x, y)\) is determined largely by \(\|x - y\|\) and find a similar thing when \(x\) and \(y\) are based on separate test points. For larger \(t\) it is clear that the quantity is both upper and lower bounded, verifying the assumptions used in this work. Across time-steps, it seems that \(\kappa(x, y)\) is largest when \(\|x - y\|\) is small with values growing as \(t\) is taken small.

5.0.0.3 Score perturbations on CIFAR-10.

Remaining in the same setting as the previous section, we consider a perturbed version of the score: \(\tilde{s}(x, t) = s(x, t) + {1}_{t_p \leqslant t \leqslant t_p + \delta} \xi \;\), where \(\xi\) is a Gaussian vector. We consider how changing \(t_p\) impacts the KL divergence and score matching losses. In Figure 6, we see that changing \(t_p\) greatly changes the character of samples, with low values of \(t_p\) producing grainy images and large values of \(t_p\) creating subtle but more global, semantically significant changes. In Figure 5 we plot how the KL divergence of samples changes as \(t_p\) is increased and we compare this with how different weightings of score matching losses change. We find that our weighting, with a crudely estimated regularity constant, better predicts the sensitivity of perturbations than the ELBO loss, whilst faithfully remaining an upper bound.

6 Conclusion↩︎

In this paper, we explored the score matching gap for DMs, i.e. the fundamental discrepancy between the learning goal used in practice and the theoretical quantity we aim to control. We first provided theoretical support for the need of additional structural assumptions to improve the classical score matching bound. We then revisited the ELBO and obtained a time-decaying version through the regularity properties of the score estimators. Then, we extended the score matching gap to two alternative topologies: the Wasserstein one, leading to similar qualitative conclusions, and the reverse KL, where regularity of the data distribution is needed instead of regularity of the score network. Finally, we empirically investigated our theory on a toy setting.

Limitations & future works. Several directions remain to be studied. In particular, tightening the log-Sobolev estimates of 3.2 might lead to improved bounds and even more informative decay terms. Alternatively, we believe that our technical insights might be exploited in the convergence analysis literature and tighten the time-dependence of several existing results.

Impact Statement↩︎

Our work is largely theoretical, and we do not expect it to have any particular ethical or societal impact.

Acknowledgements↩︎

U.S. is partially supported by the French government under the management of Agence Nationale de la Recherche as part of the “Investissements d’avenir” program, reference ANR-19-P3IA-0001 (PRAIRIE 3IA Institute). B.D., T.F., M.H, and U.S. are supported by the European Research Council Starting Grant DYNASTY – 101039676. A.D. is supported by the France 2030 program with the reference ANR-25-PEIA-0001 (THEOREM project). A.D. is funded by the European Union (ERC-2022-SYG-OCEAN-101071601). Views and opinions expressed are however those of the author(s) only and do not necessarily reflect those of the European Union or the European Research Council Executive Agency. Neither the European Union nor the granting authority can be held responsible for them. A.D. is supported by Hi! Paris and Agence Nationale de la Recherche (Grant 11-LABX-0047). This work received government funding administered by the National Research Agency (ANR) under the France 2030 program "Hi! PARIS", grant number ANR-23-IACL-0005. The authors would like to thank Giovanni Conforti for valuable comments and stimulating discussions.

The appendix is organized as follows.

• We present in 7 some additional technical background related to Fokker-Planck equations, estimation of log-Sobolev constants and reflection couplings.

• 8 is dedicated to the omitted proofs of 3.

• 9 presents the omitted proofs of 4.2 and exploits the technical background presented in 7.3.

• In 10, we present some additional results related to time-decaying generalization bounds for diffusion models, by combining our results with the generalization framework recently introduced by [55], [56].

7 Additional technical background↩︎

In this section, we provide some additional technical background related to diffusion processes, namely, Fokker-Planck equations, stability properties of the log-Sobolev inequalities, and the reflection couplings as presented in [71].

7.1 Fokker-Planck Equations↩︎

We give below a brief technical background on the Fokker-Planck equations satisfied by the backward processes. The notations are the same as in 2. Note that, in the whole paper, we use \(\partial_t\) for \(\partial / \partial t\).

Whenever the data distribution satisfies \(\mu \ll \gamma^d\), the Lebesgue density \((\overleftarrow{p}_t)_{t \geqslant 0}\) satisfies the Fokker-Planck (or forward Kolmogorov) equation [74]–[76], \[\begin{align} \label{eq:fokker-planck-backward} \partial_t \overleftarrow{p}_t = \Delta \overleftarrow{p}_t - \nabla \cdot (\overleftarrow{p}_t b_{T - t})\;, \quad t \in [0,T)\;, \end{align}\tag{11}\] where \(b_t (x) := 2\nabla \log \tilde{p}( t, x) - x\) for all \(x \in {\mathbb{R}^d}\).

Similar to 11 , under mild conditions, \((\overleftarrow{p}^{\theta}_t)_{t \geqslant 0}\) satisfies the Fokker-Planck equation, in the sense of distributions, \[\begin{align} \label{eq:fokker-planck-estimated-backward} \partial_t \overleftarrow{p}^{\theta}_t = \Delta \overleftarrow{p}^{\theta}_t - \nabla \cdot (\overleftarrow{p}^{\theta}_t b_{T - t}^\theta)\;, \quad \overleftarrow{p}^{\theta}_0 := \text{\normalfont d}\gamma^d / \text{\normalfont d}x \;, \end{align}\tag{12}\] with \(b_t^\theta (x) := s^{\theta}(t, x) - x\). In all the following, we assume that the score network \((t,x)\mapsto s^\theta(t,x)\) is jointly continuous for all \(\theta \in \Theta\).

These equations play a central role in our analysis, as they serve as the basis of the entropy flow computations, as it already appeard in existing works [45], [48].

7.2 Estimation of log-Sobolev constants↩︎

We recall the following classical stability properties of the logarithmic Sobolev inequality [12], [47]. We include proofs for the sake of completeness.

7.3 Reflection couplings and Eberle’s result↩︎

Consider an SDE of the form \(\text{\normalfont d}X_t = b(X_t) \text{\normalfont d}t + \sqrt{2} \text{\normalfont d}\mathrm{B}_t\). To obtain a Wasserstein contraction, [71] uses a reflection coupling of the form, \[\begin{align} &\text{\normalfont d}X_t = b(X_t) \text{\normalfont d}t + \sqrt{2} \text{\normalfont d}\mathrm{B}_t,\\ &\text{\normalfont d}Y_t = b(Y_t) \text{\normalfont d}t + \sqrt{2} (I - 2 e_t e_t^T) \text{\normalfont d}\mathrm{B}_t, \qquad \forall t < T,\quad X_t = Y_t, \forall t \geqslant T, \end{align}\] where \(e_t = (X_t - Y_t) / \|X_t - Y_t\|\) and \(T = \inf \{ t \geqslant 0: X_t = Y_t \}\). This coupling reflects the Brownian motion along the direction connecting the two processes, which helps control their distance. Using Itô’s lemma, they show that for any nice function \(f: \mathbb{R}_+ \to \mathbb{R}_+\), \[\begin{align} \text{\normalfont d}f(r_t) = f'(r_t) \frac{\langle X_t - Y_t, b(X_t) - b(Y_t) \rangle}{\|X_t - Y_t\|} \text{\normalfont d}t + 4 f''(r_t) \text{\normalfont d}t + 2 \sqrt{2} f'(r_t) \text{\normalfont d}\mathrm{B}_t, \end{align}\] where \(r_t = \|X_t - Y_t\|\). Defining \(\kappa\) as, \[\begin{align} \kappa(r) = \inf \bigg \{ - \frac{\langle x - y, b(x) - b(y) \rangle}{\left\Vert x - y \right\Vert^2}, ~ x,y \in {\mathbb{R}^d}, ~ \left\Vert x - y \right\Vert = r \bigg \}, \end{align}\] we have that, for contractions in the Wasserstein metric \(\mathscr{W}_f\) with rate \(c\), it is sufficient to have \[\label{eq:funcional-inequality-eberle} f''(r) - \frac{1}{4} r \kappa(r) f'(r) \leqslant- \frac{c}{4} f(r).\tag{13}\] They then prove the following lemma.

The construction of the function \(f\) and contraction rate \(c\) proceeds as follows. First, define the parameters: \[\label{eq:eberle-r0-r1} R_0 = \inf \{ R \geqslant 0: \kappa(r) \geqslant 0 \, \forall r \geqslant R \}, \qquad R_1 = \inf \{ R \geqslant R_0: \kappa(r) R (R - R_0) \geqslant 8 \, \forall r \leqslant R \}\\\tag{14}\] Next, define the auxiliary functions, \[\varphi(r) = \exp\left( -\frac{1}{4} \int_0^r s \kappa^-(s) \text{\normalfont d}s \right), \qquad \Phi(r) = \int_0^r \varphi(s) \text{\normalfont d}s,\] where \(\kappa^-(s) = \max\{-\kappa(s), 0\}\) denotes the negative part. The contraction rate \(c > 0\) is then given by \[\begin{align} \frac{1}{c} = \frac{1}{2} \int_0^{R_1} \frac{\Phi(s)}{\varphi(s)} \text{\normalfont d}s. \end{align}\] To construct the distance function \(f\), first define \(g: [0, \infty) \to \mathbb{R}_+\) by \[g(r) = 1 - \frac{c}{4} \int_0^{r \wedge R_1} \frac{\Phi(s)}{\varphi(s)} \text{\normalfont d}s.\] Finally, the strictly concave, increasing function \(f: [0, \infty) \to [0, \infty)\) is defined by, \[\begin{align} f(r) = \int_0^{r} \varphi(s) g(s) \text{\normalfont d}s. \end{align}\]

7.3.1 Special case: dissipativity at distance and 3↩︎

In this subsection, we explicit the constants and functions above in the case of dissipativity at distance (2), following the derivations of [71].

Assume that the drift \(b : {\mathbb{R}^d}\to {\mathbb{R}^d}\) satisfies the \((L,K,R/2)\)-dissipativity at distance condition, in the sense of 2. In particular, this implies that there are non-negative constants, \(L, K, R \geqslant 0\) such that \[\frac{\langle x - y, b(x) - b(y) \rangle}{\left\Vert x - y \right\Vert^2} \leqslant\begin{cases} L, & \left\Vert x - y \right\Vert \leqslant R\;, \\ - K, & \left\Vert x - y \right\Vert \geqslant R\;. \end{cases}\] Note that this implies, in particular, that \[\kappa(r) \geqslant\begin{cases} - L, & r \leqslant R\;, \\ K, & r \geqslant R\;. \end{cases}\] The next lemma provides the explicit rate of contraction for the distance \(\mathscr{W}_f\), under dissipativity at distance.

[71] provides a simplification of the contraction rate in the Lemma below.

7.3.1.1 Constants appearing in 3.

We can simplify the metric function by setting \[\kappa_0(r) := \begin{cases} - L, & r \leqslant R\;, \\ K, & r \geqslant R \;, \end{cases}\] so that \(\kappa \geqslant\kappa_0\). By noting that \(R_0\), \(R_1\), and \(1/c\) are decreasing functions of \(\kappa\), we can make the following choices for the functions and constants.

• \(R_0 = R\).

• \(R_1\) is the smallest value greater than \(R_0\) that satisfies \(K R_1 (R_1 - R) \geqslant 8\), i.e., \[R_1 = \frac{R}{2} + \sqrt{\frac{R^2}{4} + \frac{8}{K}}\;.\]

• The functions \(\varphi\) and \(\Phi\) are then given by, \[\label{eq:function-expressions-dissipativity-at-distance} \varphi(r) = \exp\left( -\frac{L}{8} (r \wedge R)^2 \right), \qquad \Phi(r) = \sqrt{\frac{2\pi}{L}} \operatorname{erf}\left( (r \wedge R) \sqrt{\frac{L}{8}} \right) + (r - R)_+ \exp\left( -\frac{LR^2}{8} \right).\tag{15}\]

The resulting function \(f\) and the contraction rate \(c\) are then given as in 7.3. In particular, we can take \[\begin{align} \frac{2}{c} \leqslant\int_0^{R_0} \left( r \wedge \sqrt{\frac{2\pi}{L}} \right) e^{\frac{Lr^2}{8}} \text{\normalfont d}r + \frac{4}{K} + \sqrt{\frac{8}{K}} R_0 e^{\frac{LR_0^2}{8}} \;. \end{align}\] Then, the formulas of 7 are valid after replacing \(R_0\) by \(R\) in all the expressions.

We use the notation \(f_{L,R,K}\) and \(c_{L,R,K}\) for the resulting function and contraction rate, respectively.

In the following lemmas, we prove some additional properties of the functions constructed above. These results complement the construction of [71] and will be particularly useful for our analysis. They may be seen as an additional technical contribution.

The following lemma allows to compare the Wasserstein distances \(\mathscr{W}_f\) and \(\mathscr{W}_1\).

8 Omitted proofs of 3↩︎

8.1 Proofs for the negative result↩︎

In this section, we present the proof of the negative result in Proposition 1. We begin by recalling the necessary prerequisites for the dynamical Schrödinger bridge (SB) problem in the Ornstein–Uhlenbeck setting and provide two technical lemmas.

For any \(\rho \in {\mathcal{P}}({\mathbb{R}^d})\), let \(R^\rho\) denote the path measure of the Ornstein–Uhlenbeck process on \([0, T]\) with initial distribution \(\rho\). For probability measures \(\rho_0, \rho_1 \in {\mathcal{P}}({\mathbb{R}^d})\), we define the dynamical Schrödinger bridge problem, \[\label{eq:sb-problem} \text{Minimize } \text{\normalfont{{KL}}}\left(Q | R^\rho \right)\;, \qquad \text{subject to } Q_0 = \rho_0,\; Q_T = \rho_1\;.\tag{16}\] Under the conditions \(\text{\normalfont{{KL}}}\left(\rho_0 | \gamma^d \right), \text{\normalfont{{KL}}}\left(\rho_1 | \gamma^d \right) < \infty\), it is known that it admits a unique solution [43]. Furthermore, due to the decomposition, \[\label{eq:sb95init95cond} \text{\normalfont{{KL}}}\left(Q | R^\rho \right) = \mathbb{E}_{Q_0} \left [ \text{\normalfont{{KL}}}\left(Q_{0:T|0} | R^\rho_{0:T|0} \right) \right ] + \text{\normalfont{{KL}}}\left(Q_0 | R^\rho_0 \right) = \mathbb{E}_{\rho_0} \left [ \text{\normalfont{{KL}}}\left(Q_{0:T|0} | R_{0:T|0} \right) \right ] + \text{\normalfont{{KL}}}\left(\rho_0 | \rho \right),\tag{17}\] we see that the solution \(Q\), will not depend on the initial law of the reference, \(\rho\): any two OU references \(R^\rho, R^{\rho'}\) with \(\text{\normalfont{{KL}}}\left(\rho_0 | \rho \right), \text{\normalfont{{KL}}}\left(\rho_0 | \rho' \right) < \infty\) share the same Markov kernel and hence, the two resulting SB problems share the same solution, which in this section, we denote by \(Q^*\). Taking \(\rho = \gamma^d\) for convenience and \(R := R^{\gamma^d}\), we recall the classical result that the path measure \(Q^*\) factorizes as, \[\label{eq:fg-decomp-general} \text{\normalfont d}Q^*(\omega) = f(\omega_0)\, g(\omega_T)\, \text{\normalfont d}R(\omega) \;,\tag{18}\] for measurable Schrödinger potentials \(f, g : {\mathbb{R}^d}\to \mathbb{R}_+\) [77]. Since \(R\) has marginal \(\gamma^d\) at every time, the pair of Schrödinger equations characterizing \(f, g\) reduces to \[\label{eq:schrod-equations} \frac{\text{\normalfont d}\rho_0}{\text{\normalfont d}\gamma^d}(z) = f(z) \cdot (P_T g) (z)\;, \qquad \frac{\text{\normalfont d}\rho_1}{\text{\normalfont d}\gamma^d}(x) = (P_T f)(x) \cdot g(x) \;.\tag{19}\] Moreover, the optimal \(Q^*\) is Markov [43] and inherits from \(R\) the structure of a diffusion process on \([0, T]\) with the same diffusion coefficient. Furthermore, its drift is given by the Doob \(h\)-transform, \[\label{eq:bridge-drift} b^Q(t, x) = -x + 2 \nabla \log (P_{T-t} g)(x)\;,\tag{20}\] which is well defined for almost all \(t\). For a more substantial introduction to Schrödinger bridge problems, we refer the read to the notes of [78] and [77] as well as the survey of [43].

The first lemma is a technical property showing that the conditional distribution \(Q^*_{T|0}\) is \(1\)-strongly log-concave.

The second lemma bounds the gradient of the Schrödinger potential and follows from Proposition 2.4 of [80].

8.1.1 Proof of 1↩︎

We instantiate the Schrödinger bridge problem 16 with \(\rho_0 = \nu\) and \(\rho_1 = \gamma^d\). Since \(\text{\normalfont{{KL}}}\left(\nu | \gamma^d \right) < \infty\), it follows that a unique solution exists and we denote the resulting bridge by \(Q^*\). We also denote the Schrödinger potentials by \(f, g\) which satisfy 18 and, as a result of 19 , satisfy \(P_T f \cdot g \equiv 1\). It is this optimal path measure \(Q^*\) which we use to choose the score function that closes the score matching gap.

By the bridging property [43], the path-space KL divergence reduces to the joint endpoints, \[\begin{align} \text{\normalfont{{KL}}}\left(R | Q^* \right) &= \text{\normalfont{{KL}}}\left(R_{0,T} | Q^*_{0,T} \right) + \int \text{\normalfont{{KL}}}\left(R_{:|0,T} | Q^*_{:|0,T} \right) \text{\normalfont d}R_{0,T}\\ &= \text{\normalfont{{KL}}}\left(R_{0,T} | Q^*_{0,T} \right)\;, \end{align}\] or equivalently, \(\text{\normalfont{{KL}}}\left(R_{:|0,T} | Q^*_{:|0,T} \right) = 0\) almost everywhere. Similarly, we can apply this property to the reference measure \(R^\mu\) starting at initial distribution \(\mu\) to obtain, \[\begin{align} \text{\normalfont{{KL}}}\left(R^\mu | Q^* \right) &= \text{\normalfont{{KL}}}\left(R^\mu_{0,T} | Q^*_{0,T} \right) + \int \text{\normalfont{{KL}}}\left(R_{:|0,T} | Q^*_{:|0,T} \right) \text{\normalfont d}R^\mu_{0,T}\nonumber\\ &= \text{\normalfont{{KL}}}\left(R^\mu_{0,T} | Q^*_{0,T} \right)\;,\label{eq:negative951} \end{align}\tag{26}\] where we use that \(R^\mu_{:|0,T} = R_{:|0,T}\). Furthermore, the chain rule for the KL divergence yields, \[\label{eq:negative952} \text{\normalfont{{KL}}}\left(R^{\mu}_{0, T} | Q^*_{0,T} \right) = \text{\normalfont{{KL}}}\left(\mu | \nu \right) + \Delta_T, \qquad \Delta_T := \int \text{\normalfont{{KL}}}\left(R^{\mu}_{T|0}(\cdot \mid x) | Q^*_{T|0}(\cdot \mid x) \right) \text{\normalfont d}\mu(x)\;,\tag{27}\] where we use that \(Q^*_0 = \nu\). Combining 26 and 27 , we obtain a correspondence between the marginal and path measure KL divergence: \[\label{eq:negative95295h} \text{\normalfont{{KL}}}\left(\mu | \nu \right) = \text{\normalfont{{KL}}}\left(R^\mu | Q^* \right) - \Delta_T.\tag{28}\]

Next, we control the quantity \(\Delta_T\). From Lemma 11, we have that for any \(z \in \mathbb{R}^d\), \(Q^*_{T|0}(\cdot|z)\) is \(1\)-strongly log-concave. Hence by the Bakry–Émery criterion [12], it must satisfy a logarithmic Sobolev inequality [def:log-sobolev-inequality], with constant \(1\). By 21 , we have \[\frac{\text{\normalfont d}R^{\mu}_{T|0}}{\text{\normalfont d}Q^*_{T|0}} (x \mid z) = P_T f(x) \, P_T g(z) \;.\] Thus, the logarithmic Sobolev inequality produces the bound, \[\label{eq:negative953} \Delta_T \leqslant\frac{1}{2} \int \int \norm{\nabla \log P_T f}^2 \text{\normalfont d}R^\mu_{T|0}(\cdot|z) \, \text{\normalfont d}\mu(z) = \frac{1}{2} \int \norm{\nabla \log P_T f}^2 \text{\normalfont d}\overrightarrow{\mu}_T\;.\tag{29}\] To apply the bound in Lemma 12, we must replace \(\overrightarrow{\mu}_T\) by \(\gamma^d\). For this, we observe that for any \(y \in {\mathbb{R}^d}\), \[\tilde{p}_{T|0}(x|y) = (1 - e^{-2T})^{-d/2}\, \exp\!\left( \frac{\Vert y \Vert^2}{2} - \frac{e^{-2T} \Vert x - e^T y \Vert^2}{2(1 - e^{-2T})} \right) \leqslant(1 - e^{-2T})^{-d/2}\, e^{\Vert y \Vert^2/2}\;.\] Since \(\text{\normalfont d}\overrightarrow{\mu}_T/\text{\normalfont d}\gamma^d(x) = \int \tilde{p}_{T|0}(x|y)\, \mu(\text{\normalfont d}y)\), we thus obtain, \[\label{eq:negative954} \frac{\text{\normalfont d}\overrightarrow{\mu}_T}{\text{\normalfont d}\gamma^d}(x) \leqslant(1 - e^{-2T})^{-d/2}\, e^{R^2/2}\;,\tag{30}\] where we use that \(\mathrm{supp}(\mu) \subseteq B_R(\mathbf{0})\). Thus, we conclude that combining 29 and 30 and applying Lemma 12, we obtain, \[\begin{align} \Delta_T &\leqslant\frac{1}{2} (1 - e^{-2T})^{-d/2}\, e^{R^2/2} \int \norm{\nabla \log P_T f}^2 \text{\normalfont d}\gamma^d\\ &\leqslant\frac{2(1 - e^{-2T})^{-d/2}\, e^{R^2/2}}{e^{2T} - 1}\, \left( \text{\normalfont{{KL}}}\left(\nu | \gamma^d \right) + e^{-2T} d \right)\;. \end{align}\] Thus, whenever \(T \geqslant\log{d}\), we have that \[\Delta_T \leqslant 2 e^{\frac{1+R^2}{2}-2T} (1 + \text{\normalfont{{KL}}}\left(\nu | \gamma^d \right)),\] and in particular, \(\text{\normalfont{{KL}}}\left(R^\mu | Q^* \right) < \infty\).

By the argument of 20 , the path measure \(Q^*\) coincides with the forward SDE, \[\text{\normalfont d}X_t = \bigl( -X_t + 2 \nabla \log (P_{T-t} g)(X_t) \bigr) \text{\normalfont d}t + \sqrt{2}\, \text{\normalfont d}\mathrm{B}_t\;, \qquad X_0 \sim \nu\;.\] We now reverse \(Q^*\) in time. Since \(\text{\normalfont{{KL}}}\left(R^\mu | Q^* \right) < \infty\), it follows from [81] that the reversed path measure \(\overleftarrow{Q^*}\) is the law of the diffusion started at \(\overleftarrow{Q^*}_0 = Q^*_T = \gamma^d\) with drift, \[\overleftarrow{b}^{Q^*}(t, y) = y - 2 \nabla \log (P_t g)(y) + 2 \nabla \log q^*_{T-t}(y)\;.\] In fact, we can state the drift more directly in terms of the \(fg\)-decomposition. From 18 , it follows that the marginals have density, \[\label{eq:Qstar-marginal} q^*_t(x) = (P_t f)(x)\, (P_{T-t} g)(x)\, \gamma^d(x)\;,\tag{31}\] where we use that, by reversibility, \(\mathbb{E}_R[f(X_0)\mid X_t = x] = P_t f(x)\). Consequently, using \(\nabla \log \gamma^d(x) = -x\), \[\begin{align} \overleftarrow{b}^{Q^*}(t, y) &= y - 2 \nabla \log (P_t g)(y) + 2 \bigl( \nabla \log (P_{T-t} f)(y) + \nabla \log (P_t g)(y) - y \bigr)\\ &= -y + 2 \nabla \log (P_{T-t} f)(y)\;. \end{align}\] Thus, by defining the score as, \[\label{eq:score-def} s^\theta(\tau, x) := 2 \nabla \log (P_\tau f)(x)\;,\tag{32}\] the reverse process \({\overleftarrow{X}}^\theta_t\) coincides with the reverse path measure \(\overleftarrow{Q^*}\).

Furthemore, the time-reversal \(\overleftarrow{R}^\mu\) of the path measure \(R^\mu\) also coincides with the true reverse process \({\overleftarrow{X}}_t\) defined in 4 . Since the KL divergence is invariant under time reversal, we have that the reversal \(\overleftarrow{R}^\mu\) of the path measure \(R^\mu\) satisfies, \(\text{\normalfont{{KL}}}\left(\overleftarrow{R}^\mu | \overleftarrow{Q}^* \right) = \text{\normalfont{{KL}}}\left(R^\mu | Q^* \right) < \infty\). Therefore, it follows from [82] that a Girsanov-type identity holds: \[\label{eq:negative955} \text{\normalfont{{KL}}}\left(\overleftarrow{R}^\mu | \overleftarrow{Q}^* \right) \geqslant\frac{1}{4} \int_0^T \mathbb{E}_{} \left[ \varepsilon_t^\theta({\overrightarrow{X}}_t)^2 \right] \text{\normalfont d}t + \text{\normalfont{{KL}}}\left(\overleftarrow{\mu}_0 | Q_0 \right) = \frac{1}{4} \int_0^T \mathbb{E}_{} \left[ \varepsilon_t^\theta({\overrightarrow{X}}_t)^2 \right] \text{\normalfont d}t + \text{\normalfont{{KL}}}\left(\overrightarrow{\mu}_T | \gamma^d \right)\;.\tag{33}\] Combining the two inequalities 28 and 33 yields the bound in the statement. 0◻

8.2 Technical lemmas↩︎

The next lemma regards the regularity of the true backward density \(\overleftarrow{p}^{\theta}\) and is inspired from [26].

8.3 Proof of 1↩︎

In this section, we prove 1. First, we recall that we use the notation \(\partial_t\) for \(\partial / \partial_t\), and \[\begin{align} b_t (x) := 2\nabla \log \tilde{p}( t, x) - x \;, \qquad b_t^\theta (x) := s^\theta(t,x) - x \;. \end{align}\]

8.4 Proofs of 3.2↩︎

8.4.1 Proof of 3↩︎

8.4.2 Proof of 5↩︎

Before presenting the proof of 5, we give the following technical lemma, which provides an estimate of the log-Sobolev constant of \(\overleftarrow{\mu}^{\theta}_t\) under a one-sided Lipschitz condition on the score network. This is a direct corollary of [52] (cited by [53]).

We can now present the proof of 5.

8.4.3 Time-dependent one-sided Lipschitz constants↩︎

In this subsection, we explore the case where the one-sided Lipschitz constant of the network is allowed to depend on time. This requires to adapt the log-Sobolev estimates of [52], [53], which we first present.

For the purpose of this operator, we introduce the operators \[\begin{align} P_{s,t} f(x) := \mathbb{E}_{} \left[ f(X_t) | X_s = x \right] \;. \end{align}\]

In the next theorem, we present the associated score-matching bound.

The above corrolary also concludes the proof of 5 by setting \(\lambda_T(t) := \eta_T(T - t)\).

8.5 Proof of 7↩︎

8.6 Proofs of 4.1↩︎

The proof of 2 is formally similar to the proof of 1, as explained below.

We give below the proof of 2.

We now present the following simple proposition, which is an estimate of the time-dependent log-Sobolev constant along the Ornstein-Uhlenbeck semigroup.

8.6.1 Proof of 8↩︎

8.6.2 Proof of 1↩︎

9 Omitted proofs of 4.2↩︎

In this section, we present the proofs of 11 and 2, which we split in several lemmas and propositions.

Let \(b_t(x) = s(t, x) - x\) and \(b^{\theta}_t(x) = s^{\theta}(t, x) - x\). Since we wish to derive contractions in the case where the bias terms are not equal, we must make a technical but important modification to the proof. In particular, the handling of what happens when the processes collide. The point of collision is the precise point where the reflection coupling is not well-defined and the application of Itô’s lemma in Euclidean norm is not valid. In the case of [71], the authors handle this using an optional stopping time, as well as the fact that the processes stick together after collision.

In the present case, we deal with it differently, borrowing the methodology of [73]. We couple the processes using a reflection coupling at distance and with a synchronous coupling when close. We define the coupling as follows: \[\begin{align} \label{eq:diffusion-coupling} &\text{\normalfont d}\overleftarrow{X}_t = b_{T-t}(\overleftarrow{X}_t) \text{\normalfont d}t + \sqrt{2} \phi_r^\delta(Z_t) \text{\normalfont d}\mathrm{B}_t + \sqrt{2} \phi_s^\delta(Z_t) \text{\normalfont d}\tilde{\mathrm{B}}_t\;,\\ &\text{\normalfont d}\overleftarrow{X}^{\theta}_t = b_{T-t}^\theta(\overleftarrow{X}^{\theta}_t) \text{\normalfont d}t + \sqrt{2} \phi_r^\delta(Z_t) (I - 2 e_t e_t^T) \text{\normalfont d}\mathrm{B}_t + \sqrt{2} \phi_s^\delta(Z_t) \text{\normalfont d}\tilde{\mathrm{B}}_t\;,\\ &Z_t = \overleftarrow{X}_t - \overleftarrow{X}^{\theta}_t\;,\nonumber \end{align}\tag{37}\] where \(e_t = (\overleftarrow{X}_t - \overleftarrow{X}^{\theta}_t) / \|\overleftarrow{X}_t - \overleftarrow{X}^{\theta}_t\|\) and for any \(\delta > 0\), we define the functions \(\phi_r^\delta: \mathbb{R}^d \to \mathbb{R}_+\) and \(\phi_s^\delta: \mathbb{R}^d \to \mathbb{R}_+\) as any Lipschitz continuous functions satisfying, \[\phi_r^\delta(z)^2 + \phi_s^\delta(z)^2 = 1, \qquad \phi_r^\delta(z) = \begin{cases} 1, & \|z\| \geqslant\delta\;, \\ 0, & \|z\| \leqslant\delta/2\;. \end{cases}\]

The following proposition completes the proof of 11.

The next corollary is a bound with respect to the Wasserstein’s distance \(\mathscr{W}_1\), that we can deduce from 10. It completes the proof of 2.

The argument is completely symmetric, but for consistency with previous sections, we state the result for the forward case.

10 Additional Results: Time-Uniform Generalization bounds for Diffusion Models↩︎

Recently, [55], [56] proposed an algorithm- and data-dependent analysis of \(\mathbb{E}_{} [ \varepsilon_t^\theta ({\overleftarrow{X}}_t) ]\) in 1 . We briefly recall their setup here. DMs are trained to minimize over \([0,T]\) the denoising score matching loss defined for \(z \in {\mathbb{R}^d}\) and \(t\in[0,T]\) as \[\begin{align} \ell_t (\theta, z) := \mathbb{E}_{} \left[ \left\Vert s^\theta(t, {\overrightarrow{X}}_t^z) - 2 \nabla \log \tilde{p}_{t | 0} ({\overrightarrow{X}}_t^z | z) \right\Vert^2 \right] \;. \end{align}\] In practice, as \(\mu\) is unknown, we use a finite dataset \(\mathbf{Z}^{(n)} := (Z_1,\dots,Z_n) \sim \mu^{\otimes n}\) sampled from \(\mu\). Then, what is minimized in practice by the learning algorithm (e.g., SGD, ADAM, ...) is the empirical denoising score matching loss, \[\begin{align} \mathscr{L}_{\mathrm{DSM}}^{(n)} (\theta,t) := \frac{1}{n} \sum_{1 \leqslant i \leqslant n} \ell_t (\theta, Z_i) \end{align}\] These authors define the generalization error as the difference between \(\mathscr{L}_{\mathrm{DSM}}^{(n)}\) and its population version, i.e., \[\begin{align} \mathscr{G}^{(n)} (\mathbf{Z}^{(n)},\theta,t) := \int \ell_t (\theta, z) \text{\normalfont d}\mu(z) - \mathscr{L}_{\mathrm{DSM}}^{(n)} (\theta,t)\;, \quad \theta \in \Theta \;. \end{align}\] Let \(\require{upgreek} {\widehat{\mu}_n}:= n^{-1} \sum_{1 \leqslant i \leqslant n} \updelta_{Z_i}\) be the empirical data distribution and \({\overrightarrow{X}}_t^{(n)}\) the process given by 3 initialized at \({\overrightarrow{X}}_0 \sim {\widehat{\mu}_n}\). Finally, we denote by \(\theta^{(n)}\) the (random) parameter learned by a learning algorithm (e.g., SGD, ADAM) optimizing \(\mathscr{L}_{\mathrm{DSM}}^{(n)}\). The following lemma is particular case of the results of [55], [56].

The proof of this lemma can be found for instance in [56] (take \(\lambda := \delta_t\) in their proof).

This corollary relates the performance of DMs to the generalization error \(\mathscr{G}^{(n)} (\theta^{(n)}, t)\) of the DSM loss with a better dependence on \(T\) compared to [55], [56]. It shows that the overall generation performance is mainly impacted by the generalization error associated with small noise levels (i.e., values of \(t\) close to \(T\) in the above equation). As we have seen above, this generalization of our results is only based on 16 and a direct application of our bounds. Therefore, it is clear that some of our other results may be generalized similarly. In particular, we could extend our Wasserstein bounds of 4.2 to the generalization error framework, which is new to the best of our knowledge.

11 Experiments details↩︎

11.1 Toy setting↩︎

We consider a simple toy distribution consisting of a uniform distribution on a circle embedded in 2-dimensional Euclidean space. We use a simple training set consisting of 8 points spread equidistant on the circle and we train a diffusion model based on a minimal feed-forward neural network. The implementation is based on a blog that can be found here and is similar to what is considered in [86]. Each tick on the x-axis of the plots in 7 represents \(1000\) epochs.

11.2 CIFAR-10↩︎

We also consider CIFAR-10, an implementation of the DDPM model from [1]. We use the configuration titledvp.ddpm.cifar10_continuous which implements the DDPM model of [4] but for the continuous-time variance preserving setting. The architecture is a U-Net [87] with the encoder and decoder each consisting of four resolution levels \((32\times32, 16\times16, 8\times8, 4\times4)\), with two residual blocks per level and utilizes self-attention at a resolution of \(16\times16\). The model is conditioned on time, with the timestep encoded via a sinusoidal embedding and injected into each residual block.

For Figure 4, we choose \(x, y\) by randomly choosing a test data point \(z_0\), taking a sample \(z \sim \tilde{p}_{t|0}(dz|z_0)\) and then setting \(x = z + r_x \xi_x, y = z + r_y \xi_y\), where \(\xi_x, \xi_y\) are standard multivariate Gaussians, and \(r_x, r_y\) are Weibull distributed scalars.

For Figures 6 we use the implementation of the DDIM sampler with the same initial Gaussian noise across values of \(t_p\). For Figure 5, we calculate the change in KL divergence by instead computing the change in log-likelihood using the bpd implementation in the codebase. Since this is a stochastic estimate, we compute 5 batches of log-likelihood values and plot the mean and the standard deviation.

References↩︎

---

1. \(\mathrm{D}\) typically denotes a divergence or a Wasserstein distance.↩︎