Score-based generative models are provably robust: an uncertainty quantification perspective

Authors: nikiforos  mimikos-stamatopoulos
benjamin j zhang
markos a katsoulakis
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-03 18:18:29

Abstract: Through an uncertainty quantification (UQ) perspective, we show that score-based generative models (SGMs) are provably robust to the multiple sources of error in practical implementation. Our primary tool is the Wasserstein uncertainty propagation (WUP) theorem, a *model-form UQ* bound that describes how the $L^2$ error from learning the score function propagates to a Wasserstein-1 ($\mathbf{d}_1$) ball around the true data distribution under the evolution of the Fokker-Planck equation. We show how errors due to (a) finite sample approximation, (b) early stopping, (c) score-matching objective choice, (d) score function parametrization expressiveness, and (e) reference distribution choice, impact the quality of the generative model in terms of a $\mathbf{d}_1$ bound of computable quantities. The WUP theorem relies on Bernstein estimates for Hamilton-Jacobi-Bellman partial differential equations (PDE) and the regularizing properties of diffusion processes. Specifically, *PDE regularity theory* shows that *stochasticity* is the key mechanism ensuring SGM algorithms are provably robust. The WUP theorem applies to integral probability metrics beyond $\mathbf{d}_1$, such as the total variation distance and the maximum mean discrepancy. Sample complexity and generalization bounds in $\mathbf{d}_1$ follow directly from the WUP theorem. Our approach requires minimal assumptions, is agnostic to the manifold hypothesis and avoids absolute continuity assumptions for the target distribution. Additionally, our results clarify the *trade-offs* among multiple error sources in SGMs.

Main paper: 
 introduction: 
 Score-based generative models (SGMs) [1, 2] are highly effective [3] , producing high quality samples, with more stable and less computationally intensive training methods than generative adversarial nets and normalizing flows [4] . The models are empirically robust to approximations and errors in learning the score function. While SGM generalization properties have been studied for idealized conditions [5, 6, 7, 8] , analyses of their robustness in practical settings remains underexplored. This paper analyzes SGMs through the regularity theory of nonlinear partial differential equations (PDEs) [9] , specifically Hamilton-Jacobi-Bellman (HJB) equations [10] . Our main result is the Wasserstein uncertainty propagation (WUP) theorem (Theorem 3.1), a versatile model-form uncertainty quantification (UQ) bound which we use to theoretically explain the robustness of SGMs to approximation errors in practical implementation. Generalization bounds for integral probability metrics (IPMs), such as the Wasserstein-1 (d 1 ) and the total variation (TV) distance follow directly from the WUP theorem.In the context of SGMs, the WUP theorem shows how an L 2 neighborhood around the true score function propagates to an IPM neighborhood around the true data distribution. By relating how various approximations in SGMs contributes to L 2 error with respect to the true score function, we establish how well the resulting SGM generalizes. Theorem 3.2 shows how the error in the learned score-function with respect to the explicit score-matching objective, or uniform-in-time L 2 error, and the choice of initial condition define the radii of d 1 and TV neighborhoods, respectively. Additionally, Theorem 3.3 addresses the case where the score function is learned from finite data using the denoising score-matching (DSM) objective and incorporates early stopping.Our bounds capture trade-offs among the errors. The ability of SGMs to generalize depends on choices such as the early stopping time and how overtraining to the DSM objective and neglecting early stopping lead to SGMs that memorize and overfit to the data. Notably, our approach relies on minimal assumptions on the data distribution, being agnostic to the manifold hypothesis. This contrasts with existing convergence results which assume the hypothesis [8] or not [5, 7] a priori. Furthermore, unlike previous work, we obtain our generalization bounds with respect to the d 1 distance directly, without appealing to the Girsanov theorem, the Kullback-Leibler divergence, the χ 2 divergence, or Pinsker's inequality [5, 11, 7] . This suggests our bounds may be sharper than those which bound stronger norms and divergences. A notable feature of our DSM generalization bound is that the error bound is a computable function of the DSM objective, contrasting with prior results typically assume the learned score is close to the truth with respect to the ESM objective [5] .The WUP theorem also enables robust uncertainty quantification (UQ) for score-based generative models, a rare capability in generative models. Robust UQ recognizes that learning complex models involves multiple sources of uncertainty due to modeling choices and imperfect data. The distributionally robust perspective [12] quantifies a uncertainty set based on divergences and probability metrics [13, 14, 15, 16] to measure the impact of model uncertainty around a baseline model. The drawback is that these sets are typically difficult to find computationally. The WUP theorem is an example of a robust UQ bound for IPMs that is computable.• We introduce the use of regularity theory of nonlinear PDEs for analyzing generative flows [9, 10] . Our key idea is using the Kolmogorov backward equation to study the evolution of IPMs under a generative flow, yielding generalization bounds. WUP is one particular application of this idea, and while we use it to study SGMs, this analysis is not limited to only SGMs. We show that the regularizing properties of the underlying Fokker-Planck equation imply SGM's robustness to errors with few assumptions on the data distribution. An intuitive explanation of our main technical contribution is provided in Section 4.• The WUP theorem (Theorem 3.1), which we state and prove, describes how an L 2 ball centered around the true score function maps to a neighborhood of IPMs (such as d 1 , TV, and MMD [17] ). WUP is a model-form uncertainty quantification bound, which maps uncertainties introduced by modeling choices when practically implementing SGM. The resulting generalization bounds WUP produces demonstrate that with properly chosen model parameters, SGM is robust to errors due to score-matching, finite sample approximation, early stopping time, and choice of the reference distribution. (Theorems 3.2, 3.3, and C.1).• When applied to SGM, WUP produces generalization bounds under minimal hypotheses on the data distribution and score function, and explains the impact of implementation errors has on generalization. Our approach is agnostic to the manifold hypothesis, applying whether or not the data distribution has a density. Moreover, it is adaptable to additional assumptions about the target distribution to yield improved generalization bounds.Convergence and generalization of SGMs have been well-studied. Many approaches [5, 7, 18] assume the target distribution is absolutely continuous with respect to a Gaussian, and obtain generalization bounds for TV, χ 2 , and d 1 by bounding the KL divergence, a stronger divergence, via the Girsanov theorem [11, 6, 18 ], Pinsker's inequality [5] , functional inequalities [7] , or other means [11, 19] .While [20] shows that kernel estimators for the score function are minimax optimal distribution estimators, they make no connections to deep learning-based SGMs, and they assume the target distribution is sub-Gaussian and is in a Sobolev space. Meanwhile, [6] has comprehensive sample complexity results that considers neural net approximations in the finite sample regime. Their results, however, are specialized to distributions in a Besov space, and, similar to other work [5, 11, 7, 19] , rely on bounding the KL divergence, which breakdown under the manifold hypothesis.Our results are similar to [8] which derives d 1 bounds directly, proving SGM convergence under the manifold hypothesis. However, the analysis assumes a particular discretization of the SGM. In [21] , an uncertainty propagation bound for the Wasserstein-2 distance is proven using a similar approach to this work. Their work, however, relies on strong, uncheckable assumptions on the score function and target measure, does not address the errors we study here, and is not extendable to IPMs. 
