Consistency Trajectory Models: Learning Probability Flow ODE Trajectory of Diffusion

Introduction

Deep generative models encounter distinct training and sampling challenges. Variational Autoencoder (VAE) (Kingma & Welling, 2013) can be trained easily but may suffer from posterior collapse, resulting in blurry samples, while Generative Adversarial Network (GAN) (Goodfellow et al., 2014) generates high-quality samples but faces training instability. Conversely, Diffusion Model (DM) (Sohl-Dickstein et al., 2015; Ho et al., 2020; Song et al., 2020b) addresses these issues by learning the score (i.e., gradient of log-density) (Song & Ermon, 2019), which can generate high quality samples. However, compared to VAE and GAN excelling at fast sampling, DM involves a gradual denoising process that slows down sampling, requiring numerous model evaluations.

Score-based diffusion models synthesize data by solving the reverse-time (stochastic or deterministic) process corresponding to a prescribed forward process that adds noise to the data (Song & Ermon, 2019; Song et al., 2020b). Although advanced numerical solvers (Lu et al., 2022b; Zhang & Chen, 2022) of Stochastic Differential Equations (SDE) or Ordinary Differential Equations (ODE) substantially reduce the required Number of Function Evaluations (NFE), further improvements are challenging due to the intrinsic discretization error present in all solvers (De Bortoli et al., 2021). Recent developments in sample efficiency thus focus on directly estimation of the integral along the sample trajectory, amortizing the computational cost of numerical solvers. Distillation models (Salimans & Ho, 2021) in Figure 1, exemplified by the Consistency Model (CM) (Song et al., 2023), presents a promising approach for estimating the integration with a single NFE (Figure 2). However, their generation quality does not improve as NFE increase, and there is no straightforward mechanism for balancing computational resources (NFE) with quality.

This paper introduces the Consistency Trajectory Model (CTM) as a unified framework simultaneously assessing both the integrand (score function) and the integral (sample) of the Probability Flow (PF) ODE, thus bridging score-based and distillation models (Figure 1). CTM estimates both infinitesimal steps (score function) and long steps (integral over any time horizon) of the PF ODE from any initial condition, providing increased flexibility at inference time. Its score evaluation capability accommodates a range of score-based sampling algorithms based on solving differential equations (Song et al., 2020b), expanding its applicability across various domains (Saharia et al., 2022a). In particular, CTM enables exact likelihood computation, setting it apart from previous distillation models. Additionally, its integral approximation capability facilitates the incorporation of distillation sampling methods (Salimans & Ho, 2021; Song et al., 2023) that involve long “jumps” along the solution trajectory. This unique feature enables a novel sampling method called $\gamma$-sampling, which alternates forward and backward jumps along the solution trajectory, with $\gamma$ governing the level of stochasticity.

CTM’s dual modeling capability for both infinitesimal and long steps of the PF ODE greatly enhances its training flexibility as well. It allows concurrent training with reconstruction loss, denoising diffusion loss, and adversarial loss within a unified framework. Notably, by incorporating CTM’s training approach with $\gamma$-sampling, we achieve the new State-Of-The-Art (SOTA) performance in both density estimation and image generation for CIFAR-10 (Krizhevsky et al., 2009) (Figure 3) and ImageNet (Russakovsky et al., 2015) at a resolution of $64\times 64$ (Table 5).

Preliminary

In DM (Sohl-Dickstein et al., 2015; Song et al., 2020b), the encoder structure is formulated using a set of continuous-time random variables defined by a fixed forward diffusion processThis paper can be extended to VPSDE encoding (Song et al., 2020b) with re-scaling (Kim et al., 2022a).,

Sampling from DM involves solving the PF ODE, equivalent to computing the integral

where $\mathbf{x}_{T}$ is sampled from a prior distribution $\pi$ approximating $p_{T}$. Decoding strategies of DM primarily fall into two categories: score-based sampling with time-discretized numerical integral solvers, and distillation sampling where a neural network directly estimates the integral.

Any off-the-shelf ODE solver, denoted as $\texttt{Solver}(\mathbf{x}_{T},T,0;\bm{\phi})$ (with an initial value of $\mathbf{x}_{T}$ at time $T$ and ending at time ), can be directly applied to solve Eq. (1) (Song et al., 2020b). For instance, DDIM (Song et al., 2020a) corresponds to a 1st-order Euler solver, while EDM (Karras et al., 2022) introduces a 2nd-order Heun solver. Despite recent advancements in numerical solvers (Lu et al., 2022b; Zhang & Chen, 2022), further improvements may be challenging due to the inherent discretization error present in all solvers (De Bortoli et al., 2021), ultimately limiting the sample quality obtained with few NFEs.

Distillation Sampling

Distillation models (Salimans & Ho, 2021; Meng et al., 2023) successfully amortize the sampling cost by directly estimating the integral of Eq. (1) with a single neural network evaluation. However, their multistep sampling approach (Song et al., 2023) exhibits degrading sample quality with increasing NFE, lacking a clear trade-off between computational budget (NFE) and sample fidelity. Furthermore, multistep sampling is not deterministic, leading to uncontrollable sample variance. We refer to Appendix A for a thorough literature review.

CTM: An Unification of Score-based and Distillation Models

To address the challenges in both score-based and distillation samplings, we introduce the Consistency Trajectory Model (CTM), which seamlessly integrates both decoding strategies. Consequently, our model is versatile and can perform sampling through either SDE/ODE solving or direct prediction of intermediate points along the PF ODE trajectory.

CTM predicts both infinitesimal changes and intermediate points of the PF ODE trajectory. Specifically, we define $G(\mathbf{x}_{t},t,s)$ as the solution of the PF ODE from initial time $t$ with an initial condition $\mathbf{x}_{t}$ to final time $s\leq t$:

$G$ can access any intermediate point along the trajectory by varying final time $s$. However, with the current expression of $G$, the infinitesimal change needed to recover the denoiser information (the integrand) can only be obtained by evaluating the $s$-derivative at time $t$, $\frac{\partial}{\partial s}G(\mathbf{x}_{t},t,s)|_{s=t}$. Therefore, we introduce a dedicated expression for $G$ using an auxiliary function $g$ to enable easy access to both the integral via $G$ and the integrand via $g$ with Lemma 1.

When $s=0$, $G(\mathbf{x}_{t},t,0)=g(\mathbf{x}_{t},t,0)$ is the solution of PF ODE at $s=0$, initialized at $\mathbf{x}_{t}$.

