Improved Analysis of Score-based Generative Modeling: User-Friendly Bounds under Minimal Smoothness Assumptions

Introduction

Generative modeling is one of the central tasks in machine learning, which aims to learn a probability distribution from data and generate data from the learned distribution. Score-based generative modeling (SGM) has achieved state-of-art performance in data generation tasks [SE19, SSK+20, SDME21, DN21], surpassing other models like generative adversarial networks (GAN) [GPAM+14], normalizing flows [RM15], variational autoencoders [KW14], and energy-based models [ZML16]. Due to the impressive sample quality, SGM has great potential in various applications, including computer vision [DN21, RBL+21], natural language processing [AJH+21], inverse problems [SSXE22, CSY21], molecular graph modeling [SLXT21, GRG+22], reinforcement learning [WHZ22], and solving high-dimensional PDEs [BVE22].

The key idea of SGM is to use a forward process to diffuse the data distribution to some prior (often the standard Gaussian), and learn a backward process to transform the prior to the data distribution by estimating the score functions of the forward diffusion process. Such a procedure provides an expressive and efficient way to model high-dimensional distributions for two reasons: 1) It is easy to construct a forward process that converges fast to the Gaussian, no matter how complex the data distribution is. For example, the Ornstein-Uhlenbeck (OU) process has stationary distribution equal to the standard Gaussian and converges rapidly. 2) Several scalable score matching methods such as denoising score matching [Vin11] and sliced score matching [SGSE19] allow us to learn the score function for use by the backward process.

While SGM has achieved great success in practice, theoretical understanding of the power of SGM is far from complete. Recent works [LLT22b, CCL+22] established that when an accurate score estimator is given, SGM can sample from general distributions with polynomial complexity and without requiring structural assumptions such as log-concavity or functional inequalities. (By polynomial complexity we mean that the running time is polynomial and the final error depends polynomially on the score estimation error and other parameters.) This is surprising in the sampling context, as it implies a sharp contrast between SGM and sampling dynamics with gradient flow structure (such as Langevin dynamics), where convergence rates depend crucially on the structure of the data distribution. In this paper, we further establish the effectiveness of SGM by showing that convergence with reasonable rates requires very weak smoothness conditions. Indeed, we obtain a logarithmic dependence on the smoothness, or no dependence when comparing against a slightly perturbed data distribution.

We use $p_{t}(x)$ to denote the density of $x_{t}$. In particular, $p_{T}$ is close to the white noise distribution $\mathcal{N}(0,I_{d})$. Then $x_{t}$ also satisfies the reverse SDE

However, we cannot directly simulate (3) since the score function $\nabla\log p_{t}$ is not available. Thus we learn the score function $\nabla\log p_{t}$ from the noisy data. First, we parameterize the score function within a function class such as that of neural networks, $s_{\theta}(x,t)$. Then we optimize one of the score-matching objectives (denosing score matching [Vin11] is often used; see appendix A for details), from which we obtain a score estimator $s_{\theta}$ such that the $L^{2}$ score estimation error

is small. Using the estimated score, we can generate samples from an approximation of the reverse SDE starting from the prior distribution:

The Choice of Forward Process.

We focus on the case $f(x,t)=-\frac{1}{2}x,\,g(t)\equiv 1$. The choice of $f(x,t)$ matches the choice in the original paper [SSK+20], though our analysis may be adapted for some other choices of drift terms; the choice of constant variance function does not cause any loss of generality since the changing the variance function is equivalent to rescaling time (when $f$ does not depend on $t$). In this case, the forward process becomes the Ornstein-Uhlenbeck process, which has an explicit conditional density:

Moreover, the Ornstein-Uhlenbeck process converges exponentially to the standard Gaussian distribution:

Time Discretization.

In practice, we need to use a discrete-time approximation for the sampling dynamics (4). Let $\delta=t_{0}\leq t_{1}\leq\cdots\leq t_{N}=T$ be the discretization points, where $\delta=0$ for the normal setting and $\delta>0$ for the early-stopping setting. For the $k$-th discretization step ($1\leq k\leq N$), we denote $h_{k}:=t_{k}-t_{k-1}$ as the step size. We will compare different choices of discretization points and identify the optimal choice in different settings.

Let $t_{k}^{\prime}=T-t_{N-k}$ be the corresponding discretization points in the reverse SDE. We consider two types of discretization schemes, which are widely used in existing work.

The expontential integrator scheme [SME21, ZC22]: by using the semi-linear structure of (2), we discretize only in the nonlinear term and retain the continuous dynamics arising from the linear term:

for $k=0,\ldots,N-1$, which is solved explicitly by

where $\eta_{k}\sim\mathcal{N}(0,I_{d})$.

2 Related Work

