Source-linked AI summary

Near-Linear Accuracy Bounds for Moreau--Yosida Unadjusted Langevin Sampling

Yuchen Xin, Zhihua Zhang

arXiv:2609.40193v1cs.LG

TL;DR

For composite log-concave targets, smoothing improves approximation but makes discretization accuracy depend on the smoothing parameter. The paper analyzes MYULA with a Poisson-based estimate and proves near-linear accuracy complexity under an explicit step-size condition.

  • Problem

    Moreau smoothing improves approximation of the target, but increases the global smoothness bound, making control of discretization accuracy as smoothing decreases a central issue.

    Source: §1 · p. 2
  • Method

    The analysis tests a second-order stationary residual against a Poisson solution to convert it into a Wasserstein bound, then combines the estimate with smoothing bias and contraction.

    Source: §1.1 · p. 3, §4 · p. 11
  • Results

    The invariant-measure bias is eO(h), and for fixed model parameters and initialization MYULA requires eO(ε^-1) iterations to achieve √m W2(μN, π) ≤ ε.

    Source: §6 · p. 24
  • Takeaways & Limitations

    The bounds apply to classical MYULA without additional higher-order smoothness assumptions, with one gradient evaluation and one exact proximal evaluation per iteration.

    Source: §1 · p. 2, §3.2 · p. 10, §6 · p. 24

    The sufficient step-size restriction is stronger than hLλ ≤ c.

    Source: §3.2 · p. 10

Abstract

from arXiv · show

We establish near-linear accuracy bounds for the classical Moreau--Yosida unadjusted Langevin algorithm (MYULA). The target is $\pi\propto e^{-f-g}$, where $f\in C^2(\mathbb{R}^d)$ is $m$-strongly convex with Lipschitz gradient and $g$ is convex and globally Lipschitz. Under an explicit parameter-dependent step-size condition, we bound the invariant-measure bias relative to the Moreau-smoothed target by $\widetilde O(h)$, with only logarithmic dependence on the inverse smoothing parameter in the error coefficient. Combining this estimate with the Moreau approximation bias and Wasserstein contraction gives $\widetilde O(\varepsilon^{-1})$ iterations to make the $N$th-iterate law $\mu_N$ satisfy $\sqrt m\,W_2(\mu_N,\pi)\le\varepsilon$, for fixed model parameters and initialization. We bound the stationary error directly, without assuming third derivatives or a Lipschitz Hessian. Each iteration uses one gradient evaluation and one exact proximal evaluation. The key idea in our analysis is to convert a second-order stationary residual into a Wasserstein bound using a Poisson-based estimate.

1 Introduction

MYULA smooths a nonsmooth convex penalty, but its accuracy for the original target depends on controlling how the smoothing parameter affects discretization. The paper develops bounds for this trade-off, a Poisson-based proof strategy, and comparisons with related sampling guarantees.

1 Introduction

  • Reducing the Moreau smoothing parameter improves approximation of the target but increases the global smoothness bound, so discretization error must be controlled jointly with smoothing.MYULA replaces g with its Moreau envelope gλ; the resulting smoothness bound is Lλ = Lf + λ^-1.Source: §1 · p. 2

1.1 Contribution

  • Under its step-size condition, the invariant-measure bias has near-linear order in h, with only logarithmic dependence on λ^-1 in the error coefficient.The sufficient step-size restriction is stronger than hLλ ≤ c.Source: §1.1 · p. 3
  • For fixed model parameters and initialization, the sufficient iteration complexity is eO(ε^-1) for √m W2(μN, π) ≤ ε, improving the cited eO(ε^-4/3) guarantee.The comparison is with Xin and Zhang [9]; it concerns the fixed-parameter accuracy exponent, not uniformly better dependence on every model parameter.Source: §1.1 · p. 3, §1.2 · p. 3
  • The proof tests a second-order stationary residual against a Poisson solution, shifting derivatives onto the test function to obtain a Wasserstein bound without higher-order smoothness assumptions.Weighted energy estimates retain the Hessian as a weight, while drift–heat identities and flux estimates control the residual.Source: §1.1 · p. 3

1.2 Related Work

  • Xin and Zhang [8] obtain eO(ε^-2) complexity for structured penalties, while Xin and Zhang [9] study the same assumptions and transition but obtain eO(ε^-4/3) complexity.The latter analysis gives invariant-measure bias O(h) + eO(h^3/4); this paper replaces it with eO(h) under a stronger step-size restriction.Source: §1.2 · p. 3
  • First-order ULA bias for a fixed smooth target does not directly yield near-linear accuracy complexity after Moreau smoothing; the gain comes from joint smoothing–discretization analysis.The related-work passage applies a first-order estimate to f + gλ and chooses λ proportional to ε/G².Source: §1.2 · p. 4
  • Fan et al. [4] obtain polylogarithmic W2 accuracy dependence under their initialization and subroutine assumptions, using approximate rejection sampling rather than a deterministic proximal update.That guarantee concerns a different algorithm and computational model; the result here is for classical MYULA with one exact proximal evaluation per step.Source: §1.2 · p. 4