Indeed, the $G$’s expression in Lemma 1 is naturally linked to the Taylor approximation to the integral:

for any $s\leq t$. Here, it is evident that $g$ includes all residual terms in Taylor expansion, which turns to be the discretization error in sampling. The goal of CTM is to approximate this $g$-function using a neural network $g_{\bm{\theta}}$ and estimate the solution trajectory with the parametrization inspired by Lemma 1 as follows:

We remark that this parametrization satisfies the initial condition $G_{\bm{\theta}}(\mathbf{x}_{t},t,t)=\mathbf{x}_{t}$ for free, leading to improved training stabilityEnsuring the initial condition’s satisfaction is crucial for stable training. Directly estimating $G$ with a network leads to rapid divergence, causing instability as the network output may deviate arbitrarily from the initial condition.. Appendix C.1 offers further insights into this parametrization.

2 CTM Training

To achieve trajectory learning, CTM should match the model prediction to the ground truth $G$ by

for any $s\leq t$. We opt to approximate $G$ by solving the empirical PF ODE with a pre-trained score model $D_{\bm{\phi}}$. Our neural network is then trained to align with the reconstruction:

In a scenario with no ODE discretization and no score approximation errors, Solver perfectly reconstructs the PF ODE trajectory, and comparing the prediction and reconstruction in Eq. (3) leads $G_{\bm{\theta}^{*}}(\mathbf{x}_{t},t,s)=G(\mathbf{x}_{t},t,s)$ at optimal $\bm{\theta}^{*}$, given sufficient network flexibility. The same conclusion holds by matching the local consistency:

where $\texttt{sg}(\cdot)$ is stop-gradient. With the initial condition ($\mathbf{G}_{\bm{\theta}}(\mathbf{x}_{t},t,t)=\mathbf{x}_{t}$) satisfied, matching Eq. (4) avoids collapsing to the trivial solution (Proposition 4 in Appendix B.1).

To estimate the entire solution trajectory with higher precision, we introduce soft matching, illustrated in Figure 5, ensuring consistency between prediction from $\mathbf{x}_{t}$ and the prediction from $\texttt{Solver}(\mathbf{x}_{t},t,u;\bm{\phi})$ for any $u\in[s,t)$:

As $u=s$, Eq. (3) enforces global consistency matching, i.e., a reconstruction loss.

As $u=t-\Delta t$, Eq. (4) is local consistency matching. Additionally, if $s=0$, it recovers CM’s distillation loss.

To quantify the dissimilarity between $G_{\bm{\theta}}(\mathbf{x}_{t},t,s)$ and $G_{\texttt{sg}(\bm{\theta})}(\texttt{Solver}(\mathbf{x}_{t},t,u;\bm{\phi}),u,s)$ and enforce Eq. (5), we could use either distance in pixel space or in feature space. However, the pixel distance may overemphasize the distance at large $s$ due to the diffusion scale, requiring time-weighting adjustments. Furthermore, feature distance requires a time-conditional feature extractor, which can be expensive to train. Hence, we propose to use a feature distance $d$ in clean data space by comparing

which leads the model’s prediction, at optimum, to match with the empirical PF ODE’s solution trajectory, defined by the pre-trained DM (teacher), see Appendix B (Propositions 3 and 5) for details.

3 Training Consistency Trajectory Models

Empirically, regularizing $\mathcal{L}_{\text{CTM}}$ with $\mathcal{L}_{\text{DSM}}$ improves score accuracy, which is especially important in large NFE sampling regimes.

On the other hand, CTM, distilling from the teacher model, is constrained by the teacher’s $D_{\bm{\phi}}$ performance. This challenge can be mitigated with adversarial training to improve trajectory estimation. The one-step generation of CTM enables us to calculate the adversarial loss efficiently, in the similar way of conventional GAN training:

where $d_{\bm{\eta}}$ is a discriminator. This adversarial training allows the student model (CTM) to beat the teacher model (DM). To summarize, CTM allows the integration of reconstruction-based CTM loss, diffusion loss, and adversarial loss

in a single training framework, by optimizing $\min_{\bm{\theta}}\max_{\bm{\eta}}\mathcal{L}(\bm{\theta},\bm{\eta})$. Here, $\lambda_{\text{DSM}}\geq 0$ and $\lambda_{\text{GAN}}\geq 0$ are the weighting functions, see Algorithm 1.

Sampling with CTM

CTM enables score evaluation through $g_{\bm{\theta}}(\bm{x}_{t},t,t)$, supporting standard score-based sampling with ODE/SDE solvers. In high-dimensional image synthesis, as shown in Figure 6’s left two columns, CTM performs comparably to EDM using Heun’s method as a PF ODE solver.

CTM additionally enables time traversal along the solution trajectory, allowing for the newly introduced $\gamma$-sampling method, refer to Algorithm 3 and Figure 7. Suppose the sampling timesteps are $T=t_{0}>\cdots>t_{N}=0$. With $\mathbf{x}_{t_{0}}\sim\pi$, where $\pi$ is the prior distribution, $\gamma$-sampling denoises $\mathbf{x}_{t_{0}}$ to time $\sqrt{1-\gamma^{2}}t_{1}$ with $G_{\bm{\theta}}(\mathbf{x}_{t_{0}},t_{0},\sqrt{1-\gamma^{2}}t_{1})$, and perturb this denoised sample with forward diffusion to the noise level at time $t_{1}$. It iterates this back-and-forth traversal until reaching to time $t_{N}=0$.

Our $\gamma$-sampling is a new distillation sampler that unifies previously proposed sampling techniques, including distillation sampling and score-based sampling.

Figure 7-(c): When $\gamma=0$, it becomes the deterministic distillation sampling that estimates the solution of the PF ODE. A key distinction between the $\gamma$-sampling and score-based sampling is that CTM avoids sampling errors by directly estimating Eq. (2). However, score-based samplers like DDIM (1st-order Euler solver) or EDM (2nd-order Heun solver) are susceptible to discretization errors from Taylor approximation, especially with small NFE. (the leftmost column of Figure 6). Deterministic nature as $\gamma=0$ ensures the sample semantic preserved across NFE changes, visualized in the rightmost column of Figure 6.

Figure 7-(b): When $0<\gamma<1$, it generalizes the EDM’s stochastic sampler (Algorithm 2). Appendix B.3 shows that $\gamma$-sampling’s sample variances scale proportionally with $\gamma^{2}$.

