Fast samplers for Inverse Problems in Iterative Refinement models

Authors: kushagra  pandey
ruihan  yang
stephan  mandt
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-15 09:40:47

Abstract: Constructing fast samplers for unconditional diffusion and flow-matching models has received much attention recently; however, existing methods for solving *inverse problems*, such as super-resolution, inpainting, or deblurring, still require hundreds to thousands of iterative steps to obtain high-quality results. We propose a plug-and-play framework for constructing efficient samplers for inverse problems, requiring only *pre-trained* diffusion or flow-matching models. We present *Conditional Conjugate Integrators*, which leverage the specific form of the inverse problem to project the respective conditional diffusion/flow dynamics into a more amenable space for sampling. Our method complements popular posterior approximation methods for solving inverse problems using diffusion/flow models. We evaluate the proposed method's performance on various linear image restoration tasks across multiple datasets, employing diffusion and flow-matching models. Notably, on challenging inverse problems like 4x super-resolution on the ImageNet dataset, our method can generate high-quality samples in as few as *5* conditional sampling steps and outperforms competing baselines requiring 20-1000 steps. Our code will be publicly available at https://github.com/mandt-lab/c-pigdm.

Main paper: 
 introduction: 
 Iterative refinement models, such as diffusion generative models and flow matching methods [Sohl-Dickstein et al., 2015 , Ho et al., 2020 , Song et al., 2020 , Lipman et al., 2023 , Albergo and Vanden-Eijnden, 2023] , have seen increasing popularity in recent months, and much effort has been invested in accelerating unconditional sampling in these models [Pandey et al., 2024 , Shaul et al., 2024 , Sauer et al., 2024 , Karras et al., 2022 , Salimans and Ho, 2022 , Zhang and Chen, 2023 , Lu et al., 2022 , Song et al., 2021] . However, while most efficient samplers have been designed in the unconditional setup, current methods for solving inverse problems, such as deblurring, inpainting, or super-resolution, still require hundreds to thousands of neural network evaluations to achieve the highest perceptual quality. Moreover, in addition to a score function evaluation, a class of existing methods for solving inverse problems using pre-trained unconditional iterative refinement models often involves expensive Jacobian-vector products [Song et al., 2022 , Chung et al., 2022a] , making a single sampling step quite expensive and therefore, intolerably slow for most practical applications. Given an initial sampling latent x ts at time t s , our sampler projects the diffusion/flow dynamics to a more amenable space for sampling using a projector operator Φ which is conditioned on the degradation operator H and the sampling guidance scale w. The diffusion/flow sampling is then performed in the projected space. Post completion, the generated sample in the projected space is transformed back into the original space using the inverse of the projection operator, yielding the final generated sample. We define the form of the operator Φ in Section 2.2. Conditional Conjugate Integrators can significantly speed up sampling in challenging inverse problems and can generate high-quality samples in as few as 5 NFEs as compared to existing baselines, which require from 20-1000 NFEs (see Section 3). This paper presents a principled framework for designing efficient samplers for guided sampling in iterative refinement models, accelerating existing samplers like ΠGDM by an order of magnitude. We present our framework for inverse problems where the degradation operator is known and might be corrupted with additional noise. Crucially, our transformations do not require any re-training and merely rely on some algebraic manipulations of the equations to be simulated.Intuitively, we expand on the concept of Conjugate Integrators [Pandey et al., 2024] by projecting the conditional generation process in inverse problems to another space that might be better conditioned for faster sampling (See Figure 1 ). To this end, we separate the linear and non-linear components in the generation process and parameterize the transformation by analytically solving the linear coefficients. By the end of the sampling procedure, we map back to the original sampling space, leading to the concept of Conditional Conjugate Integrators that apply to various iterative refinement models such as diffusion models, flows, and interpolants.In more detail, our main contributions are as follows.• Conditional Conjugate Integrators: We repurpose the recently proposed Conjugate Integrator framework [Pandey et al., 2024] for fast guided sampling in iterative refinement models (diffusions and flows) for linear inverse problems and refer to it as Conditional Conjugate Integrators. Next, we design a specific parameterization of the proposed framework, which encodes the structure of the linear inverse problem in the sampler design itself.• Theoretical Analysis: Our parameterization exhibits theoretical properties that help us identify key parameters for sampler design. More specifically, we show that our parameterization (by design) enables recovering high-frequency details early on during sampling. This further enables fast-guided sampling while maintaining good sample quality in the context of inverse problems.• Empirical Results. Empirically, we show that our proposed sampler significantly improves over baselines in terms of sampling efficiency on challenging benchmarks across inverse problems like super-resolution, inpainting, and Gaussian deblurring. For instance, on a challenging 4x superresolution task on the ImageNet dataset, our proposed sampler achieves better sample quality at 5 steps, compared to 20-1000 steps required by competing baselines.Additionally, we extend the proposed framework for noisy and non-linear inverse problems with qualitative demonstrations.2 Fast Samplers for Inverse Problems using Diffusions/Flows. 
