Maximum Likelihood Training of Score-Based Diffusion Models
Year: 2021

Abstract:
Score-based diffusion models synthesize samples by reversing a stochastic process that diffuses data to noise, and are trained by minimizing a weighted combination of score matching losses. The log-likelihood of score-based diffusion models can be tractably computed through a connection to continuous normalizing flows, but log-likelihood is not directly optimized by the weighted combination of score matching losses. We show that for a specific weighting scheme, the objective upper bounds the negative log-likelihood, thus enabling approximate maximum likelihood training of score-based diffusion models. We empirically observe that maximum likelihood training consistently improves the likelihood of score-based diffusion models across multiple datasets, stochastic processes, and model architectures. Our best models achieve negative log-likelihoods of 2.83 and 3.76 bits/dim on CIFAR-10 and ImageNet 32ˆ 32 without any data augmentation, on a par with state-of-the-art autoregressive models on these tasks.

heading: 1 Introduction
text: Score-based generative models [44, 45, 48] and diffusion probabilistic models [43, 19] have recently achieved state-of-the-art sample quality in a number of tasks, including image generation [48, 11], audio synthesis [5, 27, 37], and shape generation [3]. Both families of models perturb data with a sequence of noise distributions, and generate samples by learning to reverse this path from noise to data. Through stochastic calculus, these approaches can be unified into a single framework [48] which we refer to as score-based diffusion models in this paper.
The framework of score-based diffusion models [48] involves gradually diffusing the data distribution towards a given noise distribution using a stochastic differential equation (SDE), and learning the time reversal of this SDE for sample generation. Crucially, the reverse-time SDE has a closed-form expression which depends solely on the time-dependent gradient field (a.k.a., score) of the perturbed data distribution. This gradient field can be efficiently estimated by training a neural network (called a score-based model [44, 45]) with a weighted combination of score matching losses [23, 56, 46] as the objective. A key advantage of score-based diffusion models is that they can be transformed into continuous normalizing flows (CNFs) [6, 15], thus allowing tractable likelihood computation with numerical ODE solvers. ˚Equal contribution.
35th Conference on Neural Information Processing Systems (NeurIPS 2021).
Compared to vanilla CNFs, score-based diffusion models are much more efficient to train. This is because the maximum likelihood objective for training CNFs requires running an expensive ODE solver for every optimization step, while the weighted combination of score matching losses for training score-based models does not. However, unlike maximum likelihood training, minimizing a combination of score matching losses does not necessarily lead to better likelihood values. Since better likelihoods are useful for applications including compression [21, 20, 51], semi-supervised learning [10], adversarial purification [47], and comparing against likelihood-based generative models, we seek a training objective for score-based diffusion models that is as efficient as score matching but also promotes higher likelihoods.
We show that such an objective can be readily obtained through slight modification of the weighted combination of score matching losses. Our theory reveals that with a specific choice of weighting, which we term the likelihood weighting, the combination of score matching losses actually upper bounds the negative log-likelihood. We further prove that this upper bound becomes tight when our score-based model corresponds to the true time-dependent gradient field of a certain reverse-time SDE. Using likelihood weighting increases the variance of our objective, which we counteract by introducing a variance reduction technique based on importance sampling. Our bound is analogous to the classic evidence lower bound used for training latent-variable models in the variational autoencoding framework [26, 39], and can be viewed as a continuous-time generalization of [43].
With our likelihood weighting, we can minimize the weighted combination of score matching losses for approximate maximum likelihood training of score-based diffusion models. Compared to weightings in previous work [48], we consistently improve likelihood values across multiple datasets, model architectures, and SDEs, with only slight degradation of Fréchet Inception distances [17]. Moreover, our upper bound on negative log-likelihood allows training with variational dequantization [18], with which we reach negative log-likelihood of 2.83 bits/dim on CIFAR-10 [28] and 3.76 bits/dim on ImageNet 32ˆ32 [55] with no data augmentation. Our models present the first instances of normalizing flows which achieve comparable likelihood to cutting-edge autoregressive models.

heading: 2 Score-based diffusion models
text: Score-based diffusion models are deep generative models that smoothly transform data to noise using a diffusion process, and synthesize samples by learning and simulating the time-reversal of this diffusion. The overall approach is illustrated in Fig. 1.