We highlight two recent papers [CCL+22, LLT22b]. Both papers provide convergence guarantees with polynomial complexity without relying on any structural assumptions on the data distribution such as log-concavity or a functional inequality. In particular, the analysis of [CCL+22] is based on the Girsanov change of measure framework and the authors consider the following two settings: 1) The score functions in the whole trajectory of the forward process satisfy the Lipschitz condition with a uniform Lipschitz constant. 2) The data distribution has bounded support. Although the smoothness condition on the forward process seems mild, it may be hard to check whether the uniform bound for the Lipschitz constants scales polynomially w.r.t. the dimension $d$. In fact, this is a property of the whole process, related to tail bounds of the data distribution. The work [LLT22b] alternatively uses the idea of excluding bad sets in order to reduce to the setting of an $L^{\infty}$-accurate score estimator. This results in a worse dependence on the problem parameters; however, they do relax the smoothness condition on the whole trajectory to one on only the data distribution, and the bounded support assumption to sufficient tail decay.

Many other works have provided convergence analyses, but do not achieve polynomial complexity except in restricted settings, for example relying on functional inequalities (thus precluding multi-modal distributions) [BMR20, LLT22a, WY22], manifold hypotheses [DeB22], or $L^{\infty}$-accurate score estimates [DTHD21]. In the setting where only an $L^{2}$-accurate score estimate of the data distribution is given, [KHR22] give a statistical lower bound which shows it is in general impossible to accurately sample the distribution. This highlights the fact that having score estimates for multiple distributions—e.g., the data distribution with different amounts of noise added—is necessary for efficient sampling; this is done in practice and in our analysis. In a different direction, SGM is also related to recent work on algorithmic stochastic localization [AMS22], in which for the spin glass models under consideration, the score function (i.e., the posterior mean) can be accurately estimated using approximate message passing.

3 Our Contributions

In this paper, we quantitatively show that an $L^{2}$-accurate score estimator is enough to guarantee that the sampling dynamics (5), (6) result in a distribution close to the data distribution in various regimes. Our results combine the advantages of [CCL+22, LLT22b]: under weak assumptions on the data distribution and the score estimator, we provide a concise analysis and refined guarantees for the convergence of SGM under several settings, described below and summarized in Table 1.

Revisiting the setting where the Lipshitz constant of $\nabla\log p_{t},\,0\leq t\leq T$ is uniformly bounded (the trajectory-smooth setting), we provide three refinements compared to [CCL+22]: 1) We sidestep the technical issue of checking Novikov’s condition and provide a reverse KL divergence guarantee, which is stronger than a TV guarantee. 2) For the exponential integrator scheme, the number of steps dependends logarithmically rather than polynomially on the second moment. 3) We do not assume the data distribution has finite KL divergence wrt the standard Gaussian.

Non-smooth setting.

By adding an extra truncation step on the algorithm, we also obtain a pure Wasserstein bound depending on the tail decay of the data distribution, significantly improving the prior result [LLT22a, Theorem 2.2].

Finally, we consider the intermediate assumption of smoothness of $\nabla\log p_{0}$, rather than the whole forward process as in [CCL+22]. In this case, we can bound discretization error in the low-noise regime so that early stopping is not required. We combine the smooth and non-smooth analyses to bound the number of steps logarithmically in $L$, the Lipschitz constant of $\nabla\log p_{0}$.

Furthermore, we analyze difference choices of discretization schemes and step-size schedules (equivalently, different variance functions). This may help guide the practical implementation of SGM.

4 Notations

For random vectors, we denote $\|\cdot\|_{\psi_{k}}:=\|\|\cdot\|\|_{\psi_{k}}$. We use $x\asymp y$ if there exist absolute constants $C_{1},C_{2}>0$ such that $C_{1}y\leq x\leq C_{2}y$. Write $x\lesssim y$ to mean $x\leq Cy$ for an absolute constant $C>0$, and define $x\gtrsim y$ analogously.

Notations for the Forward Process.

Let $P$ be the data distribution and $p_{0}$ be its density (if it exists). For $0<t\leq T$, let $p_{t}$ be the density of $x_{t}$ defined in the forward process (1) with $f(t,x)=\frac{1}{2}g(t)^{2}x_{t}$. Define $\sigma_{t}$ as the conditional variance of $x_{t}$ given $x_{0}$, i.e.,

Notations for Reverse Processes.

Let $s(x,t)$ be the estimated score function. The reverse processes arising in our setting are defined as follows:

Let $\hat{y}_{t}$ be the discrete approximation of $y_{t}$ defined in (5) or (6) starting from $\hat{y}_{0}\sim\mathcal{N}(0,I_{d})$. We use $\hat{q}_{t}$ to denote the density of $\hat{y}_{t}$.

Main Results

We first consider the trajectory smoothness assumption, where we strengthen the result of [CCL+22]. Then, we state our results for more general settings in various regimes.

All the results rely on $L^{2}$-accuracy of the score estimator:

The learned score function $s(x,t)$ satisfies for any $1\leq k\leq N$,

First, we improve result of [CCL+22] for the trajectory-smooth setting, weakening the assumptions and strengthening the conclusion.

Suppose that Assumptions 1,2,3 hold. If $L\geq 1$, $h_{k}\leq 1$ for $k=1,\ldots,N$ and $T\geq 1$, using uniform discretization points yields the followings

Using exponential integrator scheme (6), we have