2 Preliminaries

The preliminaries set notation and weak-derivative conventions, then establish tools for Moreau regularization, drift-based density evolution, Poisson equations, transport, and density-ratio moments.

2 Preliminaries

  • P2(Rd) denotes probability measures with finite second moment, and the section fixes conventions for norms, matrix order, and eO notation.eO suppresses fixed powers of logarithms, but no polynomial parameter factors.Source: §2 · p. 4, §2.1 · p. 5

2.1 Weak derivatives

  • Weak derivatives are defined distributionally; when represented by a locally integrable function, they agree almost everywhere without implying classical differentiability everywhere.The divergence formulas require only local integrability, allowing density equations to be used by integration by parts.Source: §2.1 · p. 5

2.2 Moreau regularization

  • Moreau regularization of convex globally Lipschitz g yields gλ with globally Lipschitz gradient, while its Hessian bound is understood weakly.The potential Uλ used below is m-strongly convex with Lλ-Lipschitz gradient.Source: §2.2 · p. 5, §2.3 · p. 5

2.3 Drift maps, pushforwards, and density evolution

  • For 0 ≤ tLλ < 1, the drift map Tt is a global bi-Lipschitz bijection, and its pushforward preserves the L1 norm of integrable densities.Probability pushforwards give the law of TtX; the free heat semigroup adds an independent N(0, 2tId) increment.Source: §2.3 · p. 6

2.4 Poisson equations and transport

  • A Poisson equation converts stationary density residuals into differences of expectations, while strong convexity supplies Poincaré control and solution regularity without a continuous Hessian or third derivative.The weak solution is unique and mean-zero; the approximation result permits its use in residual identities initially stated for smooth compactly supported tests.Source: §2.4 · pp. 7–8

2.5 Entropy and density-ratio moments

  • Strong convexity also yields a log-Sobolev inequality, and Rényi divergence is equivalently expressed through the Lα norm of the density ratio.For α > 1, ∥ν/πλ∥Lα(πλ) = exp{(α − 1)Rα(ν∥πλ)/α}.Source: §2.5 · p. 9

3 The Moreau–Yosida Unadjusted Langevin Sampling

Under the stated convexity and smoothness assumptions, the analysis gives fixed-smoothing invariant-bias control and near-linear iteration bounds for MYULA. For fixed model parameters and initialization, the iteration bound is eO(ε^-1).

3.1 Problem formulation and assumptions

  • The analysis assumes f is strongly convex with Lipschitz gradient and g is convex and globally Lipschitz; each iteration uses one gradient and one exact proximal evaluation.Source: §3.1 · pp. 9–10

3.2 Main results

  • The first theorem controls invariant-measure bias at fixed smoothing under a sufficient step-size restriction and gives a unique invariant law.Source: §3.2 · p. 10
  • For fixed model parameters and initialization, N(ε)=eO(ε^-1); for fixed Lf/m and G/√m, the bound is eO(d/ε).Here eO suppresses logarithmic factors.Source: §3.2 · p. 11

4 Proof of the Main Results

The proof combines analytic and stationary-flux estimates with a Poisson-based residual-to-transport bound, then uses smoothing bias and contraction to establish end-to-end accuracy. [Actually no markers]

4.1 The analytic estimates

  • Drift-moment bounds and finite-order density-ratio estimates supply the integrability needed for the later stationary and transport arguments.The drift-moment estimates are uniform in the smoothing parameter, and the Euler interpolation’s equal endpoint laws support a Rényi-dissipation estimate using the log-Sobolev inequality for πλ.Source: §4.1 · p. 12

4.2 The stationary residual

  • Stationarity yields weak flux identities for the residual, with the distributional Laplacian of the averaged density represented in L1.The stationary residual is expressed through vector and matrix fluxes; the argument does not assume a classical Laplacian of the averaged density.Source: §4.2 · pp. 13–15

4.3 A residual-to-transport estimate

  • A Poisson-based estimate converts second-order residual identities into a Wasserstein bound without requiring classical derivatives of the fluxes or weak Hessian.The residual estimate applies when the stated density and weighted-norm conditions hold.Source: §4.3 · p. 16

4.4 Proof of the fixed-smoothing theorem

  • For fixed smoothing, the proof verifies the residual estimate’s density and flux hypotheses and invokes it to bound the invariant-law error.Under the step-size condition, the chain has a unique invariant law with a strictly positive smooth density.Source: §4.4 · pp. 18–19