heading: 2.1 Diffusing data to noise with an SDE
text: Let ppxq denote the unknown distribution of a dataset consisting of D-dimensional i.i.d. samples. Score-based diffusion models [48] employ a stochastic differential equation (SDE) to diffuse ppxq towards a noise distribution. The SDEs are of the form
dx “ fpx, tqdt` gptqdw, (1) where fp¨, tq : RD Ñ RD is the drift coefficient, gptq P R is the diffusion coefficient, and w P RD denotes a standard Wiener process (a.k.a., Brownian motion). Intuitively, we can interpret dw as infinitesimal Gaussian noise. The solution of an SDE is a diffusion process txptqutPr0,T s, where r0, T s is a fixed time horizon. We let ptpxq denote the marginal distribution of xptq, and p0tpx1 | xq denote the transition distribution from xp0q to xptq. Note that by definition we always have p0 “ p when using an SDE to perturb the data distribution.
The role of the SDE is to smooth the data distribution by adding noise, gradually removing structure until little of the original signal remains. In the framework of score-based diffusion models, we choose fpx, tq, gptq, and T such that the diffusion process txptqutPr0,T s approaches some analytically tractable prior distribution πpxq at t “ T , meaning pT pxq « πpxq. Three families of SDEs suitable for this task are outlined in [48], namely Variance Exploding (VE) SDEs, Variance Preserving (VP) SDEs, and subVP SDEs.

heading: 2.2 Generating samples with the reverse SDE
text: Sample generation in score-based diffusion models relies on time-reversal of the diffusion process. For well-behaved drift and diffusion coefficients, the forward diffusion described in Eq. (1) has an
associated reverse-time diffusion process [1, 16] given by the following SDE
dx “ “ fpx, tq ´ gptq2∇x log ptpxq ‰ dt` gptqdw̄, (2) where w̄ is now a standard Wiener process in the reverse-time direction. Here dt represents an infinitesimal negative time step, meaning that the above SDE must be solved from t “ T to t “ 0. This reverse-time SDE results in exactly the same diffusion process txptqutPr0,T s as Eq. (1), assuming it is initialized with xpT q „ pT pxq. This result allows for the construction of diffusion-based generative models, and its functional form reveals the key target for learning: the time-dependent score function ∇x log ptpxq. Again, see Fig. 1 for a helpful visualization of this two-part formulation. In order to estimate ∇x log ptpxq from a given dataset, we fit the parameters of a neural network sθpx, tq, termed a score-based model, such that sθpx, tq « ∇x log ptpxq for almost all x P RD and t P r0, T s. Unlike many likelihood-based generative models, a score-based model does not need to satisfy the integral constraints of a density function, and is therefore much easier to parameterize. Good score-based models should keep the following least squares loss small
JSMpθ;λp¨qq :“ 1
2
ż T
0
Eptpxqrλptq ‖∇x log ptpxq ´ sθpx, tq‖ 2 2sdt, (3)
where λ : r0, T s Ñ Rą0 is a positive weighting function. The integrand features the well-known score matching [23] objective Eptpxqr‖∇x log ptpxq ´ sθpx, tq‖ 2 2s. We therefore refer to Eq. (3) as a weighted combination of score matching losses.
With score matching techniques [56, 46], we can compute Eq. (3) up to an additive constant and minimize it for training score-based models. For example, we can use denoising score matching [56] to transform JSMpθ;λp¨qq into the following, which is equivalent up to a constant independent of θ:
JDSMpθ;λp¨qq :“ 1
2
ż T
0
Eppxqp0tpx1|xqrλptq ∥∥∇x1 log p0tpx1 | xq ´ sθpx1, tq∥∥22sdt. (4)
Whenever the drift coefficient fθpx, tq is linear in x (which is true for all SDEs in [48]), the transition density p0tpx1 | xq is a tractable Gaussian distribution. We can form a Monte Carlo estimate of both the time integral and expectation in JDSMpθ;λp¨qq with a sample pt,x,x1q, where t is uniformly drawn from r0, T s, x „ ppxq is a sample from the dataset, and x1 „ p0tpx1 | xq. The gradient ∇x1 log p0tpx1 | xq can also be computed in closed form since p0tpx1 | xq is Gaussian. After training a score-based model sθpx, tq with JDSMpθ;λp¨qq, we can plug it into the reverse-time SDE in Eq. (2). Samples are then generated by solving this reverse-time SDE with numerical SDE solvers, given an initial sample from πpxq at t “ T . Since the forward SDE Eq. (1) is designed such that pT pxq « πpxq, the reverse-time SDE will closely trace the diffusion process given by Eq. (1) in the reverse time direction, and yield an approximate data sample at t “ 0 (as visualized in Fig. 1).