In particular, choosing $T=\log\left(\frac{M_{2}+d}{\epsilon_{0}^{2}}\right)$ and $N=\Theta\left(\frac{dT^{2}L^{2}}{\epsilon_{0}^{2}}\right)$ makes this $\widetilde{O}\left(\epsilon_{0}^{2}\right)$.

Using the Euler-Maruyama scheme (5), we have

For the exponential integrator, the error consists of three parts: the error of the forward process, the score matching error, and the discretization error, detailed in Section 3.

The extra conditions on $L,h_{k},T$ in the above theorem are introduced to present the result more concisely, and are not a limitation of the analysis.

Comparing to the exponential integrator scheme, the Euler-Maruyama scheme causes an additional high-order discretization error term related to the second-order moment of the data distribution. This implies a separation between the exponential integrator scheme and the Euler-Maruyama scheme: the error of the exponential integrator scheme scales logarithmically in the second moment of the data distribution (as it suffices for $T$ to increase by $O(\log M_{2})$), while the error of the Euler-Maruyama scheme scales linearly.

2 Results for General Distributions with Early Stopping

We now consider the most general setting: we provide convergence guarantees for any distribution that has a bounded second-order moment, without introducing any structural assumptions or smoothness conditions. Hence, our results are applicable to the case that the score function is non-smooth or even not well defined, like distributions supported on a low-dimensional manifold.

Due to our weak assumptions, the backward process (2) may have very bad properties when $t$ is close to , so we need to employ early stopping. For any small constant $\delta>0$, we show that running the sampling dynamics (6) for time $T-\delta$ will result in a distribution close to $p_{\delta}$ in KL divergence. Note that in general, it is impossible to obtain KL or TV closeness to $P$ as this requires matching exactly the support of $P$.

We provide the convergence bound for general discretization and further quantify the bound for several specific choices.

There is a universal constant $K$ such that the following hold. Suppose that Assumptions 1 and 2 hold and the step sizes satisfy

Define $\Pi:=\sum_{k=1}^{N}\frac{h_{k}^{2}}{\sigma_{t_{k-1}}^{4}}$. For $T\geq 2,\delta\leq\frac{1}{2}$, the exponential integrator scheme (6) with early stopping result in a distribution $\hat{q}_{T-\delta}$ such that

In particular, for exponentially decreasing step size $h_{k}=c\min\{t_{k},1\}$, where $c\leq\frac{1}{Kd}$ (or, equivalently $\frac{{\log\left({\frac{1}{\delta}}\right)+T}}{N}\leq\frac{1}{Kd}$), then (8) holds and

Choosing $T=\log\left(\frac{M_{2}+d}{\epsilon_{0}^{2}}\right),N=\Theta\left(\frac{\left(\log\left(\frac{1}{\delta}\right)+T\right)^{2}d^{2}}{\epsilon_{0}^{2}}\right)$ makes this $\widetilde{O}\left(\epsilon_{0}^{2}\right)$.

In addition, for Euler-Maruyama scheme (5), the same bounds hold with an additional term $M_{2}\sum_{k=1}^{N}h_{k}^{3}$ term in the right hand side of (9).

The technical condition (8) is required for the change-of-measure argument in Lemma 13.

By rescaling time, choosing constant variance function $g\equiv 1$ and exponentially decreasing step size is equivalent to choosing exponential $g$ and constant step size. We state the theorem with constant $g$ for convenience (with an exponential choice of $g$, we would only reach the data distribution $P$ at time $t=-\infty$).

The key difficulty in analyzing general distributions is that the discretization error is hard to control without the Lipschitz condition on $\nabla\log p_{t}$. Our approach is to use a high-probability bound for the Hessian matrix $\nabla^{2}\log p_{t}$ with a change of measure. This approach works well for constant-order $t$, while in the low-noise regime the bound will explode as $t$ tends to 0. We overcome the blow-up of discretization error by early stopping.

When $t$ goes to 0, the regularity of $\nabla\log p_{t}$ becomes worse so slowing down the SDE leads to a smaller discretization error. In the result of Theorem 2, the term $\Pi=\sum_{k=1}^{N}\frac{h_{k}^{2}}{\sigma_{t_{k-1}}^{4}}$ in the upper bound (9) depends on the choice of discretization points. In particular,

If we choose uniform discretization $h_{k}=c$, the dependence on $\frac{1}{\delta}$ becomes linear.

[SSK+20] considers variance function $g(t)=\sqrt{t}$ with uniform discreitzation. This is equivalent to using constant variance function with quadratic discretization points $t_{k}=(\delta+kh)^{2}$ for appropriate $h$. This choice of discretization points induces a linear step size and our Theorem results in a square-root dependence on $\frac{1}{\delta}$.

In Theorem 2, by using exponentially decaying (and then constant) step size, we reduce this error to a logarithmic dependence. Indeed, the term $\Pi$ achieves its minimum (up to a constant) under our choice of discretization points.

See Appendix B for details. Although under our assumptions, the theory suggests that exponentially decreasing step sizes are optimal, other issues may arise in practice. We leave an experimental comparison of different $g$’s or step sizes to future work.