The optimal choice of $\gamma$ depends on practical usage and empirical configuration (Karras et al., 2022; Xu et al., 2023). Figure 8 demonstrates $\gamma$-sampling in stroke-based generation (Meng et al., 2021), revealing that the sampler with $\gamma=1$ leads to significant semantic deviations from the reference stroke, while smaller $\gamma$ values yield closer semantic alignment and maintain high fidelity. In contrast, Figure 9 showcases $\gamma$’s impact on generation performance. In Figure 9-(a), $\gamma$ has less influence with small NFE, but the setup with $\gamma\approx 0$ is the only one that resembles the performance of the Heun’s solver as NFE increases. Additionally, CM’s multistep sampler ($\gamma=1$) significantly degrades sample quality as NFE increases. This quality deterioration concerning $\gamma$ becomes more pronounced with higher NFEs, shown in Figure 9-(b), potentially attributed to error accumulation during the iterative long “jumps” for denoising. We explain this phenomenon using a 2-step $\gamma$-sampling example in the following theorem, see Theorem 8 for a generalized result for $N$-steps.

Let $t\in(0,T)$ and $\gamma\in$. Denote $p_{\bm{\theta}^{*},2}$ as the density obtained from the $\gamma$-sampler with the optimal CTM, following the transition sequence $T\rightarrow\sqrt{1-\gamma^{2}}t\rightarrow t\rightarrow 0$, starting from $p_{T}$. Then D_{TV}\left(p_{\text{data}},p_{\bm{\theta}^{*},2}\right)=\mathcal{O}\big{(}\sqrt{T-\sqrt{1-\gamma^{2}}t+t}\big{)}.

When it becomes $N$-steps, $\gamma=1$-sampling iteratively conducts long jumps from $t_{n}$ to for each step $n$, which aggregates the error to be $\mathcal{O}(\sqrt{T+t_{1}+\cdots+t_{N}})$. In contrast, such time overlap between jumps does not occur in $\gamma=0$-sampling, eliminating the error accumulation, resulting in $\mathcal{O}(\sqrt{T})$ error, see Appendix C.2. In summary, CTM addresses challenges associated with large NFE in distillation models with $\gamma=0$ and removes the discretization error in score-based models.

Experiments

We evaluate CTM on CIFAR-10 and ImageNet $64\times 64$, using the pre-trained diffusion checkpoints from EDM for CIFAR-10 and CM for ImageNet as the teacher models. We adopt EDM’s training configuration for $\mathcal{L}_{\text{DSM}}(\bm{\theta})$ and employ StyleGAN-XL’s (Sauer et al., 2022) discriminator for $\mathcal{L}_{\text{GAN}}(\bm{\theta},\bm{\eta})$. During training, we employ adaptive weights $\lambda_{\text{DSM}}$ and $\lambda_{\text{GAN}}$, inspired by VQGAN (Esser et al., 2021) to balance DSM and GAN losses with the CTM loss. For both datasets, we utilize the DDPM architecture. For CIFAR-10, we take EDM’s implementation; and for ImageNet, CM’s implementation is used. On top of these architectures, we incorporate $s$-information via auxiliary temporal embedding with positional embedding (Vaswani et al., 2017), and add this embedding to the $t$-embedding. This training setup (Appendix D), along with the deterministic sampling ($\gamma=0$), allows CTM’s generation to outperform teacher models with NFE $1$ and achieve SOTA FIDs with NFE $2$.

CIFAR-10 CTM’s NFE $1$ generation excels both EDM and StyleGAN-XL with FID of $1.73$ on conditional CIFAR-10, and CTM achieves the SOTA FID of $1.63$ with $2$ NFEs, surpassing all generative models. These results are obtained with the implementation based on the official PyTorch code of CM. However, retraining CM with this official PyTorch code yields FID of $10.53$ (unconditional), higher than the reported FID of $3.55$. Additionally, CTM’s ability to approximate scores using $g_{\bm{\theta}}(\mathbf{x}_{t},t,t)$ enables evaluating Negative Log-Likelihood (NLL) (Song et al., 2021; Kim et al., 2022b), establishing a new SOTA NLL. This improvement can be attributed, in part, to CTM’s reconstruction loss when $u=s$, and improved alignment with the oracle process (Lai et al., 2023a).

ImageNet CTM’s generation surpasses both teacher EDM and StyleGAN-XL with NFE 1, outperforming previous models with no guidance (Dhariwal & Nichol, 2021), see Figure 11 for the comparison of CTM with the teacher model. Notably, all results in Tables 5 and 5 are achieved within $30$K-$100$K training iterations, requiring only $10\%$ of the iterations compared to CM and EDM.

Classifier-Rejection Sampling CTM’s fast sampling enables classifier-rejection sampling. In the evaluation, for each class, we select the top 50 samples out of $\frac{50}{1-r}$ samples based on predicted class probability, where $r$ is the rejection ratio. This sampler, combined with NFE $1$ sampling, consumes an average of NFE $\frac{1}{1-r}$. In Figure 10, CTM, employing cost-effective classifier-rejection sampling, shows a FID-IS trade-off comparable to classifier-guided results (Ho & Salimans, 2021) achieved with high NFEs of 250. Additionally, Figure 12 confirms that samples rejected by the classifier exhibit superior quality and maintain class consistency, in agreement with the findings of Ho & Salimans (2021). We employ the classifier at resolution of $64\times 64$ provided by Dhariwal & Nichol (2021).

2 Qualitative Analysis

CTM Loss Figure 13 highlights the advantages of employing the proposed soft consistency matching in Eq. (5) during CTM training. It outperforms the local consistency matching (Eq. (4)). Additionally, it demonstrates comparable performance to the global consistency matching (Eq. 3) with NFE $1$, superior performance with large NFE. Furthermore, soft matching is computationally efficient, enhancing the scalability of CTM.

DSM Loss Figure 14 illustrates two benefits of incorporating $\mathcal{L}_{\text{DSM}}$ with $\mathcal{L}_{\text{CTM}}$. It preserves sample quality for small NFE unless DSM scale outweighs CTM. For large NFE sampling, it significantly improves sample quality due to accurate score estimation. Throughout the paper, we maintain $\lambda_{\text{DSM}}=1$ based on insights from Figure 14, unless otherwise specified.