heading: 3 Likelihood of score-based diffusion models
text: The forward and backward diffusion processes in score-based diffusion models induce two probabilistic models for which we can define a likelihood. The first probabilistic model, denoted as pSDEθ pxq, is given by the approximate reverse-time SDE constructed from our score-based model sθpx, tq. In particular, suppose tx̂θptqutPr0,T s is a stochastic process given by
dx̂ “ “ fpx̂, tq ´ gptq2sθpx̂, tq ‰ dt` gptqdw̄, x̂θpT q „ π. (5)
We define pSDEθ as the marginal distribution of x̂θp0q. The probabilistic model pSDEθ is jointly defined by the score-based model sθpx, tq, the prior π, plus the drift and diffusion coefficients of the forward SDE in Eq. (1). We can obtain a sample x̂θp0q „ pSDEθ by numerically solving the reverse-time SDE in Eq. (5) with an initial noise vector x̂θpT q „ π. The other probabilistic model, denoted pODEθ pxq, is derived from the SDE’s associated probability flow ODE [32, 48]. Every SDE has a corresponding probability flow ODE whose marginal distribution at each time t matches that of the SDE, so that they share the same ptpxq for all time. In particular, the ODE corresponding to the SDE in Eq. (1) is given by
dx dt “ fpx, tq ´ 1 2 gptq2∇x log ptpxq. (6)
Unlike the SDEs in Eq. (1) and Eq. (2), this ODE describes fully deterministic dynamics for the process. Notably, it still features the same time-dependent score function ∇x log ptpxq. By approximating this score function with our model sθpx, tq, the probability flow ODE becomes
dx̃ dt “ fpx̃, tq ´ 1 2 gptq2sθpx̃, tq. (7)
In fact, this ODE is an instance of a continuous normalizing flow (CNF) [15], and we can quantify how the ODE dynamics transform volumes across time in exactly the same way as these traditional flow-based models [6]. Given a prior distribution πpxq, and a trajectory function x̃θ : r0, T s Ñ RD satisfying the ODE in Eq. (7), we define pODEθ as the marginal distribution of x̃θp0q when x̃θpT q „ π. Similarly to pSDEθ , the model p ODE θ is jointly defined by the score-based model sθpx, tq, the prior π, and the forward SDE in Eq. (1). Leveraging the instantaneous change-of-variables formula [6], we can evaluate log pODEθ pxq exactly with numerical ODE solvers. Since pODEθ is a CNF, we can generate a sample x̃θp0q „ pODEθ by numerically solving the ODE in Eq. (7) with an initial value x̃θpT q „ π.
Although computing log pODEθ pxq is tractable, training pODEθ with maximum likelihood will require calling an ODE solver for every optimization step [6, 15], which can be prohibitively expensive for large-scale score-based models. Unlike pODEθ , we cannot evaluate log p SDE θ pxq exactly for an arbitrary data point x. However, we have a lower bound on log pSDEθ pxq which allows both efficient evaluation and optimization, as will be shown in Section 4.2.

heading: 4 Bounding the likelihood of score-based diffusion models
text: Many applications benefit from models which achieve high likelihood. One example is lossless compression, where log-likelihood directly corresponds to the minimum expected number of bits needed to encode a message. Indeed, popular likelihood-based models such as variational autoencoders and normalizing flows have already found success in image compression [51, 20, 21]. Despite some known drawbacks [50], likelihood is still one of the most popular metrics for evaluating and comparing generative models.
Maximizing the likelihood of score-based diffusion models can be accomplished by either maximizing the likelihood of pSDEθ or p ODE θ . Although p ODE θ is a continuous normalizing flow (CNF) and its loglikelihood is tractable, training with maximum likelihood is expensive. As mentioned already, it requires solving an ODE at every optimization step in order to evaluate the log-likelihood on a batch of training data. In contrast, training with the weighted combination of score matching losses is much more efficient, yet in general it does not directly promote high likelihood of either pSDEθ or p ODE θ .
In what follows, we show that with a specific choice of the weighting function λptq, the combination of score matching losses JSMpθ;λp¨qq actually becomes an upper bound on DKLpp } pSDEθ q, and can therefore serve as an efficient proxy for maximum likelihood training. In addition, we provide a related lower bound on log pSDEθ pxq that can be evaluated efficiently on any individual datapoint x.

heading: 4.1 Bounding the KL divergence with likelihood weighting
text: It is well-known that maximizing the log-likelihood of a probabilistic model is equivalent to minimizing the KL divergence from the data distribution to the model distribution. We show in the following theorem that for the model pSDEθ , this KL divergence can be upper bounded by JSMpθ;λp¨qq when using the weighting function λptq “ gptq2, where gptq is the diffusion coefficient of SDE in Eq. (1).
Theorem 1. Let ppxq be the data distribution, πpxq be a known prior distribution, and pSDEθ be defined as in Section 3. Suppose txptqutPr0,T s is a stochastic process defined by the SDE in Eq. (1) with xp0q „ p, where the marginal distribution of xptq is denoted as pt. Under some regularity conditions detailed in Appendix A, we have
DKLpp } pSDEθ q ď JSMpθ; gp¨q2q `DKLppT } πq. (8)
Sketch of proof. Let µ and ν denote the path measures of SDEs in Eq. (1) and Eq. (5) respectively. Intuitively, µ is the joint distribution of the diffusion process txptqutPr0,T s given in Section 2.1, and ν represents the joint distribution of the process tx̂θptqutPr0,T s defined in Section 3. Since we can marginalize µ and ν to obtain distributions p and pSDEθ , the data processing inequality gives DKLpp } pSDEθ q ď DKLpµ } νq. From the chain rule for the KL divergence, we also have DKLpµ } νq “ DKLppT } πq ` EpT pzqrDKLpµp¨ | xpT q “ zq } νp¨ | x̂θpT q “ zqqs, where the KL divergence in the final term can be computed by applying the Girsanov theorem [34] to Eq. (5) and the reverse-time SDE of Eq. (1).
When the prior distribution π is fixed, Theorem 1 guarantees that optimizing the weighted combination of score matching losses JSMpθ; gp¨q2q is equivalent to minimizing an upper bound on the KL divergence from the data distribution p to the model distribution pSDEθ . Due to well-known equivalence between minimizing KL divergence and maximizing likelihood, we have the following corollary. Corollary 1. Consider the same conditions and notations in Theorem 1. When π is a fixed prior distribution that does not depend on θ, we have
´Eppxqrlog pSDEθ pxqs ď JSMpθ; gp¨q2q ` C1 “ JDSMpθ; gp¨q2q ` C2, where C1 and C2 are constants independent of θ.
In light of the result in Corollary 1, we henceforth term λptq “ gptq2 the likelihood weighting. The original weighting functions in [48] are inspired from earlier work such as [44, 45] and [19], which are motivated by balancing different score matching losses in the combination, and justified by empirical performance. In contrast, likelihood weighting is motivated from maximizing the likelihood of a probabilistic model induced by the diffusion process, and derived by theoretical analysis. There are three types of SDEs considered in [48]: the Variance Exploding (VE) SDE, the Variance Preserving (VP) SDE, and the subVP SDE. In Table 1, we summarize all these SDEs and contrast their original weighting functions with our likelihood weighting. For VE SDE, our likelihood weighting incidentally coincides with the original weighting used in [48], whereas for VP and subVP SDEs they differ from one another.
Theorem 1 leaves two questions unanswered. First, what are the conditions for the bound to be tight (become an equality)? Second, is there any connection between pSDEθ and p ODE θ under some conditions? We provide both answers in the following theorem. Theorem 2. Suppose ppxq and qpxq have continuous second-order derivatives and finite second moments. Let txptqutPr0,T s be the diffusion process defined by the SDE in Eq. (1). We use pt and qt to denote the distributions of xptq when xp0q „ p and xp0q „ q, and assume they satisfy the same assumptions in Appendix A. Under the conditions qT “ π and sθpx, tq ” ∇x log qtpxq for all t P r0, T s, we have the following equivalence in distributions
pSDEθ “ pODEθ “ q. (9) Moreover, we have
DKLpp } pSDEθ q “ JSMpθ; gp¨q2q `DKLppT } πq. (10)
Sketch of proof. When sθpx, tq matches ∇x log qtpxq, they both represent the time-dependent score of the same stochastic process so we immediately have pSDEθ “ q. According to the theory of probability flow ODEs, we also have pODEθ “ q “ pSDEθ . To prove Eq. (10), we note that DKLpp } pSDEθ q “ DKLpp}qq “ DKLppT }qT q´ şT 0 d dtDKLppt } qtqdt “ DKLppT }πq´ şT 0 d dtDKLppt } qtqdt. We can now complete the proof by simplifying the integrand using the Fokker–Planck equation of pt and qt followed by integration by parts.
In practice, the conditions of Theorem 2 are hard to satisfy since our score-based model sθpx, tq will not exactly match the score function ∇x log qtpxq of some reverse-time diffusion process with the initial distribution qT “ π. In other words, our score model may not be a valid time-dependent score function of a stochastic process with an appropriate initial distribution. Therefore, although score matching with likelihood weighting performs approximate maximum likelihood training for pSDEθ , we emphasize that it is not theoretically guaranteed to make the likelihood of p ODE θ better. That said, pODEθ will closely match p SDE θ if our score-based model well-approximates the true score such that sθpx, tq « ∇x log ptpxq for all x and t P r0, T s. Moreover, we empirically observe in our experiments (see Table 2) that training with the likelihood weighting is actually able to consistently improve the likelihood of pODEθ across multiple datasets, SDEs, and model architectures.

