A New Perspective on Reverse Diffusion for Monte Carlo Sampling

Abstract: This paper introduces a novel perspective on the use of reverse diffusion processes for sampling from unnormalized densities. The central idea is to embed the target density as the marginal at the initial time of a suitably constructed diffusion process evolving over a finite horizon. In contrast to existing approaches, the proposed methodology involves neither time discretization error nor score function estimation, so that Monte Carlo variability is the only source of approximation. A key theoretical result characterizes the Radon-Nikodym derivative of the reverse diffusion transition distribution with respect to that of an Ornstein-Uhlenbeck (OU) process. This representation provides a tractable change-of-measure formulation and serves as the foundation for two distinct classes of Monte Carlo algorithms. The first class approximates the reverse transition distribution via a sequence of pseudo-marginal Metropolis-Hastings MCMC algorithms. The resulting scheme produces an approximate i.i.d. sample from the target distribution and is fully parallelizable, as trajectories can be generated independently. The second class consists of MCMC algorithms targeting the joint law of the whole diffusion path in $[0,T]$, for a suitably chosen horizon $T$. The proposed samplers combine three types of updates. One update simulates the diffusion forward in time according to an OU dynamics, conditional on its initial value. The remaining two update the backward component via Metropolis-type steps: one conditions on the terminal value at time $T$ and the other one does not. In both cases, acceptance probabilities are implemented using Barker-type Bernoulli factory constructions. The proposed methods perform well for targets with multimodality and complex dependence structures, providing a scalable and efficient alternative to the widely used random-walk Metropolis algorithm.

• The paper introduces an exact reverse diffusion approach that eliminates both discretization error and score function approximation.
• The SPARK and Ping-Pong algorithms are developed to robustly sample from complex, multimodal, and unnormalized densities.
• Empirical results show accurate mode allocation and scalable performance in high-dimensional, heavy-tailed targets.

Reverse Diffusion as Exact Monte Carlo: New Algorithms for Sampling Unnormalized Densities

Introduction and Motivations

The paper "A New Perspective on Reverse Diffusion for Monte Carlo Sampling" (2606.04357) addresses the intrinsic challenges associated with sampling from complex unnormalized densities—specifically, highly multimodal distributions and those with intricate dependence structures. Traditional MCMC methods such as random-walk Metropolis-Hastings and component-wise Gibbs samplers often become inefficient under such conditions due to slow mixing and metastability, particularly in high-dimensional or correlated spaces. While gradient-based methods (Langevin samplers, HMC) and advanced population-based and tempering approaches mitigate these issues to varying degrees, they have limitations: heavy reliance on discretization, gradient evaluations, or expensive tuning, and they typically encounter difficulties with separated modes.

Recent advances in generative modeling have leveraged denoising diffusion and score-based techniques, often combining SDEs with deep neural network-based score estimation. Despite their empirical success, these samplers are still approximate, incurring bias due to discretization or score approximation and entailing computational burdens from neural training and inference.

This paper proposes an alternative framework, reformulating sampling from an unnormalized density as exact simulation of reverse-time diffusions. The central innovation is the elimination of both time-discretization error and score function estimation. The only approximation comes from Monte Carlo variability. The methodology leverages a change-of-measure identity for reverse diffusions, yielding a tractable Radon-Nikodym derivative with respect to an Ornstein-Uhlenbeck (OU) process. This theoretical structure supports two classes of Monte Carlo algorithms: (1) a sequential pseudo-marginal strategy termed SPARK, and (2) a global extended-space MCMC strategy called Ping-Pong.

Theoretical Foundations

The methodological cornerstone is the exact representation of the law of a reverse-time SDE, starting from the target distribution as the initial condition. Notably, the Radon-Nikodym derivative of the reverse law with respect to a reference OU process is tractable, bypassing the need for a learned or approximate score. Concretely, the forward process follows:

$$dX_s = -\frac{1}{2}X_s\, ds + dB_s, \quad X_0 \sim p_0,$$

with $p_0$ the unnormalized target distribution. The reverse-time SDE admits an explicit change-of-measure to an OU process, governed by the marginal density $p(t,x)$ at intermediate times, for which the paper provides an unbiased Monte Carlo estimator.

In particular, for $p_0(x) = c\,\pi_u(x)$ (with known unnormalized $\pi_u$), the following unbiased estimator is introduced:

$$p(t,x) = c\,e^{\frac{d}{2}t} \mathbb{E}_{Z_{t,x}}[\pi_u(Z_{t,x})]$$

where $Z_{t,x} \sim \mathcal{N}(e^{t/2}x, (e^t-1)I_d)$. Importance sampling refinements and theoretical results on estimator variance under different tail conditions are given, covering both sub-Gaussian and heavy-tailed targets.

Change-of-Measure Identity

A general change-of-measure theorem (Theorem 1) relates the reverse diffusion law to the reference process, permitting exact construction of transition and marginal densities. The corollary for the OU setting yields the following proportionality (for a path $\omega_{0:T}$):

$$\frac{d\bar P}{d\bar Q}(\omega_{0:T}) \propto \frac{p_0(\omega_T)}{q(\omega_0)}$$

