MAP Estimation with Denoisers: Convergence Rates and Guarantees
Authors: Scott Pesme
Giacomo Meanti
Michael Arbel
Julien Mairal
Venue: NeurIPS 2025 poster ()
Keywords: MAP estimation
proximal operator
denoiser model
plug and play methods
Primary area: optimization

Abstract: Denoiser models have become powerful tools for inverse problems, enabling the use of pretrained networks to approximate the score of a smoothed prior distribution. These models are often used in heuristic iterative schemes aimed at solving Maximum a Posteriori (MAP) optimisation problems, where the proximal operator of the negative log-prior plays a central role. In practice, this operator is intractable, and practitioners plug in a pretrained denoiser as a surrogate—despite the lack of general theoretical justification for this substitution. In this work, we show that a simple algorithm, closely related to several used in practice, provably converges to the proximal operator under a log-concavity assumption on the prior $p$. We show that this algorithm can be interpreted as a gradient descent on smoothed proximal objectives. Our analysis thus provides a theoretical foundation for a class of empirically successful but previously heuristic methods

MAP Estimation with Denoisers: Convergence Rates and Guarantees
Scott Pesme
Giacomo Meanti
Michael Arbel
Julien Mairal
Univ. Grenoble Alpes, Inria, CNRS, Grenoble INP, LJK firstname. Lastname@inri.....
Abstract
Denoiser models have become powerful tools for inverse problems, enabling the use of pretrained networks to approximate the score of smoothed prior distribution. These models are often used in heuristic iterative schemes aimed solving Maximum a Posteriori (MAP) optimisation problems, where the proximal operator of the negative log-prior plays central role. practice, this operator is intractable, and practitioners plug in a pretrained denoiser as a surrogate-despinn the lack of general theoretical justification for this substitution. In this work, we show that a simple algorithm, closely related to several used in practice, provably converges the proximal operator under a log-concyvity assumption on the prior p. We show that this algorithm can be interpreted as a gradient descent on smoothed proximal objectives. Our analysis thus provides a theoretical foundation for class of empirically successful but previously heuristic methods.
Introduction
Inverse problems are ubiquitous in scientific and engineering fields involving image acquisition. In many such problems, the object of interest is not directly observed but instead undergoes a degradation process-suue as blurring, downsampling, or noise corruption. The goal is to reverse this degradation and recover the original image.
classical approach formulates this task as an optimisation problem balancing two terms: a data fidelity term, modelling the observation process, and a regularisation term, encoding prior knowledge about the solution. Historically, regularisers such as total variation or wavelet sparsity were hand- crafted [Mallat, 1999]. While effective some extent, recent approaches often rely on data-driven priors, using pretrained denoisers and generative models. In particular diffusion and flow-based models offer powerful ways learn the true image distribution p from large datasets.
This opens the door to principled formulations like Maximum a Posteriori (MAP) estimation:
arg min Xf(x) In p(x), ZERd
(MAP)
which corresponds the posterior mode under the likelihood p(y x) xx exp(-A((x)) and prior p(x). In practice, however, this optimisation problem is extremely challenging to solve: evaluating the score -V p(x) is often intractable, the term In p(x)) be severely ill-conditioned, and the data fidelity term f(x) is frequently not strongly convex. A wide range of methods have been proposed address these problems, and many of them perform remarkably well empirically. Yet. these methods do not come with the guarantee of actually minimising the MAP objective, making their success difficult to interpret.
natural class of algorithms for addressing the MAP optimisation problem are proximal splitting methods [see, e.g., Beck and Teboulle, 2009, Figueiredo al., 2007, Combettes and Pesquet, 201 1], which are particularly effective when dealing with objectives that combine smooth and non-smooth
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

components. These methods alternate between two steps: one that follows the gradient of the data fidelity term, and another that incorporates prior knowledge through what is known as "proximal update" a correction step informed by the prior distribution.
However, for prior models relying on unknown data distribution, this proximal update is extremely difficult to compute exactly. To circumvent this, a popular line of work introduced by Venkatakrishnan al. [2013] known as Plug-and-Play (PnP) replaces the intractable proximal step with a pretrained denoising neural network. PnP methods have shown excellent empirical performance in a wide range of inverse problems. But despite their success, they come with significant caveat: the denoiser is not designed to match the proximal operator it replaces. As result, the overall algorithm no longer corresponds solving the original MAP estimation problem, which limits its interpretability and makes it hard to analyse theoretically unless strong constraints are imposed on the denoiser [Hurault et al., 2022, Sun et al., 2021, Hertrich et al., 2021, Cohen et al., 2021].
More recently, a new wave of approaches has emerged which view inverse problems as a sampling task, see [Delbracio and Milanfar, 2023, Chung al., 2023, Boys al., 2024] among others, moving further away from traditional optimisation frameworks. One example is the Cold Diffusion [Bansal al., 2023] algorithm, which combines denoising steps with corruption steps towards the observed data, with decreasing intensity. While this method often produces high-quality results in practice, especially with a small number of steps. it also lacks strong convergence guarantees and may diverge during extended runs with default parameters [Delbracio and Milanfar, 2023].
In this work, revisit denoising-basee iterative schemes from a theoretical perspective, focusing on the case where the negative log-density p is log-concave and potentially ill-conditioned. Specifically, show that a simple algorithm originally proposed by Bansal et al. [2023] with appropriate step- sizes converges to the proximal operator of the negative log-density, and we establish corresponding convergence rates. Having a reliable approximation of the proximal operator enables its integration into broader MAP estimation frameworks, akin to Plug-and-Plan methods, but now supported by a rigorous theoretical foundation.
Our contribution:
establishing convergence rates for MAP estimation. In this work, for
suitable choice of sequences of noise levels o 2 0 and weights as (0, 1), we consider the following recursion to compute the proximal operator of - In p point y E Rdd (MMSE Averaging)
26+1 (1 ax) MMSE,, (Ix) a.yy
with MMSE, = E[X X tE -====:00
=], =2,
where the expectation taken over X p and practice, the theoretical minimum mean square error denoiser MMSE, can be approximated by a neural network which has been trained to match the MMSE denoiser.
/N La) conditionally X + or 2. In
Each iterate in the recursion is a weighted average between a denoised version of the current point and the original input y. echoing the structure of methods like Cold Diffusion [Bansal et al., 2023]. What makes this recursion striking is that, for appropriate choices of weights ar and vanishing noise levels og 0, it can be rewritten (see Proposition I)-via the Tweedie formula [Efron, 2011]
Ik+1 T 0k Fo, (IK),, with Fzx -
| 2 - In Poz (z).
where Po denotes the convolution of the prior with a Gaussian of variance σ². Under this reinterpre- tation, the recursion corresponds to gradient descent on a sequence of smoothed objectives (Fau), )x converging to the true proximal objective F(x) == _ Ily x1|2 In p(z) whose minimiser the proximal point prox_, In (w). This perspective enables a rigorous convergence analysis: as + O, each update more closely resembles step on F. and the iterates can be shown to converge to its minimiser.
We show that, under a log-concavity assumption on and bound on the third derivative of In p. the iterates Fx of the MMSE Averaging recursion converge to the true proximal point the following rate (see Theorem 1):
ITe prox_rx (v))) Ö(1/))
where Ö(-) hides logarithmic factors. Importantly, our convergence bound does not rely on the L-smoothness constant of the negative log-prior In p. which could be arbitrarily large.

---

This result provides theoretical grounding for algorithms that previously lacked variational in- terpretation, establishing a direct connection between heuristic denoising schemes and principled optimisation algorithms. Crucially, it yields an explicit method to approximate the proximal operator of the negative log-prior-a a central building block in many optimisation frameworks for inverse problems [ /Venkatakrishnan et al., 2013, Romano et al., 2017, Hurault et al., 2021]. Once available, this proximal operator can be readily integrated into broader algorithms, such as proximal gradient descent and its accelerated variants [Beck and Teboulle, 2009]. In Theorem 2, we demonstrate exactly this by plugging our approximation into a proximal gradient method to solve the MAP problem.
The proof of convergence with explicit rates of the MMSE Averaging iterates towards the proximal operator, while conceptually intuitive, requires a careful blend of inexact optimisation analysis and tools from partial differential equations most notably the heat equation- to control how the minimiser of the smoothed objectives F. evolves with the noise level.
2 Related Works
Our work shares similar motivations with much of the literature on Plug-and-Play (PnP) methods for inverse problems [Venkatakrishnan et al., 2013]. The PnP literature is vast, and for particularly clear and comprehensive overview, we refer the reader to the PhD thesis of Samuel Hurault [Hurault, 2023]. PnP methods replace the proximal operator prox-y -y In p (v) with a generic denoiser D.. typically parameterised by the noise level Ø. A wide variety of denoisers have been used, including classical approaches [Dabov et al., 2007, Zoran and Weiss, 2011], CNN-based denoisers [Zhang et al., 2021, Kamilov et al., 2023, Zhang et al., 2017] and, more recently, diffusion models [Graikos et al., 2022. Zhu et al., 2023]. These methods are often combined with different optimisation schemes (e.g.. PGD [Terris et al., 2020], ADMM [Romano et al., 2017], HQS [Zhang al., 2017I) and adapted to different specific inverse problems. Several works [Sreehari et al., 2016, Gavaskar and Chaudhury, 2020, Nair et al., 2021, Xu et al., 2020] show that a variety of PnP algorithms converge, however they cannot guarantee that the denoiser is a proximal operator, let alone the proximal operator of the correct functional. Furthermore the convergence proofs often rely on restrictive assumptions on the denoising model [Reehorst and Schniter, 2018). Indeed, the denoiser is usually trained [Zhang et al., 2021, Meinhardt et al., 2017] to minimize the MSE and converges to the MMSE estimator which can be very different from the MAP [Gribonval, 2011].
Gaussian noise assumptions
Gradient step (GS) denoisers [Cohen al., 2021, Hurault et al., 2021] parameterize D. where 9- is neural network. It is then possible to show that D, is indeed the proximal operator of an explicit functional [Hurault et al., 2022], but this function is unfortunately not the negative log prior as desired. Similarly, Hauptmann et al. [2024] link linear denoisers to the proximal operator of a regularization functional, which is however again not Iny.
Vg
Two recent theoretical works share our concerns about existing PnP methods and strive to learn the correct proximal operator: Fang et al. [2023] replace the usual MSE loss by a proximal matching loss which is guaranteed in the limit to yield pn Though elegant, they do not establish any convergence rate, and their training procedure only approximates the desired limit without providing a bound on the approximation error. Using an approach somewhat close to ours, Laumont et al. [2023] introduce PnP-SGD, which performs stochastic gradient descent on smoothed version of the proximal objective F,. However, by keeping the smoothing parameter fixed (o = σ), their method only converges to the proximal operator of the smoothed density and the convergence rate depends on the smoothness constant of Fa, which can be arbitrarily large and lead to slow convergence as explained in this work.
The second class of approaches which are receiving more and more attention in the context of solving inverse problems are conditional diffusion methods. These algorithms are typically based on modifying the smoothed prior score V In (=.,)) through a pretrained diffusion model- into the posterior score V In Pn (a=- v). Coupled with sampling along the reverse diffusion SDE this allows to generate samples from the desired probability distribution. Dhariwal and Nichol [2021] propose to use a classifier to estimate In p(y z)) Jalal et al. [2021] approximate (y xx)= R p(v .)) obtained through the explicit likelihood term under Gaussian noise, the DPS algorithm [Chung et al., 2023] approximates the mean of the smoothed log prior with the Tweedie formula and Boys et al. [2024] additionally approximates the standard deviation. All such methods aim to sample from the posterior distribution rather than identify its maximum. Moreover, they rely on approximations that are difficult to control, offering no guarantees of sampling from the true posterior. Although
3

---

asymptotic guarantees can be achieved with more sophisticated algorithms [Wu et al., 2023], these methods are not designed recover the MAP estimate.
Using flow matching, Zhang et al. [2024] approximate the MAP solution directly, without relying on the proximal operator. Instead, they construct trajectory that trades off between the prior and data fidelity terms, but no convergence rates are given. Finally, Ben-Hamu et al. [2024] solve a similar problem, but additionally need an expensive backpropagation step through an ODE every step.
3 Main Result: Convergence Towards the Proximal Operator
We begin by showing that the MMSE Averaging recursion corresponds to gradient descent on sequence of smoothed approximations the proximal objective F. We then show that these smoothed objectives are significantly better conditioned than the original unsmoothed problem. Finally, we prove convergence of the iterates and provide explicit convergence rates.
3.1 From MMSE Averaging to Gradient Descent on Smoothed Proximal Objectives
We can connect the recursion in MMSE Averaging to the negative log-prior - In p by leveraging the celebrated Tweedie identity (see for example Efron [201 I]), which links the MMSE denoiser the gradient of the log-density of smoothed version of the prior. Specifically, if pa denotes the Gaussian convolution of p with centred Gaussian of variance o² (i.e. the density of X or), then:
MMSE, (2) + o (2).
Plugging the above identity into the MMSE Averaging recursion allows expressing the iterate update in terms of the score of the smoothed density Pas. which already resembles gradient descent update:
(1 - ark o
(Ix y)
2 In In (or((())) (1x))
Fk+1
Okk
a
Rearranging the terms in the above expression naturally leads to the following simple observation: Proposition The MMSE Averaging recursion with choice of weights ag 1/(k: + 2) and noise sequence of = - + 1) can be rewritten:
Øk+1 Fx - og Fo, (#k), with Fo, (x) = [1y sr12 In P.r. (x)).
In this form, the recursion is naturally interpreted sequence of smoothed proximal objectives (Foz), )x and with stepsizes (ax)k- This reformulation not only enables clean convergence analysis but also offers new perspective on the MMSE Averaging recursion: as 0j -> 0, one can hope that the iterates approach the minimiser of the original (unsmoothed) proximal objective:
gradient descent algorithm applied to
F(x) == 2 ly x1| p(ar).
(Proximal Objective)
Moreover, we argue that this smoothed approach leads faster convergence than applying standard gradient descent directly to the original, potentially badly conditioned Proximal Objective.
3.2 Good Conditioning Properties of Fo
Compared to the original objective F, the function F,, enjoys much better properties. In particular, the next result shows that F, -moothh, with smoothness controlled by the noise level a. Proposition 2. For any a > 0, the function Fg is Lm-shhoth, with
T
Lg 1
02
The proof can be found in Section A and a simple consequence of known results on the Hessian of - In Pr. This smoothing effect introduces a natural trade-off: for large the objective F, becomes easier to minimise thanks to an improved smoothness, but the minimiser F, may deviate significantly from that of the original problem. While this smoothness property holds for any density

---