heading: 4.2 Bounding the log-likelihood on individual datapoints
text: The bound in Theorem 1 is for the entire distributions of p and pSDEθ , but we often seek to bound the log-likelihood for an individual data point x. In addition, JSMpθ;λp¨qq in the bound is not directly computable due to the unknown quantity ∇x log ptpxq, and can only be evaluated up to an additive constant through JDSMpθ;λp¨qq (as we already discussed in Section 2.2). Therefore, the bound in Theorem 1 is only suitable for training purposes. To address these issues, we provide the following bounds for individual data points.
Theorem 3. Let p0tpx1 | xq denote the transition distribution from p0pxq to ptpxq for the SDE in Eq. (1). With the same notations and conditions in Theorem 1, we have
´ log pSDEθ pxq ď LSMθ pxq “ LDSMθ pxq, (11)
where LSMθ pxq is defined as
´Ep0T px1|xqrlog πpx 1qs ` 1
2
ż T
0
Ep0tpx1|xq ” 2gptq2∇x1 ¨ sθpx1, tq ` gptq2 ∥∥sθpx1, tq∥∥22 ´ 2∇x1 ¨ fpx1, tqıdt,
and LDSMθ pxq is given by
´ Ep0T px1|xqrlog πpx 1qs ` 1
2
ż T
0
Ep0tpx1|xq ” gptq2 ∥∥sθpx1, tq ´∇x1 log p0tpx1 | xq∥∥22ıdt
´ 1 2
ż T
0
Ep0tpx1|xq ” gptq2 ∥∥∇x1 log p0tpx1 | xq∥∥22 ` 2∇x1 ¨ fpx1, tqıdt.
Sketch of proof. For any continuous data distribution p, we have ´Eppxqrlog pSDEθ pxqs “ DKLpp } pSDEθ q ` Hppq, where Hppq denotes the differential entropy of p. The KL term can be bounded according to Theorem 1, while the differential entropy has an identity similar to Theorem 2 (see Theorem 4 in Appendix A). Combining the bound of DKLpp } pSDEθ q and the identity of Hppq, we obtain a bound on ´Eppxqrlog pSDEθ pxqs that holds for all continuous distribution p. Removing the expectation over p on both sides then gives us a bound on ´ log pSDEθ pxq for an individual datapoint x. We can simplify this bound to LSMθ pxq and LDSMθ pxq with similar techniques to [23] and [56].
We provide two equivalent bounds LSMθ pxq and LDSMθ pxq. The former bears resemblance to score matching while the second resembles denoising score matching. Both admit efficient unbiased estimators when fp¨, tq is linear, as the time integrals and expectations in LSMθ pxq and LDSMθ pxq can be estimated by samples of the form pt,x1q, where t is uniformly sampled over r0, T s, and x1 „ p0tpx1 | xq. Since the transition distribution p0tpx1 | xq is a tractable Gaussian when fp¨, tq is linear, we can easily sample from it as well as evaluating ∇x1 log p0tpx1 | xq for computing LDSMθ pxq.
Moreover, the divergences ∇x ¨ sθpx, tq and ∇x ¨ fpx, tq in LSMθ pxq and LDSMθ pxq have efficient unbiased estimators via the Skilling–Hutchinson trick [42, 22].
We can view LDSMθ pxq as a continuous-time generalization of the evidence lower bound (ELBO) in diffusion probabilistic models [43, 19]. Our bounds in Theorem 3 are not only useful for optimizing and estimating log pSDEθ pxq, but also for training the drift and diffusion coefficients fpx, tq and gptq jointly with the score-based model sθpx, tq; we leave this avenue of research for future work. In addition, we can plug the bounds in Theorem 3 into any objective that involves minimizing ´ log pSDEθ pxq to obtain an efficient surrogate. Section 5.2 provides an example, where we perform variational dequantization to further improve the likelihood of score-based diffusion models.
Similar to the observation in Section 4.1, LSMθ pxq and LDSMθ pxq are not guaranteed to upper bound ´ log pODEθ pxq. However, they should become approximate upper bounds when sθpx, tq is trained sufficiently close to the ground truth. In fact, we empirically observe that´ log pODEθ pxq ď LSMθ pxq “ LDSMθ pxq holds true for x sampled from the dataset in all experiments.

heading: 4.3 Numerical stability
text: So far we have assumed that the SDEs are defined in the time horizon r0, T s in all theoretical analysis. In practice, however, we often face numerical instabilities when t Ñ 0. To avoid them, we choose a small non-zero starting time ą 0, and train/evaluate score-based diffusion models in the time horizon r , T s instead of r0, T s. Since is small, training score-based diffusion models with likelihood weighting still approximately maximizes their model likelihood. Yet at test time, the likelihood bound as computed in Theorem 3 is slightly biased, rendering the values not directly comparable to results reported in other works. We use Jensen’s inequality to correct for this bias in our experiments, for which we provide a detailed explanation in Appendix B.

heading: 4.4 Related work
text: Our result in Theorem 2 can be viewed as a generalization of De Bruijin’s identity ([49], Eq. 2.12) from its original differential form to an integral form. De Bruijn’s identity relates the rate of change of the Shannon entropy under an additive Gaussian noise channel to the Fisher information, a result which can be interpreted geometrically as relating the rate of change of the volume of a distribution’s typical set to its surface area. Ref. [2] (Lemma 1) builds on this result and presents an integral and relative form of de Bruijn’s identity which relates the KL divergence to the integral of the relative Fisher information for a distribution of interest and a reference standard normal. More generally, various identities and inequalities involving the (relative) Shannon entropy and (relative) Fisher information have found use in proofs of the central limit theorem [24]. Ref. [31] (Theorem 1) covers similar ground to the relative form of de Bruijn’s identity, but is perhaps the first to consider its implications for learning in probabilistic models by framing the discussion in terms of the score matching objective ([23], Eq. 2).

heading: 5 Improving the likelihood of score-based diffusion models
text: Our theoretical analysis implies that training with the likelihood weighting should improve the likelihood of score-based diffusion models. To verify this empirically, we test likelihood weighting with different model architectures, SDEs, and datasets. We observe that switching to likelihood weighting increases the variance of the training objective and propose to counteract it with importance sampling. We additionally combine our bound with variational dequantization [18] which narrows the gap between the likelihood of continuous and discrete probability models. All combined, we observe consistent improvement of likelihoods for both pSDEθ and p ODE θ across all settings. We term the model pODEθ trained in this way ScoreFlow, and show that it achieves excellent likelihoods on CIFAR-10 [28] and ImageNet 32ˆ32 [55], on a par with cutting-edge autoregressive models.

heading: 5.1 Variance reduction via importance sampling
text: As mentioned in Section 2.2, we typically use Monte Carlo sampling to approximate the time integral in JDSMpθ;λp¨qq during training. In particular, we first uniformly sample a time step t „ Ur0, T s, and then use the denoising score matching loss at t as an estimate for the whole time integral. This
Monte Carlo approximation is much faster than computing the time integral accurately, but introduces additional variance to the training loss.
We empirically observe that this Monte Carlo approximation suffers from a larger variance when using our likelihood weighting instead of the original weightings in [48]. Leveraging importance sampling, we propose a new Monte Carlo approximation that significantly reduces the variance of learning curves under likelihood weighting, as demonstrated in Fig. 2. In fact, with importance sampling, the loss variance (after convergence) decreases from 98.48 to 0.068 on CIFAR-10, and decreases from 0.51 to 0.043 on ImageNet.
Let λptq “ αptq2 denote the weightings in [48] (reproduced in Table 1), and recall that our likelihood weighting is λptq “ gptq2. Since αptq2 empirically leads to lower variance, we can use a proposal distribution pptq :“ gptq2{αptq2Z to change the weighting in JDSMpθ; gp¨q2q from gptq2 to αptq2 with importance sampling, where Z is a normalizing constant that ensures ş
pptqdt “ 1. Specifically, for any function hptq, we estimate the time integral şT
0 gptq2hptqdt with
ż T
0
gptq2hptqdt “ Z ż T
0
pptqαptq2hptqdt « TZαpt̃q2hpt̃q, (12)
where t̃ is a sample from pptq. When training score-based models with likelihood weighting, hptq corresponds to the denoising score matching loss at time t.
Ref. [33] also observes that optimizing the ELBO for diffusion probabilistic models has large variance, and proposes to reduce it with importance sampling. They build their proposal distribution based on historical loss values stored at thousands of discrete time steps. Despite this similarity, our method is easier to implement without needing to maintain history, can be used for evaluation, and is particularly suited to the continuous-time setting.

heading: 5.2 Variational dequantization
text: Digital images are discrete data, and must be dequantized when training continuous density models like normalizing flows [12, 13] and score-based diffusion models. One popular approach to this is uniform dequantization [53, 50], where we add small uniform noise over r0, 1q to images taking values in t0, 1, ¨ ¨ ¨ , 255u. As shown in [50], training a continuous model pθpxq on uniformly dequantized data implicitly maximizes a lower bound on the log-likelihood of a certain discrete model Pθpxq. Due to the gap between pθpxq and Pθpxq, comparing the likelihood of continuous density models to models which fit discrete data directly, such as autoregressive models [55] or variational autoencoders, naturally puts the former at a disadvantage.
To minimize the gap between pθpxq and Pθpxq, ref. [18] proposes variational dequantization, where a separate normalizing flow model qφpu | xq is trained to produce the dequantization noise by optimizing the following objective
max φ Ex„ppxqEu„qφp¨|xqrlog pθpx` uq ´ log qφpu | xqs. (13)
Plugging in the lower bound on log pθpxq from Theorem 3, we can optimize Eq. (13) to improve the likelihood of score-based diffusion models.