background and notation: 
 Let dimension d ∈ N, and terminal time T > 0. Denote T d ⊂ R d to be the unit torus and RT d ⊂ R d to be the torus of radius R > 0. Let P(Ω) be the space of probability distributions on Ω, where Ω is R d or RT d . Let π ∈ P(Ω) be the target data distribution, known only through a finite number of samples {z j } N j=1 ∼ π. Let f : [0, T ] × R d → R be a vector field and σ : R → R be a positive function. Score-based generative modeling [1, 2] is based on considering two Stochastic Differential Equations (SDEs for short) whose flow maps are inverses of each other. Define Y (s), X(t) to be the forward and reverse diffusion processes over time interval s, t ∈ [0, T ] that evolve according to dYEQUATIONwhere Y (s) ∼ η π (s, •), and W (t) is a standard Brownian motion in R d . Heres, η π (s, •) is the density of Y (s) at time s where the initial distribution was π. From [22] , it is known that ifm 0 = η π (T, •), then X(t) ∼ η π (T -t, •).SGMs generated new samples from π simulating trajectories of the reverse process (2) using an approximate score function s θ ≈ s = ∇ log η π . Typically f and σ are chosen such that η π (T, •) is well-approximated by a normal distribution, which is then used as the initial distribution.The score function is learned via samples {z j } N j=1 from π by minimizing one of several scorematching objectives. Let s θ be some function where the parameters θ are learned via one of three score-matching objectives. For a path ρ : [0, T ] → P(Ω), we define the functionalsEQUATIONwhere the subscript L highlights the linear dependence of J L on the underlying measure ρ. The explicit score matching objective (ESM) is J ESM (η π , θ) = J (η π , θ). As evaluation of the true score function ∇ log η π is typically inaccessible, the implicit (ISM) J ISM (η π , θ) = J L (η π , θ) [23] or the denoising (DSM) score-matching objectives [24, 1] are computed in practiceJ DSM (η π , θ) = T 0 Ω Ω s θ -∇ log η x ′ 2 dη x ′ (s)dπ(x ′ )ds.Here η x ′ (s) denotes the probability transition kernel from x ′ to x of (1) at time s. The DSM objective is most frequently used in practice as it does not require computing derivatives of s θ . It does, however, require the evaluation of η x ′ (s) in closed form, which is typically only accessible for linear SDEs. Remark 2.1 (Choice of domain Ω). While our approach will also produce bounds when Ω = R d , we primarily focus on the torus RT d , which is equivalent to a bounded domain with periodic boundary conditions. This choice ensures that the long-time behavior of the noising process (1) converges to the uniform distribution on RT d . Therefore, we set f = 0 and σ(t) = √ 2 instead of using the Ornstein-Uhlenbeck process. We make this choice for mathematical clarity, however our results generally apply, with minor modifications, to the entire space R d and for other noising processes.Finite sample approximations of ESM, ISM, and DSM are not generally equivalent. However, [24, 25] , show that DSM becomes equivalent to ESM and ISM in the finite sample regime when η π (s) is replaced with its kernel density estimate η N (s), where η N (0) = π N =1 N N i=1 δ zi . The kernel estimate, however, does not have a well-defined score at s = 0, so the DSM objective is often integrated only for s ∈ [ϵ, T ], an example of early-stopping in SGM [1, 26] . In continuous time, this has the effect of score-matching for the mollified distribution π N,ϵ = π N ⋆ ρ ϵ , where ρ ϵ is the transition probability kernel for s = ϵ 1 . Then it can be shown thatEQUATION)3 An uncertainty quantification approach to generalization in SGMsWe describe our UQ approach to generalization in SGMs and overview our main results. Discussion of the proof methods is deferred to Sections 4 and 5. Our primary goal is to study how practical and approximation errors made in implementing SGMs translate into errors in the generative distribution relative to the true distribution.√ ρ∥ 2 2 = Ω m 0 log(m 0 ) -ρ(T ) log(ρ(T ))dx.Proof. Testing (41) against f ′ (ρ) for a smooth function f , we obtainf (ρ(T, x))dx + T 0 f ′′ (ρ)|∇ρ| 2 dxdt = f (m 0 (x))dxand the choice off (x) = x log(x) yields ρ(T, x) log(ρ(T, x))dx + T 0 |∇ρ| 2 ρ dxdt = m 0 (x) log(m 0 (x))dx,which shows the second claim. Finally, if m 0 log(m 0 ) ∈ L 1 this yields estimates on supt∈[0,T ] ∥ρ(t) log(ρ(t))∥ 1 , ∥∇ √ ρ∥ 2 2 since |ρ(T, x) log(ρ(T, x))|dx = ρ(T, x) log(ρ(T, x))dx + 2 (ρ(T, x) log(ρ(T, x))) -dx ≤ ρ(T, x) log(ρ(T, x))dx + 2e -1 |Ω|,where we used the fact that x log(x) ≥ -1 e .The rest of the subsection justifies the steps outlined in Subsection 5.2. Namely starting from a estimate of the form J (η N,ϵ , θ) ≤ e nn our goal is to obtain an estimate of the formJ (η ϵ , θ) ≤ e ′ nn .For this we use formula (42). In particular we first use the linear dependence of J L with respect to the underlying measure to bound|J L (η N,ϵ , θ) -J L (η ϵ )|.Then, instead of comparing the gradients on the difference∥∇ log(η N,ϵ )∥ 2 2 -∥∇ log(η ϵ )∥ 2 2we use Proposition B.1 from which it is enough to bound the entropy terms. The next two Propositions provide these results along with some other technical facts.Proposition B.2. Let T > 0 and π i for i = 1, 2 denote two probability measures in Ω such that ∥π i log(π i )∥ L 1 < ∞ and ρ i the corresponding solutions to (41). Then, there exists a dimensional constant C > 0 such that|J L (ρ 2 , θ) -J L (ρ 1 , θ)| ≤ CT sup t∈[0,T ] d 1 (ρ 2 (t), ρ 1 (t))∥s θ ∥ 2 C 2 ([0,T ]×Ω) ≤ CT d(π 1 , π 2 )∥s θ ∥ 2 C 2 ([0,T ]×Ω) .In the above s θ indicates the score function approximation.Proof. We recall thatJ L (ρ i , θ) = T 0 E ρ i 1 2 |s θ | 2 + div(s θ ) ds = T 0 1 2 |s θ | 2 + div(s θ )dρ i (s)ds.Therefore, ifg θ = 1 2 |s θ | 2 + div(s θ ) |J L (ρ 2 , s θ ) -J L (ρ 1 , s θ )| = T 0 1 2 |s θ | 2 + div(s θ )d(ρ 2 (s) -ρ 1 (s))ds = T 0 g θ d(ρ 2 -ρ 1 )(s)ds ≤ sup 0≤t≤T ∥∇g θ ∥ ∞ T 0 d 1 (ρ 2 (s), ρ 1 (s))ds ≤ T sup 0≤t≤T ∥∇g θ ∥ ∞ sup 0≤t≤T d 1 (ρ 2 (t), ρ 1 (t)). Since ∥∇g θ ∥ ∞ ≤ ∥s θ ∥ 2 C 2to conclude we need to show thatsup 0≤t≤T d 1 (ρ 2 (t), ρ 1 (t)) ≤ d 1 (π 2 , π 1 ).However this follows from our Lemma 3.1 applied for b 1 = b 2 = ⃗ 0 and m i = π i for i = 1, 2.Lemma B.3. Let π ϵ , πN,ϵ , 0 < δ < ϵ be as in Theorem 3.3 with ϵ < 1. There exists a dimensional constant C = C(d) > 0 such thatEQUATIONEQUATIONandEQUATIONMoreover, let C > 0 denote the constant from Proposition D.3 and T ≥ R 2 large enough such that2CR -d e -ωT R 2 ≤ 1 2 1 vol(RT d ) ,then up to a dimensional constant C > 0 we haveEQUATIONProof. Let ρ N,ϵ , ρ ϵ : [0, T ] × Ω → R, be solutions toEQUATIONfor ρ 0 = π, ρ 0 = π N respectively. Since π ϵ = Γ(ϵ) ⋆ π, πN,ϵ = Γ(ϵ) ⋆ π N we note that π ϵ = ρ ϵ (ϵ, x) and πN,ϵ = ρ ϵ,N (ϵ, x). Using the same adjoint method as in the proof of Theorem 3.1 for a function ψ : Ω → R we consider the function ϕ : [0, ϵ] × Ω → R given byEQUATIONTesting against ϕ in the equation for λ = ρ ϵ -ρ ϵ,N yields(π N,ϵ -π ϵ )ψ(x)dx = ϕ(0, x)d(π -π N )(x).First we consider ∥ψ∥ ∞ ≤ 1 which from (31) , yields∥∇ϕ(0)∥ ∞ ≤ ∥ψ∥ ∞ √ ϵ ≤ C √ ϵ .Thus, for any ∥ψ∥ ∞ ≤ 1 we obtain(π N,ϵ -π ϵ )ψ(x)dx ≤ Cd 1 (π, π N ) √ ϵ , which implies ∥π N,ϵ -π ϵ ∥ L 1 (ω) ≤ Cd 1 (π, π N ) √ ϵ .Next, by considering Lip(ψ) ≤ 1 which implies∥∇ϕ(0)∥ ∞ ≤ 1 we obtain d 1 (π N,ϵ , π ϵ ) ≤ d 1 (π, π N ). For estimate (45) we note that πN,ϵ log(π N,ϵ ) -π ϵ log(π ϵ ) = 1 0 (1 + log(sπ N,ϵ + (1 -s)π ϵ ))ds(π N,ϵ -π) ≤ (1 + | log(δ)|)|π N,ϵ -π ϵ |. Therefore, ∥π N,ϵ log(π N,ϵ ) -π ϵ log(π ϵ )∥ L 1 (Ω) ≤ (1 + | log(δ)|)∥π N,ϵ -π ϵ ∥ L 1 ≤ C(1 + | log(δ)|) d 1 (π N , π) √ ϵ .For estimate (46) we note that for any convex function hh(x) -h(y) ≥ ∇h(y) • (x -y) =⇒ h(y) -h(x) ≤ ∇h(y) • (y -x).Hence by applying the above to h(x) = x log(x) we obtain log(η N,ϵ (T ))η N,ϵ (T ) -log η(T ))η(T )dx ≤ (1 + log(η N,ϵ (T )))d(η N,ϵ (T ) -η ϵ (T ))≤ ∥1 + log(η N,ϵ (T ))∥ ∞ ∥η N,ϵ (T ) -η ϵ (T )∥ L 1 .We now note that for any g :RT d → R, ∥g∥ ∞ ≤ 1 g(x)d(η N,ϵ -η ϵ )(T ) = Γ(T ) ⋆ gd(π N,ϵ -π ϵ ) ≤ ∥∇Γ(T ) ⋆ g∥ ∞ d 1 (π N,ϵ , π ϵ ) ≤ C √ T d 1 (π N,ϵ , π ϵ ).Taking the supremum yields∥η N,ϵ (T ) -η ϵ (T )∥ L 1 ≤ C √ T d 1 (π N,ϵ , π ϵ ) ≤ d 1 (π N,ϵ , π ϵ ).Finally, from Proposition D.3 above applied to ρ = η N,ϵ -1 vol(RT d ) we have η N,ϵ (T, x) ≥ -∥ρ(T )∥ ∞ + 1 vol(RT d ) ≥ 1 2 1 vol(RT d ) as well as η N,ϵ (T, x) ≤ 1 vol(RT d ) + ∥ρ(T )∥ ∞ ≤ 2 1 vol(RT d )by our assumptions on T > 0. Therefore,∥1 + log(η N,ϵ (T ))∥ ∞ ≤ C(1 + d log(R))and the result follows.Remark B.4. It is worth noting that the empirical sample π N could have been regularized by any smooth kernel ρ. However the choice of the heat kernel which provides the same result as early stopping behaves quite nicely with respect to the metric d 1 since it is a contraction. This is evident in Lemma B.3 above where we see that we do not pay any additional cost due to the mollification in (43).In fact the above techniques show the following general result regarding ESM to DSM bounds. Combining the above we can state the following general result, about ESM to DSM bounds. Theorem B.5. Let π 1 , π 2 ∈ P(Ω) denote two probability densities, such thatEQUATIONand δ > 0 such that π i (x) ≥ δ. Define ρ i : [0, T ] × Ω → R as the solutions toEQUATION)Then, up to a dimensional constant C = C(d) > 0EQUATIONB.1 Proof of Theorem 3.3Now we can prove Theorem 3.3.Proof. From Proposition B.1 we obtainJ (η ϵ , s θ * ) = J L (η ϵ , s θ * ) + 2∥∇ √ η ϵ ∥ 2 2 = J (η N,ϵ , s θ * ) + 2 ∥∇ √ η ϵ ∥ 2 2 -∥∇ η N,ϵ ∥ 2 2 ∥ + J L (η ϵ , s θ * ) -J L (η N,ϵ , s θ * ). By assumption we have the boundJ (η N,ϵ , s θ * ) ≤ e nn . By Proposition B.1 2 ∥∇ √ η ϵ ∥ 2 2 -∥∇ η N,ϵ ∥ 2 2 = Ω π ϵ (x) log(π ϵ (x))-π N,ϵ (x) log(π N,ϵ (x))dx- Ω η ϵ (T ) log(η ϵ (T ))-η N,ϵ (T ) log(η N,ϵ (T ))dx.From Lemma B.3 we have that∥π ϵ log(π ϵ ) -πN,ϵ log(π N,ϵ )∥ L 1 (Ω) ≤ C(1 + | log(δ)|) d 1 (π N , π) √ ϵ .Using the same argument and the fact that from maximum principle η ϵ (T ), η N,ϵ ≥ δ > 0, we can show∥η(T ) log(η(T )) -η N,ϵ log(η N,ϵ (T ))∥ L 1 (Ω) ≤ (1 + | log(δ)|)∥η ϵ (T ) -η N,ϵ (T )∥ L 1 (Ω) ≤ C √ T d 1 (π ϵ , πN,ϵ ) ≤ C √ T d 1 (π, π N ).Furthermore, by Proposition B.2 we haveJ L (η ϵ , s θ * )-J L (η N,ϵ , s θ * ) ≤ CT d 1 (π N,ϵ , π ϵ )∥s θ ∥ 2 C 2 ([0,T ]×Ω) ≤ CT d 1 (π N , π)∥s θ ∥ 2 C 2 ([0,T ]×Ω) .Putting everything together we haveJ (η, s θ * ) ≤ e nn + C 1 + | log(δ)| √ ϵ + 1 √ T + T ∥s θ ∥ 2 C 2 ([0,T ]×Ω) d 1 (π N , π).As an immediate Corollary we have the following. Corollary B.6. Using the same notation as in Section 3.5, if e nn > 0 is such that 0 ≤ J (η N,ϵ , θ * ) < e nn (52)EQUATIONwhere for a dimensional constantEQUATIONRemark B.7 (On regularity implications of NN approximation). It is important to note that the norm of our NN approximation ∥∇s θ ∥ ∞ is in fact dependent on ϵ > 0 as we can see from the following∥s θ ∥ C 1 ≳ ∥s θ ∥ L 2 (η N,ϵ ) ≳ ∥∇ log(η N,ϵ )∥ L 2 (η N,ϵ ) - √ e nnand ∥∇ log(η N,ϵ )∥ L 2 (η N,ϵ ) blows up as ϵ → 0.We finish this subsection of technical results with an observation about the implications on the existence of a smooth NN-approximation. Proposition B.8 (ESM implies regularity). Let π, s θ and e nn > 0 be as above for Ω = RT d . Assume moreover, that ∥∇s θ ∥ ∞ < ∞. Then, π = π(x) admits a density and in fact ∥π log(π)∥L 1 (Ω) < ∞. Proof. Since ∥∇s θ ∥ ∞ < ∞, for some constant C = C(R, d, e nn ) > 0 it holds 2∥∇ √ η π ∥ 2 L 2 = T 0 |∇ log(η π )| 2 dη π (s, x)ds ≤ C.Moreover, by the regularizing properties of the diffusion we have Using the previous results we can in fact show a result about the average error in Theorem B.5. The main observation here is that by Jensen's inequality when we take expectation with respect to the sample π N we no longer require a lower bound on our densities. However, as we will see we require a very restrictive assumption on the norm of our score function approximation Theorem C.1. (Average DSM generalization) Let e nn , A > 0 and assume that for each sample πN from π there exists a s θ * such thatEQUATION)Let m g (T ) be the generated distribution from πN (which is also random since πN is random). Let C > 0 be the dimensional constant appearing in Proposition D.3 andT ≥ R 2 large enough such that 2CR -d e -ωT R 2 ≤ 1 2 1 vol(RT d ) , then E d 1 (π, m g (T )) ≤ CR 3 2 (1 + √ A)(R 2 e -ωT R 2 + T e ′ nn ), wheree ′ nn = ϵ + e nn + C N 1 2d T RA 2 + 1 √ T (1 + d log(R)) .Remark C.2. It is important to emphasize in the above that we are assuming that we may pick both an error e nn > 0 and a bound A > 0 on the gradient of s θ globally for all random samples πN . Moreover, s θ is a weak approximation of ∇ log(η N,ϵ ) which although smooth for all ϵ > 0 it does exhibit a blow-up as ϵ > 0 tends to zero. Therefore, since by triangular inequality up to a constantC = C(R, d) > 0 CA ≥ ∥s θ ∥ C 2 ≥ ∥s θ ∥ L 2 (η N,ϵ ) ≥ ∥D η N,ϵ ∥ L 2 -e nnthe constant A > 0 will also exhibit a blow-up which can make the approximation worse. Therefore there could be a potential trade-off between NN approximation e nn becoming small and the value of the norm A.First we require a preliminary lemma. Lemma C.3. Under the same assumptions as in Theorem C.1 we have thatEQUATIONProof. We note thatJ (η ϵ , s θ * ) = J L (η ϵ , s θ * ) + 4∥D √ η ϵ ∥ 2 2 = J (η N,ϵ , s θ * ) + 2 ∥D √ η ϵ ∥ 2 2 -∥D η N,ϵ ∥ 2 2 ∥ + J L (η, s θ * ) -J L (η N,ϵ , s θ * ) .By assumption we have the boundJ (η N,ϵ , s θ * ) ≤ e nn .By Proposition B.12 ∥D √ η ϵ ∥ 2 2 -∥D η N,ϵ ∥ 2 2 = Ω π ϵ (x) log(π ϵ (x)) -πN,ϵ (x) log(π N,ϵ (x))dx - Ω η ϵ (T ) log(η ϵ (T )) -η N,ϵ log(η N,ϵ (T ))dx.Taking expectations with respect to the sample πN we getE Ω π ϵ (x) log(π ϵ (x))-π N,ϵ log(π N,ϵ )dx = Ω π ϵ (x) log(π ϵ (x))-E πN,ϵ (x) log(π N,ϵ (x)) dx ≤ Ω π ϵ (x) log(π ϵ (x)) -E πN,ϵ (x) log(E πN,ϵ )dx = 0.In the above we used Jensen's inequality along with the fact thatE πN,ϵ = π ϵ .Moreover, by Lemma B.3 we have thatE η N,ϵ log(η N,ϵ (T )) - Ω η ϵ (T ) log(η ϵ (T ))dx ≤ C √ T (1 + d log(R))E d 1 (π, πN ) ≤ C √ T N 1 d (1 + d log(R)).Finally, from Proposition B.2 we haveE J L (η, s θ * ) -J L (η N,ϵ , s θ * ) ≤ CT E ∥s N θ ∥ 2 C 2 d 1 (π, πN ) ≤ CT N 1 d A 2Thus,E J (η ϵ , θ * ) ≤ e nn + C N 1 d T A 2 + 1 √ T (1 + d log(R)) ,and the result follows.• Case 1: Assume that 0 ≤ s 0 < T and we will show that for C > 0 large enough this leads to a contradiction. At (s 0 , x 0 ) we have the optimality conditions ∇z(s 0 , x 0 ) = C∇u(s 0 , x 0 ) ∆z(s 0 , x 0 ) ≤ C∆u(s 0 , x 0 ) ∂ t z(s 0 , x 0 ) ≤ C∂ t u(s 0 , x 0 ). Hence,C(-∂ t u(s 0 , x 0 ) -α(s 0 )∆u(s 0 , x 0 ) + b • ∇u(s 0 , x 0 )) ≤ (-∂ t z(s 0 , x 0 ) -α(s 0 )∆z(s 0 , x 0 ) + b • ∇z(s 0 , x 0 ))and using the respective equations for u, z we obtain that at (s 0 , x 0 )-C α(s 0 ) 2 |∇u| 2 ≤ -α(s 0 )|∇ 2 u| 2 -⟨∇u • ∇b, ∇u⟩ -α(s 0 )∇u • ∇z.Using the optimality conditions yieldsα(s 0 ) C 2 |∇u| 2 + α(s 0 )|∇ 2 u| 2 ≤ -⟨∇u • ∇b, ∇u⟩ ≤ ∥∇b∥ ∞ |∇u| 2 , which is a contradiction if C ≥ 2∥∇b∥ ∞ M . Thus for C 0 = 2∥∇b∥ ∞ Mthe function z -C 0 u achieves its maximum at some point (T, x 0 ).• Case 2: Assume that s 0 = T . Then, for every(t, x) ∈ [0, T ] × Ω we have z(t, x) ≤ w(T, x 0 ) + Cu(t, x) ≤ 1 2 |∇v| 2 + C∥v∥ ∞ + C∥u + ∥ ∞and by maximum principle we also have that ∥u∥ ∞ ≤ ∥v∥ ∞ .Thus, up to a universal constant|∇u(t, x)| 2 ≲ (1 + ∥∇b∥M )∥v∥ C 1 ,which proves estimate (58). Now we prove estimate (59). The proof follows the previous strategy only this time we consider the function w(t, x) = (T -t)z -Cu, where again C > 0 will be determined later. Assume that the maximum occurs at some point (s 0 , x 0 ).1. Case 1: Assume that s 0 < T . Looking again the corresponding optimality conditions and using the equations for u and z we have at (s 0 , x 0 )Cα(s 0 ) 2 -1 |∇u| 2 ≤ -(T -s 0 )⟨∇u • ∇b, ∇u⟩ ≤ T ∥∇b∥ ∞ |∇u| 2 ,which is a contradiction if C ≥ 2T M ∥∇b∥ ∞ + 2M . Assume moreover that for some M > 0 0 < 1 M ≤ α(t).Then up to a dimensional constant C > 0 we have the following estimatesEQUATIONProof. Since ψ ≥ 1 by the Maximum Principle we have that1 ≤ ϕ(t, x) ≤ ∥ψ∥ ∞ .Thus we can define u : [0, T ] × Ω → R by the formula ϕ(t, x) = e -u(t,x) 2 ⇐⇒ u(t, x) = -2 log(ϕ(t, x)),We then have The following result is classical and establishes exponential rate of convergence in d 1 for solutions to the heat equation in the uniform distribution on the torus. We provide a proof for the readers convenience.EQUATIONProposition D.3. Let Ω = RT d and m i ∈ P(Ω), i = 1, 2 be two probability measures. Define ρ i : [0, ∞) × Ω → [0, ∞) by ∂ t ρ i -∆ρ i = 0 in (0, ∞) × Ω, ρ i (0) = m i in Ω.Then, there exists constants C = C(d) > 0, ω = ω(d) > 0 such that d 1 (ρ 2 (t), ρ 1 (t)) ≤ CRe -ωt R 2 d 1 (m 2 , m 1 ) for t ≥ 0.In particular by choosing m 2 =EQUATIONProof. First we assume that R = 1 and recall that if Γ(t, x) denotes the heat kernel on the unit torus T d then a constant C = C(d) > 0 we have the normalizing effect|∇Γ(t) ⋆ g| ≤ C ∥g∥ ∞ √ t .Fix ψ : T d → R such that Lip(ψ) ≤ 1 and ψ(0) = 0. Such a function ψ will also satisfy ∥ψ∥ ∞ ≤ 1. Repeating the same argument as above with initial conditions m i = ρ i (T ), we obtaind 1 (ρ 2 (2T ), ρ 1 (2T )) ≤ C √ T d 1 (ρ 2 (T ), ρ 1 (T )) ≤ C √ T 2 d 1 (m 2 , m 1 )and by induction for all k ∈ N we haved 1 (ρ 2 (kT ), ρ 1 (kT )) ≤ C √ T k d 1 (m 2 , m 1 ) = C k T -k 2 d 1 (m 2 , m 1 ).Let θ ≥ 0 and fix a k ∈ N such that kT ≤ θ < (k + 1)T.Since the heat operator defines a contraction in d 1 we have that the result follows.Here, ⋆ denotes convolution.(1 + ∥∇s θ ∥ ∞ ) Re -ωT R 2 d 1 π, 