GAN Loss Analogous to the DSM loss, Figure 15 illustrates the advantages of incorporating the GAN loss for both small and large NFE sample quality. Figure 11 demonstrates that CTM can produce samples resembling those of EDM (teacher), with GAN refining local details. Throughout the paper, we adopt the warm-up strategy for GAN training: deactivate GAN training with $\lambda_{\text{GAN}}=0$ for warm-up iterations and then activate GAN training with $\lambda_{\text{GAN}}=1$, in line with the recommendation from VQGAN (Esser et al., 2021). This warm-up strategy is applied by default unless otherwise specified.

Conclusion

CTM, a novel generative model, addresses issues in established models. With a unique training approach accessing intermediate PF ODE solutions, it enables unrestricted time traversal and seamless integration with prior models’ training advantages. A universal framework for Consistency and Diffusion Models, CTM excels in both training and sampling. Remarkably, it surpasses its teacher model, achieving SOTA results in FID and likelihood for few-steps diffusion model sampling on CIFAR-10 and ImageNet $64\times 64$, highlighting its versatility and process.

Acknowledgement

We sincerely acknowledge the support of everyone who made this research possible. Our heartfelt thanks go to Koichi Saito, Woosung Choi, Kin Wai Cheuk, and Yukara Ikemiya for their assistance.

References

Appendix A Related Works

DMs excel in high-fidelity synthetic image and audio generation (Dhariwal & Nichol, 2021; Saharia et al., 2022b; Rombach et al., 2022), as well as in applications like media editing, restoration (Meng et al., 2021; Cheuk et al., 2023; Kawar et al., 2022; Saito et al., 2023; Hernandez-Olivan et al., 2023; Murata et al., 2023). Recent research aims to enhance DMs in sample quality (Kim et al., 2022b; a), density estimation (Song et al., 2021; Lu et al., 2022a), and especially, sampling speed (Song et al., 2020a).

Fast Sampling of DMs

The SDE framework underlying DMs (Song et al., 2020b) has driven research into various numerical methods for accelerating DM sampling, exemplified by works such as (Song et al., 2020a; Zhang & Chen, 2022; Lu et al., 2022b). Notably, (Lu et al., 2022b) reduced the ODE solver steps to as few as $10$-$15$. Other approaches involve learning the solution operator of ODEs (Zheng et al., 2023), discovering optimal transport paths for sampling (Liu et al., 2022), or employing distillation techniques (Luhman & Luhman, 2021; Salimans & Ho, 2021; Berthelot et al., 2023; Shao et al., 2023). However, previous distillation models may experience slow convergence or extended runtime. Gu et al. (2023) introduced a bootstrapping approach for data-free distillation. Furthermore, Song et al. (2023) introduced CM which extracts DMs’ PF ODE to establish a direct mapping from noise to clean predictions, achieving one-step sampling while maintaining good sample quality. CM has been adapted to enhance the training stability of GANs, as (Lu et al., 2023). However, it’s important to note that their focus does not revolve around achieving sampling acceleration for DMs, nor are the results restricted to simple datasets.

Consistency of DMs

Score-based generative models rely on a differential equation framework, employing neural networks trained on data to model the conversion between data and noise. These networks must satisfy specific consistency requirements due to the mathematical nature of the underlying equation. Early investigations, such as (Kim et al., 2022c), identified discrepancies between learned scores and ground truth scores. Recent developments have introduced various consistency concepts, showing their ability to enhance sample quality (Daras et al., 2023; Li et al., 2023), accelerate sampling speed (Song et al., 2023), and improve density estimation in diffusion modeling (Lai et al., 2023a). Notably, Lai et al. (2023b) established the theoretical equivalence of these consistency concepts, suggesting the potential for a unified framework that can empirically leverage their advantages. CTM can be viewed as the first framework which achieves all the desired properties.

Appendix B Theoretical Insights on CTM

In this section, we explore several theoretical aspects of CTM, encompassing convergence analysis (Section B.1), properties of well-trained CTM, variance bounds for $\gamma$-sampling, and a more general form of accumulated errors induced by $\gamma$-sampling (cf. Theorem 2).

We first introduce and review some notions. Starting at time $t$ with an initial value of $\mathbf{x}_{t}$ and ending at time $s$, recall that $G(\mathbf{x}_{t},t,s)$ represents the true solution of the PF ODE, and $G(\mathbf{x}_{t},t,s;\bm{\phi})$ is the solution function of the following empirical PF ODE.

Here $\bm{\phi}$ denotes the teacher model’s weights learned from DSM. Thus, $G(\mathbf{x}_{t},t,s;\bm{\phi})$ can be expressed as

CTM’s practical implementation follows CM’s one, utilizing discrete timesteps $t_{0}=0<t_{1}<\cdots<t_{N}=T$ for training. Initially, we assume local consistency matching for simplicity, but this can be extended to soft matching. This transforms the CTM loss in Eq. (7) to the discrete time counterpart:

In the following theorem, we demonstrate that irrespective of the initial time $t_{n}$ and end time $t_{m}$, CTM $G_{\bm{\theta}}(\cdot,t_{n},t_{m};\bm{\phi})$, will eventually converge to its teacher model, $G(\cdot,t_{n},t_{m};\bm{\phi})$.

Define $\Delta_{N}t:=\underset{n\in[\![1,N]\!]}{\max}\left\{\left|t_{n+1}-t_{n}\right|\right\}$. Assume that $G_{\bm{\theta}}$ is uniform Lipschitz in $\bm{x}$ and that the ODE solver admits local truncation error bounded uniformly by $\mathcal{O}((\Delta_{N}t)^{p+1})$ with $p\geq 1$. If there is a $\bm{\theta}_{N}$ so that $\mathcal{L}_{\text{CTM}}^{N}(\bm{\theta}_{N};\bm{\phi})=0$, then for any $n\in[\![1,N]\!]$ and $m\in[\![1,n]\!]$

Similar argument applies, confirming convergence along the PF ODE trajectory, ensuring Eq. (4) with $\bm{\theta}$ replacing $\texttt{sg}(\bm{\theta})$:

Convergence of Densities.

In Proposition 3, we demonstrated point-wise trajectory convergence, from which we infer that CTM may converge to its training target in terms of density. More precisely, in Proposition 5, we establish that if CTM’s target $\mathbf{x}_{\text{target}}$ is derived from the teacher model (as defined above), then the data density induced by CTM will converge to that of the teacher model. Specifically, if the target $\mathbf{x}_{\text{target}}$ perfectly approximates the true $G$-function:

Then the data density generated by CTM will ultimately learn the data distribution $p_{\text{data}}$.

The uniform Lipschitzness of $G_{\bm{\theta}}$ (and $G$),