4.5 Proof of the end-to-end complexity theorem

  • The end-to-end proof allocates the accuracy budget among Moreau bias, invariant-measure bias, and convergence from the initial law, using contraction for the last term.The proof combines the invariant-measure estimate with the Moreau approximation bound to obtain the iteration guarantee.Source: §4.5 · p. 20

5 Stationary Flux Estimates

This section establishes stationary-flux bounds by controlling the averaged drift, expanding the heat evolution finitely, and deriving bounds for the remaining fluxes.

5.1 The averaged drift flux

  • The averaged drift flux is bounded by C√rR using moment and density-ratio estimates.The estimate follows by combining weighted norms of the drift and invariant density with the stated moment bounds.Source: §5.1 · p. 21

5.2 A finite heat expansion

  • A finite heat expansion, derived from the stationary drift balance, provides a route to estimating the gradient of the invariant density.The expansion is an exact finite identity; its smooth representative is obtained by Gaussian convolution.Source: §5.2 · pp. 21–22

5.3 Proof of the flux bounds

  • The remaining flux bounds are obtained by expressing the heat-path gradient through the invariant-density gradient and averaged drift, then applying the established estimates.The section uses this route to bound H; the bound on E follows from drift-moment and pushforward estimates.Source: §5 · p. 20, §5.3 · pp. 23–24

6 Conclusion

The conclusion states an invariant-measure bias guarantee for MYULA and combines it with Moreau approximation and Wasserstein contraction to obtain an end-to-end iteration guarantee.

6 Conclusion

  • The invariant-measure bias relative to the smoothed target is eO(h) under an explicit step-size condition, with only logarithmic dependence on λ^-1 in the error coefficient.The resulting end-to-end guarantee is eO(ε^-1) iterations to achieve √m W2(µN, π) ≤ ε for fixed model parameters and initialization.Source: §6 · p. 24
  • The guarantees leave classical MYULA unchanged and require no additional higher-order smoothness assumptions.The conclusion also states that each iteration uses one gradient evaluation and one exact proximal evaluation.Source: §6 · p. 24

A.1 Target drift moments

This appendix section develops drift-moment estimates using integration by parts and the curvature and Lipschitz assumptions.

A.1 Target drift moments

  • Higher drift moments are controlled using integration by parts, the curvature bound on f, and the bound ∥∇gλ∥≤G.The argument applies to ∥∇f∥^(s−2)∇f and then uses Hölder’s inequality to obtain the moment bound.Source: §A.1 · pp. 24–25
  • The resulting even-order moment growth yields an exponential estimate for the smoothed-target drift.The passage states the bound ∥∇Uλ∥^(2k) integrated against πλ is at most (18kR²)^k.Source: §A.1 · p. 25

A.2 A Gibbs change of variables for the heat operator

The appendix develops heat-kernel moment, derivative, and pushforward estimates under explicit time restrictions, then uses them in the stationary density-ratio proof.

A.2 A Gibbs change of variables for the heat operator

  • A.2 uses a Gibbs change of variables and Jensen’s inequality to bound heat-averaged moments by 2 when the time and moment order satisfy the stated restriction.The argument covers scalar, vector, and matrix densities with their corresponding norms.Source: §A.2 · pp. 25–26

A.3 Gaussian derivatives without an extra dimension factor

  • A.3 derives scalar gradient and vector-field divergence bounds without an extra dimension factor, using Gaussian integration by parts and heat-kernel estimates.The bounds apply under time restrictions tied to the moment order and R².Source: §A.3 · pp. 26–27

A.4 Pushforward estimates by an integrated Jacobian

  • A.4 bounds heat-kernel pushforwards using integrated Jacobians of the drift maps, under a step condition involving Lλ and R².The maps are controlled through their Jacobians and the Gibbs density ratio.Source: §A.4 · p. 28

B Finite-Order Stationary Density Ratios

  • The stationary density-ratio proof combines weighted drift and curvature moments, a drift-increment bound, and Rényi dissipation with the stationary endpoint condition.The appendix follows this three-step structure under the restrictions of (4.4).Source: B Finite-Order Stationary Density Ratios · p. 29

B.1 Weighted energy for the density ratio

The proof regularizes the potentials to justify the density-ratio calculation, then establishes weighted moments and drift-increment estimates for use in the dissipation argument.