source of errors in score-based generative modeling: 
 We attribute errors of SGM to the following six sources:• Data distribution is only accessible via samples e 1 : The target distribution is only known through a finite set of samples. In practice, the regularity of the score function and data distribution, such as Lipschitzianity or whether the distribution has a density, is unknowable.• Choice of score-matching objective e 2 : In practice, score-matching objectives are approximated via samples. The DSM objective is most frequently used as it avoids computing derivatives of the score function. While DSM and ISM are equivalent given the exact density evolution η(x, s), they differ when approximated with finite samples. Previous analysis typically assumes that the learned score function is close to the true score function within some L 2 distance, i.e., a ball determined by J ESM . In contrast to previous work [11, 7, 5] , we show in Theorem B.5 how training through DSM translates into a bound for ESM. Moreover we prove Proposition B.8 which states that if there exists a neural net that well-approximates the true score then that the target density is necessarily regular.• Expressiveness of the score function e 3 . The expressivity of neural net approximations to the score function will depend on the particular expressivity of the parametrization.• Choice of reference distribution e 4 . With access to the exact score function, the denoising process (2) produces trajectories that sample from π only if the initial condition starts at η π (•, T ). In practice, however, the initial distribution is usually a Gaussian approximation of η π or, in the case of the OU noising process, its corresponding stationary distribution.• Early stopping of the denoising process e 5 . In practice, the score-matching objective is integrated from s ∈ [ϵ, T ], for small ϵ, instead of s = 0 [1, 26] . This prevents the SGM from memorizing the training data [27, 28] , and is equivalent to running the denoising process for t ∈ [0, T -ϵ]. Early stopping is crucial when the data distribution is supported on a low-dimensional manifold, where it has no density with respect to Lebesgue measure. Previous analyses of SGM often adjust ϵ to optimize generalization bounds [6, 20] .• Discretization error e 6 . The denoising SDE must be solved via a numerical scheme. Previous work [5, 11, 7, 8] have considered the impact of discretization error on the generalization abilities of SGM. While our analysis does not consider the discretization error, our approach can be extended to study discretization error through the use of modified equations [29] .Let d be an integral probability metric (IPM) between measures ν 1 , ν 2 ∈ P(Ω)EQUATIONfor some function space X . We assume X = {ψ : Ω → R, ∥∇ψ∥ ∞ ≤ 1} or {ψ : Ω → R, ∥ψ∥ ∞ ≤ 1}, corresponding to the Wasserstein-1 (d 1 ) and total variation distances, respectively. Recall that if ν 1 , ν 2 have densities, then their TV distance is equivalent to the L 1 distance between their densities.Let π be the target data distribution and take the stationary distribution of the noising process on the torus m g (0) = 1 vol(RT d ) as the initial condition. Given a learned score function s θ , let m g (T ) be the generative distribution produced by SGM. We study how IPMs between π and m g (T ) relate to errors from the first five sources of errors discussed above, i.e., we seek bounds of the form d(m g (T ), π) ≤ F(e 1 , e 2 , e 3 , e 4 , e 5 ).Our key insight is that for score-based generative modeling, we can derive bounds of the form (6) via a model-form uncertainty quantification bound for the generative Fokker-Planck equation. The Wasserstein Uncertainty Propagation theorem formalizes this insight.The WUP theorem is a general statement about how L2 neighborhoods in the space of drift functions map to neighborhoods in P(Ω) defined by IPMs. Theorem 3.1 (Wasserstein Uncertainty Propagation).Let Ω = RT d or R d . Let b 1 , b 2 : [0, T ] × Ω → R d be given with ∥∇b 1 ∥ ∞ < ∞, and m 1 , m 2 ∈ P(Ω). If m i for i = 1, 2 are given by ∂ t m i -∆m i -div(m i b i ) = 0, m i (0) = m i (7)then, up to a universal constant C > 0, we have the following:EQUATIONandEQUATION•EQUATIONEquation ( 8) is a particularly notable result as we bound the TV distance between m 1 (T ) and m 2 (T ) in terms of a weaker d 1 metric between m 1 and m 2 . This is due to the regularizing effects of diffusion processes, which is showcased in the proof. See Section 4 and Section A.1 for full details.We use WUP to produce generalization bounds when the only errors are due to the choice of the initial condition and the approximation of the score function with respect to the ESM objective. Theorem 3.2 (ESM bounds). Let π ∈ P(Ω) where Ω = RT d for some R > 0. Moreover let e nn > 0.Assume that the learned score function s θ is such that J (η π , θ) ≤ e nn . Then for b θ = s θ (T -t, x) the generated distribution m g (T ) ≈ π whereEQUATIONsatisfiesd 1 (π, m g (T )) ≤ CR 1 vol(RT d )Error from choice of reference measure (e4)+ √ e nn.(EQUATIONIf the stronger estimate supEQUATION)Here, the constants C, ω depend only on the dimension d.Applying WUP for b 1 and b 2 to be the true and learned score function, respectively, with Proposition D.3 on the convergence of the noising process to the stationary measure immediately yields the estimates above. Notice that the TV estimate ( 13) is comparable to Theorem 2 of [5] , which also assumes a uniform-in-time L 2 -accurate score function. Again, notice that we are able to bound the TV distance between m g (T ) and π using the weaker d 1 distance between π and 1 vol(RT d ) .In practice, the score is learned through samples using the DSM objective with an early stopping parameter [1, 26] . SGM is effective at producing samples from distributions supported on (or near) low dimensional manifolds. Our d 1 generalization bound describes how early stopping aids in generalization.Theorem 3.3. (Pointwise DSM generalization) Let b θ = s θ (T -t, x) and m g : [0, T ]×RT d → R begiven by (11) . Assume the learned score function is such thatJ DSM (η N,ϵ , θ) = J (η N,ϵ , θ) < e nn . Let 0 < δ < ϵ be such that δ ≤ πN,ϵ , δ < π ϵ . Then up to a dimensional constant C = C(d) > 0, d 1 (π, m g (T )) ≲ √ ϵEarly stopping (e5)EQUATIONwheree ′ nn ≲ √ e nnDSM score function error (e3)+ 1 + | log(δ)| √ ϵ + T ∥s θ ∥ 2 C 2 ([0,T ]×Ω) d 1 (π N , π)Finite sample error (e1)Choice of SM objective (e2).The pointwise DSM generalization bound applies to every finite training sample of size N . A crucial part of this theorem relates the error in the practical DSM objective function to the ESM error needed to apply the WUP theorem. We connect the DSM objective early stopping to the ESM error with the mollified distribution π ϵ in Lemma B.5. This result is agnostic to the manifold hypothesis for π, as long as ϵ > 0. Such bounds will blow up if the KL divergence were used instead. In fact, previous results that use the KL divergence to bound the TV distance [5, 19, 11, 6] will produce vacuous bounds for the d 1 distance under the manifold hypothesis as their d 1 generalization bounds are derived by bounding the KL divergence. Remark 3.4 (Density lower bound). Similar to [6] , our DSM generalization bound relies on a density lower bound assumption. We note however, that our DSM generalization bound holds for any random sample {z i } N i=1 ∼ π. This density lower bound assumption can be removed via Jensen's inequality if we consider the expected d 1 distance between π and m g (T ) over random empirical measures. See Theorem C.1 and its proof in Section C. Remark 3.5 (Trade-offs and memorization). Theorem 3.3 captures trade-offs when training SGMs and memorization effects. To minimize the error from early stopping, we can let ϵ → 0. However, empirically, without early stopping, SGMs overfit to the kernel approximation and memorize the training data [27, 28, 25] . The bound is vacuous when ϵ = 0 regardless of whether the distribution lies on a low-dimensional manifold. As ϵ → 0, training the DSM to be small implies that the score function must approximate a rough function with large Lipschitz constant, which will increase the bound. This shows that overtraining the DSM objective may not necessarily yield a better generative model. Moreover, suppose that π N = π, then as ϵ → 0 and e nn → 0, we have that d 1 (π, m g (T )) → 0, indicating the model memorizes the training data. Our results corroborate those of [25] . 