background and problem statement: 
 Diffusion models define a continuous-time forward process (usually with an affine drift) to convert data x 0 ∈ R d into noise. A learnable reverse process is trained to generate data from noise. In this work, we only consider deterministic reverse processes specified as an ODE [Song et al., 2020] ,EQUATION)The score is usually intractable and is approximated using a parametric estimator s θ (x t , t), trained using denoising score matching [Vincent, 2011 , Song and Ermon, 2019 , Song et al., 2020] . Analogously, one-sided stochastic interpolants [Albergo and Vanden-Eijnden, 2023 ] define an interpolant2 x t = α t x 1 + γ t z, where x 1 ∼ p data ,and z ∼ N (0, I) to define a transport map between the generative prior (typically an isotropic Gaussian) and the data distribution. Interestingly, the one-sided interpolant induces a vector field b(x t , t) = E[ αt x 1 + γt z|x t ], where αt , γt represent the time derivatives of α t and γ t , respectively. The vector field b(.) is typically learned using a neural network approximation b θ (x t , t). The deterministic interpolant process can then be specified as dx t = b θ (x t , t) dt. Numerically solving these deterministic generative processes with a sufficient sampling budget can generate plausible samples from noise.Problem Statement. Given a noisy linear degradation process (we will consider non-linear processes later) with a degradation operator H specified over an unobserved data point x 0 ,EQUATIONthe goal is to recover the original signal x 0 . Additionally, given an unconditional pre-trained diffusion or flow matching model, one approach for solving inverse problems is to infer the posterior distribution over the data given the degraded observation, i.e., p(x 0 |y) ∝ p(y|x 0 )p(x 0 ) by simulating the conditional reverse process dynamics i.e.Diffusion:EQUATIONFlows:dx t = b(x t , y, t)dt,where ∇ xt log p(x t |y) and b(x t , y, t) are the conditional score and velocity estimates, respectively. One approach could be to directly model the conditional score or velocity estimates using a conditional iterative refinement model [Saharia et al., 2022a,b] . However, such approaches are problem-dependent, requiring expensive training pipelines to account for the lack of generalization across inverse problems. Additionally, such methods rely on the availability of paired (x t , y) measurements, which can be expensive to acquire. Alternatively, problem-agnostic methods leverage pre-trained unconditional iterative refinement models to estimate the conditional score or velocity fields and can generalize to different inverse problems without extra training. In this work, we restrict our discussion to the latter and discuss estimating conditional score/velocity fields next. We include a formal proof for the result in Eqn. 4 from an interpolant perspective in Appendix A.1. Since the conditional score and velocity estimates require approximating the term ∇ xt log p(y|x t ), we discuss its estimation next.Estimation of the Noise Conditional Score ∇ xt log p(y|x t ): The noise conditional distribution p(y|x t ) can be represented as p(y|x t ) = p(y|x 0 )p(x 0 |x t )dx 0 . For problem-agnostic models, it is common to approximate the posterior p(x 0 |x t ) using an unimodal Gaussian distribution [Chung et al., 2022a , Song et al., 2022] . In this work, we restrict our discussion to the posterior approximation in ΠGDM [Song et al., 2022] and its flow variant [Pokle et al., 2023] (named as ΠGFM in our work), p(x 0 |x t ) ≈ N (x 0 , r 2 t I d ), which yields the following estimate of the conditional score:EQUATION)where x0 is the first-order Tweedie's moment estimate [Stein, 1981] . Our choice of using the ΠGDM approximation is motivated by its expressive posterior approximation p(x 0 |x t ) compared to other methods such as DPS or MCG. This makes it an excellent starting point for low-budget sampling.Conjugate Integrators The main idea in conjugate integrators [Pandey et al., 2024] is to project the diffusion dynamics in Eqn. 1 into another space where sampling might be more efficient. The projected diffusion dynamics can then be solved using any numerical ODE solver. On completion, the dynamics can be projected back to the original space to generate samples from the data distribution. To this end, Pandey et al. [2024] introduce an invertible time-dependent affine transformation xt = A t x t . Interestingly, conjugate samplers have theoretical connections to prior work in fast sampling for unconditional diffusion models [Song et al., 2021 , Zhang and Chen, 2023 , Lu et al., 2022] . We refer the readers to Pandey et al. [2024] for exact details.Next, we design conjugate integrators for linear inverse problems. For simplicity, we discuss noiseless inverse problems, σ y = 0, and defer the discussion of noisy inverse problems to Section 2.4. Furthermore, due to space constraints, we present our analysis for diffusion models and defer the discussion of flows to Appendix B. Lastly, without loss of generality, we assume the standard score network parameterization, s θ (x t , t) = C out (t)ϵ θ (x t , t) where C out (t) is the notation from the score precondition defined in Karras et al. [2022] .A straightforward way to define conditional conjugate integrators is to treat the score estimate ∇ xt log p(y|x t ) as a black-box i.e., ignore the structure of the inverse problem. For this case, we formally specify the conjugate integrator formulation as, Proposition 1. (Extended ΠGDM) For the conditional diffusion dynamics defined in Eqn. 3, introducing a diffeomorphism, xt = A t x t , where,EQUATIONinduces the following projected diffusion dynamics,EQUATIONwhere H † = H ⊤ (HH ⊤ ) -1 and P = H ⊤ (HH ⊤ ) -1 H represent the pseudoinverse and the orthogonal projector operators for the degradation operator H. (Proof in Appendix A.2)Similar to Pandey et al. [2024] , the matrix B t is a design choice. We refer to the formulation in Eqn. 7 as Extended ΠGDM since for B t = 0, the ODE in Eqn. 7 becomes equivalent to the ΠGDM formulation proposed in Song et al. [2022] . This is because, for B t = 0, Conjugate Integrators are equivalent to DDIM [Song et al., 2021] (See Pandey et al. [2024] for proof). Therefore, the projected diffusion dynamics in Eqn. 7 already present a more generic framework for designing samplers for inverse problems over ΠGDM. In this work, we only explore the parameterization in Eqn. 7 for B t = 0 and hence refer to it simply as ΠGDM (analogously ΠGFM for flows; see Appendix B).One characteristic of the formulation in Eqn. 7 is the black-box nature of the conditional score ∇ xt log p(y|x t ). However, the inherent linearity in the conditional score can be used to design better conditioned (more on this in Section 2.3) conjugate integrators, which we illustrate formally in the form of the following result.Proposition 2. (Conjugate ΠGDM) Given a noiseless linear inverse problem with σ y = 0, a design matrix B : [0, 1] → R d×d , and the conditional score ∇ xt log p(y|x t ) approximated using Eqn. 5, introducing the transformation xt = A t x t , whereEQUATIONinduces the following projected diffusion dynamics:EQUATIONwhere exp(.) denotes the matrix exponential, H † , and P are the pseudoinverse and projector operators (as defined previously). Proof in Appendix A.3.In this case, the coefficients Φ y , Φ j , and Φ s depend on time t and the degradation operator H (See Appendix A.3 for full definitions). Intuitively, by including information about the degradation operator H and the guidance scale in the transformation A t in Eqn. 8, we incorporate the specific structure of the inverse problem in the sampler design, which can have several advantages (more on this in Section 2.3). Moreover, the matrix B t is a design choice of our parameterization (we will discuss exact choices in Section 2.2.2). We refer to this parameterization as C-ΠGDM (analogously C-ΠGFM for flows; see Appendix B). In this work, we restrict our discussion to this parameterization and discuss some practical and theoretical aspects next.Choice of Diffusions and Flows: While our proposed integrators are applicable to generic diffusion processes [Dockhorn et al., 2022, Pandey and Mandt, 2023] and flows [Ma et al., 2024] , we restrict follow-up discussion to VP-SDE [Song et al., 2020] diffusion for which [Liu et al., 20223, Lipman et al., 2023] for which α t = t, γ t = 1 -t. For our score network parameterization, we set C out (t) = -1/σ t , corresponding to the standard ϵ-prediction [Ho et al., 2020 , Song et al., 2020] parameterization in diffusion models.F t = -1 2 β t I d , G t = √ β t I d and OT-flowsChoice of B t : Similar to Pandey et al. [2024] , we set B t = λI d , where λ is a time-invariant scalar hyperparameter tuned during inference for optimal sample quality.Choice of w t : Similar to prior work [Song et al., 2022 , Pokle et al., 2023 Having an extra multiplicative factor of µ 2 t (for VP-SDE) or α 2 t (for flows) stabilizes the numerical computation of coefficients in Eqn. 9 before sampling. We tune the static guidance weight w during inference for optimal sample quality.Choice of Start Time: Given a degradation output y, it is common to start diffusion or flow sampling at τ < T or τ > 0, respectively [Chung et al., 2022b , Song et al., 2022 , Pokle et al., 2023] . Consequently, we initialize the diffusion sampling process as x τ = µ τ H † y + σ τ z. Analogously for flows, we initialize sampling atx τ = α τ H † y + γ τ z.Choice of the ODE Solver: Unless specified otherwise, we use the Euler discretization scheme for C-ΠG(D/F)M samplers.We illustrate a generic C-ΠGDM sampling routine in Algorithm 1 and include additional implementation details in Appendix D. Next, we present some theoretical aspects of our proposed method.With the simplifications in Section 2.2.2, the transformation A t in Eqn. 8 simplifies to:EQUATIONwhere P = H ⊤ (HH ⊤ ) -1 H is an orthogonal projection operator. v l = hAt n Bt n A -1 tn x + (Φ n+1 y -Φ n y )y ▷ Linear drift 13: v nl = (Φ n+1 s -Φ n s )ϵ θ (x, tn) + (Φ n+1 j -Φ n j ) ∂xϵ θ (x, tn)(H † y -P x0)▷ Non-Linear drift 14:x = x + v l + v nl ▷ Euler Update 15: end for return x = A -1 t N x ▷ Project back to original space when doneComputing A t : While computing the matrix exponential in Eqn. 10 might seem non-trivial, it has several interesting properties that make it tractable to compute. More specifically, the matrix exponential in Eqn. 10 can be simplified as (Proof in Appendix A.4),EQUATIONwhere exp(.) in Eqn. 11 represents the scalar exponential. Furthermore, the integrals in Eqn. 11 are trivial to compute analytically or numerically, making A t easier to compute. Moreover, A -1 t can also be compactly represented as,EQUATIONand is also tractable to compute. Due to the tractability of A t and A -1 t , the projected diffusion dynamics in C-ΠGDM are straightforward to simulate numerically.Intuition behind A t : Next, we analyze several theoretical properties of the transformation matrix A t in Eqn. 11. More specifically,EQUATIONSince P = H ⊤ (HH ⊤ ) -1 H is an orthogonal projector, the matrix I d -P is also an orthogonal projector which projects any vector v in the nullspace of P . Therefore, we can decompose the state x t into two orthogonal components x t = P x t + (I d -P )x t . Plugging this form in Eqn. 13,EQUATIONIntuitively, near t = T (i.e., at the start of reverse diffusion sampling), for a large static guidance weight w, exp(κ 2 (t)) → 0. In this limit, from eqn. 14, xt ≈ (I d -P )x t . This implies that for a large guidance weight w, the diffusion dynamics are projected into the nullspace of the projection operator P . Intuitively, for an inverse problem like superresolution, this implies that near the start of the diffusion process, the projected diffusion dynamics correspond to the denoising of the highfrequency details missing in P x t . This is because the projector operation, P x t = H † Hx t can be interpreted as the pseudoinverse of the noisy degraded state x t , and, therefore, (I d -P )x t represents the high-frequency details missing from the signal component in P x t .Moreover, near t = 0 (i.e., near the end of reverse diffusion sampling), assuming the guidance weight w is not too large, both coefficients exp(κ 1 (t)) and exp(κ 2 (t)) → 1, which implies xt ≈ x t . This implies that near t = 0, diffusion happens in the original space, which can prevent over-sharpening artifacts towards the end of sampling. Therefore, we hypothesize that a large w can also lead to over-sharpened results near the end of sampling, resulting in artifacts in the generated samples. Therefore, introducing the projection A t as defined in Eqn. 10, introduces a tradeoff in the choice of w to control for sample quality. Lastly, since the parameter λ controls the magnitude of xt , it exhibits a similar tradeoff. Indeed, we will empirically demonstrate these tradeoffs in Section 3.3. While our discussion has been limited to diffusion models, a similar theoretical intuition also holds for flows (See Appendix B for proof).While our discussion has been primarily in the context of noiseless linear inverse problems, the conditional Conjugate Integrator framework can also be extended to develop samplers for noisy linear and non-linear inverse problems. We provide a more detailed explanation for the same in App. C. 