This structure allows the correct target measure to be enforced across the entire diffusion path.

Algorithmic Frameworks

SPARK: Sequential Pseudo-marginal Sampling

SPARK (Sequential Pseudo-marginal Algorithm with Reverse-time Kernels) constructs approximate i.i.d. samples from the target via a sequence of conditional MCMC or PMMH updates along a carefully designed time grid, reverse-propagating from the Gaussian reference at time $T$ to the target at time $p_0$0. Each update conditions on the result of the subsequent step, and only requires local MCMC mixing due to the pathwise construction, analogous in spirit to sequential SMC methods but with crucial differences—namely the absence of discretization error and the reliance only on unbiased estimators.

The SPARK procedure can be visualized through the progressive refinement of the particle distribution as it evolves from the well-behaved terminal marginal back to the complex target, as illustrated below.

Figure 1: Samples generated by SPARK at each time in the schedule.

Strong empirical performance is noted even on sharply multimodal and highly dependent targets. Notably, in a 10D Gaussian mixture, SPARK accurately reproduces the component weights and achieves correct allocation across modes.

Figure 2: Sample generated by SPARK projected onto the first two coordinates. Empirical component weights are $p_0$1, $p_0$2, and $p_0$3. Red crosshairs indicate the true component means.

Ping-Pong: Extended Space MCMC with SPIDER Proposal

The Ping-Pong methodology constructs a single Markov chain on the entire diffusion path, whose marginal at time 0 is the desired target. Proposals are generated using a multi-level sequential importance-resampling mechanism (SPIDER), and a Barker-type Metropolis correction is applied, using divide-and-conquer Bernoulli factory techniques to handle unbiased acceptance probabilities involving intractable proposal densities.

The implemented update strategies include:

1. Forward simulation: Exact propagation of the forward OU process.
2. Backward update (Reverse SIR proposal): Propose the reverse path and accept/reject using Barker's rule.
3. Global path update: Propose entire paths and accept via Barker.

Ping-Pong offers exact sampling in the extended space and is amenable to embedding as a Gibbs step within more complex samplers, as it does not require explicit gradients of the target density. Empirically, on challenging radial distributions, Ping-Pong produces diverse trajectories that reflect accurate exploration and mixing across modes.

Figure 3: Sample (left) and trace plot of the radius (right) generated by Ping-Pong.

Numerical Results and Practicalities

Benchmarks demonstrate that SPARK can generate thousands of approximately independent samples within seconds to minutes, whereas Ping-Pong MCMC, while more expensive per sample, provides an exact extended-space Markov chain with robust mixing. Empirical component weights are correctly recovered in high-dimensional mixtures. Both approaches remain efficient and scalable with increasing complexity and multimodality.

Key tuning considerations addressed include:

• Rescaling/standardization of the target prior to diffusion
• Choice of schedule (denser at early times for better regularization)
• Importance proposal selection based on the tail behavior of $p_0$4
• Variance and acceptance diagnostics for pseudo-marginal steps
• Optimization of the time horizon $p_0$5, balancing prior smoothing with computational cost.

Both frameworks exhibit superior robustness for multimodal, heavy-tailed, or non-Gaussian dependence structures compared to traditional random-walk-based MCMC.

Implications, Theoretical and Practical

The presented approach establishes, for the first time, an unbiased, discretization-free, and zeroth-order Monte Carlo sampling framework based on reverse diffusions for arbitrary unnormalized densities. This is a significant theoretical extension of the denoising diffusion paradigm, shifting the focus from approximate generative modeling (with neural approximations and discretization bias) to exact inference with only Monte Carlo variability.

Key implications include:

• Broader applicability: No need for differentiable or tractable gradients, enabling use on combinatorial or black-box targets.
• Theoretical tractability: The pathwise change-of-measure and unbiased density estimation enable rigorous analysis of method error and convergence.
• Parallelization: Algorithmic independence of replications (SPARK) and pathwise proposals (Ping-Pong) naturally scale in parallel and distributed environments.
• Embedding into hierarchies: The extended-space approach of Ping-Pong can be directly integrated as a Metropolis-within-Gibbs update for complex Bayesian models.

Some limitations include potentially higher computational cost, especially when importance estimators have high variance or when the acceptance probabilities become small for aggressive proposals. However, the versatility in the absence of the gradient and the removal of discretization bias provide a compelling alternative for targets where those costs justify the methodological improvements.

Conclusion

This work introduces a principled, zeroth-order, and unbiased framework for Monte Carlo sampling based on reverse-time diffusions, circumventing both time-discretization error and the need for gradient or score function estimation. Through two algorithmic schemes—SPARK and Ping-Pong—the approach demonstrates compelling empirical results on difficult targets, emphasizing robustness to multimodality and complex dependencies. The methods enrich the algorithmic toolkit for Bayesian computation, highlighting the value of rigorous change-of-measure representations and unbiased density estimation in modern scalable simulation. Future work is poised to explore further efficiency enhancements, improved path proposal mechanisms, and broader applicability to non-diffusive reference processes.