function p. obtaining convergence guarantees requires stronger assumptions. In particular, we will focus on the case where p log-concave and satisfies regularity conditions. Although this assumption is clearly idealised and does not hold for many practical distributions, it offers a manageable setting for theoretical analysis.
Assumption 1. The density p is log-concave, and stricily positive IR..
In particular, this ensures that - Inp p is well-defined and convex over R so that the Prox- imal Objective function F is 1-strongly convex and admits a unique minimiser, denoted by prox - Inp(() = arg min F. Furthermore, the stability of log-concavity under convolution (a prox_r
special case of the Prékopa-Leindler inequality, see [Saumard and Wellner, 2014, Proposition 3.5.J) ensures that In Po is convex for all > O, and hence that Fa is 1-strongly convex. Along with Proposition this allows quantify how much the smoothing improves the conditioning of the objective in the following proposition.
Proposition 3. Under Assumption 1, for a O, the function Fg Lo -mmoothh and P -o-stonngly convex with L. +/O0² and pla = 1. The condition number of Fa is therefore most
La - (1+
Kig
Po
This result highlights key benefit of the smoothed proximal objective: as increases the function Fg becomes significantly better conditioned, with the condition number Kg decreasing toward 1 as + 0o. For example, setting V already yields condition number of K.7 2.
Next, we impose an assumption on the third derivative the log-prion which is crucial in our analysis for controlling the Lipschitz continuity of the map o² + arg min Fg. Without such control, it would be difficult to establish any meaningful convergence guarantees for the iterates of MMSE Averaging.
Assumption 2. The prior p three times differentiable and the third derivative of in p is bounded. We denote by M 2 0 the quantity:
V |>'n M,
sup 2E3d
where for A E Rdxdxd, ||||| = (Eis,k A? A2th) 1/2 corresponds to the Frobenius norm.
This assumption controls how skewed and "non-quadratic" the to control the stability of the minimisers prox- ln (3) == arg min Fs varies. Also note that an upper bound on the third derivative does not imply an upper bound on the second one: indeed for Gaussian prior p. its negative log likelihood is simple quadratic, which have arbitrarily large L-smoothness, while its third derivative trivially
is, and we make it in order
3.3 Convergence of the MMSE Averaging Iterates Towards the Proximal Operator
Leveraging the upper bound on the condition number of the objectives (Fo)o20. we obtain the following convergence result on the iterates 2ą of the MMSE Averaging recursion:
Theorem 1 (Convergence to the Proximal operator). Under Assumptions I and 2, let prox_ In (v) denote the unique solution of the Proximal Objective problem. Then, the MMSE Averaging iterates with parameters ag 1/(k + 2), of T/(k
(k- + 1) and initialised at zo
satisfy:
(In k)) +7
|| prrxx top (v)Il
[lw 1pwxxprxx... p +
Comparison with naive GD: illustration with a Gaussian prior. The most important part of our result that the convergence bound does depend on the L-smoothness of In p. which could be arbitrarily large. The convergence rate depends only on a bound on the third derivative of In which may remain moderate even when the second derivative is large. This is unlike gradient descent (GD) applied directly to the proximal objective F, whose rate scales poorly with the L-smoothness of - In We illustrate this with toy yet instructive case Gaussian prior, for which the third derivative of the log likelihood is trivially zero, yet the second derivative can be arbitrarily large. Let p the density of a d-dimensional Gaussian (O, H-1), with H positive definite matrix whose smallest eigenvalue we arbitrarily consider to be and whose largest eigenvalue L > can be
5

---

arbitrarily large. In this setting the negative log-prior In is a quadratic with Hessian H and F a quadratic too with Hessian equal to (I+ + -H). The corresponding smoothness constant of F is therefore Lr = 1 + 7 I., and the strong convexity constant is μr + T. Since [ can be arbitrarily large, gradient descent on requires an arbitrarily small (and con-practical) step size a < 1/ L.p. For = F,, the iterates satisfy the standard convergence bound:
| prox- (y)(l)l - (1 (1-1s11111
FF
k/2
1y
prox In In((v)]l
LF
leading to an iteration complexity of L- log(L/))) to reach e-accuracy, From Theorem 1, since M = the MMSE Averaging iteration converges much faster, with rate (1/kk) (i.e. iteration complexity O(1/E)), which is tight up the log term (see Section A.3).
Parameter-free algorithm. A key practical advantage of our result is that it guarantees convergence for parameter-free choice of weights ag and noise levels o.e Specifically, these sequences depend only on the chosen regularisation parameter I and do not require any knowledge of smoothness or Lipschitz constants, condition number, or other This makes the algorithm particularly simple to use and eliminates the need for costly hyperparameter tuning.
properties of the prior distribution p.
Sketch of proof. The proof (given in Section A.3) combines techniques for approximate gradient optimization and a priori estimates on the solution to a partial differential equation. We begin by applying the standard descent lemma to the smoothed objective Fo,, which yields a contraction towards its minimiser at a rate determined by the condition number Kay which is controlled through Proposition 3. guaranteeing consistent progress. However, because the minimiser of Fr,, changes with Tb. we must control how much it drifts over the iterations. To do this, we study the evolution of the minimiser of F, as a function of a by analysing the differential equation it satisfies. This made possible by the fact that Po satisfies the heat equation. The resulting ODE for arg min F, involves the quantity VA In Pa. which we are able to bound uniformly in a by M M by carefully analysing the parabolic inequality satisfied by ||V3 In (x())l p and using the bound from Assumption 2 for = O. Summing the incremental drift contributions and combining them with the contraction bound yields the final convergence result toward the true proximal point.
Link with cold diffusion. There is a notable similarity between our algorithm and a heuristic approach introduced in Bansal et al. [2023], which generates images by inverting known degradation. When the degradation operator is defined as a linear interpolation between the degraded image and the clean image (as explained in Section 6.2 in Delbracio and Milanfar (2023]), cold diffusion initialises at Fo = y and applies the following iteration for a fixed number of steps N:
k N
Ik+1 oz Do k) + agys with ag
where Do is trained denoiser, as for our recursion MMSE Averaging. However, note that the choice Or := k/N differs from the schedule used in our theoretical analysis. While this empirical scheme yields strong results for very small N. it lacks convergence guarantees and tends to diverge the number of iterations increases. We suspect that this instability may be due to the fact that the fixed ratio K/N does not necessarily correspond to a well-behaved weighting policy.
Comparison with standard random smoothing techniques. The smoothing that appears through - Im differs significantly from classical random smoothing approaches (c.g.. Nesterov and Spokoiny [2017]). In standard random smoothing, the goal is regularise a possibly non- smooth function h by convolving it with a Gaussian, yielding a smooth approximation h=(2)= Eg-A(0.02 ri [= + =)]. This smoothed function h inherits favourable differentiability properties that are well understood and can be leveraged in zeroth-order or gradient- based optimisation. In In Pr. where
contrast, our approach considers the logarithm of a smoothed Po is the Gaussian convolution of a density p. This subtle change has major impact: the logarithm does not commute with convolution, and the resulting function exhibits different analytic properties. As result, existing results from the standard random smoothing literature cannot be directly applied.

---

Extension to priors supported on an affine subspace. Our analysis naturally extends to the case where the prior distribution μ is supported on an affine subspace S R of dimension < representing first step toward modelling the assumption that clean images lie manifold within the ambient space. Indeed, assuming that the restriction ja to S admits a positive log-concave density p with respect to the r-dimensional Lehesgue measure on S, the smoothed density Pa then defined over id and can naturally be decomposed into a Gaussian term orthogonal to S and convolution restricted to S. Specifically, for any point E IRd the smoothed density Po (=) factorizes into a Gaussian penalty for the distance Importantly, this decomposition allows derivatives intrinsic to As a result, Theorem still holds but with ambient dimension replaced by the effective dimension r' d. We formally prove this in Section A.4.
low-dimensional
to S, and the intrinsic smoothing of p along express the third-order derivatives of In Pa in terms of
Extension when using approximate scores. In practice we do not have access to the exact V In Po but only an approximation the score, often provided by trained neural network 9o R In P.. In this more realistic case, the MMSE Averaging recursion becomes 2x+11 2g O&k (V (Ix) + where & == In Poz (2x) gr, (arz) denotes the approximation error at step k. Assuming that these errors are uniformly bounded along the trajectory, i.e. ||Ex|| < &. we can show that the iterates converge to point at distance O(g) from the true proximal point, with the same rate as in Theorem I.
We refer to Section A.3 for the proof.
From Approximate Proximal Operators to MAP Estimation
We now return the original MAP optimisation problem, recalled here:
arg min Af(x) In p(x).
=E3d
We denote the objective by J(x) = Xf(x) - -Ip((r) and work under the following assumption on the data fidelity term
Assumption 3. The data fidelity term f is convex, lower-bounded, and L
This a mild assumption that holds for many common data fidelity terms, such as f(x) = 211 | Ar- yl which is L -smooth with L = 1/Amax (AT A). Note that we do not require f to be strongly convex. Under this assumption, we denote aMA E arg min J any minimiser of J.
Algorithm. When the proximal operator is accessible, minimising can be achieved using proximal gradient descent, starting from (0) = y:
(n+1) prox (2(()(() (x x
(PGD)
Under Assumptions and 3, the classical result of Beck and Teboulle [2009] (see their Theorem 3.1) guarantees that for step size < 1/(ALy), the following convergence rate holds:
||y - TMAP 1
(xx - J(MMAP) S
2rn
In our setting, however, we not have direct access to the exact proximal operator Ir Instead, we compute approximate version using the MMSE Averaging recursion. Given sequence (km))21 specifying the number of internal iterations used to approximate each proximal step, this leads naturally to an inexact proximal gradient descent algorithm.
7

---

Algorithm Approximate Proximal Gradient Descent (Approx PGD)
Require: Noisy image y, Step size parameter number of inner steps (kn)n21
Initialise: z:00) +-1 y
for 2, do
Data fidelity gradient descent step (x+1) - [()) ~
20
Approximate proximal step flm+1) -oxx-1lp(( (=0 (
for k 0 bx+1 do
E+I
Og f
k+2
[n+1)
f- (11-a))MMB MMSE,, (2k (((att)) + 2 (n+1)
2x+1
end for
((++1)
(n+1)
f- 2k
end for
We prove the following convergence result for the approximate proximal gradient descent iterates from Algorithm 1.
Theorem 2 (Convergence towards the MAP estimator with explicit bounds). For T S Xe, and number of steps in the inner loop which increases kn Le nlta for C, > O, under Assumptions I to 3 the approximate proximal gradient descent iterates (((")))rrrn Algorithm I satisfy:
((x - (AMAP) S SO { and [| x(10)||
1 1+7
S
1=1
where ((())_ prox_ (AV denotes the exact proximal update at iteration n.
The constants hidden in the ((1/n) and ((1/n) terms depend explicitly on the problem parameters and are given in detail Section A.5.
Comment on the convergence bound. This result provides a meaningful convergence guarantee in the context of MAP estimation. Since we do not assume strong convexity of f. is more natural to measure progress through convergence in function value rather than in the iterates themselves. However, a direct bound on (()(n)) J cannot be expected in general: because the iterates ((() are only approximate updates of the true proximal points (n)), even a small error between ((n) and
((()) can result in a large discrepancy objective value due to the potentially poor conditioning of
Instead, our analysis shows that the iterates (n) are close to exact proximal iterates (x))), whose average MAP error provably small, As result, even though we cannot directly control J/ (á(").
we ensure that the iterates are close the iterates (n) which provably converge (in average) towards the optimum.
Sketch of proof. We start from the classical descent inequality for proximal gradient updates. Since we use approximate proximal steps ((m), we quantify the error ((r)) = (()) - (x())) using Theorem 1 and bound its impact on the objective. Summing over iterations and controlling the errors yiclds the O(1nn) rate for the objective. The second bound follows directly from the convergence of the inner loop to the true proximal operator thanks to Theorem 1. Note that although our proof follows a similar strategy to that of Schmidt al. [2011], which analyses inexact proximal gradient methods, their results do not directly apply becaus our approximation guarantee from Theorem concern the iterates and not the objective function values.
Finally, note that while we consider an approximate version of proximal gradient descent, one could also analyse its accelerated counterpart, in the spirit of FISTA Beck and Teboulle [2009], which would yield faster convergence rates under the same assumptions. We leave this direction for future work.

---

Numerical Visualisations
To better understand the effect of smoothing on the proximal objective- and how influences the gradient descent trajectory- consider simple two-dimensional example where the prior is
Gaussian distribution with a highly anisotropic covariance Σ 1/L. for L > 1. In this
setting, the density p( Ø2) is sharply concentrated around the Ii -i-axis and rapidly decays as soon 22 + O. The corresponding proximal objective F(x) 2 ly =|| In p(x) is then quadratic function with Hessian equal to
V F=)))=1++22 (100- 11+
1+TL)
As illustrated in Figure this severe ill-conditioning leads gradient descent on F to stagnate, making very little progress toward the true proximal point prox In (v).
However, smoothing the prior leads significant change in behaviour. Since Pa corresponds the convolution of p with a Gaussian of variance o², it remains Gaussian with covariance E, E+ 12. The smoothed proximal objective Fg (x) 1 fly x| r In 1l (I)) with Hessian
is then also quadratic, but now
(z) (1 7 02)
7/(1
1 1 + +.o2)
As increases, this Hessian interpolates between the poorly conditioned V2F and the well- conditioned identity matrix Iz. This transition is clearly visualised in Figure I, which shows how the level curves of F. become more isotropic as a increases. However, while smoothing improves conditioning, it also causes the minimiser prox Io (v) = arg min F, to drift away from the solu- tion prox In (v) = arg min F which we ultimately aim to recover (the red triangle in Figure 1). This highlights the need for decreasing schedule of d within the recursion: benefit from better conditioning at early stages while still converging to the correct solution. This strategy leads to significantly improved optimisation performance.
shown Figure 2, gradient descent applied to the sequence of smoothed objectives (Fr, )å, using the step size and noise schedule specified in Proposition 1, converges rapidly to the desired solution.
T/2
y
422 axis
X prox Io (N)
223 2x15
22 Miis
X
X
axis
T1 axis
T1 axis
Figure 1: Visualisation of the level curves of the smoothed proximal objective (x(x)
x 1
7 In Po (z) for different values of The unsmoothed objective F is poorly conditioned (left plot).
but the conditioning improves significantly as a increases.

---

2'k
prox-; Inp(y)
GD on F
GD (F
T22 axis
10-1
GD on F
GD on (Foz)e
10a
10a
10a
T1 axis
Number of iterations
Figure 2:
Illustration of the iterate trajectories (left plor) and convergence rates (right plot) of naive gradient descent on F (which has condition number x 500) versus gradient descent on the smoothed objectives (Fa.)a. using toy 2D Gaussian prior. Gradient descent on F, using stepsize a 0.8/ F (chosen for better visualisation), suffers from poor conditioning and makes little progress toward the optimal solution prox-r In (1). In contrast, gradient descent the smoothed objectives (Frr)) ) converges rapidly, clearly exhibiting O(1/k) rate.
Conclusion
In this work, prove that the iterative denoising- based scheme MMSE Averaging converges to the proximal operator of the negative log-prior In a central component in MAP estimation for inverse problems. We show that, under suitable choices of averaging weights or and noise levels ds. the algorithm can be interpreted as gradient descent on a sequence smoothed proximal objectives. Leveraging this perspective, prove that the iterates converge the true proximal point at rate of Ö(1/k), under the assumption that the prior p log-concave and has bounded third derivatives.
This result offers a principled foundation for a class of denoising- schemes and connects them to classical optimisation theory. Importantly, it provides explicit way approximate the proximal operator of In p. enabling the use of standard proximal methods to solve the MAP problem. We demonstrate this by incorporating our approximation into proximal gradient descent and deriving convergence guarantees for the resulting algorithm.
Despite these advances, our theoretical guarantees rely strong assumptions most notably that the prior log-concave, sufficiently smooth, and supported on all of RR.. Extending the analysis to more realistic settings, such as non-convex priors or those supported on low-dimensional manifolds, an exciting direction future work.
Acknowledgements
S. Pesme would like to warmly thank Filippo Santambrogio, who kindly and spontaneously responded email asking for help with Lemma 5. The proof present is entirely based the email exchange we had. S. Pesme is also grateful to Nikita Simonov, who generously replied to a similar email with several insightful suggestions for approaching the same lemma. Finally, S. Pesme would like to thank Loucas Pillaud-Vivis for the many valuable discussions they had regarding the proof' of Lemmas and 6.
This work was supported by ERC grant number 101087696 (APHELEIA project) as well as by the ANR project BONSAI (grant ANR-23-CER 23-CC23-0011-011
10

---

References
Arpit Bansal, Eitan Borgnia, Hong-Min Chu, Jie Li, Hamid Kazemi, Furong Huang, Micah Goldblum, Jonas Geiping, and Tom Goldstein. Cold diffusion: Inverting arbitrary image transforms without noise. Advances in Neural Information Processing Systems, 08841259-412888 2023.
Amir Beck and Mare Teboulle. A fast iterative shrinkagg-threshooldnng algorithm for linear inverse problems. SIAM journal imaging sciences, 2(1):183-22) 202, 2009.
Heli Ben- -emmH Omri Puny, Itai Gat, Brian Karrer, Uriel Singer, and Yaron Lipman. D-flow: Differentiating through flows for controlled generation. In Proceedings of the 41st International Conference on Machine Learning. volume 235 of Proceedings of Machine Learning Research, PMLR, 21-27 Jul 2024.
pages
Benjamin Boys, Mark Girolami, Jakiw Pidstrigach, Sebastian Reich, Alan Mosca, and O. Deniz Akyildiz. Tweedie moment projected diffusions for inverse problems, 2024.
Hyungjin Chung, Jeongsol Kim, Michael Thompson Mecann, Mare Louis Klasky, and Jong Chul Ye. Diffusion posterior sampling for general noisy inverse problems. In International Conference Learning Representations, 2023.
Regev Cohen, Yochai Blau, Daniel Freedman, and Ehud Rivlin. It has potential: Gradient-driven
denoisers for convergent solutions inverse problems. Advances in Neural Information Processing Systems, 34:18152-188111 1816,, 2021.
Patrick L Combettes and Jean-Chriss Christophe Pesquet.
Proximal splitting methods in signal processing. Fixed-point algorithmss for inverse problems in science and engineering, pages 185-212, 2011.
Kostadin Dabov, Alessandro Foi, Vladimir Katkovnik, and Karen Egiazarian. Image denoising by sparse 3-d transform-domait collaborative filtering. IEEE Transactions on Image Processing, 16 (8):2080-2095, 2007. doi: 10. (009/11P2000.9902233.
Mauricio Delbracio and Peyman Milanfar. Inversion by direct iteration: An alternative to denoising diffusion for image restoration. arXiv preprint arXiv:2303. 2023.
Prafulla Dhariwal and Alexander Quinn Nichol. Diffusion models beat GANs on image synthesis. In A. Beygelzimer, Y. Dauphin, P. Liang, and Wortman Vaughan, editors, Advances in Neural Information Processing Systemx, 2021.
Bradley Efron. Tweedie's formula and selection bias, Journal of the American Statistical Association, 600-1661,, 2011.
Lawrence C Evans. Partial differential equations, volume 19. American Mathematical Society, 2022.
Zhenghan Fang, Sam Buchanan, and Jeremias Sulam. What's in a prior? learned proximal networks for inverse problems. arXiv preprint arXiv:2310. 14344, 2023.
Mário AT Figueiredo, Robert Nowak, and Stephen J Wright. Gradient projection for sparse reconstruction: Application to compressed sensing and other inverse problems. IEEE Journal of 0908886--97,0000
selected topics in signal processing,
Ruturaj G Gavaskar and Kunal N Chaudhury. Plug-and-play ista converges with kernel denoisers. IEEE Signal Processing Letters, 27:610-61 614, 2020.
Alexandros Graikos, Nikolay Malkin, Nebojsa Jojic, and Dimitris Samaras. Diffusson models as plug-and-play priors. In Alice H. Oh, Alekh Agarwal, Danielle Belgrave, and Kyunghyun Cho, editors, Advances in Neural Information Processing Systems, 2022.
Rémi Gribonval. Should penalized least squares regression be interpreted maximum posteriori estimation? IEEE Transactions on Signal Processing, 59(5):2405-24 2410, 2011.
A Hauptmann, S Mukherjee, CB Schönlieb, and F Sherry. Convergent regularization in inverse problems and linear plug-and-play denoisers. Foundations of Computational Mathematics, 2024. doi: 10. N0.17833/CMM.108

---

Johannes Hertrich, Sebastian Neumayer, and Gabriele Steidl. Convolutional proximal neural networks and blug-and-pul algorithms. Linear Algebra and its Applications, 631:203-234, 2021.
Samuel Hurault. Convergent tttu--ndd-layy methods for image inverse problems with explicit and nonconvex deep regularization. PhD thesis, Université de Bordeaux, 2023.
Samuel Hurault, Arthur Leclaire, and Nicolas Papadakis. Gradient step denoiser for convergent plug-and-play In International Conference on Learning Representations, 2021.
Samuel Hurault, Arthur Leclaire, and Nicolas Papadakis. Proximal denoiser for convergent plug- and-play optimization with nonconvex regularization. In International Conference on Machine Learning, pages 9483-9505. PMIR, 2022.
Ajil Jalal, Marius Arvinte Giannis Daras, Eric Price, Alexandros G Dimakis, and Jon Tamir. Robust compressed sensing mri with deep generative priors. In M. Ranzato, A. Beygelzimer, Y. Dauphin, P.S. Liang, and Wortman Vaughan, editors, Advances in Neural Information Processing Systems, volume 34, pages 14938-14954, Curran Associates, Inc., 2021.
Ulugbek S Kamilov, Charles A Bouman, Gregery T Buzzard, and Brendt Wohlberg. Plug-and-piae play methods for integrating physical and learned models in computational imaging: Theory, algorithms, and applications. IEEE Signal Processing Magazine, 40(1):85-97, 2023.
Rémi Laumont, Valentin De Bortoli, Andrés Almansa, Julie Delon, Alain Durmus, and Marcelo Pereyra. On maximum a posteriori estimation with plug & play priors and stochastic gradient descent. Journal of Mathematical Imaging and Vision, 65(1):140-163, 2023.
Ki-Ahm Lee and Juan Luis Vázquez. Geometrical properties of solutions of the porous medium equation for large times. Indiana University Mathematics Journal, pages 991-1016, 2003.
Stéphane Mallat. wavelet tour of signal processing. Elsevier, 1999.
T. Meinhardt, M. Moeller, C. Hazirbas, and D. Cremers. Learning proximal operators: Using denoising networks for regularizing inverse imaging problems. In 2017 IEEE International
Conference on Computer Vision (ICCV), pages 1799-1808, 2017.
Pravin Nair, Ruturaj Gavaskar, and Kunal Narayan Chaudhury. Fixed-point and objective conver- gence of plug-and-play algorithms. IEEE Transactions on Computational Imaging, 7:337-348, 2021.
Yurii Nesterov. Introductory lectures on convex optimization: A basic course, volume 87. Springer Science & Business Media, 2013.
Yurii Nesterov and Vladimir Spokoiny. Random gradient-free minimization of convex functions. Foundations of Computational Mathematics, 17(2):527-566, 2017.
Bernt Oksendal. Stochastic differential equations: an introduction with applications. Springer
Science & Business Media, 2013.
Edward T Reehorst and Philip Schniter. Regularization by denoising: Clarifications and new interpretations. IEFE transactions on computational imaging, 5(1):52-67, 2018.
Yaniv Romano, Michael Elad, and Peyman Milanfar. The little engine that could: Regularization by denoising (red). SIAM Journal on Imaging Sciences, 10(4): 2017. doi: 10. 1377/ 10284.
Adrien Saumard and Jon A Wellner. Log-concavity and strong log-concavity: review. Statistics surveys, 8:45, 2014.
Mark Schmidt, Nicolas Roux, and Francis Bach. Convergence rates of inexact proximal-gradien methods for convex optimization. Advances neural information processing systems, 24, 2011.
Suhas Sreehari, S. Venkatakrishnan, Brendt Wohlberg, Gregery T. Buzzard, Lawrence F. Drummy, Jeffrey P. Simmons, and Charles A. Bouman. Plug-and-play priors for bright field electron tomography and sparse interpolation. IEEE Transactions on Computational Imaging, 2(4):408- 423, 2016. doi: 10. 10.1009TC1220112599977
12

---

Yu Sun, Zihui Wu, Xiaojian Xu, Brendt Wohlberg, and Ulugbek S Kamilov. Scalable plug-an Blug-and-pp -an---d-pp and-play admm with convergence guarantees. IEEE Transactions on Computational Imaging. 7:849-863, 2021.
Matthieu Terris, Audrey Repetti, rean-Christophe Pesquet, and Yves Wiaux. Building firmly nonex- pansive convolutional neural networks. ICASSP 2020 2020 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), pages 8658 8662, 2020.
Singanallur V Venkatakrishnan, Charles A Bouman, and Brendt Wohlberg. Plug-and-pun priors for model based reconstruction. In 2013 IEEE global conference on signal and information processing, pages 945-948. IEEE, 2013.
Luhuan Wu, Brian L. Trippe. Christian A Naesseth, John Patrick Cunningham, and David Blei. Practical and asymptotically exact conditional sampling in diffusion models. In Thirty-seventh Conference on Neural Information Processing Systems, 2023.
Xiaojian Xu, Yu Sun, Jiaming Liu, Brendt Wohlberg, and Ulugbek Kamilov. Provable convergence of plug-and-pla priors with mmse denoisers. IEEE Signal Processing Letters, 11-12800-244 2020.
Kai Zhang, Wangmeng Zuo, Shuhang Gu, and Lei Zhang. Learning deep emn denoiser prior for image restoration. In IEEE Conference Computer Vision and Pattern Recognition, pages 3929-3938. 2017.
Kai Zhang, Yawei Li, Wangmeng Zuo, Lei Zhang, Luc Van Gool, and Radu Timofte, Plug-and-play image restoration with deep denoiser prior. IEEE Transactions on Pattern Analysis and Machine Intelligence, 44(10):630-6336 2021.
Yasi Zhang, Peiyu Yu, Yaxuan Zhu, Yingshan Chang, Feng Gao. Ying Nian Wu, and Oscar Leong. Flow priors for linear inverse problems via iterative corrupted trajectory matching. arXiv preprint a405v: 188/6, 2024.
Yuanzhi Zhu, Kai Zhang, Jingyun Liang, Jiezhang Cao, Bihan Wen, Radu Timofte, and Luc Van Gool. Denoising diffusion models for plug-and ana- image restoration. In IEEE Conference on Computer Vision and Pattern Recognition Workshops (NTIRE), 2023.
Daniel Zoran and Yair Weiss. From learning models of natural image patches whole image restoration. In 201 International Conference on Computer Vision, pages 479-486, 2011. doi: 10.1109/CCv.2011..66
13

---

NeurIPS Paper Checklist
1. Claims
Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Question:
Answer: [Yes]
Justification:
An informal version of our main results is provided in the introduction.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No
NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings.
It is fine to include aspirational goals as motivation as long as is clear that these goals attained by the paper.
2. Limitations
Question:
Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification:
The limitations of each of our assumptions are discussed after stating these
assumptions. Additional limitations are also discussed in the discussion section.
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 encouraged to create separate "I "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 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 few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution is low images are taken in low lighting. Or 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, worse outcome might 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 impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory assumptions and proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
14

---

Answer: [Yes]
Justification: all assumptions are clearly stated and theoretical results are proved in the appendix.
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 cross- referenced.
All assumptions should be clearly stated or referenced the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but they appear in the supplemental material, the authors are encouraged provide 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.
4. Experimental result reproducibility
Question: Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that affects the main claims and/or conclusions of the paper (regardless whether the code and data are provided or not)?
Answer: [Yes]
Justification:
experimental details are given in Section 5.
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 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, descríbing the architecture fully might suffice, or if the contribution is specific model and empirical evaluation, it may be necessary to either make possible for others to replicate the model with the same dataset, 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.
While NeurIPS does not require releasing code, the conference does require submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) the contribution is primarily a new algorithm, the paper should make clear how reproduce that algorithm. (b) If the contribution is primarily new model architecture, the paper should describe the architecture clearly and fully.
(c) the contribution is new model (e-g., a large language model), then there should either be way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome 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 possible for other researchers to have some path to reproducing verifying the results.
5. Open access to data and code
15

---

Does the paper provide open access to the data and code, with sufficient instruc- tions faithfully reproduce the main experimental results, as described in supplemental material?
Question:
Answer: [Yes]
Justification: code provided in the supplementary material and will be publicly released after publication.
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/
cetminisssooooooooiioioooooooioi oeeiionooloooooo oolicy) for more details.
While encourage the release of code and data, understand that this might not be possible, so "No" an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurlPS code and data submission guidelines (https: / //nips. cc/ wwc//ullc///aeeeccceeeceecceelee http//uudee//eeeee Ceaees/CooSSubmieec cessinnoolccyy for more details.
The authors should provide instructions data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible in supplemental material (appended to the paper) recommended, but including URLs to data and code permitted.
6. Experimental setting/details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: see 4.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core of the paper 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, as supplemental material.
7. Experiment statistical significance
Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Question:
Answer: [No]
Justification: the impact of randomness very small 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, confi-
dence intervals, statistical significance tests, 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 with given experimental conditions).
16

---

The method for calculating the error bars should be explained (closed form formula, call 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 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 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful figures symmetric error bars that would yield results that are out of range (e.g. negative rates).
show in tables or
If error bars are reported in tables or plots, The authors should explain the text how they were calculated and reference the corresponding figures tables in the text.
8. Experiments compute resources
For each experiment, does the paper provide sufficient information on the puter resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?
Question:
Answer: [No]
Justification: the experiments are small and run on personal laptop.
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 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 failed experiments that didn't 't make it into the paper).
9. Code of ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurlPS Code of Ethics https pertps://neete http///eeuripp.eee http://uurrsssccppuiiieeeeeeeeee htcp/pulii/Ethicsuuiiiiiiiiieiei
Answer: [Yes]
Justification:
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 cireumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., if there is special consid- cration 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: [NA]
Justification: no societal impact, as it is this paper remains theoretical with no impact on broader society
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 why the paper does not address societal impact.
17

---

Examples of negative societal impacts include potential malicious unintended uses (e.g., disinformation, generating fake profiles, surveillance), faireess 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 direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that 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 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 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 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 monitor how system learns from feedback over time, improving the efficiency and accessibility of ML).
Safeguards
Question: Does the paper deseribe 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: the paper poses no such risks
Guidelines:
The answer NA means that the paper poses no such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Question:
Answer: [NA]
Justification: the paper does not use existing assets
Guidelines:
The answer NA means that the paper does use existing assets. The authors should cite the original paper that produced the code package dataset.
The authors should state which version of the asset used and, possible, include URL.
The name of the license (e.g., CC-BY 4.0) should be included for each asset.
For scraped data from a particular source (e.g., website). the copyright and terms of service of that source should be provided.
If assets are released, the license, copyright information, and terms use in package should be provided. For popular datasets, paperswi wwwode.hcoe...c. com datasets has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
18

---

For existing datasets that are re-packaged, both the original license and the license of the derived asset (if has changed) should be provided.
If this information is not available online, the authors are encouraged to reach out to the asset's s creators.
13. New assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
Justification: the paper does not release new assets
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details of the dataset/code/moddl as part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset used.
At submission time, remember to anonymize your assets (if applicable). You can either create anonymized URL or include anonymized zip file.
14. Crowdsourcing and research with human subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well as details about compensation (if any)?
Answer: [NA]
Justification: the paper does not involve crowdsourcing nor research with human subjects
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included in the main paper.
According the NeurIPS Code of Bthies, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional review board (IRB) approvals or equivalent for research with human subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or equivalent approval/review based the requirements of your country or institution) were obtained?
Answer: [NA]
Justification: the paper does not involve crowdsourcing nor research with human subjects
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Depending on the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurlPS Code of Ethics and the guidelines for their institution.
For initial submissions, not include any information that would break anonymity (if applicable), such the institution conducting the review.
19

---

16. Declaration of LLM usage
Question: Does the paper describe the usage of LLMs if it is an important, original, or non-standard component of the core methods in this research? Note that if the LLM is used only for writing, editing, or formatting purposes and does not impact the core methodology. scientific rigorousness, originality of the research, declaration is not required.
Answer: [NA]
Justification: the core method development in this research does not involve LLMs as any important, original, or non-standard components.
Guidelines:
The answer NA means that the core method development in this research does not involve LLMs as any important, original, or non-standard components.
Please refer to our LLM policy (https: / /ttp///ersp////// theurips.cc/Confeeee /eeeneee/////e/////////// 20LM/LLM for what should or should be described.
20

---

Organisation of the appendix.
1. In Section A, provide the proofs of Proposition 2 and Theorems and 2. We also show that Theorem is tight in the case Gaussian priors, and extend to priors supported on a low-dimensional affine subspace, and provide convergence guarantee in the more realistic setting where we do not have access to the true score function but only to an approximation.
2. In Section B, we provide several lemmas which enable control + zyg
A Proofs of Proposition 2 and Theorems and 2
A.1 Preliminary results
We start by the following proposition establishing that - In Pa is convex and 2 -mooth..
Proposition 4. Fix a > O. Under Assumption + In pa (xr) is conver with a Hessian satisfying:
in 0² [Ia
σ² Var(e|X de = 2)] = 02 La
Proof. The convexity of - In (x) follows directly by the classical fact that convolution of log-concave densities with Gaussian is still log-concave (see [Saumard and Wellner, 2014, Proposition 3.5)). The fact that the Hessian is upper-bounded by 3t2 La is a direct consequence of an identity which can be seen as a "second order Tweedie formula" (e.g. Lemma A. in Gribonval [2011] or Lee and Vázquez [2003] equation 5.8.):
1 La Var(s|X to = =2]
Im
Id.
where denotes standard d-dimensional Gaussian random variable (8 NN Ia)) and the matrix inequality is due to the positiveness of the covariance matrix. For completeness we give the proof of the second order Tweedie identity. From the standard Tweedie identity (see, e.g. Efron [201 1J) have that:
E oE =
-Vlnpo(2)
0²
(2 x) 2)p(z(2)dz
0²
Rd
Φ0 x: =11))(()) dz,
x))
JRa 2 1/- x''l x''l)p(z'ddd"
where de (2) exp(-
2 Notice that (2(2)= = (z) We can now compute the Hessian
207 of Inpo. letting X, X
+ σE:
(2 02 12002ppp ((p)((() 2) 02 / (z x)p(z/z2dd) 82)
Ia
σ²
Rd
0² [dd
Ia
(E((X, - x) = 22-EE[X - X X = 2)002))
0²
σ²
1
la
Var(e|Xo = =)))
which concludes the proof.
Now, we recall and prove Proposition 2, which is direct consequence of Proposition 4. Proposition For any a > 0, the function F is T Tmo-shh with
La 1+ o2
Proof. The result directly follows from Proposition which implies that In Pa is 1/00 -smoohh,,, so that F, is Lo -smooth with Lg + z.
D
21

---

A.2 Analysis of the MMSE Averaging iterates
We start by recalling our main result Theorem 1, which provides a convergence rate towards the proximal operator of the MMSE Averaging recursion.
Theorem 1 (Convergence to the Proximal operator). Under Assumptions I and 2, lei prox-r in (v) denote the unique solution of the Proximal Objective problem. Then, the MMSE Averaging iterates with parameters a = 1/(k + 2), of T/(k + 1) and initialised at 20 satisfy:
(In + [1y - - tu +
|| - prox lup (1)l -
+
Proof. From Proposition 3. are guaranteed that Fo, is strongly convex and smooth with
22 To avoid heavy notations, we denote Fop prox -n
Poz 1, Los =
No: k +
ok
(v) arg min Fo, well as x² prox_ (n) arg min F, note that these quantities are well defined and unique by the strong convexity of F... and F.
Recall that due Proposition one step of the MMSE Averaging recursion can be seen one step of gradient descent on Fr, with stepsize d ボ+2 which exactly corresponds to og = 1/ T.o.. Hence, at iteration k a standard convex optimisation result (see Theorem 2.1.15 in Nesterov [2013])
guarantees the contraction:
|| I0e S
1
Hov
1/2
|| - xta
+
Loas
Kop -1
1/2
1|
Kox +
1
1/2
|| x=t11 II
(1)
We now use the triangle inequality to write:
1: + 1 1/2
| 2-a
k+33
(1
-'r_, + ||=02_1 1)).
(2)
And we clearly see that we need to be able to control the regularity of a + xg This is done in Proposition 9, where we show that x; is Lip schitz in o²:
| xt,ll2 s C(of - -
for 02 Si S /T and where C = +| yll + 7M Vd. Since o% S VT, we can use this bound and insert it in inequality (2) to get:
+ 1 1/2
|
(1 -0 || (oh-1- C).
+
It remains unroll the inequality until h: and using the fact that TT
it1 1+3
(++1)+++)
(K+2)(443)*
V6
(j 1)(j +2)
|| 21,11 1
|0
xx || +
- -0)11-03)c.
/ (k + 2) + 3)
(Ir + 2)( + 3)
And from inequality (1) we have that ||
- 2
5
| xt,11 ||. Since In = y/ this leads to:
V2
(i+1)(ii 22)
| 11 V(k (k
11
xx +
Σ
(kt 2k+ +3) 103-1-00 03)C.
33)
j=1
22

---

Now since of k+I. we have that 103-1-07) = j(i+ii hence for k 2 1:
v2
(j+ 1)(j +2) TC
|0 - Foell
11y/
23ol +
V(k + 2) +3)
(k + 2)(k + +3) 505 +1)
v2
TC Σ j+2
1: 2 y x6ull
k+ j(
And we can simply bound:
k j+2
k 1Ú+1)
k
Σ
f 5 j+1
S 1 + Σ
1 + In(k),
J=1
1=1
j=1
Therefore
V2
(2 + In(()))FC
| - Towll S
F +
k ++ 2
h: +
Now using the triangular inequality Proposition 9 which bounds |
|| x²1 || v ||
Iza + ||
x11 and using
x² || S C
we get that:
v2 I1y x
(2 In()))
TC
|k++ 2011
k +
k+2
h &+1
And using the triangular incquality again:
|yy - x0
S | x + 1x - -00ll
S y x || + C
x || + TC,
where the second inequality is due to Proposition 9. Therefore:
/2 ly
x=1
(In + v2 T TC,
[ - x²1
k: +
+ 1
v2 11 x=1l ||
(Ink) + 4 TC.
k k+1
+ 1
Plugging the definition of = 1 || - y|l + TTM Vd we can finally write:
(In k) +
| - x=1| <
| - all +V
+ 1
which concludes the proof.
This next proposition proves the tightness of Theorem 1 (up to constants and the log-term) in the case of Gaussian prior. Here we assume that p the density of d-dimensional Gaussian M(,, E), with Σ a positive definite matrix. Without loss of generality, we can assume that the Gaussian is centered: ie.. 0.
Proposition 5 (Exact convergence rate for Gaussian priors.). Under the assumption that the prior p is a d-dimensional centered Gaussian (0, E), then we have that the MMSE Averaging recursion with ag 1/ + 2), of = T T/kk + 1) initialised at Dy = y satisfies the identity:
y prox -mp(y)) (v)
K prox p (mp(v)
k ++1
Proof. In this setting, the negative log-prior - Inp is a quadratic with Hessian H/ = 1-11,am and F is quadratic:
F(x) 211yy
z=|²
T
-z
23

---

Its minimiser is given by:
x = prpx-rlup()) = (T + -
And since Po ~ -0 la). the smoothed objective writes:
For(2) 2 1 - =112
- (2
and its gradient is:
V (x) = y 7(Σ oa - xx.
We now prove the result by induction. For 0, we have zo and the base case trivially holds. Inductive step: The inductive hypothesis provides that:
Th k: x*).
Using the identity x² (I + + 1ーー1ー11 1)- ys have:
y リ-x² -Σ-
Ib x²
E E-12"
k+
Then,
(E+ oTaa xxk ( K Ta) -1
(I k+1 --
E xx²
y x
so that:
+ (y - x²) = y - x
V xx y (y Now from Proposition 1, the update writes:
y +
(y - x²)
h +
k +
IK = Tb
VF (IK)
k+2
= =x" y- x²
-xx
k+1 1
(k 1)(k 22
= x + (y x")
2+2
This completes the inductive step, and hence the proof.
A.3 Extension when using approximate scores
In practice, when using trained denoiser, we do not have access to the exact score V In Po, but only an approximation 9o V In Po. In this more realistic case, the update rule becomes:
Fk+1 Fg as (r* y - -To. (Ix))
where we use the approximation go instead of the true score V Inp This recursion rewrites
Fk+1 26 (u( (#x) + T
(Noisy recursion)
where Ex In In (#x) For (Ixx) denotes the approximation error at step k. Assuming that these errors are uniformly bounded along the trajectory, i.e. ||Ex]] S &. we can show that the iterates converge to point at distance (E) from the true proximal point, with the same rate as in Theorem 1.
Proposition (Convergence with approximate scores). Under Assumptions and let prox In (v) denote the unique solution of the Proximal Objective problem. If the score approximation errors satisfy [|Ex]] < for all k. then the Noisy recursion iterates with parameters a 1/(k + 2). of = (/ 1) and initialised at Do y satisfy:
|| prox in (v)]
(In (Ink) k) + 7
[lly prox mnx +2 MV + V
k 1+1
24

---

Proof. To avoid heavy notations, we denote 2y = prox- In Pog x² prox-
(v) = arg min Fox well as arg min F, note that these quantities well defined and unique by the strong convexity of Fax and F.
Let ēx+1 z4 0g F (Ix) be the noiseless step. Following the exact same arguments as in the proof of Theorem 1, we get the single step contraction
k
1/2
| 2
10 x ll.
k 3
Since 2k+1 ēk+1 T& the triangle inequality leads to:
| 2t,
h: + 1/2 k k+3/ 3
I Izaull
+ k+ 11&xl1.
T
(3)
Next, in the noiscless case, decompose as:
|x - x |1 1 |:00 - - || + 1
For S [ Soa 1 || C(o2-1--00 -
with C= = * ||2 yll + MM v and 1 03 k[k+1]
Plugging this into the previous inequality
and unrolling from j = 1 to kt gives
v6
(j + +)
| - -2
||1 - xza +
(k + 2)(k -
V(k + 2)(k: +3)
j=1
( 2)(+
||E/1|-
(k + 2)( + 3) j 2
From inequality (3) with ks have that | Hop get:
v3 IT 2t,1l ll + 111l. Since Tp y.
V2
(j+ 1)(j +2)
|0 11
y 20 Σ
(k+ + ()(k+)) -3) - o3)
(k + 2)(k +3)
(j (j+ +3)
|6||--
(k: 2) +3) j+2
The second sum is bounded exactly in the noiseless case:
(++ +)(++++) -C
(2 (2+In
(k+ 2) + 3) 10++1)
k+2
For the noise sum, using ||ر S & and V
(++2)++ [K+2)6++] J
13/2
we obtain
(kK+2)(k+3)
k
(k ( 2) 3) j+2 |l <
2)j 3) т
/3/2(k + 1)TE
$ / Vē...
V(k + 2)
Putting things together, exactly as in the proof of Theorem 1, we obtain for all 2 1:
k) +
(Ia** vll + ++++ Vžrs.
| x11
1: +
Thus, the iterates converge to an O neighbourhood of prox (1) with the same ((1kk) rate as in the noiseless case.
25

---

A.4 Extension to distributions supported on affine subspaces of R
In this subsection prove that Theorem I can naturally be extended the case where the prior distribution is supported on an affine subspace of dimension r < d. in which case the dimension d which appears in the upperbound reduces the effective dimension r. Formally, we assume that the clean images r are drawn from a probability distribution f on F satisfying the following:
Assumption 4. There exists affine subspace S R of dimension S d such that the probability distribution il E P (RR)) satisfies:
μ is supported
S: p(R \S = 0. Moreover, the restriction of μ to S admits density p + IR+ with respect to the r-dimensional Lebesgue measure on S. By abuse of notation, extend to Rd by setting p(z) 0 for E Rd S.
p(x) > for all E S.
p is log-concave.
denote the Gaussian kernel on Rd of variance o², now let Co = (2mr2)1/2 such that Se (x(x) Cg. The smoothed density function Pa IRd RR,, then writes, for all 2 Rd:
Let (x:) = exp ||=]|(2)
20
(=(2) 1 / 0,(2 (z(2=2)))0 =z)
Cd
= Ch. p(x) lz x) dæ.
For E IRd, let denote the orthogonal projection of = on S. Using orthogonality, notice that:
Po (z) (0(2-21) 1
[
C J p(z))o r)dz.
Cgt
Therefore, for z S, letting Po (2)
d Co
Is (( x)dr denote the convolution of p with the
Gaussian kernel over S. we get that
po 12- 21112
(o InC.
20²
And importantly:
(s) - -s
where the V and As denote the intrinsic gradients and Laplacians on S.
Therefore using Lemma 5 for Po we have the following upper bound: In
sup zERd
| (n))
sup 2LES
$ Vr sup
V Inp(2)))]-
fEES carries through, with the ambient dimension
From this point onward, the proof of Theorem replaced by the effective dimension r.
A.5 Analysis of the approximate PGD Algorithm
We now restate and prove the convergence of the approximate PGD algorithm towards the MAP estimator. The following is a restatement of Theorem with explicit constants.
Theorem 3 (Convergence towards the MAP estimator with explicit bounds). For S XLy and number of steps in the inner loop which increases as kn [c- [e-n!+ty n! for e, n > 0, the approximate proximal gradient descent iterates (i(n)) nn from Algorithm satisfy:
Σ [2)(1) - -((MAAp) S
1
L(11ee
[Iv I + Σ || + 2R,
Σ |E:
2rn
1=1
1=1
[m((n)) - x(()|| < (1 n) (n()))(())))) + In(c) +
Ry,er
c. nl
26

---

where x(n) = prox_r In tn (8(1-1) TXV fV corresponds to the true proximal mapping, and where the quantities Bm Li=1 ||Eil| and Ei=1 || are explicitly upper bounded in Lemma 1. For, e.g., 1) = 1 and c 10, the bounds become:
1
[xx()) - J 2 =k
300 - 11 -MMP 12 + 600 (FAIIV(((xxx))))+++++
2 (n)(n) + 10 (6. 11 =MAPll + 12- (F2AI//VV((AM)))) + ++
| - |1 x
n²
Proof. For S Xt.'' the classic inequality after step of the true proximal descent x r(
prox-rlap(E" (E(n) -TVV provides that (see, e.g. equation 3.6 in Beck and Teboulle [2009]):
((x(n+2)) J - (1 MMPP 1 - | (n= M
1
(4)
27
Now for n 2 1, let En = =((n) ((() correspond to approximation error which can be quantified
using Theorem 1. Letting J : (dMM)) for 2 1, inequality (4) can be expanded as:
I(zlat1)) - J - 2r
(a(ln) TMAPI 2| || xMall + ||
1|2 + 2(z(()) ( aMAP!)
( (la(m1) *MAP!! 112 - ||r(( xMAPl1 ||en 112 + 2/|en || -ahAAPl)
-
27
27 (Ix(n) ZMAPP 21 - [|x 1|=((tt) 2MaPl1 || 11 + 22,,cllEnnl)
1
where the second inequality is due to the Cauchy-Sehrve Schwarz inequality, and the bound [(x(m)) - -MPP 1 B. is due to Lemma 1. It remains to sum this inequality from i = 1 to and add inequality
with n = 0 to get:
Σ (J(( JJ) Jn) S
1 (1
-MMAP!' -1)x0 - TMAP +
Σ
||::11 + 2.. Σ leill)
2r
2m1
$ 2r
1 [1 MMAP
112 +
1|0:112 2 2
Σ
20
where the second inequality is due to Lemma Diving by n leads to the first result. The second for which the upper bound is given in Lemma 1. O
comes from the fact that ||=,||| ||
x(0)| ||
The following lemma provides bound on this approximation error at each step, along with bounds on other useful quantities.
Lemma 1. For S xt., and a number of steps the inner loop which increases as kn [c-n n= for c, > O. let (f(m)), denote the approximate proximal gradient descent iterates from Algorithm and let En = {(() I(n) denote the approximation error iteration n, where r(a) = proxo_ In (8(m TXV -XV is the true proximal point. Then it holds that:
(1 (1+1) 1) In(n) + In(()) +
[|x(10) zMArl & Ryg.
||
Ryges
C. 1-
|enl| S Sn.e* Rynes
||ell2 S Tye R2
n=1
n=1
where
Ry By, +
Byy = exp(2S, 11 xMAPl + Sy. (F\I/V((rrAAP)) + + P2Mv))
1+1 (1 17 (In(c) + 7))
Sn
2(In(c) + 7)2 4(1 + 7)2
Th = e² (2n + 1)3
1+
2n +
27

---

For, e.g., 17 c= = 10, these quantities simply become:
R Bij, = 60 lly MiPll 120 (F>IVV(=MAP))) +
Sne R Tn
Proof. From inequality (4), for n 2 we have that:
||x( MMPll - || 1) zMAPll
(5)
|| - 11(-1))1 + ||| - MPPll
1)1
| + ||, #MsrPl-
Furthermore, from Theorem I, since e. 2 kn S c. n2, we get for n 2 1:
(In kn) + [|
||s/||| -1 (1) s
-XVS(3((-1)) xz (0) +
kn + 1
(1 + n) In(n) + In(e) +7 [| - -
- - )) V.d.
e.n 2+
(6)
Now, we use the triangle inequality to write:
|| - (Ia - TAVf) +
<|1x:(00 xAAPP + IEMA (Ia (VV
(7)
+ 11(aa- - (ZM) - (la 1l4-7A (=
Now, since "MAP satisfies the fixed point property TMAP prox from the definition of ((1), we can write:
and
(n) IaAll ||prox_ In ([V))
(La
S [I( - -V- (AV=
where the inequality due the non-expansiveness of the proximal operator. Inequality 7 then becomes
[xx( -La- ((
21a - (1f]=-)A (La - (1--AA((( FAIlV
$ 2 s(n-1)/| + +>//V
where the second inequality is because La - -Af is Lipschitz for - 1/(AL)). Therefore, injecting this bound in the inequality 6, we get for 2 1:
(1 77) In(n) + In(c) [2| (n-))
|| 1
xMAPP +
+
c.n
(1 + In(n) + In()) +7 [2]len-//| |l + 211x(1-11) zMAll + = V )l + -2M
c.1.+o
(8)
where the second inequality still holds for n 1 with the convention Eo and xo f0 = y. Now adding the inequality [|x(")) *MArll S ||en-1]l i|l (R-))= *MAPll from inequality (5) to the above inequality 8, and letting w,, ||e,.]| ||x(m) #MMpl for 2 0, we get the following recursive inequality for n 2 1:
Wy, S (1 + 20.)wn-1 CnA,
where
Cn (1+ In(n) + In(c) +7 A= +rMva, wo Ily IMAPl|-
It now remains to unroll the recursive inequality on w.. which is done in the auxiliary Lemma 2 to obtain:
Wn S exp(2S,, (wo + ASne),
28

---

where
Sn -=1 (1 (1+1) 17- (In(c) + 7)),
cr2
Putting things together we get the following uniform bound on w.n
Wn I Byg = exp(2S,.e) [|v MAPll + Sq,e" (MAIV +
From the definition of wi, = || || + ||_,()) - MMalll, we trivially get that ||x((n) xmMAPl I By.c. and now from inequality (8) get, for n 2 1:
|| + 17) In(n) In(c) +
[2Br.c + T T +
C. c.n
Letting Ry == 2B n. + A [[ + + Bye we prove the two first inequalities of the statement.
Now to bound Σn=1 ||ا ll we simply reuse the bound obtained on Σi C S S in the proof of
Lemma 2 obtain:
| ll S Sn.c Ry...
n=1
Finally for Σn= || upperbound:
(1+ 1) In(n) + In(() +
2
2(1+
00
In (n)
2(In(c + 7) 0c
S
c. 1
c²
n2 (1(+))
c²
2((+))))))))
n=1
n=1
We now bound the two series using integrals:
In (k)
In (x) dr 2
[++))
(1()))))
(2m + 1)31
2014m dr = 1 +
(1++g
2n +
==1
Putting everything together, we obtain the bound:
Σ ||
4(1 1/) (2 + 1)3
2(In(c) + +7)
(1
2n + ++)))
1
1
RGy.:
6²
which concludes the proof.
Lemma 2. The recursive inequality
w/ (1, wn-1 C.A, where Cn (1+ 7) [n(n))) + In(e) +7
c-nltn
unrolls as:
Wn S exp(22,) (wo ASp,),
where
Sn. = 1 1+ (1+1- 7)).
ca2
Proof. We iteratively apply the inequality to obtain:
S II (1 + 2C2) + C. II (1 2C,), i=1
1=1
j=i+1
with the convention that empty products are equal to 1.
We now bound the product IT (1 2C2) ) by using the inequality log(1 + x)) S x to get:
"
n
log (1 + 202) Σ
log(1 (1 +0,) < Σ 20,,
j=1
29

---

hence,
a II 20,) (0,) - exp
2 2Σ C3
==1
To bound the sum Ej=1 C2. we split the numerator:
Σ
C
+++ 17
In(e)
+ 80
Σ
1
jztn
j=1
==1
We use the known bounds:
Σ L
s1+
000
1
dt
1+ .
x
Int t 11+
dt
1 1
11+1
j=2
which gives:
Σ j=1
Cj S
1+1
(In(c) (l) (1+ ;)
cr2
1 1+ on?
(1 + n (In(c) + +7) =: Sn
Then we have:
a
n Σc. II (1+20;) S Sn.c exp(2S,.)))
(1 S exp (2S,.o),
j=1
i=1
j=++1
Plugging these into the expression for wn yields the final bound:
Wn S cxp(22,,c) (wo + ASn.)
30

---

B Controlling a + F
The goal of this appendix is show that the minimiser z is Lipschitz-continncn with respect to o². To establish this, we need to control how the objective function F, evolves as a changes. natural way to approach this is through PDE perspective, since the smoothed density P, satisfies the heat equation. This connection allows us to describe how Pr. its logarithm, and its gradient (i.e., the score function) evolve with respect to o².
Throughout this appendix, we use the following notation for differential operators acting on functions f: I R:
f f denotes the gradient of f. a vector RR,
V f denotes the Hessian of f. a x d matrix of second-order partial derivatives,
V denotes the third-order derivative tensor f, a rank-3 tensor in IR [ddxdxd
Af tr(V2f) f) denotes the Laplacian of
The first lemma provides several PDEs satisfied by Pa. In Po. and the score function V Pa
Lemma 3. Let p(x) be probability density and denote by Po (1:) its convolution with isoropic centered Gaussian of variance o². For >> 0, it holds that (x) > 0 for all x E Bd and follows the heat equation:
Opr 1
Apo
0g²
Moreover, In follows the following partial differential equation:
0ln Po
1 (A + I[Vinpoll')
Og
2
Taking the gradient in the previous equation we get that the score functions follow:
0V (Vin)a(r
( + 2 In (z))
Og²
Proof. Standard results (see, e.g., [Evans, 2022, Chapter 2]) guarantee that (a, + po (x()) is O on R Fd and satisfies the heat equation:
Op. Ap
dō² 2
By differentiating lnp w.r.t. o² and using the above, we directly have:
Ol Po
1 Ape
dG² 2 Po
To get the PDE satisfied by In Po notice that:
Apo
A Anp
Po
Using both equation above directly yields:
Oln Pa
Anpp +
a0²
Taking the gradient in the above identity leads to the last partial differential equation of the Lemma and concludes the proof.
This next lemma justifies the use of smoothed gradient descent by confirming that, as the smoothing parameter - 0, the minimisers of the smoothed objectives original (non-smoothed) objective F. In other words, the limit of the smoothed minimisers coincides with the proximal point ultimately aim to recover.
F, converge to the minimiser of the
31

---

Lemma 4. Recall that we define
1 F(a)= = l1lly x||² - and (x(x) 2 by x|| Po (2))
Recall that prox-rlapy)) = arg min F(x) and that prox-- lu (0) arg min Fo(x). It holds =ER2
sERd
that
pupp- (v)
+ 11ox--11p10 o-+0
Proof. Let K be compact set, since is continuous and p(x) > 0 on K (Assumption 1), we have that there exists such that re p(x) Z a. Now since p Lipschitz continuous on K,
Theorem in Nesterov and Spokoiny [2017] ensures that supzek (P=(r) - p(z)l - 0. Therefore
o-30
for small enough inf=e (x) 2 a/2 and from standard inequalities on the logarithm:
[P((x) - p(x)) min(po(2)) (I), S (x()) p(z)).
[l(p=(<)) - In(p(x))] <
Therefore suPeck In(po(x)) In(p(z)))
- 0 on all compact sets K, and trivially:
o-0
sup (x) F(x) o-2 0.
=eK
To ease notations, let let be the minimiser of F& and x" that of F. Note that such minimisers exist and are unique since F6 and F' are strongly convex by Proposition 4. Consider the values Fx(x). By optimality of r; we know that Fg (x²) S (x²).. Moreover, since F, + F uniformly on compact sets, have Fx(x") + F(**), in particular, the sequence (Fo (xt)) is uniformly bounded above:
(x S Fo(x)) F(2")+
for a small enough. Now, assume that || + oo along some sequence. Since the functions Fa are all 1-strongly convex, they can all be lower bounded by the same quadratic and we would have Fo + oo, contradicting the bound above. Therefore, the sequence (x"))=E(Urr) is bounded, and thus contained in a fixed compact set K C TR...
Since F. + F' uniformly K, any cluster point Too of (xz)) satisfies
F(*x)) lim Fo Fx())) a-
< lim Fo (x*) F(x").
a=
Therefore, by uniqueness of the minimiser of F, must be that dio = x= so that x _ x²
a
The next proposition establishes the existence and smoothness of the solution path x as a function of a.
Proposition 7 (Existence of the smooth solution path). Recall that
1 2 1yy x1| In (x))..).
(x)
Denote by x; the minimiser of Fo for any o > Then 0² + na is continuously differentiable on (0, r] and satisfies the following ordinary differential equation:
dx d - - -
Proof. By smoothness of the solution the heat equation (see, e.g., [Evans, 2022, Chapter 2)), we have that ++ F, (x) is differentiable for any 0 and (o 2) + V F.F, (x) is jointly differentiable on R x Rd, Then, by Proposition 4, we have that the Hessian 2 F. (x) is invertible and satisfies: (x) = Ia- We can then apply the implicit function theorem, which guarantees the existence of unique solution path o² + =; to the implicit equation: VF, = that differentiable on (), By strong convexity of F.. this solution path coincides with the minimisers of F, for all > 0. The
ODE followed by a² - x is obtained by taking the derivative with respect o² of the identity
F., (x)) = 0.
32

---

Proposition 8 (Bound on the solutions). Let x, == arg zin (o(x), then for o² T it holds that
1 ||=1 yll 1ly prox-rlnp(y) (v)ll + 2 M
Proof. Let write 2 Proposition 7, we have that Hy follows the differential equation:
de do (note that the derivative is with respect to o² and not σ). From
f - - = = 0,2 Vmppo
0,>
(9)
where the last equality follows from Lemma Furthermore, recalling the optimality condition satisfied by 24 i... Vin (xt) = } (x² - y). if follows that:
1 Qa - 8) + B,,
(10)
2T
where the matrix Qo and vector B, are given by:
Q. == In (Ix)) La] V Impo 20
(11)
Bg = Inp. + = - Invip (3)))
(12)
Here, the matrix Qo is positive semi-definite since 2 In Po (z) is positive by Proposition 4. Now from eq. (10), we get:
b v|1²
(4 - -y)
2 do²
27 |1=;
Uli. + (Bo, y)
(Ba x2 y)
- BBllll> y|l.
From the upperbound lo s M /d which follows from Lemma 5, we directly have
that |B,l| || I F. vd. Injecting this bound in the above inequality and dividing both sides by | yll yields:
dll - vll S I Mvd. 2
2
Integrating of the above inequality from 0 o², using that lim,- x3 from Lemma 4, we get:
prox (1)
| yll lly prox_rlnp (v)ll 1
0 -MV
2
1
Il pno (v)ll + MVd.
2
where the last inequality is since we consider 0² ST. Proposition 9 (Lipschitz continuity of o² x²)), Let x = arg aig (x(x), then for of < of < 5, holds that:
| xx,, (of 03) [
prox + TMV
And taking 02 + 0 in the above inequality:
||=; prox r (x))l S 0² [ll prox (m)
33

---

Proof. Recall from Equation (9):
2/ -- (mpo + + IV (x) Vilp
Now, by Proposition 4, have that (npo = and spectral norm bound the inverse yields:
v (xpo -
and:
|||~v po (aap (z)p
Putting things together obtain that:
| <
V Iapo
+
(13)
V Impp (x
Mvd.
(14)
where the second inequality is due Lemma 5. Now recall that the optimality condition which define I is V Vl = (2 - y). Plugging this equality in the upperbound we get that:
| S yy x
5 M Mvd
2
prox in (v)l TMdd,
where the last inequality is due to Proposition 8. From here suffices to notice that, for σ 2 02 0:
| - x=,1l
4 do
| do
of
| Ill
(o]
o
prox-inupvv) + STMVil
which proves the first statement. The second follows from the fact that wt prox_.. (n) by Lemma 4.
This last result is the technical lemma this work. It establishes that the third derivative of the smoothed log-density In can be uniformly controleed-indeeendetty of a. This regularity bound essential for tracking how the minimisers x evolve as varies.
Lemma 5. For all 0, it holds that supzERd VAln S VdM.
We would like to emphasise again that the following proof is entirely based on the computations and insights that Filippo Santambrogio generously shared with us response to an email we sent asking for ideas on how to approach this result. The proof is technical and relies several surprising simplifications that Filippo identified.
Proof. To simplify notations, throughout the proof we let == o² and V(t, == (Ip,, = In pa (2)) correspond the convex potential associated to p.. The proof first relies on showing that V (t, æ)]l must be maximal for t = 0.
34

---

Establishing a parabolic inequality for V x)|l.
From Lemma 3, we have that the potential
V follows the following PDE:
1 (AV
0 V =
2 - |VV])
2
For i, [d], we let wigh Desu V. which therefore follows:
wewikk = 2 (Amigh Oiyk \V
Now let Wijk wiek. multiplying the previous equation by Wigh we get:
Oruigh = Wijk wiyh 1ip [
1
(Auigt (Awiyk) WisAOInxIVVIIII
2
1
$ Auijk (
2
Summing over i,j, and letting S(t, x) = [tt
Eigr Mijk, we have that:
a,s - 2 AS
Σ Gjk
It remains to control the last term in the inequality. Since /VI2 (OV), taking the third derivative with respect to j. le we get that:
Oijk |VVI2= =2 Σ a.V. OigiaV + +d.neV daV + d.a V 0y. + Oegi V daaV
2(VV, Vwigk) +2
Σ WyM daV + was 0.V + wast V.
Multiplying the equality by Wijk and summing over i, k we get:
Σ
V = (VV, VS) + 2 tk
Σ Waykelgaa VaV + wighWas VnV + Wiyxeee 2u V.
Ljke
However notice that from the convexity of V(a,)) :) for all 0, get that:
Σ
WaykWya VzV) 2 0,
jk
20
which implies that the function S(t,x) = (t, 12
satisfies the following parabolic incquality
0,5 < 7 As (VV.
(15)
Proving that S is maximal for - To prove that S must attain its maximum for 0, let us fix t and for & t.) we let stt, S(t, f, and (t, V(t t, correspond to the "reversed time" counterparts of S and V. Adapting Equation (15), the parabolic inequality satisfied by $ is:
A5 (VV. V 0,5 -
(16)
For t E (0, we now consider the following stochastic differential equation:
dX, VV X,)dt+ dB,,
(17)
initialised at X,=0 = 2p for some 2g E Rdd From Lemma 6, we are guaranteed the existence and uniqueness of strong solution to this stochastic differential equation over [0, t1]. We can then apply the Itô formula to Š(t, Xt):
1 A50, X,)dl as(1, X,) = 0,5(1, x,) di + (v5(1, X), + 1
0,,((, X,) dt (vš(t. Xx), (V,,, + A X.)dt+ + (3)(,x.)),
2 (V)(t,X,), dB,)
35

---

where the last inequality is due to the parabolic inequality on S from eq. (16). Now integrating from t 0 ti we obtain:
/ (1)(
S(4, S(0, X.-0) +
=(00,z0) + / (.),
Since the expectation of the stochastic integral is 0, and recalling that S(t,x) = S(t1 - t, we obtain:
B[S(0, X,, E xx,)) 2 s(0, ((1, vo)-
remains that sup, a00 < from Assumption 2 obtain that:
sup S(0, 2) (0, Xx,) 20), =ER*
Since this inequality holds for all zo E Rd and t 0, get that:
sup
S(t, = sup 5(0, Vt 0. =E
=E3"
Therefore, recalling that S(1,x) = V all 2 0:
| lu (z)))()))1 we finally have that for
sup =ERd
| S sup ||
SER
From V to /AA/. From the Shuchy-Schwartz inequality, one gets:
d
d (ugff) <d d Σ =
|
Ouisf
S d
1,=1
i.j.==1
which concludes the proof.
Lemma For a horizon time l > O, let -n denote the backward-time log-density defined over [0, Rd. Then for all initialisation X Rd, the stochastic differential equation defined in Equation (17) which we recall here:
dX. -VV X.)dt dd.,
has unique strong solution over [0,
Proof. From Proposition 4 we have for all x E RRd:
> -V 3 =
Therefore 21, = V(h - I, x) satisfies:
1
0 x v2V(1, 1 h カ- la-
This entails that for all e > 0, V is globally Lipschitz for E [0,4 (1 s):
[V(t,)) - vi x')l || x' ||.
which ensures the existence of unique strong solution over [0, t1 E] (see e.g. Theorem 5.2.1 in Oksendal [2013]) and hence over [0, t1). It remains to show that X, does not blow up t
36

---

Proving that X, is bounded over (0,t To do so, we consider the Lyapunov 2 || roll?, for which the Itō formula provides that:
1 d||x; zoll² = (dX,, -o) d dt
2
= ()) - - + dt (dB,, X₁ zo).
2
Now recall that for all t, the function æ + ((t, x) is convex (Proposition 4) and hence we have the inequality (VV (t, x'), x - x²) _ ((t, V x'), which leads to:
d =oll2 - (V
V1 xx) dt + dt (dB,, X; zp)-
Recalling the integral definition of Po Po (z) = f2 p(2)0. density function of variance x² inf Po (Io) Puii (zo) > 0 (since p is assumed strictly positive over R from Assumption 1).
2)d2, where do denotes gaussian have that sup, Po (x) S Pmax sup,, p(x) well as
Therefore
d zoll² I Cdt 2(X, Io, (Bt)
with C = 2 12(maax PPiiiaaa///iiiiii (rr)) + d. Now integrating from 0 to < t we obtain:
| xo - Ct + 2 / (Xe- Io. dBp)
S Ct M,
(18)
where 2x xo, dBp) a continuous-time martingale. Bounding M, over (0,
Taking the expectation in the last inequality get:
BlllX, xoll²] Ct Ct..
Now notice that due to the Itō isometry, we have that:
E[M 4E J. XX zoll ''t]'tt'] = 4 L. E[IIXr - zall 't s 4Ct.
We now apply Doob's 's martingale inequality to the process M²:
E[M?
4Ct;
P sup M² 2 A2)
<
r'st
A2
A
And since
sup M3 2 1²) = U sup M > A? ecc
n21 #ec
where the sequence of events are monotonically increasing, we obtain that:
4C
TP sup Mf A²)
lim P sup
M 2 12
s
n-000
A
Therefore lim P Mu 2 A²) = 0 which translates into:
4-000 tKh
P (sup < x xx)=1 =I tch
Due to inequality 18, this means that the trajectories (X, (X,(w))ee(ooa)) are bounded for almost all w. Therefore, due to the continuity of VV x) over R x Rd the path t + (t,, Xx(W)) bounded on [0, f1). Hence, for almost all w,
Xi(w) = 1.' ( + B.(w)
must admit limit when t + t Hence X; extends continuously to and X,, (w) still satisfies the integral form of the SDE. Hence a strong solution exists on the whole interval [0, f.. Unicity over [0, t1] follows from unicity over [0, t2) and taking the limit in t7.
37


Reviews:
[{"review_id": "TNk3qaarDo", "reviewer": "Reviewer_wDZG", "summary": "The authors examine a simple but practical algorithm aimed at solving the MAP optimization objective, $\\arg\\min_{x \\in \\mathbb{R}^d} \\bigl(\\lambda f(x)-\\ln p(x)\\bigr)$. In general, this optimization objective can be solved using proximal gradient descent, but typically, the proximal operator for $\\ln p(x)$ is intractable or unavailable. The authors first propose a simple recursion, called MMSE averaging, that estimates this proximal operator by running gradient descent on a sequence of proximal operators of decreasingly noised distributions. Given some assumptions on log concavity and higher-order derivative boundedness, the authors show this recursion tracks the true proximal operator value. The authors next substitute this recursive proximal calculation for MAP optimization and show $O\\left(\\frac{1}{k}\\right)$ convergence.", "questions": "- For learned denoisers (e.g. in diffusion models), could the authors comment on the advantage of using the MMSE averaging algorithm as opposed to directly using access to the parametrized score at a small noise level ($\\nabla \\log p(x, \\sigma=\\text{\\rm small})$)?\n\n- What is the sensitivity of the convergence guarantees to the specified noise schedule $\\sigma_k$ in Theorems 1/2?\n\n- Could the authors comment on the benefit of MMSE averaging relative to using a preconditioner with gradient descent when $\\log p(x)$ is tractable?\n\n- What heuristics are there for setting the iteration number in Algorithm 1 for practical implementation?\n\n- Small point, but should $\\alpha_k$ be set to $\\frac{k+1}{k+2}$ in Proposition 1 and elsewhere?", "limitations": "Limitations are clearly addressed.", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 4, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strengths:\n- The paper presents definitive theoretical contributions towards a highly relevant \"tilted\" MAP optimization problem and identifies a simple algorithm to approximate the optimization. \n\n- The explanation of the algorithm as using the noised proximal objectives for an \"annealed\" or \"smoothed\" series of gradient descent steps is quite interesting and provides insight towards the compatibility of traditional optimization techniques with pretrained diffusion models.\n\n- Proofs and presentation are very clear and easy to follow.\n\n\nWeaknesses:\n- If $\\log p(x)$ is intractable, typically convolving with Gaussian noise will also be intractable, so the algorithm seems to be practically relevant when trained denoisers/score functions are learned, but learned well at the noise levels required.", "was_revised": true, "final_justification": "- The paper is a strong theoretical contribution to a practically relevant problem, and rebuttals from the authors towards points/questions raised confirmed the original positive impression.", "rebuttal": "We wish to thank the reviewer for their comments and thoughtful engagement with our work. The questions raised in the review in are addressed below.\n - **Intractable $\\\\log p$.** We agree and would like to clarify that, in practice, we do not compute $\\\\log p\\_\\\\sigma$ directly. Instead, modern denoising models (*e.g.*, trained MMSE denoisers, score-based diffusion models, flow-based models) are explicitly trained to approximate the score $\\\\nabla \\\\log p\\_\\\\sigma$ for a range of noise levels. These models can be trained on randomly sampled noise levels, which enables them to interpolate effectively for intermediate values, including those required by our method. Hence, our analysis applies to practical implementations where this score is approximated via learned denoisers.\n - **Comparison with using $\\\\nabla \\\\log p(x, \\\\sigma = \\\\mathrm{small})$ directly.** Using a fixed small noise level to approximate the score of the original prior is a valid approach, but as discussed in the paper, it requires performing very small step sizes that can lead to arbitrarily slow convergence to the MAP of a *smoothed* density. In contrast, the scheme we analyze leverages large $\\\\sigma$ in early iterations and gradually anneals to smaller $\\\\sigma$ to reach the *true* proximal point with $O(1/k)$ rate (cf. Theorem 1 and paragraph right after, lines 202-218). Importantly, this scheme does not rely on any smoothness parameter about $\\\\log p$ (*e.g.*, Lipschitz constants), which are typically unknown in practice. There are of course limitations due to the assumptions made about $\\\\log p$, but, as far as we know, this is the first convergence result with known rate obtained for such optimization schemes with non-trivial (non-Gaussian) data distributions.\n\n - **Sensitivity to the noise schedule $\\\\sigma\\_k$ in Theorems 1/2.** Our analysis can indeed easily be adapted to a more general class of decreasing noise schedules. In particular, for sequences $\\\\sigma\\_k$ such that $\\\\sum\\_k \\\\sigma\\_k = \\\\infty$ and corresponding step sizes $\\\\gamma_k = 1 / (1 + \\\\tau / \\\\sigma_k^2)$, we obtain a general convergence rate of the form: $O \\\\Big( \\\\min(\\\\sigma\\_k, \\\\exp( - \\\\sum\\_{j = 1}^k \\\\sigma_j / \\\\tau ) \\\\Big)$ for the smoothed gradient descent recursion of Proposition 1. Thus, any $\\\\sigma\\_k$ decreasing at most as $\\\\tau / k$ ensures a $O(1/k)$ rate. The specific choice $\\\\sigma_k^2 = \\\\tau / (k+1)$ is made to emphasize the connection with the MMSE Averaging scheme $x\\_{k+1} = \\\\alpha\\_k \\mathrm{MMSE}\\_{\\\\sigma\\_k}(x\\_k) + (1 - \\\\alpha\\_k) y$, while other choices of noise schedules and averaging weights will lose this equivalent formulation. We will clarify this in the paper.\n\n - **Comparison to preconditioned GD when $\\\\log p$ is tractable.** If $\\\\log p$ and its gradient were tractable, it is very likely that better optimization schemes could be used, but it would fall beyond the typical use case where $p$ represents an unknown data distribution. That being said, theoretical convergence with a preconditioner $P$ would depend on the condition number $\\\\kappa_P = \\\\sup\\_x \\\\lambda\\_{\\\\mathrm{max}}(P \\\\nabla^2 \\\\log p(x)) / \\\\inf\\_x \\ \\\\lambda\\_{\\\\mathrm{min}}(P \\\\nabla^2 \\\\log p(x))$ with a rate of $O(\\\\kappa\\_P \\\\log(1/\\\\varepsilon))$. In the best cases this can lead to a considerable speed-up with $\\\\kappa_P$ close to $1$ but in practice for very badly conditioned functions $\\\\log p$ the matrix inversion (the optimal $P$ being $(\\\\nabla^2 \\\\log p(x))^{-1}$) could lead to severe numerical errors. In contrast, the MMSE averaging scheme does not require any knowledge of $\\\\log p$ and does not explicitly require any sort of matrix inversion, as long as the MMSE estimator can be computed/well approximated.\n\n - **Heuristics for the number of iterations in Algorithm 1.** A full empirical study is beyond the scope of our paper, but it is indeed an important future direction to investigate. Existing works such as *Delbracio & Milanfar (2023)* suggest that 5–20 inner denoising steps are often sufficient for high-quality results. We expect our method to behave similarly: the inner loop (MMSE Averaging) benefits from rapid convergence due to improved conditioning in early iterations. We will add a comment to the paper with guidance and references.\n\n - **$\\\\alpha\\_k$ Correction.** Thank you for catching this. Indeed, there is a typo: the correct choice used throughout the theoretical analysis is $\\alpha_k = \\\\frac{k+1}{k+2}$, not $\\frac{1}{k+2}$. We will correct this to be consistent in the paper."}, {"review_id": "SPcsPGXCPy", "reviewer": "Reviewer_LqTm", "summary": "This paper deals with an algorithm to compute efficiently the proximal operator of some minus log-density $- \\ln p$, by using a kind of multi-scale gradient descent that operates on some blurred versions $-\\ln p_\\sigma$ where $p_\\sigma$ is the convolution of $p$ with the density of the Gaussian $\\mathcal{N}(0,\\sigma^2 I)$.\nThe proposed method is formulated as an \"MMSE averaging\" algorithm, which iterates relaxed MMSE denoising with prior $p$.\nThe authors prove the convergence of this algorithm for a particular choice of step sizes, under the hypothesis that $p$ is positive log-concave and that the third derivative of $\\ln p$ is bounded.\nThey also give an example illustrating the benefit of this algorithm with a quadratic $-\\ln p$ corresponding to an ill-conditioned matrix $A$ (in this case, isotropic smoothing improves the conditioning of the matrix and thus the convergence of gradient descent).\nThe authors also integrates this algorithm in a modified loop (called Approx PGD) in order to target the maximum a posteriori (MAP) of an inverse problem with data-fidelity $f$ and prior $p$. The convergence of this second algorithm is proven under similar hypothese.", "questions": "- Could the authors clarify their contribution on Approx PGD compared to existing results about inexact proximal algorithms?\n\n- The proposed method looks quite similar to other algorithms based on simulated annealing. Can the authors compare the proposed algorithm to other ones based on simulated annealing?\n\n- Can the authors confirm that it is really mandatory for their contribution to consider a log-prior $\\ln p$, and not just any concave function? \n\n- Are there some interesting non-Gaussian priors on which the proposed results apply? For such a non-Gaussian case, I think that we do not expect the third derivative to be much less than the second one.", "limitations": "yes", "rating": 5, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 2, "originality": 3, "strengths_and_weaknesses": "# STRENGTHS\n\n- This paper contains interesting thoughts on the way we use MAP or MMSE denoisers in iterative algorithms (like PnP) to solve inverse problems. It makes also a nice connection with recent algorithms based on Cold Diffusion.\n\n- The proposed algorithm is well detailed and relatively easy to implement (except maybe that the number of iterations must increase super-linearly, and that the true MMSE denoiser is rarely tractable).\n\n- The paper is well written with a nice exposition of the main results.\n\n- This paper tackles a very important issue in the literature about PnP and Diffusion, which is to get convergence guarantees of algorithms that operates with a decreasing schedule of noise levels $\\sigma_k$. \n\n\n# WEAKNESSES\n\n- The idea of using an inner loop to approximate the proximal operator in PGD is not new, as shown by the introductive section of the cited paper (Schmidt, Roux,  Bach, 2011) (who quoted several other older references on this specific topic). \n\n- The convergence guarantee applies to an algorithm based on the true MMSE denoiser (for several noise levels $\\sigma_k$). In practice, we cannot access exactly these MMSE denoisers, so that another source of error comes at stake. Therefore, the practical benefit of the given convergence bounds seems unclear.\n\n- The convergence guarantees only concern the convex case, but in many applications the most interesting priors are non convex. It appears that the proposed algorithm helps to cope with ill-conditioning of the log prior. In its current form, the paper does not explain how it would be helpful to cope with issues related to non-convexity.\n\n- The only experiment proposed in Section 5 concerns a toy example, which relies on a very well known counter-example of badly conditioned fixed step size gradient descent on a quadratic functional. \n Since a large part of the introduction is devoted to solving inverse problems with plug-and-play or diffusion algorithms, one could have expected to see an experiment showing the benefit of the proposed algorithm in such a setting.\n If the main contribution of the paper is about convergence results on proximal computations, then the paper should give more context about other optimization approaches, and could give less details on the potential applications in imaging.\n\n- The paragraph of the introduction related to PnP methods omits the fact that some of the PnP methods (like the ones of (Cohen et al., 2021a) or (Hurault et al., 2022)) are already formulated with a denoiser that is designed to approximate the MMSE (and not the MAP). And these algorithms already provide convergence guarantees. It would be interesting to know how the proposed algorithm compares with these more direct methods (which do not have an inner loop). Also, it would be interesting to see how the proposed convergence analysis relates to the one given in\n\n> \"Provable Convergence of Plug-and-Play Priors with MMSE Denoisers\". Xu et al. IEEE Signal Processing Letters, 2020.\n\n\n## Minor Remarks and Typos\n\n- Saying that $- \\ln p$ is ill conditioned in the introduction seems unclear. (In my opinion, it is commonly accepted to talk about an \"ill-conditioned matrix\" but less common to talk of an \"ill-conditioned log-prior\").\n\n- It is a bit exaggerated to say that \"a new wave of approaches\" focus on viewing \"inverse problems as a sampling task\". This line of reasoning was already present in older papers, see for example \n\n> Louchet, Moisan. \"Posterior Expectation of the Total Variation model: Properties and Experiments\". SIAM IS, 2013.\n\n- Some of the theoretical results are stated without recalling the hypothese, which make their reading and understanding more difficult. Similarly, when they recall the \"classical result by Beck and Teboulle [2009]\", they could recall the required hypothesis.\n\n- About the comparison of convergence speed: when $\\epsilon \\to 0$, the complexity $O(\\log(1 / \\epsilon))$ is better than $O(1 / \\epsilon)$ so that the authors' comment is unclear. \n\n- The step size $\\gamma_k$ is not defined when it first appears.\n\n- $J^*$ is not introduced.\n\n- A few sentences are unclear:\n\n> they cannot guarantee that the denoiser is a proximal operator, let alone the proximal operator of the correct functional.\n\n> This assumption controls how skewed and “non-quadratic” the log-prior is", "was_revised": true, "final_justification": "After considering the authors' reponses to the reviewers' comments, I increase my score from 4 Borderline Accept to 5 Accept.\nThis paper brings an interesting contribution to the field of convergence analysis of signal/image restoration algorithms.\nThe proposed convergence result applies only with a log-concave prior, but seems nevertheless interesting compared to existing results.\nUnfortunately, due to this constraining hypothesis, this theoretical result does not apply to recent non log-concave priors that are widely used in the plug-and-play or diffusion literature. But it is likely that the proposed analysis could be extended to a more general setting in the future.", "rebuttal": "We wish to thank the reviewer for their comments and thoughtful engagement with our work. The questions raised in the review in are addressed below.\n\n - **Compare to existing inexact proximal algorithms.** We of course do not claim novelty in using an inner loop for approximate proximal computations, but we will make sure this is clear in the revised manuscript. As noted by the reviewer, this approach is classical, going back to *Schmidt et al. (2011)* and earlier cited in Section 4. As discussed in lines 299–302, our analysis differs in that the approximation guarantees for the underlying proximal solver are on the *iterates* (from Th. 1), not the objective values. This distinction slightly alters the analysis, particularly in controlling error propagation. Nevertheless, we do not wish to claim this as a major contribution but more as a natural consequence of how our approximation of the proximal operator can be integrated into standard PGD frameworks.\n\n - **Denoiser approximation error.** Indeed in practice we only have access to approximations of the MMSE denoiser. Our analysis can be extended to take into account this additional source of error: if each time the denoiser is applied, the score estimate is off by $\\\\varepsilon_k$ (i.e. we get $\\\\nabla \\\\log p\\_{\\\\sigma\\_k}(x\\_k) + \\\\varepsilon\\_k$), and if $|\\\\varepsilon\\_k| \\\\leq \\\\varepsilon$ is uniformly bounded, then we can prove convergence to a neighborhood of radius $O(\\\\varepsilon)$ around the true proximal point at a $1/k$ rate. We will include this extension as a comment in the main text, referring to the appendix for a formal statement.\n\n - **Dealing with non-convexity.** Our theoretical results rely on the convexity of the negative log prior. Handling the general non-convex setting is outside the scope of our work as it would certainly lead to several additional technical hurdles. One possible way forward is to note that convolution with a Gaussian intuitively \"*convexifies*\" a probability distribution to some extent. This could significantly help during optimization in practice, and could potentially be useful for theoretical analysis: there exists a minimum $k$ (dependent on the prior $p$) for which $p\\_{\\\\sigma\\_k}$ is log-concave. For example, *Stochastic Localization via Iterative Posterior Sampling, L. Grenioux et al., 2024* show in a sampling setting that this *convexification* effect can be used to help sample from non-convex distributions belonging to a restricted class. Another possible direction will be to exploit existing tools allowing to leverage convergence results on convex auxiliary functions, as done by Paquette et al. in *Catalyst Acceleration for Gradient-Based Non-Convex Optimization*. These directions are perfectly relevant, but beyond the scope of our paper. As noted in the discussion with Rev. wDZG, there are of course limitations to our work, but, as far as we know, this is the first convergence result with known rate obtained for such optimization schemes with non-trivial (non-Gaussian) data distributions.\n\n - **Experiments.** We agree that the toy example in Section 5 is minimal and serves only illustrative purposes, to convey intuition about the conditioning improvements and the trajectory of the iterates. Our primary goal was indeed to provide convergence results and the connections to imaging arised mainly from noting the similarity between MMSE Averaging and PnP algorithms in imaging settings (as well as notably *Cold Diffusion*), where they have shown good empirical performance. We will expand the introduction by adding more context about optimization aspects of the problem. In particular, note that the gradient of (Proximal Objective) is not Lipschitz, hence gradient descent converges as $O(\\\\frac{1}{\\\\sqrt{k}})$, and acceleration will not improve this rate. Splitting methods on the other hand will require the prox of $\\\\log p$, which results in a chicken and egg problem. Finally, higher order methods are also intractable since they will require for example the Hessian of $\\\\log p$ which is not generally known. Hence the rate of $O(\\\\frac{1}{k})$ obtained in Theorem 1 is a strict improvement for the specific problem we're dealing with.\n\n - **Convergence guarantees in PnP.** In the introduction and in the related works section, we noted that there exist many PnP methods with convergence guarantees. However, we also noted that they do not converge towards the MAP solution, but rather to the fixed point of some other operator (different from the proximal operator of $-\\\\log p$). Our goal is to bridge this gap by showing that one can compute the true MAP minimiser via MMSE denoisers.\n  The analysis in *Provable Convergence of Plug-and-Play Priors with MMSE Denoisers, Xu et al., 2020*, instead follows a different path: the assumptions are much weaker (no convexity required), but convergence is only shown to a stationary point of a function which is not the proximal operator (note the regularizer explicitly stated in Eq. 8 is not equal to $-\\\\log p(x)$). We will add comparison with this paper to the related work section, since it is indeed relevant.\n\n - **Is it necessary to work with $\\log p$ (i.e., a log-prior), rather than a general concave function?** Yes. Our work requires access to the score function of $-\\\\log p\\_\\\\sigma$. If instead one starts from a general function $h$, the analogous quantity would be $\\\\log (\\\\exp(-h) * \\\\mathcal{N})$, which is generally intractable as it doesn't have a natural Tweedie-based reformulation. Indeed, as is often the case when dealing with score functions, being able to apply the Tweedie formula is fundamental.\n\n - **Non-Gaussian priors.** While strong, the assumptions needed for our results to hold do not rule out several interesting non-Gaussian priors. Two examples are \n   1. the smoothed TV prior (where the smoothing can be a simple approximation of the absolute value $|x| \\\\approx \\\\sqrt{x^2 + \\\\delta^2} - \\\\delta$, or the Pseudo-Huber loss) which can be used as a prior for images penalizing discontinuities in their gradient,\n    2. a uniform distribution over an interval, convolved with a Gaussian to smooth out the discontinuities at the edges.\n\n   Both satisfy Assumptions 1 and 2, and the third derivative is bounded by a quantity which depends on $\\\\delta$ for 1. and on the standard deviation of the Gaussian in 2. Unlike in the Gaussian case the 3rd derivative is not uniformly lower than the 2nd on the whole of $\\\\mathbb{R}$. However, note that this does not make the naive GD on the *unsmoothed* proximal objective a good algorithm, because we cannot compute $\\\\nabla \\\\log p(x)$, but only $\\\\nabla \\\\log p\\_\\\\sigma(x)$ for $\\\\sigma > 0$.\n\n - **Comparison to simulated annealing.** Both simulated annealing (SA) and the MMSE Averaging algorithm involve decreasing a \"temperature\" or noise level which decreases over time to ensure convergence. In our case this serves to deal with non-smoothness of the objective, while in the case of SA it helps deal with non-convexity (or even non-continuity issues). While the principle of reducing the noise is the same, the context of the two algorithms remains quite different: MMSE Averaging is based on gradient descent while SA performs Gibbs sampling at each step to decide the next iterate. Furthermore, while $\\\\sigma\\_k$ in MMSE Averaging does not introduce additional noise, and thus the algorithm is deterministic, SA is inherently a stochastic algorithm. We will add a short comment to clarify the connections with SA.\n\n### Concerning the minor remarks:\n\n - We agree that \"ill-conditioned log-prior\" is not standard and we will clarify this. It is taken from the standard use of \"ill-conditioned objective function\" from the optimization literature.\n\n - We will revise the sentence about “a new wave of approaches” and cite earlier works such as *Louchet & Moisan (2013)*, clarifying that this perspective has an older history which has recently gained renewed attention. Thanks for pointing this out.\n\n - We will restate all theoretical results with full hypotheses for clarity, and similarly recall the assumptions behind the *Beck and Teboulle (2009)* result.\n\n - Regarding the comparison between $O(\\\\log(1/\\\\varepsilon))$ and $O(1/\\\\varepsilon)$: our point concerns the impact of **conditioning**. We will clarify this in the discussion.\n\n - We will fix missing notations like $\\\\gamma_k$ and $J^\\\\star$, and rewrite unclear sentences accordingly."}, {"review_id": "pKSsz0485V", "reviewer": "Reviewer_vJDD", "summary": "The submitted paper follows a plug-and-play (PnP) strategy [VBW 2013] to solve the maximum a posteriori (MAP) problem (MAP). In the strategy, the proximal gradient of the negative-log prior is replaced with the recursion (MMSE Averaging). The main results are Theorems 1 and 2. Theorem 1 claims that the MMSE averaging converges to the minimizer of the proximal objective (Proximal Objective) at the rate $O(1/k)$, with $k$ denoting the number of iterations in the MMSE averaging. Theorem 1 requires two technical assumptions: In Assumption 1, the prior density is assumed log-concave and strictly positive on $\\mathbb{R}^{d}$, i.e. its support is equal to $\\mathbb{R}^{d}$. Assumption 2 postulates the existence and boundedness of the third derivates of the log prior. Theorem 2 claims that approximate proximal gradient descent (PGD) obtained by plugging in the MMSE averaging converges to the true MAP solution at the rate $O(1/k)$ in a special sense presented in Theorem 2. A numerical visualization of the MMSE averaging is also presented with a two-dimensional Gaussian toy model.", "questions": "--Line 172:\n$\\ln p$ should be $-\\ln p$ since $\\ln p$ is not convex but concave.  \n\n--Proposition 1:\nIs the weight not $\\alpha_k = 1/(k+2)$ but $1 - \\alpha_k = 1/(k+2)$? Please check it. If it is true, this typo appears in Theorem 1 and Algorithm 1 while it does not provide any changes in the theorems. \n\n--Proposition 3 and Theorem 2:\nAdd the phrase “under Assumption 1” to Proposition 3. Theorem 2 is also in the same situation. Theorem 2 requires Assumptions 1 and 2 since Theorem 1 is used in its proof.", "limitations": "yes", "rating": 4, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 1, "originality": 2, "strengths_and_weaknesses": "Originality:\nThe MMSE averaging itself is not novel. For instance, a similar algorithm can be found in [Feng etl al. 2024, Algorithm 1]. As pointed out in the submitted paper, the MMSE averaging is inspired by a diffusion model like [Bansal et al. 2023]. Algorithm 1 in the submitted paper is also in the same situation: It is PnP itself. Thus, the novel results of the submitted paper are Theorems 1 and 2.\n \n[Feng etl al. 2024] T. Li, H. Feng, L. Wang, L. Zhu, Z. Xiong and H. Huang, \"Stimulating Diffusion Model for Image Denoising via Adaptive Embedding and Ensembling,\" in IEEE Transactions on Pattern Analysis and Machine Intelligence, vol. 46, no. 12, pp. 8240-8257, Dec. 2024.\n\nSignificance:\nThe convergence results themselves in Theorems 1 and 2 are good. However, they require strong assumptions in Assumptions 1 and 2, such as the log-concavity, three times differentiability, and its boundedness for the prior density. For instance, the $\\ell_1$ regularization (or the Laplace prior) is not differentiable at the origin. In my opinion, these strong assumptions exclude almost all priors requiring PnP. ReLU cannot be used in the denoising network. Thus, the practical significance of Theorems 1 and 2 is questionable. \n\nTheorem 2 is weak from a different point of view. It only claims the average $k^{-1}\\sum_{i=1}^{k}J(x_i)$ converges to the minimum of the MAP problem at the rate $O(1/k)$. While the convergence to the minimum is guaranteed, $J(x_i)$ is not guaranteed to converge at the rate $O(1/k)$. \n\nQuality & Clarity:\nThe manuscript itself is well written and easy to follow. For a further improvement, see my comments in Questions", "was_revised": true, "final_justification": "Authors' rebuttal has resolved my concern about the significance of the assumptions. Thus, I increase my score.", "rebuttal": "We thank the reviewer for their time and thoughtful feedback. Below we respond to the points raised, highlighting the changes we intend to introduce in the manuscript in response to the comments above.\n\n - **Originality of MMSE Averaging.** We fully agree that the algorithm is not novel, and we try to clearly state this in the paper, also explicitly linking it to Cold Diffusion (Bansal et al., 2023) and Indi (Delbracio & Milanfar, 2023). We would like to thank the reviewer for pointing out the paper of Feng et al., 2024 which also provides interesting insights connecting iterative denoising to diffusion models. However, our *main contribution is to establish rigorous convergence guarantees* under specific choices of noise levels and averaging weights of the iterative denoising framework. It can be seen as a theoretical justification for many iterative algorithms which have been proposed with empirical and heuristic bases, but with no rigorous justification.\n  \n   While the algorithm itself has been widely used before, the fact that the averaging recursion converges towards the proximal operator is not obvious, and requires noting its equivalence with gradient descent on a sequence of smoothed objectives (Prop. 1 in the paper).\n   Likewise, Algorithm 1 is indeed a PnP method and it highlights one of the potential applications which is enabled by having access to a converging approximation to the proximal operator of $\\\\log p$ (MMSE Averaging). Contrary to many successful PnP methods which approximate the proximal operator without guarantees on its accuracy (e.g. directly replacing it with a denoiser), we provide a theoretical proof that — under our assumptions — the PnP algorithm we propose converges to the MAP estimate.\n - **Significance under assumptions.** Concerns about the strength of the assumptions are indeed valid. Nevertheless, our goal was to provide a clear and rigorous theoretical foundation in a setting where the analysis is quite delicate, even under strong assumptions. Establishing convergence to the *true* proximal operator and MAP estimator is highly non-trivial and not common in the literature. We believe this is an important first step towards understanding and extending such results to more realistic, potentially non smooth, settings (see also response to Rev. LqTm on the potential extensions to non-convex settings). \n  \n   Note that the constraints on the prior distribution do not directly influence the constraints on the network used for MMSE denoising: while ReLU networks are not smooth, they can be used to approximate smooth functions effectively and hence used in our framework; we will clarify this distinction in the paper as it can indeed be confusing.\n   \n    To give an example of a more realistic prior which would satisfy our assumptions, we can start with the TV prior, commonly used for images as it penalizes discontinuities in the data. The TV prior for a discrete 1D signal would be $p(x) \\\\propto \\\\exp(-\\\\sum\\_{i=1}^{n-1} |x\\_{i+1} - x\\_i|)$. Since we need the prior to be differentiable, we can use a smooth approximation to the absolute value $|x| \\\\approx \\\\sqrt{x^2+\\\\delta^2}-\\\\delta$ for small $\\\\delta$. Then we can calculate the second and third derivatives to verify that the distribution is indeed log concave, and that the bound on the third derivative depends on $\\\\delta$. This prior is certainly more realistic than the Gaussian prior for image data and is (or was, before deep learning priors emerged) widely used in inverse problems.\n - **Theorem 2 convergence to the average.** Although last-iterate convergence would be great to have, it typically requires a noticeably more complex analysis. For a good discussion on the difficulties of last-iterate analysis see *\"Last iterate convergence of SGD for Least-Squares in the Interpolation regime\"*, A. Varre, L. Pillaud-Vivien, N. Flammarion, 2021. In general average-iterate analysis is quite standard in the context of analyzing optimization algorithms, for example see *\"Convergence Rates of Inexact Proximal-Gradient Methods for Convex Optimization\", M. Schmidt, N. Le Roux, F. Bach, 2011* who prove a similar theorem in a slightly different setting using the average iterate. Note that *best-iterate* convergence can be obtained easily from average iterate convergence. We will highlight this point in the revised manuscript.\n - **Typo in $\\\\alpha\\_k$.** Thanks for catching this, as you said the correct value is $\\\\alpha\\_k = \\\\frac{k+1}{k+2}$. Algorithm 1 contained a few inaccuracies, we have restated it in the rebuttal to reviewer 3djj.\n - **Statement of assumptions.** We will modify the theorem statements to be explicit about the required assumptions in the revised manuscript."}, {"review_id": "aO9QAQpLJw", "reviewer": "Reviewer_3djj", "summary": "Denoisers are often used as priors to solve inverse problems, typically by formulating the problem as a maximum a posteriori (MAP) estimation with a negative log-prior. While pre-trained denoisers commonly approximate this prior, there's generally no theoretical basis for this practice. This paper addresses that gap, demonstrating that a simple, practical algorithm (similar to the commonly used) provably converges to the proximal operator under a log-concavity assumption of the prior. The work also offers novel interpretations of the optimization problem as a gradient descent on a smoothed proximal objective.", "questions": "I have one questions that could help me clarify my understanding of the work.\n\nMy main question is regarding Algorithm 1. The inner-loop that does the MMSE averaging is with fixed weights and noise levels. This means this loop is solving a fixed point. I don't see the connection of this fixed point with the MMSE averaging equation where $\\alpha_k$ evolves at each iteration. Same with the connection with Cold Diffusion (where the weights evolve at every step).\n\nI would also like the authors to comment on the significance. This is theoretic paper without any real numerical experiment. The theoretical result is interesting but still makes some assumptions that are not generally hold. This limits the impact. \n\nBesides this question, I was wondering if there's a way to push the discretization and the numerical Algorithm 1, into an ODE formulation  (e.g, $\\tau \\to 0$). This could help with the understanding/intuition and also to develop connections to other works (diffusion, PF-ODE). See for example Appendix B in Delbracio and Milanfar 2023 where they analyze the case of a Gaussian prior. \n\nNotation:\n In Algorithm 1. $k$ represents the iteration on the outer loop, but within the text $k$ represents a different thing (see e.g., Theorem 2). This complicates the understanding.", "limitations": "Limitations are not explicitly stated but the paper discusses future work about extending the analysis to more realistic cases ( L343-L346).", "rating": 4, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 2, "originality": 3, "strengths_and_weaknesses": "Strengths:\n- The paper is generally well written, the connections with relevant work as correctly presented in a nice manner.\n- The paper introduces connections to other works (e.g., Cold Diffusion) which is something quite interesting.\n\nWeaknesses:\n - The main weaknesses is that the paper doesn't show any real application of the theory. The theoretical result is interesting though, so it might be that there's no need for showing applications, but having these will make the paper much stronger. \n- There are some minor presentation issues (repeated notation, for example \"k\" seems to mean something different in the main body of text than on the Algorithm 1.", "was_revised": true, "final_justification": "This is an interesting theory paper, well written, and which drives connections to relevant work. The main limitations is that it doesn't show any real application of the theory, but the theoretical result is still interesting (though based on some strong assumptions). The paper had minor typos that are going to be corrected in the revised version. My recommendation thus is borderline accept.", "rebuttal": "We thank the reviewer for their comments and thoughtful engagement with our work.\nThe concerns expressed in the review are addressed below.\n\n - **Algorithm 1.** Thanks for catching this, there were indeed some typos in Algorithm 1 (note that the analysis was done on the correct algorithm). The corrected algorithm uses a schedule for noise levels $\\sigma_k$ and weights $\\alpha_k$. We will also fix the overloading of variable $k$ for the iteration index in the main text which is confusing: we will use $i$ for the PGD algorithm loop and $k$ for the MMSE Averaging loop (it is fixed in the algorithm below). The correct algorithm looks like this:\n\n    $\\\\hat{x}^{(0)} \\\\gets y$\n\n    for $i = 0, 1, \\\\dots$ do\n\n    $\\\\quad \\\\hat{x}\\_0^{(i+1)} = \\\\hat{x}^{(i)} - \\\\tau \\\\lambda \\\\nabla f(\\\\hat{x}^{(i)})$\n\n    $\\\\quad \\\\text{for}~k = 0, 1, \\\\dots, n_{i+1} - 1$ do\n\n    $\\\\qquad\\\\sigma\\_k \\\\gets \\\\sqrt{\\\\frac{\\\\tau}{k + 1}}$\n\n    $\\\\qquad\\\\alpha_k \\\\gets \\\\frac{k+1}{k+2}$\n\n    $\\\\qquad\\\\hat{x}\\_{k+1}^{(i+1)} \\\\gets \\\\alpha_k \\\\mathrm{MMSE}\\_{\\\\sigma\\_k}(\\\\hat{x}\\_k^{(i+1)} ) + (1 - \\\\alpha\\_k) \\\\hat{x}\\_{0}^{(i+1)}$\n\n    $\\\\quad\\\\hat{x}^{(i+1)} \\\\gets \\\\hat{x}^{(i+1)}\\_{n\\_{i+1}}$\n\n - **Continuous-time version.** It is indeed possible to derive a continuous-time analogue of our gradient descent recursion on the smoothed proximal objectives. Specifically, one can write \n    \\\\[ \\\\dot{x}_t = - \\\\gamma\\_t \\nabla F\\_{\\\\sigma\\_t}(x\\_t) = - \\gamma\\_t (x\\_t - y) + \\\\gamma\\_t \\\\tau \\\\nabla \\\\ln p\\_{\\\\sigma\\_t}(x\\_t) \\\\]\n  \n    where $\\\\sigma_t$ and $\\\\gamma_t$ are time-dependent parameters. \n    Now assuming that $\\\\gamma_t = (1 + \\\\tau / \\\\sigma_t^2)^{-1}$, this can be equivalently rewritten as \n    \\\\[\\\\dot{x}\\_t = - x\\_t + \\\\gamma\\_t y + (1 - \\\\gamma\\_t) \\\\mathrm{MMSE}\\_{\\\\sigma\\_t}(x\\_t) \\\\] \n    which would correspond to the continuous-time equivalent of our MMSE averaging sequence.\n    Assuming that $p(x)$ is a Gaussian (hence the MMSE has a closed form), and that $\\sigma_t$ decreases as $\\tau / t$, we can solve the ODE in closed form and show that it converges to the proximal solution. \n    However, while the ODEs have a similar shape to the one derived in Appendix B of *Delbracio & Milanfar (2023)*, the two algorithms have different objectives and are not entirely comparable: the proximal solution for us, samples from $p(x)$ for InDI. Of course the type of convergence proved in InDi is weaker than what we prove in that it has a) restrictive assumptions requiring Gaussian $p(x)$ and b) is based on a continuous limit which can't be achieved in practice, while ours is based on a discrete recursion.\n\n - **Significance.** We agree that the assumptions (log-concavity, differentiability) do limit direct applicability. However, as with many theoretical papers, we believe this work serves as a first step toward understanding the behaviour of practical denoising-based schemes from a variational optimisation perspective, especially for a problem where no convergence rates have been proven so far. Besides, even under strong assumptions, establishing convergence to the *true* proximal operator and MAP estimator is highly non-trivial. We hope this theoretical foundation can guide future extensions to more general priors.\n\n- **Lack of experiments.** Our primary goal in this paper is to provide theoretical foundations for a class of heuristics already used in practice—such as in *Cold Diffusion*—which have demonstrated strong empirical performance. Nevertheless, we agree that a thorough empirical study is very valuable, and we believe it would hold as a separate work by itself."}]

Author rebuttal:

Meta-review: This paper concerns MAP optimization in the case that the proximal operator for the log-prior is unavailable but an exact MMSE denoiser is available.  The authors propose an MMSE averaging recursion that denoises at a sequence of increasingly smaller noise levels, and they prove convergence to the true proximal operator under certain assumption on log-concavity and boundedness of higher-order derivatives.

Four experts reviewed the paper and joined the authors in many rounds of rebuttal and discussion.  Although reviewers LqTm, vJDD, and 3djj all noted that the required assumptions (e.g., log-concavity, thrice differentiability, exact MMSE denoising) limit the applicability of the results in practical settings, the four reviewers all felt that the paper makes a strong theoretical contribution and is well written.  Furthermore, they appreciated the connections to cold diffusion.
Consequently, the paper can be accepted on the condition that the authors incorporate the reviewers' suggestions into the camera-ready version.
Decision: Accept (poster)