experiments: 
 Next, we empirically demonstrate that our proposed samplers C-ΠGDM/GFM outperform recent baselines on linear image restoration tasks regarding sampling speed vs. quality tradeoff. We then present ablation experiments highlighting the key parameters of our samplers. Lastly, we present design choices for solving noisy and non-linear inverse problems using our proposed framework.Models and Dataset: For diffusion models, we utilize an unconditional pre-trained ImageNet [Deng et al., 2009] checkpoint at 256×256 resolution from OpenAI [Dhariwal and Nichol, 2021] 3 . For evaluations on the FFHQ dataset Karras et al. [2019] , we use a pre-trained checkpoint from Choi et al. [2021] also at 256×256 resolution. For flow model comparisons, we utilize three publicly available model checkpoints from Liu et al. [20223] 4 , trained on the AFHQ-Cat [Choi et al., 2020] , LSUN-Bedroom Yu et al. [2015] , and CelebA-HQ [Karras et al., 2018] datasets. Each flow model was trained at a pixel resolution of 256 × 256. For diffusion models, we conduct evaluations on a 1k subset of the evaluation set. For flows, we conduct evaluations on the entire validation set.Tasks and Metrics: We evaluate our samplers qualitatively (see Figure 2 ) and quantitatively on three challenging linear inverse problems under the noiseless setting. Firstly, we test Image Super-Resolution, enhancing images from bicubic-downsampled 64 × 64 pixels to 256 × 256 pixels. Secondly, we assess Image Inpainting performance on images with a fixed free-form center mask. Lastly, we evaluate our samplers on Gaussian Deblurring, applying a Gaussian kernel with σ = 3.0 across a 61 × 61 window. We evaluate the performance of each task based on three perceptual metrics: FID [Heusel et al., 2017] , KID [Bińkowski et al., 2018] and LPIPS [Zhang et al., 2018] .We assess the sample quality of our proposed C-ΠGDM and C-ΠGFM samplers using 5, 10, and 20 sampling steps (denoted as Number of Function Evaluations (NFE)).We conduct an extensive search to optimize the parameters w, λ and τ to identify the best-performing configuration based on sample quality. For diffusion baselines, we include DDRM [Kawar et al., 2022] , DPS [Chung et al., 2022a] , and ΠGDM [Song et al., 2022] . As recommended for DPS [Chung et al., 2022a] , we use NFE=1000 for all tasks. For DDRM, we adhere to the original implementation and run it with η b = 1.0 and η = 0.85 at NFE=20. We test our implementation of ΠGDM (see Section 2.2), with NFE values of 5, 10, and 20 and use the recommended guidance schedule of w t = r 2 t across all tasks. For flow models, we consider the recently proposed method inspired by ΠGDM running on OT-ODE path by Pokle et al. [2023] (which we refer to as ΠGFM; see Appendix B), and similarly run it with NFE values of 5, 10, and 20. We optimize all baselines by conducting an extensive grid search over w and τ for the best performance (in terms of sample quality).We present the results of our method applied to inverse problems in Table 1 , specifically using the CelebA-HQ dataset for flow-based models and the ImageNet dataset for diffusion-based models. For a comprehensive review of additional results across different datasets, please refer to Appendix E. Our method consistently surpasses other approaches across all sampling budgets (indicated by NFE) for the inpainting task. Similarly, our flow-based sampler (C-ΠGFM) exhibits superior perceptual quality for image super-resolution at NFEs of 5 and 10. The ΠGFM model only reaches comparable performance at higher NFEs. Remarkably, our diffusion-based sampler C-ΠGDM outperforms all baselines across the entire range of NFEs. Notably, C-ΠGDM outperforms competing baselines requiring 20-1000 NFEs in just 5 sampling steps on the challenging ImageNet dataset, demonstrating a significant speedup in sampling speed while preserving sample quality. A similar pattern is observed in the image deblurring task, where the performance of ΠGDM/ΠGFM approaches that of our method only when the NFE is increased to 20 steps.LPIPS↓ KID×10 -3 ↓ FID↓ C-ΠGFM ΠGFM C-ΠGFM ΠGFM C-Interestingly, we observe a plateau in performance improvements at NFE=20 for both super-resolution and deblurring tasks using our method. This suggests that while our method efficiently utilizes the iterative model under a deterministic path with an Euler solver, further enhancements in performance, particularly at higher NFEs, might require integrating stochastic sampling techniques or more advanced solvers. This potential next step could unlock further gains from our approach in complex image processing tasks.Figure 2 presents a qualitative comparison between our proposed method and the ΠG(D/F)M baseline.The inpainting results in the first column reveal that ΠGFM tends to introduce gray artifacts within the inpainted areas. This issue may stem from the initialization of the parameter τ ; optimal performance is achieved when τ ≥ 0.2, as established during our parameter tuning phase and corroborated by Pokle et al. [2023] . Consequently, insufficient NFE means ΠGFM cannot effectively eliminate the artifacts associated with the inpainting mask in our experiments. For image super-resolution, our method excels in restoring fine details, particularly evident in high-frequency image components such as human hair and wheat ears. Similarly, for the deblurring task, our method qualitatively outperforms the baseline, especially in mitigating the over-smoothing artifacts (Figure 2 , last column). Additional examples are provided in Appendix E.4.In this section, we further explore the impact of the hyperparameters w, λ, and τ , which were identified during our tuning phase and link to the theoretical insights discussed in Section 2.3. We recognize that τ is particularly task-specific and relatively straightforward to adjust. For instance, tasks such as inpainting require a smaller τ to prevent masking artifacts, whereas tasks like super- From Figure 3 we make the following observations. Firstly, for both C-ΠGDM and C-ΠGFM samplers, we observe that the optimal value of λ can differ from λ = 0. This illustrates the usefulness of parameterizing B t in our sampler design. On the contrary, ΠGDM or ΠGFM samplers do not have this flexibility and, therefore, yield sub-optimal sample quality at different sampling budgets. Secondly, we observe that deviating from the optimal λ can lead to degradation in sample quality. More specifically, we observed that deviating from our tuned value of λ leads to either over-sharpening artifacts or blurry samples (See Figs. 7, 11) . This is intuitive since λ controls the scale of the transformation xt = A t x t (see Eqn. 14) and thus plays a significant role in conditioning the projected diffusion dynamics. We observe a similar tradeoff on varying the static guidance weight w where a large magnitude of w can lead to over-sharpened artifacts while a very small guidance weight can lead to blurry samples (See Figs. 6, 10) . These empirical observations are consistent with our theoretical analysis in Section 2.3, confirming our theoretical intuition on the role of the sampler parameters w and λ.Question: For each theoretical result, does the paper provide the full set of assumptions and a complete (and correct) proof?Answer: [Yes] Justification: See Section 2 and Appendix A Guidelines:• The answer NA means that the paper does not include theoretical results.• All the theorems, formulas, and proofs in the paper should be numbered and crossreferenced. • All assumptions should be clearly stated or referenced in the statement of any theorems.• The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide a short proof sketch to provide intuition. • Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material. • Theorems and Lemmas that the proof relies upon should be properly referenced. 