Wasserstein+KL Guarantee.

Notice that when $\delta$ is small, $p_{\delta}$ is only a small perturbation (in Wasserstein distance) of the data distribution $P$. Then stopping the algorithm at appropriate $\delta$ results in a distribution that is close in KL divergence to a distribution that is close to $P$ in Wasserstein distance, and we obtain the following.

for an appropriate absolute constant $K$. Here, $M(x)=\exp(\frac{\delta}{2})x,\,Q=M_{\sharp}\hat{q}_{T-\delta}$.

Corollary 3 implies an upper bound for the bounded Lipschitz metric between the data distribution and $\hat{p}_{t_{0}}$(as mentioned in [CCL+22]):

Note our improved dependencies compared with [CCL+22, Corollary 3] and [LLT22b, Theorem 2.1].

While the smoothness assumption is relaxed, our analysis induces an additional $d$-factor in place of the Lipschitz constant of $\nabla\log p_{t}$ compared to Theorem 1. This $d$-factor comes from the high-probability bound for the Hessian matrix (see Lemma 12). However, [CCL+22, Theorem 5] suggests that the lower bound of the discretization error scales linearly on $d$. We leave open the problem of closing the gap between the dimension dependence in the upper and lower bounds.

Pure Wasserstein Guarantee.

We can also obtain a pure Wasserstein guarantee by following [LLT22a, Theorem 2.2]. For this, we need to include an extra truncation step on the algorithm output, i.e., for some choice of $R$, replacing any sample $\hat{y}_{T-\delta}\sim\hat{q}_{T-\delta}$ falling outside $B_{R}(0)$ by 0. In addition, we need to assume some concentration for $P$, so that samples from $P$ lie in $B_{R}(0)$ with high probability.

Note that the appropriate $R$ in (10) exists under mild tail conditions on the data distribution $P$. For example:

3 Result for Smooth Data Distributions

We further provide convergence analysis for smooth $p_{0}$ without using early stopping. As mentioned in Subsection 2.2, the early stopping technique is employed to bound the discretization error in the low-noise regime. We can alternatively bound this error by using the smoothness condition on $p_{0}$:

We bound the discretization error in two different time regimes: Choosing an appropriate constant $\delta_{0}>0$, when $t>\delta_{0}$, we use a high-probability Hessian bound and a change of measure argument similar to the analysis in the early stopping setting; for $t<\delta_{0}$, we alternatively derive a Lipschitz constant bound for $\nabla\log p_{t}$ (stated in Lemma 14) based on Assumption 4.

There is a universal constant $K$ such that the following holds. Under Assumptions 1, 2, and 4 hold, by using the exponentially decreasing (then constant) step size $h_{k}=c\min\{\max\{t_{k},\frac{1}{L}\},1\}$, $c=\frac{\log L+T}{N}\leq\frac{1}{Kd}$, the sampling dynamic (6) results in a distribution $\hat{q}_{T}$ such that

Choosing $T=\log\left(\frac{M_{2}+d}{\epsilon_{0}^{2}}\right)$ and $N=\Theta\left(\frac{d^{2}(T+\log L)^{2}}{\epsilon_{0}^{2}}\right)$ makes this $\widetilde{O}\left(\epsilon_{0}^{2}\right)$.

In addition, for Euler-Maruyama scheme (5), the same bounds hold with an additional $M_{2}\sum_{k=1}^{N}h_{k}^{3}$ term.

Comparing to Theorem 1, this result only depends on the Lipschitz constant of $\nabla\log p_{0}$ rather than the uniform Lipschitz constant bound for $\nabla\log p_{t},\,0\leq t\leq T$. We also ease the dependency on $L$ from $L^{2}$ to $\log^{2}L$ for optimal choice of variance function or step size, so the requirement on the smoothness of the data distribution is significantly relaxed: even if the Lipschitz constant $L$ scales exponentially on $d$, we can still obtain a polynomial complexity guarantee. Note that we do pay an extra $d$ factor compared to Theorem 1.

Proof sketches

We sketch the proofs of the main theorems using the exponential integrator discretization, and give complete proofs in Appendices C and D. We first consider the smooth setting, and then describe the modifications for the non-smooth case. Our main technical novelty lies in the arguments for the non-smooth setting, we also streamline the arguments in the smooth setting and use an interpolation rather than Girsanov approach that gives KL divergence bounds.

The first source of error arises from the mismatch between the distribution of the forward process $p_{T}$ at time $T$, and our Gaussian initialization for the reverse process, $\hat{q}_{0}=\gamma_{d}$. We can separate out this term using the chain rule for KL divergence:

The first term can be bounded using exponential mixing of the forward (Ornstein-Uhlenbeck) process towards the standard Gaussian. In conjunction with the fact that after constant time, the KL-divergence is bounded by $O(d+M_{2})$, we obtain (Lemma 9)

The remaining term can be written as a sum, again using the chain rule for KL divergence, by comparing the continuous process with the estimated, discrete process through a chain of intermediate processes where we run the continuous process until time $t_{k}$. We can interpolate the discrete processes to realize them as SDE’s. If Novikov’s conditions are satisfied, Girsanov’s Theorem then applies to bound the KL divergence in terms of the squared difference of the drift terms between the processes.