regularity theory of hamilton-jacobi-bellman pdes enables uncertainty quantification in sgms: 
 A recent result in [30] established connections between generative flows with the theory of PDEs, more specifically the theory of mean field games. We continue investigating these connections and showcase how one may study generative modeling algorithms by obtaining stability estimates for the Fokker-Planck equation. We provide a proof sketch for our Wasserstein Uncertainty Quantification theorem (Theorem 3.1). Our strategy involves (1) constructing test functions for the IPMs using the Kolmogorov backward equation, (2) choosing the suitable function space for the terminal data depending on the desired IPM, and (3) obtaining gradient estimates of the test functions via Bernstein estimates. The theorem relies on the gradient of the test function remaining bounded for t ∈ [0, T ), which is guaranteed by the regularizing properties of diffusion processes. See A.1 for full details about the proof.From (7) , we aim to compute bounds for d(mEQUATION)Let ϕ : [0, T ] × Ω → R be a test function in space and time. We integrate in space and time against the PDE ( 16) and integrate by parts to move the derivatives on to ϕ which yieldsEQUATION+ T 0 Ω m 2 ∇ϕ • (b 1 -b 2 )dxdt = 0Notice that if we choose the test function ϕ to satisfy the Kolmogorov backward equation (KBE)-∂ t ϕ -∆ϕ + b 1 • ∇ϕ = 0 in [0, T ) × Ω, ϕ(T, x) = ψ(x)in Ω (18) with terminal condition ψ ∈ X , then from (17), we haveEQUATION)The equality above is valid for any terminal condition ψ. Depending on the choice of function space X for ψ, we obtain bounds on different IPMs. Taking the supremum over X we haveEQUATION)Recall that ϕ is related to ψ via the KBE (18) To obtain bounds for d, we need to bound the two terms in (20) , which depend on the choice of function space X and assumptions on the drift terms.We split our proof sketch into two parts; the first part focuses on deriving L 1 estimates, while the second derives d 1 estimates.L 1 estimates. We first note that because of the regularizing properties of (18) we can obtain bounds on Ω λ(0, x)ϕ(0, x)dx with weaker norms on ϕ. Notice thatEQUATIONwhere Y denotes a generic space of functions, and Y ′ is its dual. In what follows, we show that regularizing effects of (18) takes functions in X to functions with more regularity Y such that Y is compactly embedded in X . This in turn implies that ∥ • ∥ Y ′ is weaker than d and so we are able to bound the stronger norm d by the weaker norm.• To prove (8), wed = ∥ • ∥ L 1 (Ω) , which corresponds with X = {ψ : Ω → R, ∥ψ∥ ∞ ≤ 1}.Observe that we can take Y = {ψ : Ω → R, ∥∇ψ∥ ∞ ≤ 1}, in which case we obtainΩ λ(0, x)ϕ(0, x)dx = ∥∇ϕ(0, x)∥ ∞ Ω (m 1 -m 2 ) ϕ(0, x) ∥∇ϕ(0, x)∥ ∞ dx (22) ≤ d 1 (m 1 , m 2 )∥∇ϕ(0, x)∥ ∞ .Notice that this bound crucially depends on showing that ϕ(0, x) is Lipschitz. • To prove (9), we can simply take X = Y, and obtainEQUATIONFor ( 8) and ( 9), Cauchy-Schwarz on the spatial integral shows the second term can be bounded asEQUATION≤ sup 0≤t≤T ∥(b 1 -b 2 )(t)∥ L 2 (m 2 (t)) T 0 ∥∇ϕ(t, •)∥ ∞ dt.Notice here that it is crucial to produce estimates for ∇ϕ.Wasserstein-1 (d 1 ) estimates. To prove (10) , we choose X = Y = {ψ : Ω → R, ∥∇ψ∥ ∞ ≤ 1}, the space of 1-Lipschitz functions. We obtain (22) again, except the terminal data ψ also has Lipschitz constant equal to 1. Applying Cauchy-Schwarz in space and time implies the second term is boundedEQUATION≤ T sup 0≤t≤T ∥∇ϕ(t, •)∥ ∞ ∥b 1 -b 2 ∥ L 2 (m 2 ) .We again see that we require estimates for ∇ϕ, and we need them to be bounded.To obtain estimates for ∇ϕ, we could differentiate (18) to derive a PDE for ∇ϕ. The resulting PDE, however, will have ∂ t (∇ϕ) grow linearly with ∇ϕ, and so if we apply (reverse) Gronwall's inequality, the resulting estimates for ∇ϕ will grow exponentially in time. To avoid this exponential time dependence, we first perform a Hopf-Cole transform u = -2 log ϕ on (18) to derive the HJB equation for uEQUATIONHere we provide an example of classical Bernstein estimates (Proposition D.1 and Corollary D.2) to obtain bounds for ∇ϕ without using Gronwall's inequality [10] . The main idea is to derive a PDE (60) for the function z = 1 2 |∇u| 2 by taking the gradient of ( 26) and then taking the inner product with ∇u. Then by showing that the function w(t, x) = z -Cu for sufficiently large C attains its maximum at (T, x 0 ), we can show thatEQUATION)Using the maximum principle for (18) , we find∥u∥ ∞ = ∥2 log ψ∥ ∞ , yielding z(t, x) ≤ C∥ log ψ∥ C 1 .Assuming that ψ is Lipschitz continuous implies boundedness of z and therefore ∇ϕ for all time. A similar result holds when ψ ∈ L ∞ only (see ( 59)). Detailed proofs and related boundes are provided in Proposition D.1. Applications of the bounds to derive the WUP is provided in Section A.1. Remark 4.1 (The regularizing role of stochasticity). In our analysis, the stochasticity in SGMs provides two types of regularizing effects. The first is early stopping, which adds a small amount of Gaussian noise of the data distribution. This, in effect, is equivalent to running the noising process for a short amount of time and immediately mollifies the initial empirical distribution so that it has a smooth density. Second, the Laplacian is the key mechanism that regularizes the test function in (18) [9] , which then allows us to bound, for example, the stronger ∥ • ∥ L 1 norm by the weaker d 1 -norm. Recall that in PDE theory, the stochasticity of a generative flow manifests as the Laplacian operator in the Fokker-Planck and Kolmogorov backward equations [9] . Without the Laplacian the regularizing the test functions in (18) would not be possible in general and, in fact, we would not have access to long time behavior results.5 Proof sketches -Score-based generative models are robust to errorsWe now provide sketches of the proofs for the main generalization bounds in Theorems 3.2 and 3. While Theorem 3.2 has no explicit assumptions on π, the assumption that the approximate score function s θ is close in L 2 to the true score implicitly implies regularity of π. Specifically, if s θ is Lipschitz (which is true for neural networks in practice) and satisfies J (π, θ) < e nn , then π must have finite entropy. This implies that π necessarily has a density. See Proposition B.8 for the formal statement and proof. 