related works: 
 Fast Unconditional Sampling: Recent research has significantly advanced the efficiency of the sampling process in unconditional diffusion/flow models [Song et al., 2020 , Lipman et al., 2023 , Manduchi et al., 2024] . One line of research involves designing efficient diffusion models to improve sampling by design [Karras et al., 2022 , Dockhorn et al., 2022 , Pandey and Mandt, 2023 , Song et al., 2023] . Since our treatment of conditional Conjugate Integrators is quite generic, our method is readily compatible with most advancements in diffusion model design. Another line of work focuses on distilling a student model from a teacher model, enabling sampling in even a single step [Salimans and Ho, 2022 , Meng et al., 2023 , Sauer et al., 2024] . However, since these methods require expensive re-training, there has been a significant interest in the development of fast samplers applicable to pretrained diffusion/flow models [Liu et al., 2022 , Pandey et al., 2024 , Shaul et al., 2024 , Zhang and Chen, 2023 , Lu et al., 2022 , Song et al., 2021 , Gonzalez et al., 2023] . Our work falls under the latter line of research, where we develop fast conditional samplers that can be applied to pretrained diffusion models.Conditional Iterative Refinement Models have become prevalent for tasks requiring controlled generation. These models often involve training specialized conditional diffusion models [Saharia et al., 2022b , Yang and Mandt, 2023 , Kong et al., 2021 , Pandey et al., 2022 , Preechakul et al., 2022 , Rombach et al., 2022 , Podell et al., 2023 , Ramesh et al., 2022 , Peebles and Xie, 2023 , Ma et al., 2024 , Esser et al., 2024 , Chen et al., 2023] and may incorporate classifier-free guidance [Ho and Salimans, 2021] or classifier guidance [Dhariwal and Nichol, 2021, Song et al., 2020] for conditional sampling. These approaches have also spurred research into solving inverse problems related to various image degradation transformations, such as inpainting and super-resolution [Kawar et al., 2022 , Chung et al., 2022a , Song et al., 2022 , Mardani et al., 2023 , Pokle et al., 2023] . Although these methods demonstrate promising outcomes, they are typically bottlenecked by a costly sampling process, emphasizing the need for a fast sampler to address inverse problems efficiently. Recent work Xu et al. [2024] employs a consistency model Song et al. [2023] to enhance posterior approximation, but incorporating an additional model may deviate from our proposal of using a single pre-trained model. DPM-Solver++ [Lu et al., 2023] also tackles the problem of accelerating guided sampling in diffusion models. However, unlike [Lu et al., 2023] , we incorporate the structure of the inverse problem in the sampler design.Question: Does the paper fully disclose all the information needed to reproduce the main experimental results of the paper to the extent that it affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?Answer: [Yes] Justification: See the beginning of Section 3 Guidelines:• The answer NA means that the paper does not include experiments.• If the paper includes experiments, a No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not. • If the contribution is a dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable. • Depending on the contribution, reproducibility can be accomplished in various ways.For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution is a specific model and empirical evaluation, it may be necessary to either make it possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g., in the case of a large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed. (d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome to describe the particular way they provide for reproducibility.In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers to have some path to reproducing or verifying the results. 
discussion: 
 We present a generic framework for designing samplers for accelerating guided sampling in iterative refinement models. In this work, we explore a specific parameterization of this framework, which incorporates the structure of the inverse problem in sampler design. We provide a theoretical intuition behind our design choices and empirically justify its effectiveness in solving linear inverse problems in as few as 5 sampling steps compared to 20-1000 NFEs required by competing baselines. While our method can serve as an important step toward designing fast-guided samplers, there are several important future directions. Firstly, our parameterization of the transform A t can be more expressive by learning it directly during the sampling stage. Secondly, in this work, we consider inverse problems with a known degradation operator. Extending our framework for solving blind inverse problems could be an important research direction. Lastly, it would be interesting to adapt our solvers to techniques for solving inverse problems in latent diffusion models [Rout et al., 2024] to enhance sampling efficiency further.Broader Impact: While our work has the potential to make synthetic data generation accessible, the techniques presented in this work should be used responsibly. Moreover, despite good sample quality in a limited sampling budget, restoration can sometimes lead to artifacts in the generated sample which can be undesirable in some domains like medical image analysis. Proof. Given a one-sided interpolant, x t = α t x 1 + γ t z, x 1 ∼ p data , z ∼ N (0, I) satisfying regularity conditions as stated in [Albergo et al., 2023] , and a degraded signal y generated using Eqn. 2, the conditional velocity field b ( x t , y, t) can be approximately estimated from the unconditional velocity field b θ (x t , t) as [Pokle et al., 2023] ,EQUATIONwhere w t represents a time-dependent scalar guidance schedule and σt , αt represent the first-order time derivatives of σ t and α t , respectively. Our proof consists of two parts: Firstly, we establish the connection between the unconditional velocity field b(x t , t) for the one-sided interpolant and the score function s(x t , t) associated with the marginal distribution p(x t ). Secondly, we use this connection to estimate the conditional velocity b(x t , y, t) in terms of b(x t , t) to establish the required result.Connection between b(x t , t) and s(x t , t): For the one-sided interpolant, by definition,EQUATIONTaking the expectation w.r.t p(x 1 , z) on both sides conditioned on the noisy state x t , we have,EQUATIONFurthermore, we have the following result from Albergo et al. [2023] ,EQUATIONwhere s(x t , t) represents the score function. From Eqns. 17, 18, it follows,EQUATIONIntuitively, the above result represents Tweedie's estimate [Stein, 1981] for estimatingx1 = E[x 1 |x t ]in the context of one-sided stochastic interpolants. Next, the one-sided interpolant also induces an unconditional velocity field specified as:EQUATIONwhere γt , αt represent the first-order time derivatives of γ t and α t , respectively. Substituting the result from Eqn. 19 into Eqn. 20, we have the following result,EQUATIONThis concludes the first part of the proof.Estimating the conditional velocity b(x t , y, t) in terms of b(x t , t). For transport conditioned on y, the conditional velocity can be expressed as (following the result in Eqn. 21):EQUATIONApproximating the unconditional velocity b(x t , t) using a parametric estimator b θ (x t , t), we get the required result.EQUATIONA.2 Proof of Proposition 1We restate Proposition 1 for convenience, Proposition. For the conditional diffusion dynamics defined in Eqn. 3, introducing the transformation xt = A t x t induces the following projected diffusion dynamics.EQUATIONEQUATIONwhere H † = H ⊤ (HH ⊤ ) -1 and P = H ⊤ (HH ⊤ ) -1 H represent the pseudoinverse and the orthogonal projector operators for the degradation operator H.Proof. We have the following form of the conditional diffusion dynamicsEQUATIONGiven an affine transformation which projects the stateEQUATIONTherefore, by the Chain Rule of calculus,EQUATIONSubstituting the ODE in Eqn. 30 in Eqn. 32EQUATIONEQUATIONSince, we have s θ (x t , t) = C out (t)ϵ θ (x t , t), the above equation can be simplified as,EQUATIONFurther parameterizing,EQUATIONwhich yields the required diffusion ODE in the projected space:EQUATIONWe have the following approximation for the conditional score ∇ xt log p(y|x t ) (for the noiseless case σ y = 0)EQUATIONwhere H † = H ⊤ (HH ⊤ ) -1 and P = H ⊤ (HH ⊤ ) -1 H represent the pseudoinverse and the orthogonal projector operators for the degradation operator H. Substituting this form of the conditional score in projected diffusion dynamics, we have,EQUATIONThis concludes the proof.We restate the theorem here for convenience.Proposition. For a noiseless linear inverse problem with σ y = 0 and the conditional score approximated using Eqn. 5, introducing the transformation xt = A t x t also induces the following projected diffusion dynamics.EQUATIONA t = exp t 0 B s -F s + w s r -2 s 2µ 2 s G s G ⊤ s P ds dΦ y = - w t r -2 t 2µ t A t G t G ⊤ t H † (44) dΦ s = - 1 2 A t G t G ⊤ t I d - w t r -2 t σ 2 t µ 2 t A t P C out (t) dΦ j = - w t r -2 t σ 2 t 2µ t A t G t G ⊤ t C out (t) (45)where exp(.) denotes the matrix exponential, H † , and P are the pseudoinverse and projector operators (as defined previously).Proof. The proof consists of two parts. Firstly, we simplify the conditional score ∇ xt log p(y|x t ).Secondly, we plug the simplified form of the conditional score into the conditional diffusion dynamics and develop conjugate integrators.Exploiting the linearity in ∇ xt log p(y|x t ): From the definition of the score ∇ xt log p(y|x t ):EQUATION)where x0 is the Tweedie's estimate of E(x 1 |x t ) given by:EQUATIONwhere µ t , σ t are the mean coefficient and standard deviation of the perturbation kernel p(x t |x 1 ) = N (µ t x 1 σ 2 t I d ), respectively, and Σ t = r 2 t (HH ⊤ ) is the variance of the ΠGDM approximation of p(y|x t ) (for the noiseless case i.e. σ y = 0). Therefore,∇ xt log p(y|x t ) = ∂ x0 ∂x t ⊤ H ⊤ Σ -1 t (y -H x0 ) (48) = 1 µ t (I d + σ 2 t ∇ xt s θ (x t , t) =S θ (xt,t) )H ⊤ Σ -1 t (y -H x0 ) (49) = 1 µ t H ⊤ Σ -1 t (y -H x0 ) + σ 2 t S θ (x t , t)H ⊤ Σ -1 t (y -H x0 ) (50) = 1 µ t H ⊤ Σ -1 t (y - 1 µ t H(x t + σ 2 t s θ (x t , t))) + σ 2 t S θ (x t , t)H ⊤ Σ -1 t (y -H x0 ) (51) = 1 µ t H ⊤ Σ -1 t y - 1 µ 2 t H ⊤ Σ -1 t Hx t Linear Terms (52) - σ 2 t µ 2 t H ⊤ Σ -1 t Hs θ (x t , t)) + σ 2 t µ t S θ (x t , t)H ⊤ Σ -1 t (y -H x0 )Non-Linear Terms (53) where S θ denotes the second-order derivative of the score function s θ (x t , t). Therefore, the conditional score ∇ xt log p(y|x t ), can be decomposed into a combination of linear and non-linear terms.Next, we use this decomposition to design conjugate integrators for noiseless linear inverse problems.Conjugate Integrator Design: From Eqn. 30, the conditional reverse diffusion dynamics can be specified as:EQUATIONPlugging in the form of the conditional score in Eqn. 53 in the above equation, we have,EQUATIONGiven an affine transformation which projects the state x t to xt ,EQUATIONthe projected diffusion dynamics can be specified as:dx t dt = dA t dt + A t (F t + w t 2µ 2 t G t G ⊤ t H ⊤ Σ -1 t H) x t - w t 2µ t A t G t G ⊤ t H ⊤ Σ -1 t y (61) - 1 2 A t G t G ⊤ t I d - w t σ 2 t µ 2 t A t H ⊤ Σ -1 t H s θ (x t , t) - w t σ 2 t 2µ t A t G t G ⊤ t S θ (x t , t)H ⊤ Σ -1 t (y -H x1 ) (62)Furthermore, the score network is parameterized as s θ (x t , t) = C out (t)ϵ θ (x t , t). Consequently, S θ (x t , t) = C out (t)∂ t ϵ θ (x t , t). Lastly, we reparameterize Σ t = r 2 t Σ where Σ = HH ⊤ . Plugging these parameterizations in the projected diffusion dynamics, we have,EQUATIONWe then parameterize,EQUATIONThis implies,EQUATIONwhere exp(.) denotes the matrix exponential. Furthermore, we parameterize,EQUATIONWith this parameterization, the projected diffusion dynamics can be compactly specified as follows:EQUATIONThis concludes the proof.We restate the result for convenience. The matrix exponential in Eqn. 10A t = exp t 0 λ + 1 2 β s ds I d - w 2 t 0 β s ds P (72)can be simplified as,EQUATION)where P is the orthogonal projector corresponding to the degradation operator H.Proof. We have,A t = exp(κ 1 t I d + κ 2 t P ) = exp(κ 1 t I d ) exp(κ 2 t P ) (74)The above result follows since I d P = P I d (commutative under multiplication). Moreover, we can further simplify the matrix exponential in P as follows,exp(κ 2 t P ) = ∞ i=0 (κ 2 t ) i P i i! (75) = I d + ∞ i=1 (κ 2 t ) i P i i! = I d + ∞ i=1 (κ 2 t ) i i! P (76)The above result follows from the property of orthogonal projectors P 2 = P . Therefore,exp(κ 2 t P ) = I d + ∞ i=1 (κ 2 t ) i i! P = I d + ∞ i=0 (κ 2 t ) i i! -1 P (77) exp(κ 2 t P ) = I d + exp(κ 2 t ) -1 P (78)which concludes the proof.We restate the result here for convenience. Proposition. For the noisy inverse problem,EQUATIONgiven the transformation A t for the noiseless case as defined in Eqn. 11, the corresponding noisy transformation can be approximated as,A σy t = A t + κ 3 (t)H † (H † ) ⊤ + O(σ 4 y ) (80) κ 3 (t) = wσ 2 y 2 t 0 β s r 2 s ds exp κ 1 (t) + κ 2 (t) -1 (81)Consequently, the inverse of the transformation A σy t can be approximated as, (AEQUATIONProof. We have,A σy t = exp t 0 λ + 1 2 β s ds - w 2 t 0 β s H ⊤ (HH ⊤ + σ 2 y r 2 t I d ) -1 Hds (83)From perturbation analysis, we introduce the following first-order approximation,EQUATIONEQUATIONEQUATIONSubstituting this approximation in the expression for A σy t (and ignoring terms in O(σ 4 y )),EQUATIONFrom the definition of the matrix exponential, it can be shown that the A σy t in Eqn. 89 can be approximated as:EQUATIONIgnoring the higher-order terms, we have,EQUATIONConsequently, we can also approximate the inverse of A σy t , as follows, (AEQUATIONwhich concludes the proof.This section discusses conditional conjugate integrators in the context of flows. For brevity, we skip deriving our results for flows since the derivations can be similar to the analysis of diffusion models with minor parameterization changes. Recall that the conditional dynamics for flows are specified as follows [Pokle et al., 2023] :EQUATIONwhere b θ (x t , t) represents the pre-trained velocity field for a flow. Moreover, we restate the form of the conditional score ∇ xt log p(y|x t ) for convenience.EQUATION)where x1 represents the Tweedie's estimate of the first moment of E(x t |x 1 ),EQUATIONwhere s(x t , t) represents the score function associated with the marginal distribution p(x t ). It can be shown that x1 can also be expressed in terms of the pre-trained velocity field b θ (x t , t) as follows,EQUATIONAnalogous to diffusion models, we can design conditional conjugate samplers for flows that treat the conditional score ∇ xt log p(y|x t ) as a black box. Similar to Proposition 1, by introducing the transformation xt = A t x t , we have the projected flow dynamics,EQUATIONwhere H † = H ⊤ (HH ⊤ ) -1 and P = H ⊤ (HH ⊤ ) -1 H represent the pseudoinverse and the orthogonal projector operators for the degradation operator H. For B t = 0, the formulation in Eqn. 98 becomes equivalent to the ΠGDM formulation proposed for OT-flows in Pokle et al. [2023] . For simplicity, since in this work, we only explore the parameterization in Eqn. 98 for B t = 0, we refer to this parameterization as ΠGFM.Analogous to the discussion of C-ΠGDM samplers in Section 2.2. More specifically, given a noiseless linear inverse problem with σ y = 0, and the conditional score ∇ xt log p(y|x t ), introducing the transformation xt = A t x t , whereA t = exp t 0 B s + w s r -2 s γ t γ2 t 2α t γ t αt -γt α t P ds (100)induces the following projected flow dynamics.dx t = A t B t A -1 t xt dt + dΦ y y + dΦ b b θ (x t , t) + dΦ j ∂ xt b θ (x t , t)(H † y -P x1 ) (101)where, Φ y = - t 0 w s r -2 t γ s γs α s A s H † ds (102)EQUATIONwhere exp(.) denotes the matrix exponential, H † , and P are the pseudoinverse and projector operators (as defined previously). Lastly, the matrix B t is a design choice of our method. We specify a recipe for C-ΠGFM sampling in Algorithm 2.Here, we discuss an extension of Conditional Conjugate Integrators to noisy and non-linear inverse problems. While our discussion is primarily in the context of diffusion models, similar theoretical arguments also apply to Flows.Noisy Linear Inverse Problems: For noisy linear inverse problems of the form,EQUATION)for VPSDE diffusion, the noisy transformation A σy t can be approximated from the transformation A t for the noiseless case (i.e., σ y = 0) as illustrated in the following result (Proof in Appendix A.5):EQUATION)EQUATION)Consequently, the inverse projection (A σy t ) -1 can be approximated from A σy t from perturbation analysis.EQUATIONTherefore, the transformation matrix A σy t and its inverse (see Appendix A.5) can also be computed tractably for the noisy case. Since, for most practical purposes, σ y is pretty small, higher order terms in σ 4 y can be safely ignored, making our approximation accurate. We include qualitative examples for 4x super-resolution with σ y = 0.05 for the ImageNet dataset in Figure 8 Non-Linear Inverse Problems: For non-linear inverse problems of the form,EQUATIONsimilar to Song et al. [2022] , we heuristically re-define linear operations like H † x t , Hx t and P x t by their non-linear equivalents h † (x t ), h(x t ) and h † (h(x t )) respectively. Consequently, analogous to Eqn. 11 the projection operator for a noiseless non-linear inverse problem, in this case, can be defined as,EQUATION) where P = h † (h(.)) is non-linear 'projector" operator. For instance, in non-linear inverse problems like compression artifact removal, h(x t ) and h † (x t ) can realized by encoders and decoders. We illustrate some qualitative examples in Figure 12 . It is worth noting that this is a purely heuristic approximation, and developing a more principled framework for non-linear inverse problems within our framework remains an interesting direction for further work.In this section, we include additional practical implementation details for both C-ΠGDM and C-ΠGFM formulations.We work with the VP-SDE diffusion [Song et al., 2020] with the forward process specified as:EQUATIONThis implies, F t = -1 2 β t and G t = √ β t . For the VP-SDE the perturbation kernel is given by,EQUATION)µ t = exp - 1 2 s 0 β s ds σ 2 t = 1 -exp - s 0 β s ds (113)The corresponding deterministic reverse process is parameterized as:EQUATIONMoreover, we adopt the standard ϵ-prediction parameterization which implies C out (t) = -1/σ t . Lastly, the Tweedies estimate x0 can be specified as:EQUATIONWe choose the parameterization B t = λI d and set the adaptive guidance weight asw t = wµ 2 t r 2 t , where r 2 t = σ 2 t σ 2 t +µ 2 t .The projected diffusion dynamics are then specified as:EQUATIONwhereA t = exp t 0 λ + 1 2 β s dsI d - w 2 t 0 β s ds P (117)which further simplifies to,EQUATIONMoreover, we have,EQUATIONEQUATIONΦ j = - t 0 w s r -2 s σ 2 s 2µ s A s G s G ⊤ s C out (s)ds = t 0 wβ s µ s σ s 2 A s ds (130) = t 0 wβ s µ s σ s 2 exp(κ 1 (s)) I d + (exp(κ 2 (s)) -1)P ds (131) = t 0 wβ s µ s σ s 2 exp(κ 1 (s))ds + t 0 wβ s µ s σ s 2 exp(κ 1 (s))(exp(κ 2 (s)) -1)ds P (132) D.2 C-ΠGFM: Practical Aspects D.2.1 OT-FlowsWe work with OT-Flows [Albergo et al., 2023 , Lipman et al., 2023, Liu et al., 20223] due to its wide adoption. More specifically, the corresponding interpolant can be specified as,EQUATIONFor this case α t = t and γ t = 1 -t. Therefore, the Tweedie's estimate of E(x t |x 1 ) can be specified as (from Eqn. 97):EQUATIONWe choose the parameterization B t = λI d and set the adaptive guidance weight asw t = wα 2 t r 2 t , where r 2 t = γ 2 t α 2 t +γ 2 t .The projected diffusion dynamics are then specified as follows:EQUATIONwhere,EQUATIONwhich further simplifies to,A t = exp(κ 1 (t)) I d + (exp(κ 2 (t)) -1)P , κ 1 (t) = t 0 λds, κ 2 (t) = w 2 t 0 s(1 -s)ds(137) Moreover, we have,EQUATIONΦ j = t 0 w s r -2 s γ 2 s α s A s ds = t 0 wα s γ 2 s A s ds (148) = t 0 ws(1 -s) 2 exp(κ 1 (s)) I d + (exp(κ 2 (s)) -1)P ds (149) = t 0 ws(1 -s) 2 exp(κ 1 (s))ds + t 0 ws(1 -s) 2 exp(κ 1 (s))(exp(κ 2 (s)) -1)ds P (150)From the above analysis, most integrals are one-dimensional and can be computed in closed form or numerically with high precision. To clarify, with a predetermined timestep schedule {t i }, the coefficients Φ can be calculated offline just once and then reused across various samples. Therefore, this computation must only be done once offline for each sampling run. For numerical approximation of these integrals, we use the odeint method from the torchdiffeq package [Chen, 2018] with parameters atol=1e-5, rtol=1e-5 and the RK45 solver [Dormand and Prince, 1980] . We set the initial value Φ init = 0 for all coefficients Φ as an initial condition for both C-ΠGDM and C-ΠGFM samplers.We use the Euler method to simulate projected diffusion/flow dynamics for simplicity. However, using higher-order numerical solvers within our framework is also possible. We leave this exploration to future work.: We use uniform spacing for timestep discretization during sampling. We hypothesize our sampler can also benefit from more advanced timestep discretization techniques Karras et al. [2022] commonly used for sampling in unconditional diffusion models in the low NFE regime.It is common to add an Euler-based denoising step from a cutoff ϵ to zero to optimize for sample quality [Song et al., 2020 , Dockhorn et al., 2022 , Jolicoeur-Martineau et al., 2021] at the expense of another sampling step. In this work, we do not use last-step denoising for our samplers.We use the network function evaluations (NFE) to assess sampling efficiency and perceptual metrics KID [Bińkowski et al., 2018] , LPIPS [Zhang et al., 2018] and FID [Heusel et al., 2017] to assess sample quality. In practice, we use the torch-fidelity [Obukhov et al., 2020] package for computing all FID and KID scores reported in this work. For LPIPS, we use the torchmetrics package with Alexnet embedding.Diffusion Baselines: For DPS [Chung et al., 2022a] , we set NFE=1000 and set the step size for each task to the value recommended in Appendix D in Chung et al. [2022a] . For DDRM [Kawar et al., 2022] , we set the number of sampling steps to NFE=20 with parameters η b = 1.0 and η = 0.85 as recommended in Kawar et al. [2022] . For both DPS and DDRM we start diffusion sampling from t = T . For our implementation of Π-GDM, we set the start time parameter τ to 0.6 for super-resolution and deblurring. We set the guidance weight w t = wr 2 t where w is tuned using grid search between 1.0 and 10.0 for best sample quality for super-resolution and deblurring. For implementation of all diffusion-based baselines, we use the official code for RED-Diff [Mardani et al., 2023] at https://github.com/NVlabs/RED-diff.Flow Baselines: In developing our flow-based baseline, we adhere to the approach outlined in ΠGFM (Pokle et al., 2024) , which advocates for a consistent guidance schedule characterized by w t = w and r t =γ 2 t γ 2 t +α 2 t .For each task, we perform a comprehensive grid search over the parameters α τ = {0.1, 0.2, . . . , 0.7} and w = {1, 2, . . . , 5} (35 combinations in total) across different datasets to identify the optimal configuration that minimizes the LPIPS score. For the implementation of Flows, we use the official implementation of Rectified Flows [Liu et al., 20223] at https://github.com/gnobitab/RectifiedFlow.We include additional comparisons between our proposed samplers and competing baselines on the AFHQ-Cat (see Table 2 ), LSUN Bedroom (see Table 3 ), and the FFHQ (see Table 4 ) datasets.A note on Inpainting evaluations for ImageNet. We find that for diffusion model evaluations, the continuous sampler for ΠGDM suffers from noisy artifacts for the inpainting task. Consequently, Conjugate ΠGDM suffers from similar artifacts. Therefore, we do not report results on this task for the ImageNet dataset. Table 3 : Quantitative evaluation on 4x superresolution, inpainting, and Gaussian deblurring on the LSUN-Bedroom dataset. We note that ΠGFM fails to generate reasonable texture in the masked region even with the maximum NFE=20, so we choose not to report the results here. (See qualitative examples in Figure 9 )Here, we highlight the robustness of C-ΠGDM in both perceptual and recovery metrics in the context of inverse problems. For completeness, we provide a comparison between DPS, ΠGDM, and C-ΠGDM in terms of PSNR, SSIM, FID, and LPIPS in Tables 5 and 6 on the ImageNet-256 and FFHQ-256 datasets on the 4x super-resolution task. It is worth noting that the PSNR and SSIM scores for all methods correspond with the best FID/LPIPS scores presented in the main text for these methods. Our method achieves competitive PSNR and SSIM scores for better perceptual quality than competing baselines like DPS/Π-GDM, even for very small sampling budgets. For instance, on the FFHQ dataset, our method achieves a PSNR of 28.97 compared to 28.49 for DPS while achieving better perceptual sample quality (LPIPS: 0.095 for ours vs 0.107 for DPS) and requiring around 200 times less sampling budget (NFE=5 for our method vs 1000 for DPS). Therefore, we argue that our perceptual quality to recovery trade-off is better than competing baselines.In addition to the guidance weight w, our method also allows tuning an additional hyperparameter λ, which controls the dynamics of the projection operator (See Sections 2.3 and 3.2 for more intuition). Therefore, tuning w and λ can help traverse the trade-off curve between perceptual quality and distortion for a fixed NFE budget. We illustrate this aspect in flexibility. Moreover, tuning the guidance weight in methods like DPS could be very expensive due to its high sampling budget requirement (around 1000 NFE).Diffusion Models:1. We include additional qualitative comparisons between Π-GDM and our proposed C-ΠGDM sampler for the ImageNet dataset in Fig. 4 . 2. We include a qualitative comparison between sample quality at different sampling budgets for the C-ΠGDM sampler in Fig. 5 . 3. We qualitatively study the impact of varying w on sample quality in Fig. 6 and the impact of varying λ on sample quality in Fig. 7 . 4. We qualitatively present the performance of the C-ΠGDM sampler for noisy inverse problems in Fig. 8 . In just 5 steps, our method can also generate good-quality samples for noisy inverse problems.Flow Models:1. We include additional qualitative comparisons between Π-GFM and our proposed C-ΠGFM sampler with different sampling budget for the all three datasets in Fig. 9 . 2. We qualitatively study the impact of varying w on sample quality in Fig. 10 Guidelines:• The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper. • The authors are encouraged to create a separate "Limitations" section in their paper.• The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be. • The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on a few datasets or with a few runs. In general, empirical results often depend on implicit assumptions, which should be articulated. • The authors should reflect on the factors that influence the performance of the approach.For example, a facial recognition algorithm may perform poorly when image resolution is low or images are taken in low lighting. Or a speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon. • The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size. • If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness. • While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an important role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.Question: Does the paper provide open access to the data and code, with sufficient instructions to faithfully reproduce the main experimental results, as described in supplemental material?Answer: [Yes] Justification: The code is not available at the time of submission, but will be published later.Guidelines:• The answer NA means that paper does not include experiments requiring code.• Please see the NeurIPS code and data submission guidelines (https://nips.cc/ public/guides/CodeSubmissionPolicy) for more details. • While we encourage the release of code and data, we understand that this might not be possible, so "No" is an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e. • Providing as much information as possible in supplemental material (appended to the paper) is recommended, but including URLs to data and code is permitted. 
experimental setting/details: 
 Question: Does the paper specify all the training and test details (e.g., data splits, hyperparameters, how they were chosen, type of optimizer, etc.) necessary to understand the results? Answer: [Yes] Justification: See the beginning of Section 3 Guidelines:• The answer NA means that the paper does not include experiments.• The experimental setting should be presented in the core of the paper to a level of detail that is necessary to appreciate the results and make sense of them. • The full details can be provided either with the code, in appendix, or as supplemental material. 7. Experiment Statistical Significance Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments? Answer: [No] Justification: We don't have any error bar data in our experiments. Guidelines:• The answer NA means that the paper does not include experiments.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions). • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors).• It should be clear whether the error bar is the standard deviation or the standard error of the mean. • It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have a 96% CI, if the hypothesis of Normality of errors is not verified. • For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative error rates). • If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text. 
experiments compute resources: 
 Question: For each experiment, does the paper provide sufficient information on the computer resources (type of compute workers, memory, time of execution) needed to reproduce the experiments? Answer: [No] Justification: Our approach utilizes established models for a range of downstream tasks, ensuring that hardware variations do not affect the outcomes. Guidelines:• The answer NA means that the paper does not include experiments.• The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage.• The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute. • The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make it into the paper). 9. Code Of Ethics Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https://neurips.cc/public/EthicsGuidelines? Answer: [Yes] Justification: the research conducted in the paper conform with the NeurIPS Code of Ethics Guidelines:• The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.• If the authors answer No, they should explain the special circumstances that require a deviation from the Code of Ethics. • The authors should make sure to preserve anonymity (e.g., if there is a special consideration due to laws or regulations in their jurisdiction). 10. Broader Impacts Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed? Answer: [Yes] Justification: See Section 5 Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact. • Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations. • The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology. • If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML). 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)? Answer: [NA] Justification: Our approach utilizes published models for a range of downstream tasks, so we do not need to add any safeguard. Guidelines:In this work, we use the terms interpolant and flows interchangeably.https://github.com/openai/guided-diffusionhttps://github.com/gnobitab/RectifiedFlow 
 
 Appendix:  