B.1 Weighted energy for the density ratio

  • The auxiliary convolution scale δ is distinct from the Moreau parameter λ, and the finite integrability bounds may depend on δ.The lemma asserts no uniform bound as δ decreases to zero.Source: §B.1 · pp. 29–30
  • The density ratio u is the quantity whose α-moment is controlled; its relative score s and the drift discrepancy e_t enter the dissipation calculation.The drift discrepancy compares the current drift at Y with the drift frozen at X.Source: §B.1 · p. 30
  • Weighted bounds on ∇f and the mean Moreau curvature provide moment estimates used to control the drift increment.The smooth-component estimate uses Δf ≤ τf and the Moreau component uses the bounded-gradient property.Source: §B.1 · pp. 30–31, §B.2 · p. 33

B.2 Gaussian integration by parts with the density-ratio weight

  • Gaussian integration by parts with the density-ratio weight yields the drift-increment bound needed by the Rényi dissipation calculation.The estimates treat the smooth and Moreau components separately under step restrictions.Source: §B.2 · pp. 31–32

B.3 Proof of the density-ratio bound

The proof controls the density ratio through a Rényi-divergence dissipation estimate, using drift bounds and stationarity at the interpolation endpoints.

B.3 Proof of the density-ratio bound

  • The proof inserts the drift-increment bound into the Rényi dissipation inequality and uses the stationary endpoint condition to bound the density ratio.It establishes the estimate by controlling the dissipation from above and bounding its associated term from below.Source: §B.3 · pp. 33, 36, 37
  • The resulting density-ratio estimate (4.6) holds for the original potential, with quantitative constants independent of the auxiliary smoothing scale.The proof first obtains the estimate for a regularized potential and then passes to the limit.Source: §B.3 · p. 38

C Poisson Regularity and Negative Sobolev Norms

  • Under mI ⪯ ∇²Uλ ⪯ LλI almost everywhere, Appendix C proves Poisson regularity and equivalence of test-function classes for the negative Sobolev norm.The stated setting has Uλ ∈ C1,1(Rd).Source: C Poisson Regularity and Negative… · p. 38

C.1 Existence and second-order energy

The appendices establish second-order Poisson energy, approximation and test-function equivalence, and the regularity arguments needed for the density-ratio calculation.

C.1 Existence and second-order energy

  • The Poisson solution has a weighted L² Hessian, yielding the first inequality in (2.20); the second follows from the Poincaré inequality (2.17).The proof obtains the Hessian bound using cutoffs and nonnegative energy terms.Source: §C.1 · p. 39

C.2 Approximation by compactly supported tests

  • Compactly supported H² approximations and mollification produce smooth compactly supported tests converging in all three norms in (2.21).The mollification scales are chosen so the total approximation error is at most 1/n.Source: §C.2 · p. 40

C.3 Equivalent test classes for the negative Sobolev norm

  • The suprema defining the negative Sobolev norm agree across H¹(πλ), smooth compactly supported tests, and C¹ tests with finite weighted gradient norm.The proof establishes equality by cutoff and mollification approximations.Source: §C.3 · pp. 40–41

D Weak Regularity and Integrability Before Dissipation

  • Appendix D independently verifies mollification properties and the integrability needed for the dissipation calculation, then passes the quantitative estimate to the original potential.These arguments support the calculations in Appendix B.Source: D Weak Regularity and Integrability… · p. 41

D.1 Integration by parts and smooth approximation

  • Convolution preserves the structural bounds and step restrictions for every δ, while higher derivative bounds of the mollified potentials may depend on δ.The estimates in Appendix B use the displayed convexity, gradient, and Hessian bounds.Source: §D.1 · p. 41

D.2 Finiteness before the R´enyi calculation

The appendix establishes the finiteness and regularity bounds needed for the R´enyi calculation, then removes an auxiliary mollification while preserving the quantitative bound.

D.2 Finiteness before the R´enyi calculation

  • The regularized Euler kernel contracts in W2 and has a unique invariant law when hLλ < c < 1.The finiteness bounds used to justify analytic operations may depend on fixed parameters and δ, but do not enter the quantitative step restrictions.Source: §D.2 · p. 42
  • An exponential moment of the invariant law, uniform density and gradient bounds across t ∈ [0, h], and joint weighted-moment bounds establish the finiteness needed for the calculation.The density and gradient bounds are uniform over the interpolation interval.Source: §D.2 · pp. 42, 46, 48
  • The interpolated density is locally C1,2 on (0, h) × Rd, with continuity also checked at both time endpoints.The proof uses dominated convergence to establish the interior regularity and endpoint continuity.Source: §D.2 · pp. 46–47

D.3 Removing the auxiliary mollification

  • Removing the auxiliary mollification gives uniform-in-time W2 convergence of the regularized interpolations and preserves the full quantitative bound.The bound follows by lower semicontinuity of Rényi divergence and taking the supremum over t ∈ [0, h].Source: §D.3 · pp. 51–52
Loading 2609.40193v1…