The uniform boundedness in $\bm{\theta}$ of $G_{\bm{\theta}}$: there is a $L(\mathbf{x})\geq 0$ so that

If for any $N$, there is a $\bm{\theta}_{N}$ such that $\mathcal{L}_{\text{CTM}}^{N}(\bm{\theta}_{N};\bm{\phi})=0$. Let $p_{\bm{\theta}_{N}}(\cdot)$ denote the pushforward distribution of $p_{T}$ induced by $G_{\bm{\theta}_{N}}(\cdot,T,0)$. Then, as $N\rightarrow\infty$, $\left\lVert p_{\bm{\theta}_{N}}(\cdot)-p_{\bm{\phi}}(\cdot)\right\rVert_{\infty}\rightarrow 0$. Particularly, if the condition in Eq. (10) is satisfied, then $\left\lVert p_{\bm{\theta}_{N}}(\cdot)-p_{\text{data}}(\cdot)\right\rVert_{{\infty}}\rightarrow 0$ as $N\rightarrow\infty$.

B.2 Non-Intersecting Trajectory of the Optimal CTM

CTM learns distinct trajectories originating from various initial points $\mathbf{x}t$ and times $t$. In the following proposition, we demonstrate that the distinct trajectories derived by the optimal CTM, which effectively distills information from its teacher model ($G_{\bm{\theta}^{*}}(\cdot,t,s)\equiv G(\cdot,t,s;\bm{\phi})$ for any $t,s\in[0,T]$), do not intersect.

This implies that $\mathbf{x}_{t}\neq\mathbf{y}_{t}$, $G_{\bm{\theta}^{*}}(\mathbf{x}_{t};t,s)\neq G_{\bm{\theta}^{*}}(\mathbf{y}_{t};t,s)$ for all $s\in[0,t]$.

Specifically, the mapping from an initial value to its corresponding solution trajectory, denoted as $\mathbf{x}_{t}\mapsto G_{\bm{\theta}^{*}}(\mathbf{x}_{t},t,\cdot)$, is injective. Conceptually, this ensures that if we use guidance at intermediate times to shift a point to another guided-target trajectory, the guidance will continue to affect the outcome at $t=0$.

B.3 Variance Bounds of γ𝛾\gamma-sampling

Suppose the sampling timesteps are $T=t_{0}>t_{1}>\cdots>t_{N}=0$. In Proposition 7, we analyze the variance of

resulting from $n$-step $\gamma$-sampling, initiated at

Here, we assume an optimal CTM which precisely distills information from the teacher model $G_{\bm{\theta}^{*}}(\cdot)=G(\cdot,t,s;\bm{\phi})$ for all $t,s\in[0,T]$, for simplicity.

where $\zeta(t_{n},t_{n+1},\gamma)=\exp{\left(2L_{\bm{\phi}}(t_{n}-\sqrt{1-\gamma^{2}}t_{n+1})\right)}$ and $L_{\bm{\phi}}$ is a Lipschitz constant of $D_{\bm{\phi}}(\cdot,t)$.

In line with our intuition, CM’s multistep sampling ($\gamma=1$) yields a broader range of $\text{Var}\left(X_{n+1}\right)$ compared to $\gamma=0$, resulting in diverging semantic meaning with increasing sampling NFE.

B.4 Accumulated Errors in the General Form of γ𝛾\gamma-sampling.

We can extend Theorem 2 for two steps $\gamma$-sampling for the case of multisteps.

We begin by clarifying the concept of “density transition by a function”. For a measurable mapping $\mathcal{T}:\Omega\rightarrow\Omega$ and a measure $\nu$ on the measurable space $\Omega$, the notation $\mathcal{T}\sharp\nu$ denotes the pushforward measure, indicating that if a random vector $X$ follows the distribution $\nu$, then $\mathcal{T}(X)$ follows the distribution $\mathcal{T}\sharp\nu$.

Given a sampling timestep $T=t_{0}>t_{1}>\cdots>t_{N}=0$. Let $p_{\bm{\theta}^{*},N}$ represent the density resulting from N-steps of $\gamma$-sampling initiated at $p_{T}$. That is,

B.5 Transition Densities with the Optimal CTM

In this section, for simplicity, we assume the optimal CTM, $G_{\bm{\theta}^{*}}\equiv G$ with a well-learned $\bm{\theta}^{*}$, which recovers the true $G$-function. We establish that the density propagated by this optimal CTM from any time $t$ to a subsequent time $s$ aligns with the predefined density determined by the fixed forward process.

We now present the proposition ensuring alignment of the transited density.

Then for any $t\in[0,T]$ and $s\in[0,t]$, $p_{s}=\mathcal{T}_{t\rightarrow s}\sharp p_{t}$.

This theorem guarantees that by learning the optimal CTM, which possesses complete trajectory information, we can retrieve all true densities at any time using CTM.

Appendix C Algorithmic Details

Our parametrization of $G_{\bm{\theta}}$ is affected from the discretized ODE solvers. For instance, the one-step Euler solver has the solution of

Again, the solver scales $\mathbf{x}_{t}$ with $\frac{s}{t}$ and multiply $1-\frac{s}{t}$ to the second term. Therefore, our $G(\mathbf{x}_{t},t,s)=\frac{s}{t}\mathbf{x}_{t}+(1-\frac{s}{t})g(\mathbf{x}_{t},t,s)$ is a natural way to represent the ODE solution.

For future research, we establish conditions enabling access to both integral and integrand expressions. Consider a continuous real-valued function $a(t,s)$. We aim to identify necessary conditions on $a(t,s)$ for the expression of $G$ as:

for a vector-value function $h(\mathbf{x}_{t},t,s)$ and that $h$ satisfies:

$\lim_{s\rightarrow t}h(\mathbf{x}_{t},t,s)$ exists;

Starting with the definition of $G$, we can obtain

Suppose that there is a continuous function $c(t)$ so that

The second equality follows from the mean value theorem (We omit the continuity argument details for Markov filtrations). Therefore, we obtain the desired property 2). We summarize the necessary conditions on $a(s,t)$ as:

We now explain the above observation with an example by considering EDM-type parametrization. Consider $c_{\text{skip}}=c_{\text{skip}}(t,s):=\sqrt{\frac{(s-\sigma_{\text{min}})^{2}+\sigma_{\text{data}}^{2}}{(t-\sigma_{\text{min}})^{2}+\sigma_{\text{data}}^{2}}}$ and c_{\text{out}}=c_{\text{out}}(t,s):=\Big{(}1-\frac{s}{t}\Big{)}. Then $G(\mathbf{x}_{t},t,s)$ can be expressed as