Decision: This paper presents a numerical sampling method for solving inverse problems using diffusion and flow models. The proposed technique is based on conjugate integrators which provide high-quality solutions to linear inverse problems in a few sampling steps. Overall the reviewers have acknowledged the efficacy of the method and its rigorous treatment of the problem and they have recommended accepting the paper.

Reviews: review: paper_summary: paper_summary: This work proposes a plug-and-play based method that leverages pretrained diffusion and flow models to solve inverse problems. The proposed method called Conditional Conjugate Integrators adapts previously proposed Conjugate Integrators framework for fast sampling of diffusion and flow models to solve linear inverse problems. The key idea is to project the diffusion (flow) dynamics into another latent space where sampling can be more efficient. Upon completion, the dynamics is projected back to the original pixel space. The paper provides the mathematical forms of the projection operator for conditional diffusion (flow) dynamics, then adapts it to get the projection operator for linear inverse problems. The paper also provides tractable forms of the projection matrix and its inverse. This derivation can be seen as a more general form of previously proposed method $\Pi$GDM ($Pi$GFM). The method can be extended to nonlinear and noisy settings. The paper provides promising results on datasets like LSUN Bedrooms, AFHQ Cats, FFHQ, ImageNet etc. on tasks like inpainting, superresolution, gaussian deblurring etc.
strength_weakness: strength_weakness: The benefits of the proposed method for the settings of noisy linear inverse problem setting as well as non-linear inverse problems remain unclear. The paper does not include any quantitative results for these two problem settings.  It only provides some limited qualitative results for super-resolution for the noisy setting with $\sigma_y$=0.05 and compression inverse problem and JPEG restoration problem for non-linear settings.
questions: questions: Questions:
1. Why does Eq. 117 hold true for VPSDE? I’m missing the simplification step here. 
2. Why does Table 1 not include results of inpainting for diffusion models on ImageNet? Similarly, Table 4 in the appendix skips results for inpainting on FFHQ.