In the last step we use the triangle inequality. However, in general Novikov’s condition may not be satisfied; [CCL+22] circumvent this using an involved truncation argument which only results in a TV bound and relies on the trajectory-smooth condition (Assumption 3). We instead use a differential inequality argument which gives the same conclusion (Lemma 6, 7, Proposition 8) and is applicable to the non-smooth setting; this step requires significant technical work (Appendix F).

Second term.

Term (2) is exactly the score estimation error, and by Assumption 1, it is bounded by $T\epsilon_{0}^{2}$.

Third term.

Term (3) is the discretization error. This discretization error bound is non-trivial since in classical numerical analysis theory, the discretization error often depends exponentially on the time $T$ due to the use of Gronwall’s inequality. Our analysis our will rely on the special structure of the Ornstein-Uhlenbeck process. We note that (3) involves both a “time” and “space” discretization error (as both the time and space arguments are different). We show in Lemma 11 that this can be bounded purely in terms of the space discretization error (which streamlines the argument of [CCL+22])

The explicit form of the OU process tells us that $\alpha_{s,t}^{-1}x_{s}=x_{t}+z$, where $z$ is a Gaussian of variance $O(s-t)$. Therefore, the second term (which dominates) can be bounded as a Lipschitz constant times the second moment of a Gaussian:

Note that we crucially use the Lipschitzness of the score in this step. Plugging this bound into the sum (3) gives the final error term.

2 Non-smooth setting (Theorem 2)

Comparing Theorem 1 (smooth setting) and Theorem 2 (non-smooth setting), we note that the discretization error changes from $\frac{T^{2}L^{2}d}{N}$ to $\frac{\left({\log\left({\frac{1}{\delta}}\right)+T}\right)^{2}d}{N}$; the intuition is that $L$ is “effectively” bounded by $\sqrt{d}$. Previously, [CCL+22] assume that $P$ is supported on a ball of radius $R$ to derive a global Lipschitzness bound $\left\|{\nabla^{2}\log p_{t}}\right\|=O\left({\frac{R^{2}}{t^{2}}}\right)$ to plug into the smooth theorem.

Our main insight is that (1) because we are averaging the error over $p_{t}$, it suffices to have a high-probability rather than uniform bound on the the Hessian, and (2) such bounds are obtainable from the smoothing properties of the forward process. In fact, to bound (11), we only need Lipschitzness in a random direction, and hence a Frobenius norm bound is sufficient (Lemma 12):

(This is the weaker analogue of an operator norm bound of $O(\sqrt{d})$, which was suggested from the $L=O(\sqrt{d})$ analogy.) This incurs significant savings over a uniform bound, and in particular does not depend on boundedness or tails of $P$. We prove this by giving a Bayesian interpretation of the Hessian as the posterior variance of the noise in the score matching objective. As a purely mathematical statement about smoothing of the OU process, this result may be of independent interest.

Finally, to use (12) in (11), we actually need to bound the Hessian not just at $x_{t}$ but along the path (in direction $z$) joining $x_{t}$ and $\alpha_{s,t}^{-1}x_{s}$: for this we need a change-of-measure argument (Lemma 13) which says that the distributions of $(x_{t},z)$ and $(x_{t}+az,z)$ are close in $\chi^{2}$-divergence, for $0\leq a\leq 1$. Finally, although the bound (12) blows up as $t\to 0$, by choosing an exponentially decreasing step size and stopping at time $\delta$, we only incur a $\log\left({\frac{1}{\delta}}\right)$ dependence, similarly to the analysis of the score estimation error (Remark 1).

If we only assume $\nabla\log p_{0}$ is $L$-Lipschitz, we can still derive Lipschitzness of $\nabla\log p_{0}$ for small time $t\leq\frac{1}{L}$ (Lemma 14). For large $t\geq L$, the argument in the non-smooth case applies (and gives a bound of $O(dL)$ in (12)). Thus, we take exponentially decreasing step size until $t=1/L$, and then constant step size, and combine the analyses of Theorems 1 and 2 to obtain Theorem 5.

Conclusion

In this paper, we analyzed the theoretical properties of SGM in various regimes. We extended existing result to the most general setting and provided refined guarantees. The current analysis provides guarantees for SGM in the framework that an $L^{2}$-accurate score estimator is available. This implies the training objective in denoising score matching is suitable for learning a generative model and partially explains why SGM is empirically successful at modeling very complex distributions, like multi-mode distributions or distributions with weak smoothness condition.

We obtain guarantees for arbitrary data distributions without smoothness assumptions, by exploiting (high-probability) smoothing properties of the forward process. Besides closing the factor-$d$ gap between our upper bound and the (suggested) lower bound, it would be interesting to carry out this kind of analysis for other choices of the forward/backward processes, such as critically damped Langevin Diffusion [DVK21], to see if improved guarantees are available. ([CCL+22] show that no improvement is available only in the setting of a uniform bound on the Lipschitz constant of the score.)