Then, we can verify that $c_{\text{skip}}$ satisfies the condition in Eq. (12) and that

The DSM loss with this $c_{\text{skip}}$ becomes

However, empirically, we find that the parametrization of $c_{\text{skip}}(t,s)$ and $c_{\text{out}}(t,s)$ other than the ODE solver-oriented one, i.e., $c_{\text{skip}}(t,s)=\frac{s}{t}$ and $c_{\text{skip}}(t,s)=1-\frac{s}{t}$, faces training instability. Therefore, we set $G(\mathbf{x}_{t},t,s)=\frac{s}{t}\mathbf{x}_{t}+(1-\frac{s}{t})g(\mathbf{x}_{t},t,s)$ as our default design and estimate $g$-function with the neural network.

C.2 Characteristics of γ𝛾\gamma-sampling

Connection with SDE When $G_{\bm{\theta}}=G$, a single step of $\gamma$-sampling is expressed as:

where $\bm{\epsilon}\sim\mathcal{N}(0,\mathbf{I})$. This formulation cannot be interpreted as a differential form (Øksendal, 2003) because it look-ahead future information (from $t_{n+1}$ to $\sqrt{1-\gamma^{2}}t_{n+1}$) to generate the sample $\mathbf{x}_{t_{n+1}}^{\gamma}$ at time $t_{n+1}$. This suggests that there is no Itô’s SDE that corresponds to our $\gamma$-sampler pathwisely, opening up new possibilities for the development of a new family of diffusion samplers.

Connection with EDM’s stochastic sampler We conduct a direct comparison between EDM’s stochastic sampler and CTM’s $\gamma$-sampling. We denote $\texttt{Heun}(\mathbf{x}_{t},t,s)$ as Heun’s solver initiated at time $t$ and point $\mathbf{x}_{t}$ and ending at time $s$. It’s worth noting that EDM’s sampler inherently experiences discretization errors stemming from the use of Heun’s solver, while CTM is immune to such errors.

C.3 Trajectory Control with Guidance

We could apply $\gamma$-sampling for application tasks, such as image inpainting or colorization, using the (straightforwardly) generalized algorithm suggested in CM. In this section, however, we propose a loss-based trajectory optimization algorithm in Algorithm 4 for potential application downstream tasks.

Appendix D Implementation Details

Following Karras et al. (2022), we utilize the EDM’s skip scale and output scale for $g_{\bm{\theta}}$ modeling as

where $\text{NN}_{\bm{\theta}}$ refers to a neural network that takes the same input arguments as $g_{\bm{\theta}}$. The advantage of this EDM-style skip and output scaling is that if we copy the teacher model’s parameters to the student model’s parameters, except student model’s $s$-embedding structure, $g_{\bm{\theta}}(\mathbf{x}_{t},t,t)$ initialized with $\bm{\phi}$ would be close to the teacher denoiser $D_{\bm{\phi}}(\mathbf{x}_{t},t)$. This good initialization partially explains the fast convergence speed.

We use 4$\times$V100 (16G) GPUs for CIFAR-10 experiments and 8$\times$A100 (40G) GPUs for ImageNet experiments. We use the warm-up for $\lambda_{\text{GAN}}$ hyperparameter. On CIFAR-10, we deactivate GAN training with $\lambda_{\text{GAN}}=0$ until 50k training iterations and activate the generator training with the adversarial loss (added to CTM and DSM losses) by increasing $\lambda_{\text{GAN}}$ to one. The minibatch per GPU is 16 in the CTM+DSM training phase, and 11 in the CTM+DSM+GAN training phase. On ImageNet, due to the excessive training budget, we deactivate GAN only for 10k iterations and activate GAN training afterwards. We fix the minibatch to be 11 throughout the CTM+DSM or the CTM+DSM+GAN training in ImageNet.

We follow the training configuration mainly from CM, but for the discriminator training, we follow that of StyleGAN-XL (Sauer et al., 2022). For $\mathcal{L}_{\text{CTM}}$ calculation, we use LPIPS (Zhang et al., 2018) as a feature extractor. We choose $t$ and $s$ from the $N$-discretized timesteps to calculate $\mathcal{L}_{\text{CTM}}$, following CM. Across the training, we choose the maximum number of ODE steps to prevent a single iteration takes too long time. For CIFAR-10, we choose $N=18$ and the maximum number of ODE steps to be 17. For ImageNet, we choose $N=40$ and the maximum number of ODE steps to be 20. We find the tendency that the training performance is improved by the number of ODE steps, so one could possibly improve our ImageNet result by choosing larger maximum ODE steps.

For $\mathcal{L}_{\text{DSM}}$ calculation, we select $50\%$ of time sampling from EDM’s original scheme of $t\sim\mathcal{N}(-1.2,1.2^{2})$. For the other half time, we first draw sample from $\xi\sim[0,0.7]$ and transform it using $(\sigma_{\text{max}}^{1/\rho}+\xi(\sigma_{\text{min}}^{1/\rho}-\sigma_{\text{max}}^{1/\rho}))^{\rho}$. This specific time sampling blocks the neural network to forget the denoiser information for large time. For $\mathcal{L}_{\text{GAN}}$ calculation, we use two feature extractors to transform GAN input to the feature space: the EfficientNet (Tan & Le, 2019) and DeiT-base (Touvron et al., 2021). Before obtaining an input’s feature, we upscale the image to 224x224 resolution with bilinear interpolation. After transforming to the feature space, we apply the cross-channel mixing and cross-scale mixing to represent the input with abundant and non-overlapping features. The output of the cross-scale mixing is a feature pyramid consisting of four feature maps at different resolutions (Sauer et al., 2022). In total, we use eight discriminators (four for EfficientNet features and the other four for DeiT-base features) for GAN training.

Following CM, we apply Exponential Moving Average (EMA) to update $\texttt{sg}(\bm{\theta})$ by

However, unlike CM, we find that our model bestly works with $\mu=0.999$ or $\mu=0.9999$, which largely remedy the subtle instability arise from GAN training. Except for the unconditional CIFAR-10 training with $\bm{\phi}$, we set $\mu$ to be 0.999 as default. Throughout the experiments, we use $\sigma_{\text{min}}=0.002$, $\sigma_{\text{max}}=80$, $\rho=7$, and $\sigma_{\text{data}}=0.5$.