Suggestions:
1. All the tables and figures must explicitly state all the relevant settings of the inverse problems. It is not immediately apparent that some of these quantitative results are only for noiseless linear inverse problems. 
2. Consider including an algorithmic box that summarizes C-$\Pi$GFM and  C-$\Pi$GDM. This would provide a concise overview of the method to the readers.
limitations: limitations: The paper discusses relevant limitations.
paper_summary: paper_summary: Inverse problems like super-resolution and deblurring remain computationally inefficient with current diffusion/flow models. The paper introduces a plug-and-play framework with Conditional Conjugate Integrators that use pre-trained iterative refinement models to simplify sampling. The method generates high-quality samples in as few as 5 steps for tasks like 4× super-resolution on ImageNet, outperforming existing methods.
strength_weakness: strength_weakness: 1. The reason why the paper is restricted to the $\Pi$GDM paradigm is unclear.
questions: questions: 1. Can the method be applied to nonlinear inverse problems? If not, please clarify the reason.
2. Though the Conditional Conjugate Integrators work well and can speed up the sampling process to 5 steps, a more detailed computational comparison is needed as each iteration becomes more complicated. Can you provide more detailed computational results?
3. Can you clarify the reason why you restrict yourself to the posterior approximation in $\Pi$GDM? Though the experimental results are promising, a direct comparison with DPS is preferred.
limitations: limitations: Yes, the paper mentioned its broader impact at the end, which looks good to me.
paper_summary: paper_summary: This paper introduces a inverse problem approach to reconstruct low-resolution blurred or obfuscated images, from pre-trained generative models. The method leverages conjugate integrators to project the diffusion dynamics in an easier space, solve the inverse problem and then map back.
strength_weakness: strength_weakness: The paper is extremely dense with equations and mathematical formulations. I would suggest moving sections 2.3 and 3.4, as well as the numerical details of the experiments, to the appendix. This would allow space in the main paper for additional qualitative results and a broader discussion on the intuition behind the conjugate integrators for inverse problems. Including a figure illustrating the method would also be appreciated. While the mathematical formulation is the main essence of the paper, a visual representation or descriptive insight would be appreciated for the understanding.
questions: questions: Minor details:
- Line 608, link the approximation to Eq. (5)
limitations: limitations: I don't have any concern.