Another future direction is to explore theories beyond the framework that an $L^{2}$-accurate score estimator is available and understand the learning of a score estimator, including the approximability, sample complexity, and the training dynamics of denoising score matching. This is related to the most challenging problems in deep learning theory; advances in deep learning theory may provide some new insight into SGM.

References

Appendix A Denoising Score Matching

For $0\leq t\leq T$, the goal of score matching for $p_{t}$ is to minimize

Since the score function $\nabla\log p_{t}$ is not available, we alternatively consider a denoising score matching objective [Vin11], which is derived from integrating by parts

where $p_{t|0}$ is the conditional distribution of $x_{t}$ given $x_{0}$, and $C$ is a constant independent of $\theta$.

In this case, by noting that $\sigma_{t_{k}}^{2}\asymp\min\{1,t_{k}\}$, we have

Appendix B Discussion on Choices of Discretization Points

In this section, we consider the scaling of the term $\Pi=\sum_{k=1}^{N}\frac{h_{k}^{2}}{\sigma_{t_{k-1}}^{4}}$ in (9) under different choices of discreitzation points.

For uniform discretization(inducing constant step size) $t_{k}=\delta+kh,\,h=\frac{T-\delta}{N}$, we have

Thus the upper bound for discretization error has a linear dependence on $\frac{1}{\delta}$.

The Linear Step Size

For quadratic discretization points(inducing linear step size) $t_{k}=(\delta+kh)^{2},\,h=\frac{\sqrt{T}-\delta}{N}$, by noting that $\frac{h_{k}}{h}\asymp\sqrt{t_{k}}$, we have

Optimality of Exponential Decaying Step Size

Now we will show that the discretization points used in Theorem 2 minimizes the term $\Pi$ (up to a constant). Indeed, note that

For the term $\Pi_{1}$, let $z_{k}=\log\frac{t_{k}}{t_{k-1}}>0$, we have $\Pi=\sum_{k=1}^{n}(e^{z_{k}}-1)^{2}$. Note that $z\mapsto(e^{z}-1)^{2}$ is convex for $z>0$. By Jensen’s inequality, when the summation of $z_{k}$’s are fixed, the minimum of $\Pi_{1}$ is reached when $z_{k}$’s are identical. Equivalently, $h_{k}=ct_{k}$ for $t_{k}\leq 1$. For the term $\Pi_{2}$, we have $\Pi_{2}=\sum_{t_{k}>1}h_{k}^{2}$. Similarly, since $h\mapsto h^{2}$ is convex for $h>0$, the minimum of $\Pi_{2}$ is reached when $h_{k}$’s are identical.

Appendix C Main Proof Ingredients

The key idea of the proof is motivated by the Girsanov change of measure framework used in [CCL+22]. However, in order to avoid the technical challenge of altering the process to satisfy Novikov’s condition, we use a differential inequality-based argument instead.

where $F_{1},F_{2},g$ are continuous functions and may depend on $a$. We assume the uniqueness and regularity condition:

Define the relative Fisher information between $p_{t}$ and $q_{t}$ by

While we have written the same Brownian motion for $X$ and $Y$, as we only care about distributions, the Brownian motions can be chosen independent with each other.

In addition, the above results also hold if we replace $\hat{y}_{t}$ with that corresponding to the Euler-Maruyama scheme:

The exponential integrator scheme (6) satisfies

For $0\leq k\leq N-1$, we use the chain rule of KL divergence to obtain

This completes the proof for the exponential integrator scheme. The proof for the Euler-Maruyama scheme is similar; the only difference is the differential inequality becomes

and we can obtain the result in an analogous way. ∎

The three terms in the upper bound of Proposition 8 match the claim in Theorem 1. The first term is controlled by the exponential convergence of the forward process, which is given in the following lemma.

Notice that $x\mapsto x\log x$ is a convex function for $x>0$. Let $p_{t|0}$ be the conditional density of $x_{t}$ given $x_{0}$. For any $t>0$, we can use Jensen’s inequality to bound the entropy of $p_{t}$:

Since $x_{t}|x_{0}=y\sim\mathcal{N}(\alpha_{t}x_{0},\sigma_{t}^{2}I_{d})$, we have

From the exponential convergence of Langevin dynamics with strongly log-concave stationary distribution (see, e.g., [VW19]), we obtain

The second term in the upper bound of Proposition 8 is exactly the same as the score estimation error defined in Assumption 1. So the key challenge is to bound the third term, which is caused by the discretization error.

Suppose that $h_{k}\leq 1$ for $1\leq k\leq N$. We have

From the definition of the forward process (1), we have

where the last inequality follows from the Cauchy-Schwartz inequality. From the explicit form of the conditional density

Taking summation over $k=1,\ldots,N$, we complete the proof. ∎

For any $0\leq t\leq s\leq T$, the forward process (1) satisfies

Since $x_{s}|x_{t}\sim\mathcal{N}\left(\alpha_{t,s}{x}_{t},(1-\alpha_{t,s}^{2})I_{d}\right)$, from Lemma 20, we can rewrite $\nabla\log{p}_{s}$ as