heading: 5.3 Experiments
text: We empirically test the performance of likelihood weighting, importance sampling and variational dequantization across multiple architectures of score-based models, SDEs, and datasets. In particular, we consider DDPM++ (“Baseline” in Table 2) and DDPM++ (deep) (“Deep” in Table 2) models with VP and subVP SDEs [48] on CIFAR-10 [28] and ImageNet 32ˆ32 [55] datasets. We omit experiments on the VE SDE since (i) under this SDE our likelihood weighting is the same as the original weighting in [48]; (ii) we empirically observe that the best VE SDE model achieves around 3.4 bits/dim on CIFAR-10 in our experiments, which is significantly worse than other SDEs. For each experiment, we report ´Erlog pODEθ pxqs (“Negative log-likelihood” in Table 2), and the upper bound ErLDSMθ pxqs on´Erlog pSDEθ pxqs (“Bound” in Table 2). In addition, we report FID scores [17] for samples from pODEθ , produced by solving the corresponding ODE with the Dormand–Prince RK45 [14] solver. Unless otherwise noted, we apply horizontal flipping as data augmentation for training models on CIFAR-10, so as to match the settings in [48, 19]. Detailed description of all our experiments can be found in Appendices B and C.
We summarize all results in Table 2. Our key observations are as follows:
1. Although Theorem 3 only guarantees ErLDSMθ pxqs ě ´Erlog pSDEθ pxqs, and in general we have pSDEθ ‰ pODEθ , we still find that ErLDSMθ pxqs (“Bound” in Table 2) ě ´Erlog pODEθ pxqs (“NLL” in Table 2) in all our settings.
2. When all conditions are fixed except for the weighting in the training objective, having a lower value of the bound for pSDEθ always leads to a lower negative log-likelihood for p ODE θ .
3. With only likelihood weighting, we can uniformly improve the likelihood of pODEθ and the bound of pSDEθ on CIFAR-10 across model architectures and SDEs, but it is not sufficient to guarantee likelihood improvement on ImageNet 32ˆ 32.
4. By combining importance sampling and likelihood weighting, we are able to achieve uniformly better likelihood for pODEθ and bounds for p SDE θ across all model architectures,
SDEs, and datasets, with only slight degradation of sample quality as measured by FID [17].
5. Variational dequantization uniformly improves both the bound for pSDEθ and the negative log-likelihood (NLL) of pODEθ in all settings, regardless of likelihood weighting.
Our experiments confirm that with importance sampling, likelihood weighting is not only effective for maximizing the lower bound for the log-likelihood of pSDEθ , but also improving the log-likelihood of pODEθ . In agreement with [19, 33], we observe that models achieving better likelihood tend to have worse FIDs. However, we emphasize that this degradation of FID is small, and samples actually
have no obvious difference in visual quality (see Figs. 3 and 4). To trade likelihood for FID, we can use weighting functions that interpolate between likelihood weighting and the original weighting functions in [48]. Our FID scores are still much better than most other likelihood-based models.
We term pODEθ a ScoreFlow when its corresponding score-based model sθpx, tq is trained with likelihood weighting, importance sampling, and variational dequantization combined. It can be viewed as a continuous normalizing flow, but is parameterized by a score-based model and trained in a more efficient way. With variational dequantization, we show ScoreFlows obtain competitive negative log-likelihoods (NLLs) of 2.83 bits/dim on CIFAR-10 and 3.76 bits/dim on ImageNet 32ˆ32. Here the ScoreFlow on CIFAR-10 is trained without horizontal flipping (different from the setting in Table 2). As shown in Table 3, our results are on a par with the stateof-the-art autoregressive models on these tasks, and outperform all existing normalizing flow
models. The likelihood for CIFAR-10 can be significantly improved by incorporating advanced data augmentation, as demonstrated in [25, 41]. While we do not compare against these approaches, we believe that incorporating the same data augmentation techniques could also improve the likelihood of ScoreFlows.