theorem 3.3: pointwise dsm generalization bound: 
 In practice, the score function is learned via a Monte Carlo approximation of the DSM objective (2). To avoid overfitting to the kernel estimator and memorizing the dataset [25] , the time integral is taken only over s ∈ [ϵ, T ] where ϵ is the early stopping parameter. This is equivalent to training with the ESM objective with the true score replaced with the kernel approximation at time ϵ [20, 27] . To derive generalization bounds of SGMs trained via DSM, we establish the relationships between (1) the mollified distribution π ϵ = Γ(ϵ) ⋆ π and the true distribution π, and (2) the DSM objective J DSM (η N,ϵ , θ) = J (η N,ϵ , θ) and the ESM objective with respect to the score function of the mollified distribution J (η π ϵ , θ). Formally, this strategy involves bounding the following two termsd 1 (π, m g (T )) Generalization bound w.r.t true distribution ≤ d 1 (π, π ϵ )Early stopping error (e5)+ d 1 (π ϵ , m g (T )).Generalization bound w.r.t. mollified distribution (28) For the early stopping error, we use the regularizing properties of the heat equation to show thatd 1 (π, π ϵ ) ≤ C √ ϵ, where C only depends on the dimension d. This implies that, measured in d 1 , early stopping only incurs a nominal C √ ϵ error even if the π does not admit a density. This result would not possible if we were to study generalization error in terms of KL or TV directly.A bound for the second term is obtained by comparing the ESM objective value between the learned score function and the true score function of the early stopped distribution, J (η π ϵ , θ) to the DSM objective J (η N,ϵ , θ). We present Theorem B.5 and its corollary, which under the assumption that η N,ϵ and η π ϵ have a lower bound δ, state that if J (η N,ϵ , s θ ) < e nn , then J (η π ϵ , s θ ) < e ′ nn withEQUATIONThe main idea behind Theorem B.5 is to note that the difference between the ESM and DSM objective functions can be written as a difference in ISM objective functions plus the entropy difference between η N,ϵ and η π ϵ . Details for bounding this term is provided in Proposition B.2 and Lemma B.3.To arrive at the final result, apply the WUP theorem to derive generalization bounds for d 1 (π ϵ , m g (T )) under the assumption that J (η π ϵ , θ) < e ′ nn , along with (29) . Finally, combine this result with the error due to the early stopping d 1 (π, π ϵ ). Full details of this proof is provided in Section A.3. 