D.2 Evaluation Details

For likelihood evaluation, we solve the PF ODE, following the practice suggested in Kim et al. (2022b) with the RK45 (Dormand & Prince, 1980) ODE solver of $\texttt{tol}=1e-3$ and $t_{\text{min}}=0.002$.

Throughout the paper, we choose $\gamma=0$ otherwise stated. In particular, for Tables 5 and 5, we report the sample quality metrics based on either the one-step sampling of CM or the $\gamma=0$ sampling for NFE 2 case. For CIFAR-10, we calculate the FID score based on Karras et al. (2022) statistics. For ImageNet, we compute the metrics following Dhariwal & Nichol (2021) and their pre-calculated statistics. For the StyleGAN-XL ImageNet result, we recalculated the metrics based on the statistics released by Dhariwal & Nichol (2021), using StyleGAN-XL’s official checkpoint.

For large-NFE sampling, we follow the EDM’s time discretization. Namely, if we draw $n$-NFE samples, we equi-divide $$with$n$points and transform it (say$\xi$) to the time scale by$(\sigma_{\text{max}}^{1/\rho}+(\sigma_{\text{min}}^{1/\rho}-\sigma_{\text{max}}^{1/\rho})\xi)^{\rho}$. However, we emphasize the time discretization for both training and sampling is a modeler’s choice.

Appendix E Additional Generated Samples

As the score, $\nabla\log p_{t}(\mathbf{x})$, is integrable, the Fundamental Theorem of Calculus applies, leading to

F.2 Proof of Theorem 2

Define $\mathcal{T}_{t\rightarrow s}$ as the oracle transition mapping from $t$ to $s$ via the diffusion process Eq. (2). Let $\mathcal{T}_{t\rightarrow s}^{\bm{\theta}^{*}}(\cdot)$ represent the transition mapping from the optimal CTM, and $\mathcal{T}_{t\rightarrow s}^{\bm{\phi}}(\cdot)$ represent the transition mapping from the empirical probability flow ODE. Since all processes start at point $T$ with initial probability distribution $p_{T}$ and $\mathcal{T}_{t\rightarrow s}^{\bm{\theta}^{*}}(\cdot)=\mathcal{T}_{t\rightarrow s}^{\bm{\phi}}(\cdot)$, Theorem 2 in (Chen et al., 2022) and $\mathcal{T}_{T\rightarrow t}\sharp p_{T}=p_{t}$ from Proposition 9 tell us that for $t>s$

Here (a) is obtained from the triangular inequality, (b) and (c) are due to $\mathcal{T}_{\sqrt{1-\gamma^{2}}t\rightarrow t}\mathcal{T}_{T\rightarrow\sqrt{1-\gamma^{2}}t}=\mathcal{T}_{T\rightarrow t}$ and $\mathcal{T}_{T\rightarrow t}\sharp p_{T}=p_{t}$ from Proposition 9, and (d) comes from Eq. (13).

F.3 Proof of Proposition 3

Consider a LPIPS-like metric, denoted as $d(\cdot,\cdot)$, determined by a feature extractor $\mathcal{F}$ of $p_{\text{data}}$. That is, $d(\mathbf{x},\mathbf{y})=\left\lVert\mathcal{F}(\mathbf{x})-\mathcal{F}(\mathbf{y})\right\rVert_{q}$ for $q\geq 1$. For simplicity of notation, we denote $\bm{\theta}_{N}$ as $\bm{\theta}$. Since $\mathcal{L}_{\text{CTM}}^{N}(\bm{\theta};\bm{\phi})=0$, it implies that for any $\mathbf{x}_{t_{n}}$, $n\in[\![1,N]\!]$, and $m\in[\![1,n]\!]$

Then due to Eq. (14) and $G$ is an ODE-trajectory function that $G(\mathbf{x}_{t_{n+1}},t_{n+1},t_{m};\bm{\phi})=G(\mathbf{x}_{t_{n}},t_{n},t_{m};\bm{\phi})$, we have

Notice that since $G_{\bm{\theta}}(\mathbf{x}_{t_{m}},t_{m},t_{m})=\mathbf{x}_{t_{m}}=G(\mathbf{x}_{t_{m}},t_{m},t_{m};\bm{\phi})$, $\mathbf{e}_{m,m}=\mathbf{0}$.

Indeed, an analogue of Proposition 3 holds for time-conditional feature extractors.

Let $d_{t}(\cdot,\cdot)$ be a LPIPS-like metric determined by a time-conditional feature extractor $\mathcal{F}_{t}$. That is, $d_{t}(\mathbf{x},\mathbf{y})=\left\lVert\mathcal{F}_{t}(\mathbf{x})-\mathcal{F}_{t}(\mathbf{y})\right\rVert_{q}$ for $q\geq 1$. We can similarly derive

F.4 Proof of Proposition 5

We first prove that for any $t\in[0,T]$ and $s\leq t$, as $N\rightarrow\infty$,

We may assume $\{t_{n}\}_{n=1}^{N}$ so that $t_{m}=s$, $t_{n}=t$, and $t_{m+1}\rightarrow s$, $t_{n+1}\rightarrow t$ as $\Delta_{N}t\rightarrow\infty$.

In particular, Eq. (15) implies that when $N\rightarrow\infty$

This implies that $p_{\bm{\theta}_{N}}(\cdot)$, the pushforward distribution of $p_{T}$ induced by $G_{\bm{\theta}_{N}}(\cdot,T,0)$, converges in distribution to $p_{\phi}(\cdot)$. Note that since $\{G_{\bm{\theta}_{N}}\}_{N}$ is uniform Lipschitz

$\{G_{\bm{\theta}_{N}}\}_{N}$ is asymptotically uniformly equicontinuous. Moreover, $\{G_{\bm{\theta}_{N}}\}_{N}$ is uniform bounded in $\bm{\theta}_{N}$. Therefore, the converse of Scheffé’s theorem (Boos, 1985; Sweeting, 1986) implies that $\left\lVert p_{\bm{\theta}_{N}}(\cdot)-p_{\bm{\phi}}(\cdot)\right\rVert_{\infty}\rightarrow 0$ as $N\rightarrow\infty$. Similar argument can be adapted to prove $\left\lVert p_{\bm{\theta}_{N}}(\cdot)-p_{\text{data}}(\cdot)\right\rVert_{{\infty}}\rightarrow 0$ as $N\rightarrow\infty$ if the regression target $p_{\bm{\phi}}(\cdot)$ is replaced with $p_{\text{data}}(\cdot)$. $\blacksquare$