where $p_{t|s}$ is the conditional density of $x_{t}$ given $x_{s}$. Thus the time discretization error can be bounded by

Therefore, splitting the error into the space-discretization and the time-discretization error,

If the score functions of the forward process is smooth, i.e., Assumption 3 holds, the space-discretization error can be directly bounded using the Lipschitz condition on $\nabla\log p_{t}$.

In the general setting, we choose a early stopping time $t_{0}$ and bound the space-discretization error for $t>t_{0}$ by a high-probability bound on the Hessian matrix $\nabla^{2}\log p_{t}$ and a change of measure argument, which are worked out in section C.1.

For smooth $p_{0}$, we further bound the space-discretization error for small $t$ by providing a Lipschitz constant bound for $\nabla\log p_{t}$ when $t$ is sufficient small, which is given in section C.2.

In this subsection, we establish the high-probability bound for the Hessian matrix $\nabla^{2}\log p_{t}$ and use the high-probability bound to control the space-discretization error. This is the critical part of our analysis that allows us to prove Theorem 2.

where $\|\cdot\|_{F,\psi_{1}}=\|\|\cdot\|_{F}\|_{\psi_{1}}$ denote the sub-exponential norm of the Frobenius norm of a random matrix.

For any positive integer $p$, using the fact that $\frac{y-x}{\sigma}$ is distributed as $\mathcal{N}(0,I_{d})$ and the power mean inequality,

Using the arbitrariness of $p$, we know that

There is a universal constant $K>0$ so that the following holds. For $0\leq t\leq s\leq T,\,\frac{s-t}{\sigma_{t}^{2}}\leq\frac{1}{Kd}$, we have

We bound the difference between the value of $\nabla\log{p}_{t}$ at different points with the Hessian:

where the last inequality comes from Lemma 12. Next, we bound the second term in (18). By the data processing inequality,

Notice that $x_{t}+az_{t,s}|(z_{t,s},x_{0})\sim\mathcal{N}(\alpha_{t}^{-1}x_{0}+az_{t,s},\sigma_{t}^{2}I_{d})$ and $x_{t}|x_{0}\sim\mathcal{N}(\alpha_{t}^{-1}x_{0},\sigma_{t}^{2}I_{d})$. We can compute the chi-squared divergence explicitly:

Finally, the condition $\frac{s-t}{\sigma_{t}^{2}}\leq\frac{1}{Kd}$ implies $e^{s-t}-1\lesssim s-t$ and $\frac{e^{s-t}-1}{\sigma_{t}^{2}}\lesssim\frac{1}{Kd}$. Thus for large enough $K$(actually, $K=1$ is enough),

Combining the bound for the first and the second terms of (18), we conclude that

Plugging (19) into (17), we complete the proof. ∎

C.2 Stability of the Lipschitz Constant

In this subsection, we show that if $p_{0}$ satisfies the smoothness condition, $p_{t}$ is also smooth for sufficiently small $t$. In particular, under Assumption 4, we can choose $t_{0}\asymp\frac{1}{L}$ and an absolute constant $C$ such that for any $0\leq t\leq t_{0}$, the Lipschitz constant of $\nabla\log p_{t}$ is bounded by $CL$.

Define a density $q(x)\propto p_{0}(\alpha_{t}^{-1}x)$. Then $\nabla\log q$ is $\alpha_{t}^{-1}L$-Lipschitz. Notice that $p_{t}$ is the Gaussian perturbation of $q$. Using Lemma 22, we write the second-order score function of $p_{t}$ as

Appendix D Proofs for the Main Theorems

Now we follow the discussion in Section C and combine everything together to complete the proof of our main theorems stated in Section 2.

For $t_{k-1}\leq t\leq t_{k}$, suppose that $\nabla\log p_{t}$ is $L$-Lipschitz for $t_{k-1}\leq t\leq t_{k}$. If $L\geq 1,h_{k}\leq 1$, we have

The space-discretization error is easily bounded by the Lipschitz condition:

where the last inequality is because of $t_{k}-t\lesssim 1$. Combining Lemma 11, Lemma 21, and (20), we have

As shown in Section C, the extra terms arising in the discretization error of Euler-Maruyama scheme can be bounded by Lemma 10, so we only need to consider the exponential integrator scheme. By Proposition 8, we can bound the KL divergence between $p_{0}$ and $\hat{q}_{T}$ by

The first term in (21) is bounded by Lemma 9. Then, we apply Lemma 16 to bound the discretization error:

For uniform discretization, the above quantity is $\frac{dT^{2}L^{2}}{N}$. We complete the proof. ∎

D.2 Proof of Theorem 2

There is a constant $K$ such that the following holds. In the early stopping setting, suppose that the variance function $g$ satisfies $\frac{h_{k}}{\sigma_{t_{k-1}}^{2}}\leq\frac{1}{Kd}$ for any integer $1\leq k\leq N$. Then we have

Noticing that $\frac{h_{k}}{\sigma_{t_{k-1}}^{2}}\lesssim\frac{1}{d}$ implies $\frac{(1-\alpha_{t,t_{k}}^{-1})^{2}}{\sigma_{t}^{2}}\lesssim\frac{t_{k}-t}{d}$ and combining this with (LABEL:555) and (23), we conclude that