discussion: pde regularity theory and uq for generative modeling: 
 Our main contribution is the study of generalization in score-based generative models from the perspective of uncertainty quantification. The regularity theory of nonlinear PDEs is the key technical tool that produces our results. We emphasize that the tools we use here can be used generative models beyond SGMs, and that we have not pushed our analysis to the limits of our tools in SGMs. Moreover, we also emphasize some downstream UQ applications for SGMs that may be of future interest.A surprising result of our work is deriving bounds for the stronger L 1 distance in terms of the weaker d 1 distance (see (22) and Theorem 3.2). The key insight (Equation ( 21)) is that the evolution of observables defined by the KBE (18) regularizes the test function. We have not fully exploited the regularizing effects of (18) in this paper, as we only focused on L 1 and d 1 estimates. Improved bounds in Sobolev spaces H s . To illustrate other extensions and choices of Y in (21) , observe that in the trivial case when b 1 = b 2 = 0, (19) 2 are not identically zero we still expect such estimates, though they will depend on the regularity of b 1 . To highlight the importance of regularizing effects, note that by [31] , if π ∈ P(T d ), then for the empirical measure π N , we expectsimplifies to λ(T, x)ψ(x)dx = λ(0, x)(Γ(T ) ⋆ ψ)(x)dx, where Γ is the heat kernel. If ∥ψ∥ ∞ ≤ 1, then Γ ⋆ ψ(T ) ∈ C ∞ and we have estimates of the form d 1 (m 1 (T ), m 2 (T )) ≤ C(s, T, d)∥m 1 -m 2 ∥ H -s for all s ∈ N. When b 1 , bd 1 (π, π N ) ≲ 1 N 1 d . However if s > d 2 , ∥π -π N ∥ H -s ≲ 1 √ N .This suggests that improved regularity may influence overcoming the curse of dimensionality.Computing expectations with respect to posterior distributions is a key task in Bayesian inference. Generative modeling, in particular, has a key role in future developments of likelihood-free inference [32, 33, 34, 35] . For generative models to be trustworthy for inference, they must to be shown to be robust. The WUP theorem provides error bounds for approximating expectations with respect to some true unknown distribution, and may be significant for SGMs in likelihood-free inference.For example, suppose we wish to estimate E π h for some distribution π and observable h, and an SGM m g (T ) is constructed to approximate the expectation. Bounds of the form (6) , such as the WUP theorem, translate into guarantees on the expectations:E π h -E mg(T ) h ≤ d(m g (T ), π) ≤ F(e 1 ,(30) In the context of SGM, the inequalities are a posteriori bounds, meaning they can be computed after learning the model. Additional regularity on h may yield improved guarantees.Robust UQ methods are based on the perspective that learning any complex model will typically involve multiple sources of uncertainty due to modeling choices, model reduction or learning from imperfect data. These uncertainties are not just in parameters but are inherently present in the mathematical model itself and will propagate to any predictions. There is substantial related work in recent years using a distributional robustness perspective, [12] , using divergences or probability metrics and their variational representations to quantify the impact of model uncertainty around a baseline model that may be either learned (e.g. a generative model) or could be just an empirical distribution from an unknown true distribution. The approach can generally be described as quantifying an uncertainty set around the baseline model that the worst-case distribution belongs to via some neighborhood defined in terms of a probability divergence [13, 14] , a Wasserstein distance [15] or maximum mean discrepancy [16] . There are, however, drawbacks to each of these approaches: a divergence ball contains only distributions with the same support as the baseline, while the uncertainty set may be hard to determine practically. The WUP Theorem is a related robust UQ notion where the uncertainty ball is in an IPM, e.g. 1-Wasserstein, MMD, or TV. However, it allows us to bypass the robust UQ for stochastic processes relying on restrictive, path-space probability divergence-based approaches or Girsanov's Theorem, [36, 37] . Furthermore, the WUP Theorem bounds use PDE theory to provide a computable uncertainty set, as we also demonstrate in the case of DSM, see Theorem 3.3.Throughout the proofs we will make frequent use of the fact that if Γ : [0, ∞) × Ω → R denotes the heat kernel on Ω = R d or Ω = RT d , then there exists a dimensional constant C = C(d) > 0 such thatEQUATIONFurthermore, recall that constants are subject to change from line to line but always maintain a variable dependence as in the corresponding statement.A.1 Theorem 3.1The proof of Theorem 3.1 follows the strategy outlined in Section 4 and uses the different gradient estimates from the HJB equations provided in Section D.1.Proof. Let λ = m 1 -m 2 , which satisfiesEQUATION)For a fixed ψ : Ω → R, for which we will specify bounds latter, let ϕ : [0, T ] × Ω → R be given byEQUATION)Testing against ϕ in (32) and using (33) , we obtainEQUATION• Proof of ( 8) and ( 9): Assume that ∥ψ∥ ∞ ≤ 1 which then implies-1 ≤ ϕ ≤ 1.Since dλ(t) = 0 for all t ∈ [0, T ] up to adding a constant we can assume without loss of generality in (34) , that 1 ≤ ψ ≤ 3 and 1 ≤ ϕ ≤ 3. Then from estimate (63) in Corollary D.2, we have(T -t)∥∇ϕ(t)∥ 2 ∞ ≤ C(T ∥∇b∥ ∞ + 1) =⇒ ∥∇ϕ(t)∥ ∞ ≤ C T ∥∇b 1 ∥ ∞ + 1 √ T -t .Next, we note thatλ(T, x)ψ(x)dx = λ(0)ϕ(0, x)dx - T 0 m 2 (s)∇ϕ(s) • (b 2 -b 1 )(s)dxds ≤ d 1 (m 1 , m 2 )∥∇ϕ(0)∥ ∞ + T 0 |∇ϕ(s)| 2 m 2 (s)dx 1 2 m 2 (s)|b 2 -b 1 | 2 (s)dx 1 2 ds ≤ d 1 (m 1 , m 2 )∥∇ϕ(0)∥ ∞ + sup 0≤t≤T ∥(b 2 -b 1 )(t)∥ L 2 (m 2 (t)) T 0 ∥∇ϕ(s)∥ ∞ ds, ≤ C( T ∥∇b 1 ∥ ∞ + 1) d 1 (m 1 , m 2 ) √ T + √ T sup 0≤t≤T ∥(b 2 -b 1 )(t)∥ L 2 (m 2 (t)) .Taking the supremum over ∥ψ∥ ∞ ≤ 1 yields[∥m 2 (T ) -m 1 (T )∥ L 1 (Ω) ≤ C( T ∥∇b 1 ∥ ∞ + 1)× d 1 (m 1 , m 2 ) √ T + √ T sup 0≤t≤T ∥(b 2 -b 1 )(t)∥ L 2 (m 2 (t)) .Finally, if instead we used the boundλ(0)ϕ(0, x)dx ≤ ∥λ(0)∥ L 1 (Ω) ∥ϕ∥ ∞ ≤ C∥m 1 -m 2 ∥ L 1 (Ω)we obtain (9) .• Proof of ( 10): Let Ω = RT d and consider a ψ with ∥∇ψ∥ ∞ ≤ 1. Without loss of generality we can assume ψ(0) = R + 1 and thusψ(x) = ψ(x) -ψ(0) + ψ(0) ≥ -|x -0| + R + 1 ≥ -R + R + 1 = 1 and similarly ψ(x) ≤ 2R + 1.Therefore by maximum principle we have that1 ≤ ϕ(t, x) ≤ 2R + 1,and so∥ψ∥ C 1 = ∥ψ∥ ∞ + ∥∇ψ∥ ∞ ≤ C(1 + R).From (34) we haveλ(T, x)ψ(x)dx ≤ d 1 (m 1 , m 2 )∥∇ϕ(0)∥ ∞ + sup 0≤t≤T ∥∇ϕ(t)∥ ∞ ∥b 1 -b 1 ∥ L 2 (m 2 ) .Thus by Corollary D.2, estimate (62) we obtainλ(T, x)ψ(x)dx ≤ C (1 + ∥∇b∥ ∞ )∥v∥ 3 C 1 (d 1 (m 1 , m 2 ) + ∥b 1 -b 1 ∥ L 2 (m 2 ) ) ≤ CR 3 2 (1 + ∥∇b 1 ∥ ∞ ) d 1 (m 1 , m 2 ) + ∥b 1 -b 1 ∥ L 2 (m 2 ) ,which after taking the supremum over ∥∇ψ∥ ∞ ≤ 1 yields (10) .The proof of Theorem 3.2 follows by an application of Theorem 3.1.EQUATIONthenm 1 = m g , m 2 (t) = η π (T -t). Thus, by an application of Theorem 3.1 we obtaind 1 (m g (T ), π) ≤ CR 3 2 (1 + ∥∇s θ ∥ ∞ )(d 1 ( 1 vol(RT d ) , η π ) + ∥s θ -∇ log(η π )∥ L 2 (η π ) ) ≤ CR 3 2 (1 + ∥∇s θ ∥ ∞ ) (Re -ωT R 2 d 1 π, 1 vol(RT d ) + √ e nn ,where in the last inequality we used Proposition D.3. The second claim follows again by Theorem 3.1 and the fact thatd 1 (ν 1 , ν 2 ) = sup ∥∇g∥∞≤1 RT d gd(ν 1 -ν 1 ) ≤ R∥ν 2 -ν 1 ∥ L 1 .The proof of Theorem 3.3 follows the same strategy as in Theorems 3.1 and 3.2. The main technical difference is that we only have a bound between our generated score function s θ and the score function generated by the sample ∇ log(η N,ϵ ) of the formT 0 |s θ -∇ log(η N,ϵ )| 2 dη N,ϵ (s)ds ≤ e nnand so we need to produce a bound between s θ and the true score function ∇ log(η ϵ ). This last step is shown in Section B and in particular Theorem B.5.Proof of Theorem 3.3. Consider the function η π,ϵ : [0, T ] × RT d → R given byEQUATION)Define the driftsb θ * (t, x) := s θ * (T -t, x) and b π ϵ (t, x) := ∇ log(η π ϵ )(T -t, x),and let m ϵ (t, x) = η π,ϵ (T -t, x) which satisfiesEQUATIONFinally, we consider the distribution of our generated sample m g : [0, T ] × Ω → R which is given byEQUATIONWe have the followingEQUATIONWe will bound each of the terms separately.1. We recall that π ϵ = Γ(ϵ) ⋆ π where Γ is the heat kernel. Therefore, by the dual formulation of d 1 we have thatd 1 (π, π ϵ ) = sup ∥∇g∥∞≤1 g(x)d(π ϵ -π) = sup ∥∇g∥∞≤1 (Γ(ϵ) ⋆ g(x) -g(x))dπ(x).Although estimates on (Γ(ϵ) ⋆ g(x) -g(x)) for a function ∥∇g∥ ∞ ≤ 1 are classical for the readers convenience we provide a quick proof. Let v(t) = Γ(t) ⋆ g, which solvesEQUATIONThen for all x ∈ Ω,|v(ϵ, x) -v(0, x)| ≤ ϵ 0 |∂ t v(s, x)|ds = ϵ 0 |∆v(s, x)|ds = ϵ 0 |(∇Γ(s) ⊗ ⋆∇g)(x)|ds ≤ C∥∇g∥ ∞ ϵ 0 1 √ s ds ≤ C √ ϵ.Where in the above we used estimate (31) . Thus, we haved 1 (π, π ϵ ) ≤ C √ ϵ,where C = C(d) > 0.2. First we note that π ϵ = m ϵ (T ). Thus applying Theorem 3.1 form 1 = m g , b 1 = b θ * m 2 = m ϵ , b 2 = b π ϵ yields d 1 (m g (T ), π ϵ ) = d 1 (m g (T ), m ϵ (T )) ≲ R 3 2 (1 + ∥∇b 1 ∥ ∞ )(d 1 (m ϵ (0), 1 vol(RT d ) ) + √ T ∥b π ϵ -b θ * ∥ L 2 (m ϵ ) ).By Proposition D.3, we haved 1 (m ϵ (0), 1 vol(RT d ) ) = d 1 (η π,ϵ (T ), 1 vol(RT d ) ) ≤CRe -ωt R 2 d 1 (π ϵ , 1 vol(RT d ))≤CR 2 e -ωt R 2 .Finally, from Theorem B.5∥b π ϵ -b θ * ∥ 2 L 2 (m ϵ ) = J (η π ϵ , θ * ) ≤ e ′ nn = e nn + C 1 + | log(δ)| √ ϵ + 1 √ T + T ∥s θ ∥ 2 C 2 ([0,T ]×Ω) d 1 (π N , π)Therefore, if T ≥ 1, putting everything together, we have shown that up to a dimensional constantd 1 (π, m g (T )) ≤ C √ ϵ + R 3 2 1 + ∥∇b 1 ∥ ∞ × R 2 e -ωT R 2 + T e nn + C 1 + | log(δ)| √ ϵ + T ∥s θ ∥ 2 C 2 ([0,T ]×Ω) d 1 (π N , π)In this section, we gather all the technical results regarding the connections between ISM, DSM, and ESM score matching functionals. For the readers convenience we recall some of our notations. Let Ω ⊂ R d and m 0 ∈ P(Ω). We will denote byEQUATIONand we may drop the superscipt m 0 when it is clear from context.First we state a simple result that justifies the formulaEQUATIONwhich formally follows byJ (ρ, θ) = T 0 |s θ | 2 dρ(s) -2s θ • ∇ log(ρ) + |∇ log(ρ)| 2 dρ(s)ds = J L (ρ, θ) + T 0 |∇ρ| 2 ρ dxds = J L (ρ) + 4∥∇ √ ρ∥ 2 2 .However ∥∇ √ ρ m0 ∥ 2 may not be finite, thus we record the exact statement along with some technical facts for later use in the following Proposition.Proposition B.1. Let m 0 be a probability density in Ω such that m 0 (x) log(m 0 (x)) ∈ L 1 (Ω) and ρ : [0, T ] × Ω → R be given by (41). Then, the following hold:1. There exists a universal constant C > 0, such thatsup t∈[0,T ] ∥ρ(t) log(ρ(t))∥ L 1 (Ω) + ∥∇ √ ρ∥ 2 L 2 ([0,T ]×Ω) ≤ CT ∥m 0 log(m 0 )∥ L 1 (Ω) ,and thus (42) holds. 
 
 Appendix:  


Decision: This paper gives theoretical analysis of score-based generative models from an uncertainty quantification perspective, showing the impact of errors from different sources. They use a Wasserstein uncertainty propagation theorem using PDE theory, which is potentially of broader interest beyond SGM's. Reviewers received the paper positively and thought it was generally well-written and mathematically sound, with solid technical contributions, though noted a lack of tightness in the bounds (including some unspecified constants) and discretization analysis. I recommend acceptance, with the reviewers' suggestions to improve clarity of presentation.

Reviews: review: paper_summary: paper_summary: This paper studies the generalization error of diffusion models. The major tool is the Wasserstein uncertainty propagation theorem. With such a result and the regularity analysis in PDE, the authors establish robust analysis for diffusion models with respect to various errors.
strength_weakness: strength_weakness: 1. The first major concern is the connection and comparison to the literature. The bounds in Theorems 3.2 and 3.3 are not explicit. I suggest the authors provide clear sample complexity results w.r.t. problem parameters. Furthermore, the authors should discuss how the bounds improved the existing results in the literature. 
2. The analysis in this paper heavily relies on the PDE theory and regularity analysis. I suggest the authors discuss relevant prior work in Section 1.2. Also, I suggest the authors briefly introduce UQ used in other ML problems beyond diffusion models. Furthermore, the score approximation and estimation theory should be mentioned as that is one source of the errors.
questions: questions: I have some minor comments and questions:
1. The results established in this work are usually referred to sample complexity bounds or distribution estimation error. The authors use the name robust analysis while the meaning of "robustness" seems different from the one in robust optimization. I suggest the authors clarify this notion in the context.
2. The paper leverages the pathwise characterization of probability distribution generated by the forward process. I suggest the authors add more explanations and/or refer the readers to the literature when a PDE is introduced, e.g., Eq. (11). The same applies to all the FK and HJB equations.
limitations: limitations: NA
paper_summary: paper_summary: The paper studies the robustness of score-based generative models (SGM) to different sources of error that are relevant in practice, such as the limited expressivity of the score function representation or the choice of reference distribution. Specifically, they upper bound the Wasserstein-1 and L1 distances between target and approximate distribution (given by an SGM) in terms of these different sources or error for both denoising and explicit score-matching objectives.
strength_weakness: strength_weakness: - The main weakness of the paper is its dense presentation. The authors did a good job of motivating and describing their contributions in the introduction, but the rest of the paper is not easy to follow. I am not sure, there is much room for improvement within the limited space of 9 pages, but I would suggest including a discussion on actionable insights one can get from these bounds. Throughout the paper, the authors comment on how important robustness is for the reliability of generative models in general, and I agree, but I could not get an intuition about when one can expect these bounds to be tight so as to ensure the resulting model is trustworthy.
- While I respect the authors’ decision of going for an entirely theoretical paper, I cannot help but feel that some small scale experiments illustrating the tightness and usefulness of the bounds would be enlightening.

### Minor Issues
- Line 35: “contributes” should probably be plural here.
- Line 53: “recognizes” is misspelled.
- Line 61: no need for “of” after “study”.
- Line 100: if I’m not mistaken, the acronym SDE was not yet defined at this point in the text.
- Line 102: Albeit clear from the context, $W$ and $\eta$ have not yet been defined by this point in the text.
- Line 121: “however our results are generally apply”
- Line 331: Word missing after “used”, probably “for”.
questions: questions: - In Section 6.2., the authors discuss an application of their bounds to likelihood-free inference. To that end, could the authors elaborate on how difficult it is to estimate their bounds in practice? Is accurately estimating the Lipschitz constant of the score function the main challenge?
limitations: limitations: This is mainly a theoretical work and I cannot foresee any direct societal impact stemming from this research.
paper_summary: paper_summary: This work studies the influence of different error terms for diffusion models from a continuous perspective under $W_1$ distance. They explain the reason why the early stopping parameter $\epsilon$ would lead to a memory phenomenon of diffusion models. To achieve these results, this work proposes a WUP theorem to explain the robustness of SGMs.
strength_weakness: strength_weakness: 1.	The abstract section mentions that stochasticity is the key mechanism ensuring SGM algorithms. However, this work does not discuss it in detail in the main content. It seems that the WUP theorem does not hold for the deterministic sampling process. It would be better to discuss it in detail.
questions: questions: Please see the Weakness part.

Question 1: This work does not consider the influence of discretization error. It would be better to discuss the challenge when considering this error.

Question 2: For me, one interesting point is the balance of different terms when considering $\epsilon$. As shown in Theorem 3.3, $e_5$ term has $\sqrt{\epsilon}$ dependence and $e_2$ term has $1/\sqrt{\epsilon}$. It seems that there exists a balance between these terms when $N$ is finite. It would be better to discuss this balance in detail.
limitations: limitations: This work does not discuss the limitation and societal impact in a independent paragraph.

Comments: comment: Thanks for the detailed response. It would be better to add the discussion of W1 to the main content in the next version. I am satisfied with the answers and will maintain my positive score.
Thanks for the detailed response. I am satisfied with the answers and will maintain my positive score.
Thank you for your detailed response. I have raised the score.