Comments: title: Author Response
comment: We thank the reviewer for their response. Please find our detailed response below:

**Comparison of perceptual vs Recovery metrics:** We agree that highlighting robustness in both perceptual and recovery metrics is important in the context of inverse problems. Following previous works like DPS and $\Pi$-GDM, we include only perceptual metrics in the main text since recovery metrics like PSNR and SSIM usually favor blurrier samples over perceptual quality. So, a tradeoff exists between distortion (a.k.a recovery)  and perceptual quality. However, whether perceptual quality/recovery is preferred depends on the application. Therefore, for completeness, we provide a comparison between DPS, $\Pi$GDM, and C-$\Pi$GDM (our method) in terms of PSNR, SSIM, FID, and LPIPS in Tables 1 and 2 on the ImageNet-256 and FFHQ-256 datasets on the 4x super-resolution task. It is worth noting that the PSNR and SSIM scores for all methods correspond with the best FID/LPIPS scores presented in the main text for these methods.

Tables 1 and 2 show that our method achieves competitive PSNR and SSIM scores for better perceptual quality than competing methods, even for very small sampling budgets. For instance, on the FFHQ dataset, our method achieves a PSNR of 28.97 compared to 28.49 for DPS while achieving better perceptual sample quality (LPIPS: 0.095 for ours vs 0.107 for DPS) and requiring around 200 times less sampling budget (NFE=5 for our method vs 1000 for DPS). Therefore, we argue that our perceptual quality to recovery trade-off is much better than DPS, given our method is significantly faster than DPS.