heading: 6 Conclusion
text: We propose an efficient training objective for approximate maximum likelihood training of scorebased diffusion models. Our theoretical analysis shows that the weighted combination of score matching losses upper bounds the negative log-likelihood when using a particular weighting function which we term the likelihood weighting. By minimizing this upper bound, we consistently improve the likelihood of score-based diffusion models across multiple model architectures, SDEs, and datasets. When combined with variational dequantization, we achieve competitive likelihoods on CIFAR-10 and ImageNet 32ˆ32, matching the performance of best-in-class autoregressive models. Our upper bound is analogous to the evidence lower bound commonly used for training variational autoencoders. Aside from promoting higher likelihood, the bound can be combined with other objectives that depend on the negative log-likelihood, and also enables joint training of the forward and backward SDEs, which we leave as a future research direction. Our results suggest that scorebased diffusion models are competitive alternatives to continuous normalizing flows which enjoy the same tractable likelihood computation but with more efficient maximum likelihood training.
Limitations and broader impact Despite promising experimental results, we would like to emphasize that there is no theoretical guarantee that improving the SDE likelihood will improve the ODE likelihood, and this is explicitly a limitation of our work. Score-based diffusion models also suffer from slow sampling. In our experiments, the ODE solver typically need around 550 and 450 evaluations of the score-based model for generation and likelihood computation on CIFAR-10 and ImageNet respectively, which is considerably slower than alternative generative models like VAEs and GANs. In addition, the current formulation of score-based diffusion models only supports continuous data, and cannot be naturally adapted to discrete data without resorting to dequantization. Similarly to other deep generative models, score-based diffusion models can potentially be used to generate harmful media contents such as ‘deepfakes’, and might reflect and amplify undesirable social bias that could exist in the training dataset.
Author Contributions
Yang Song wrote the code, ran the experiments, proposed and proved Theorems 1 and 3, and wrote most of the paper. Conor Durkan proposed and proved a first version of Theorem 2, and wrote the paper. Iain Murray and Stefano Ermon co-advised the project and provided helpful edits to the draft.
Acknowledgments and Disclosure of Funding
The authors would like to thank Sam Power, George Papamakarios, Adji Dieng for helpful feedback, and Duoduo for providing her photos in Fig. 1. This research was supported by NSF (#1651565, #1522054, #1733686), ONR (N000141912145), AFOSR (FA95501910024), ARO (W911NF-21-10125), Sloan Fellowship, and Google TPU Research Cloud. This research was also supported by the EPSRC Centre for Doctoral Training in Data Science, funded by the UK Engineering and Physical Sciences Research Council (grant EP/L016427/1), and the University of Edinburgh. Yang Song was supported by the Apple PhD Fellowship in AI/ML.

Review summaries:
This paper presents likelihood weighting to improve the log-likelihood performance of the score-based model. It is widely known that the score-based model is an excellent generative model to sample an instance from the high-dimensional data distribution. However, its encoding capability has remained question. This paper argues that the presented likelihood weighting can be universally applicable to the score-based models to improve such encoding capability without hurting the generation performance. There were series of theoretic derivations to induce the likelihood weighting, and a benchmark experiment was conducted.  

The authors show that for a specific weighting of score matching losses at various noise levels, the objective upper bounds the negative log-likelihood. They further bound the log-likelihood of individual datapoints to develop a continuous-time generalization of the ELBO.  Using the proposed objective the authors optimize score-based diffusion models towards likelihood. The authors further make use of 1) importance sampling of time steps for reduced variance, 2) variational dequantization for a tighter bound on the likelihood and 3) the use of the equivalent ODE during evaluation, which allows exact likelihood computation (apart from the dequantization) for an approximation to the SDE model. As a result of this, state-of-the-art likelihood performance is presented on CIFAR-10 and ImageNet 32x32, beating even the best autoregressive models.

The paper studies score-based denoising-diffusion generative models and makes one major contribution. It shows how the training objective, a mixture or score-matching losses defined along different times of the diffusion, can be adapted to enable maximum likelihood training (using a lower bound objective). This is rigorously derived and proven analytically. The result is a specific weighting of the score-matching losses for different times t. The paper validates the derived methodology and empirically demonstrates very strong likelihood performance, state-of-the-art among models without data augmentations.