If $K\geq 2$, $c\leq\frac{1}{Kd}$, $t_{0}=\delta$, $t_{N}=T$, and $h_{k}:=t_{k}-t_{k-1}=c\min\{t_{k},1\}$, then $\frac{h_{k}}{\sigma_{t_{k}}^{2}}\lesssim\frac{1}{Kd}$ for $k=1,\ldots,N$ and

Note that $\sigma_{t}^{2}\asymp\max\{1,t\}$. We consider the sum with $t_{k}\leq 1$ and $t_{k}>1$ separately. For $t_{k}\leq 1$, we have $\frac{G_{k}}{\sigma_{k-1}^{2}}\asymp\frac{ct_{k}}{t_{k-1}}\leq\frac{2}{Kd}$ (when $K\geq 2$, so $\frac{t_{k}}{t_{k-1}}\leq 2$). Noting that the number of terms in the sum is $\lesssim\log_{1-c}(\delta)$,

For $t_{k}>1$, $\frac{h_{k}}{\min\{t_{k-1},1\}}=c\leq\frac{1}{Kd}$ and

Combining (24) and (25) gives the result. Note the number of steps is

As shown in Section C, the extra terms arising in the discretization error of the Euler-Maruyama scheme can be bounded by Lemma 10, so we only need to consider the exponential integrator scheme. From Proposition 8 we obtain

By bounding the first term in (26) with Lemma 9 and the second term in (26) with Lemma 17, we obtain (9). Further more, we can further quantify the term $\Pi=\sum_{k=1}^{N}\frac{G_{k}^{2}}{\sigma_{t_{k-1}}^{4}}$ for exponentially decaying (and then constant) step size with Lemma 18. ∎

D.3 Proof of Corollary 3 and Corollary 4

[LLT22b, Lemma 6.6] Let $\mu$ be the standard Gaussian measure on $N\left(0,I_{d}\right)$. Then

By the $\chi^{2}$ tail bound in [LM00], for $t\geq 0$

so $\|X\|^{2}$ is stochastically dominated by a random variable with cdf $F(y)=1-e^{-\frac{y-2d}{3}}$. Then letting $P_{Y}$ be the measure corresponding to $F$,

Proof of Corollary 4.

We used the triangle inequality, data processing inequality, and Pinsker’s inequality in (28), (29), and (30), respectively. Express ${x}_{\delta}=\alpha_{\delta}{x}_{0}+\sigma_{\delta}\xi$, where ${x}_{0}\sim P,\,\xi\sim\mathcal{N}(0,I_{d})$. Now

where the second inequality comes from Lemma 19. Combining (27), (31), (32) and the choice of parameters in (10), we complete the proof. ∎

D.4 Proof of Theorem 5

As shown in Section C, the extra terms arising in the discretization error of the Euler-Maruyama scheme can be bounded by Lemma 10, so we only need to consider the exponential integrator scheme. Using Proposition 8, we obtain

In the right hand side of (33), the first term is directly bounded by Lemma 9. Thus we only have to consider the second term, which is the discretization error. Let $k_{0}$ be the largest index such that $t_{k_{0}}\leq\frac{1}{L}$. By Lemma 17 and Lemma 18,

The number of steps for this part is $N-k_{0}\lesssim\frac{1}{c}(\log L+T)$. Note $k_{0}\lesssim\frac{1}{c}$ so by Lemma 16 and Lemma 14,

Thus the total discretization error is bounded by

and the total number of steps is $N\lesssim\frac{1}{c}(\log L+T)$. Given the number of steps $N$, we can choose $c=\frac{\log L+T}{N}$; plugging this in gives the bound. We complete the proof.

Appendix E Lemmas for Computing Score Functions

In this section, we provide some lemmas for the score function, which will be used in our analysis.∎

[CEL+22] If $\nabla\log p$ is $L$-Lipchitz, we have

Using Lemma 20, we rewrite the score function as

We rewrite the second-order score function as

To prove the second expression, we note that

Appendix F Technical details for Proposition 8

By the Fokker-Plank equation, the evolution of $p_{t}$ and $q_{t}$ is given by

Proof of Lemma 7(1).

In order to prove Lemma 7(2), we need the following.

By Lemma 20, we write the score function of $p_{\sigma}$ as

Proof of Lemma 7(2).

for the exponential integrator scheme, or

for the Euler-Maruyama scheme. Hence, (36) is obtained by the Monotone Convergence Theorem and we conclude the proof. Now we check the Novikov condition, which is given by

In fact, by Lemma 23 we have $\|\nabla\log{p}_{t}({x}_{t})\|_{\psi_{2}}\lesssim\sqrt{\frac{d}{\sigma_{t}^{2}}}$. Thus

When $t-t_{k}^{\prime}$ is sufficient small, we have

The second term in the right hand side of (39) is a constant so we only need to consider the first term. Note that

Thus when $t-t_{k}^{\prime}$ is sufficient small we have