**Traversing the Recovery vs Perceptual trade-off in C-$\Pi$GDM**: In addition to the guidance weight $w$, our method also allows tuning an additional hyperparameter $\lambda$, which controls the dynamics of the projection operator (See Sections 2.3 and 3.2 for more intuition). Therefore, tuning $w$ and $\lambda$ can help traverse the trade-off curve between perceptual quality and distortion for a fixed NFE budget. We illustrate this aspect in Table 3 (fixed $\lambda$ with varying $w$) and Table 4 (fixed $w$ with varying $\lambda$) for the SR(x4) task on the ImageNet-256 dataset using the PSNR, LPIPS, and FID metrics. Therefore, our method offers greater flexibility to tune the sampling process towards either good perceptual quality or good recovery for a given application while maintaining the same number of sampling steps. In contrast, other methods like DPS do not offer such flexibility. Moreover, tuning the guidance weight in DPS is very expensive in the first place due to its high sampling budget requirement (around 1000 NFE).

We will extend these comparisons for other tasks and add them in the Appendix section of our revised paper. We would also request the reviewer to reconsider their evaluation.

|                     | PSNR  | SSIM  | FID   | LPIPS |
|---------------------|-------|-------|-------|-------|
| DPS (NFE=1000)      | 23.81 | 0.708 | 38.18 | 0.252 |
| $\Pi$GDM (NFE=20)   | 21.92 | 0.646 | 37.36 | 0.222 |
| C-$\Pi$GDM (NFE=5)  | 22.32 | 0.641 | 37.31 | 0.220 |
| C-$\Pi$GDM (NFE=10) | 23.00 | 0.651 | 34.22 | 0.206 |
| C-$\Pi$GDM (NFE=20) | 23.16 | 0.654 | 34.28 | 0.207 |

**Table 1**: Comparison between different methods on ImageNet-256 for the SR(x4) task

|                     | PSNR  | SSIM  | FID   | LPIPS |
|---------------------|-------|-------|-------|-------|
| DPS (NFE=1000)      | 28.49 | 0.834 | 30.86 | 0.107 |
| $\Pi$GDM (NFE=20)   | 28.26 | 0.818 | 26.17 | 0.087 |
| C-$\Pi$GDM (NFE=5)  | 28.97 | 0.832 | 32.01 | 0.095 |
| C-$\Pi$GDM (NFE=10) | 29.03 | 0.821 | 29.07 | 0.086 |
| C-$\Pi$GDM (NFE=20) | 28.79 | 0.809 | 26.37 | 0.083 |

**Table 2**: Comparison between different methods on FFHQ-256 for the SR(x4) task.

| w  | PSNR  | LPIPS | FID   |
|----|-------|-------|-------|
| 2  | 22.91 | 0.339 | 48.48 |
| 4  | 23.37 | 0.306 | 45.03 |
| 6  | 23.49 | 0.274 | 42.68 |
| 8  | 23.44 | 0.266 | 40.96 |
| 10 | 23.28 | 0.254 | 40.27 |
| 12 | 22.89 | 0.246 | 40.13 |
| 14 | 22.74 | 0.239 | 40.16 |

**Table 3**: Illustration of the impact of $w$ for a fixed $\lambda=0.0$ on the sample recovery (PSNR) vs sample perceptual quality (LPIPS, FID) at NFE=5 for our method. The task is SR(x4) on the ImageNet-256 dataset.

| $\lambda$ | PSNR  | LPIPS | FID   |
|-----------|-------|-------|-------|
| -1.0      | 20.96 | 0.291 | 42.56 |
| -0.8      | 21.33 | 0.265 | 40.97 |
| -0.6      | 21.69 | 0.240 | 39.38 |
| -0.4      | 22.04 | 0.223 | 37.83 |
| -0.2      | 22.32 | 0.220 | 37.31 |
| 0.2       | 22.73 | 0.257 | 45.27 |
| 0.4       | 22.90 | 0.275 | 48.98 |
| 0.6       | 23.03 | 0.283 | 47.47 |
| 0.8       | 23.11 | 0.285 | 46.2  |
| 1.0       | 23.15 | 0.285 | 46.41 |

**Table 4**: Illustration of the impact of $\lambda$ for a fixed $w=15$ on the sample recovery (PSNR) vs sample perceptual quality (LPIPS, FID) at NFE=5 for our method. The task is SR(x4) on the ImageNet-256 dataset.
I appreciate the authors' efforts in preparing the rebuttal and acknowledge the theoretical contributions of the paper, which have been well-received by the other reviewers. However, I still have some concerns. The paper does not report PSNR and SSIM results, making it difficult to fully assess the quality of the proposed methods. Given that the paper utilizes the $\Pi$GDM paradigm, it is crucial to provide these metrics for a comprehensive evaluation. DPS, for example, demonstrates robustness in both perceptual and recovery metrics across various tasks, highlighting the importance of such evaluations. Considering above, I keep the original rating.
I thank the authors for their response. I would like to retain my score.