F.5 Proof of Proposition 6

Fix a $t\in[0,T]$, the solution operator $\mathcal{T}$ of Eq. (16) with an initial condition $\mathbf{x}_{t}$ is defined as

Here $L:=\sup_{t\in[0,T]}L(t)<\infty$. In particular, if $\mathbf{x}_{t}\neq\hat{\mathbf{x}}_{t}$, $\mathcal{T}[\mathbf{x}_{t}](s)\neq\mathcal{T}[\hat{\mathbf{x}}_{t}](s)$ for all $s\in[t,T]$.

Assumptions (a) and (b) ensure the solution operator in Eq. (17) is well-defined by applying Carathéodory-type global existence theorem (Reid, 1971). We denote $\mathcal{T}[\mathbf{x}_{t}](s)$ as $\mathbf{x}(s;\mathbf{x}_{t})$. We need to prove that for any distinct initial values $\mathbf{x}_{t}$ and $\hat{\mathbf{x}}_{t}$ starting from $t$, $\mathcal{T}[\mathbf{x}_{t}]\not\equiv\mathcal{T}[\hat{\mathbf{x}}_{t}]$. Suppose on the contrary that there is an $s_{0}\in[t,T]$ so that $\mathcal{T}[\mathbf{x}_{t}](s_{0})=\mathcal{T}[\hat{\mathbf{x}}_{t}](s_{0})$. For $s\in[t_{0},s_{0}]$, consider $\mathbf{y}(s;\mathbf{x}_{t}):=\mathbf{x}(t+s_{0}-s;\mathbf{x}_{t})$ and $\mathbf{y}(s;\hat{\mathbf{x}}_{t}):=\mathbf{x}(t_{0}+s_{0}-s;\hat{\mathbf{x}}_{t})$. Then both $\mathbf{y}(s;\mathbf{x}_{t})$ and $\mathbf{y}(s;\hat{\mathbf{x}}_{t})$ satisfy the following ODE

Thus, the uniqueness theorem of solution to Eq. (19) leads to $\mathbf{y}(s_{0};\mathbf{x}_{t})=\mathbf{y}(s_{0};\hat{\mathbf{x}}_{t})$, which means $\mathbf{x}_{t}=\hat{\mathbf{x}}_{t}$. This contradicts to the assumption. Hence, $\mathcal{T}$ is injective.

On the other hand, consider the reverse time ODE of Eq. (16) by setting $\tau=\tau(u):=t+s-u$, $\mathbf{y}(u):=\mathbf{x}(t+s-u)$, and $h(\mathbf{y}(u),u):=-f(\mathbf{y}(u),t+s-u)$, then $\mathbf{y}$ satisfies the following equation

Similarly, we define the solution operator to Eq. (21) as

Here $\mathbf{y}_{t}$ denotes the initial value of Eq. (21) and $\mathbf{y}(u;\mathbf{y}_{t})$ is the solution starting from $\mathbf{y}_{t}$. Due to the Carathéodory-type global existence theorem, the operator $\mathcal{S}[\cdot](s)$ is well-defined and

For simplicity, let $\mathbf{y}_{t}:=\mathbf{x}(s;\mathbf{x}_{t})$ and $\hat{\mathbf{y}}_{t}:=\hat{\mathbf{x}}(s;\mathbf{x}_{t})$. Also, denote the solutions starting from initial values $\mathbf{y}_{t}$ and $\hat{\mathbf{y}}_{t}$ as $\mathbf{y}(u;\mathbf{y}_{t})$ and $\hat{\mathbf{y}}(u;\hat{\mathbf{y}}_{t})$, respectively. Therefore, using a similar argument, we obtain

Proof of Proposition 6.

With the definition of $G(\mathbf{x}_{t},t,s;\bm{\phi})$, we obtain

F.6 Proof of Proposition 7

Let $Y$ be an i.i.d. copy of $X$. Then $h(X)$ and $h(Y)$ are also independent. Thus, $\text{cov}(X,Y)=0$ and $\text{cov}(h(X),h(Y))=0$.

The final equality follows the same reasoning as in Eq. (F.6). Likewise, we can apply the argument from Eq. (24) to show that

Therefore, $L^{-2}\text{Var}\left(X\right)\leq\text{Var}\left(X\right)\leq L^{2}\text{Var}\left(X\right)$. $\blacksquare$

Proof of Proposition 7.

Proposition 6 implies that $G_{\bm{\theta}^{*}}(\cdot,t_{n},\sqrt{1-\gamma^{2}}t_{n+1})$ is bi-Lipschitz and that for any $\mathbf{x},\mathbf{y}$

where $\zeta(t_{n},t_{n+1},\gamma)=\exp{\left(2L_{\bm{\phi}}(t_{n}-\sqrt{1-\gamma^{2}}t_{n+1})\right)}$. Proposition 7 follows immediately from the inequalities (F.6) and (F.6). $\blacksquare$

F.7 Proof of Proposition 9

$\{p_{t}\}_{t=0}^{T}$ is known to satisfy the Fokker-Planck equation (Øksendal, 2003) (under some technical regularity conditions). In addition, we can rewrite the Fokker-Planck equation of $\{p_{t}\}_{t=0}^{T}$ as the following equation (see Eq. (37) in (Song et al., 2020b))

where $\mathbf{W}_{t}:=-t\nabla\log p_{t}$.

Now consider the continuity equation for $\mu_{t}$ defined by $\mathbf{W}_{t}$

Thus, Proposition 8.1.8 of (Ambrosio et al., 2005) implies that for $p_{T}$-a.e. $\mathbf{x}$, the following reverse time ODE (which is the Eq. (2)) admits a unique solution on $[0,T]$

Moreover, $\mu_{t}=X_{t}\sharp p_{T}$, for $t\in[0,T]$. By applying the uniqueness for the continuity equation (Proposition 8.1.7 of (Ambrosio et al., 2005)) and the uniqueness of Eq. (29), we have $p_{t}=\mu_{t}=X_{t}\sharp p_{T}=\mathcal{T}_{T\rightarrow t}\sharp p_{T}$ for $t\in[0,T]$. Again, since the uniqueness theorem with the given $p_{T}$, we obtain $p_{s}=\mathcal{T}_{t\rightarrow s}\sharp p_{t}$ for any $t\in[0,T]$ and $s\in[0,t]$.