Decorrelating ReSTIR Samplers via MCMC Mutations

Monte Carlo rendering algorithms often utilize correlations between pixels to improve efficiency and enhance image quality. For real-time applications in particular, repeated reservoir resampling offers a powerful framework to reuse samples both spatially in an image and temporally across multiple frames. While such techniques achieve equal-error up to 100 times faster for real-time direct lighting and global illumination, they are still far from optimal. For instance, unchecked spatiotemporal resampling often introduces noticeable correlation artifacts, while reservoirs holding more than one sample suffer from impoverishment in the form of duplicate samples. We demonstrate how interleaving Markov Chain Monte Carlo (MCMC) mutations with reservoir resampling helps alleviate these issues, especially in scenes with glossy materials and difficult-to-sample lighting. Moreover, our approach does not introduce any bias, and in practice we find considerable improvement in image quality with just a single mutation per reservoir sample in each frame.


INTRODUCTION
The efficiency of rendering algorithms often hinges on their ability to effectively evaluate similar integrals by reusing samples across pixels [Ward et al. 1988;Lafortune and Willems 1993;Jensen 1996;Veach and Guibas 1997].In real-time path tracing, sample reuse becomes more critical since tracing rays is computationally intensive even on high-end consumer GPUs [Kilgariff et al. 2018].Moreover, while existing denoisers drastically improve image quality even at low sample counts [Chaitanya et al. 2017;Schied et al. 2017Schied et al. , 2018;;Kozlowski and Cheblokov 2021;NVIDIA 2022], they are unable to reconstruct features missing from their input samples.Thus, sample reuse is often the only means to improve sampling quality given limited computational budgets.Compared to methods that generate independent samples, reuse is also at times the only practical approach available to render challenging scenes with caustics and tricky lighting [Hachisuka and Jensen 2009;Veach and Guibas 1997].
Recent sampling algorithms for real-time ray tracing achieve massive speedups in scenes with complex illumination by sharing samples spatially within an image and temporally across frames [Bitterli et al. 2020;Ouyang et al. 2021;Lin et al. 2021Lin et al. , 2022]].These so-called

Random Replay Shift
Reconnection Shift Hybrid (Random Replay + Reconnection) Shift Fig. 2. Glossy scenes with difficult-to-sample lighting rendered using Re-STIR PT often contain correlation artifacts irrespective of the selected shift mapping strategy [Lin et al. 2022, Section 7].Artifacts result from suboptimal importance sampling and over-enthusiastically sharing a few highcontribution samples between pixels.
ReSTIR1 based techniques select  high-contribution samples from a larger streamed candidate pool of size .They do so by reformulating resampled importance sampling (RIS) [Talbot et al. 2005] in terms of weighted reservoir sampling (WRS) [Chao 1982].While RIS effectively importance samples candidates in proportion to an arbitrary target function (e.g., the integrand of the rendering equation), WRS makes resampling efficient by reducing storage costs from  () to  ( ).Repeated resampling across pixels then helps distribute important samples over several frames for estimation.Though ReSTIR derives impressive efficiency gains from correlated sampling, the benefits of repeated resampling are not indefinite.When only a few high-contribution samples have been identified, iterative spatial reuse creates blotchy artifacts as several pixels reuse the same sample (Figures 2 and 3).Such undersampling artifacts eventually fade away with temporal reuse over several frames, using a parameter to balance between minimizing pixel error and correlations via greater sample reuse (Figure 4).However, emphasizing error reduction adds lag under camera movement with dynamically changing lighting and geometry (Section 2.4), and introduces distracting low-frequency artifacts akin to those in photon mapping [Hachisuka and Jensen 2009], Metropolis Light Transport (MLT) [Veach and Guibas 1997] and Virtual Point Light (VPL) methods [Dachsbacher et al. 2014].
As spatiotemporal correlations are difficult to quantify, resolving artifacts is challenging.For instance, popular denoisers that compute first-and second-order moments (e.g., Schied et al. [2017]) are less effective given imprecise variance estimates with correlated samples.For ReSTIR, trying to reduce such artifacts by increasing the reservoir size  is also ineffective, as resampling with replacement [Chao 1982] produces duplicate samples in the presence of strong correlations (see Wyman and Panteleev [2021, Figure 19]).
Inspired by work on Sequential Monte Carlo (SMC) [Doucet et al. 2001] and Population Monte Carlo (PMC) [Cappé et al. 2004], we demonstrate that interleaving MCMC mutations with reservoir resampling (Section 3) helps alleviate correlations and impoverishment, especially in scenes with glossy materials and difficult lighting.Unlike MLT where mutations drive information sharing across pixels, our mutations instead help mitigate artifacts caused by spatiotemporal reuse.Similar to blue-noise sampling [Mitchell 1987;Georgiev and Fajardo 2016;Heitz and Belcour 2019], these mutations produce images with better visual fidelity without necessarily reducing error (Figures 9 and 11).Our approach highlights the complementary strengths of resampling and mutations for realtime rendering: resampling reuses samples with large contributions proportional to a pixel's target distribution, while mutations diversify the resampled population by locally perturbing samples in proportion to the same target distribution.Furthermore, like Veach and Guibas [1997]'s bias elimination strategy for MLT, we show that resampling eliminates the need for any burn-in period with Metropolis-Hastings (MH) mutations [Metropolis et al. 1953;Hastings 1970] (Section 2.5, Appendix A).This drives considerable image quality improvements from even a single mutation per frame for each reservoir sample (Figures 7,9 and 10).
From an implementation perspective, our approach requires only simple additions to existing ReSTIR algorithms (see Algorithm 3)we mutate reservoir samples using Metropolis-Hastings and an appropriate target function every frame after temporal reuse.This is immediately followed by an adjustment to each mutated sample's contribution weight to maintain detailed balance and ensure unbiased estimation.Overall, our contributions include: • Demonstrating how to incorporate MCMC mutations within ReSTIR to address the pitfalls of unchecked spatiotemporal sample reuse with resampling.
• Showing how to correctly adjust the RIS weight of mutated samples in an unbiased fashion for further resampling.• Situating ReSTIR in the broader paradigm of techniques that jointly apply resampling and mutations to sampling problems, such as MLT, SMC and PMC (see Table 2).
We start with the key building blocks of our approach in the next section, and postpone discussion about related work to Section 6 for better context when comparing with our method.

BACKGROUND
The rendering equation [Kajiya 1986] gives the outgoing radiance  out leaving a point  in the direction .Expressed as an integral over directions, it is Here  e is the emitted radiance,  in (,   ) is the incoming radiance from the direction   ,  (, ,   ) is the BSDF and   is the angle between   and the surface normal at . Absent participating media, the incident radiance  in is defined recursively as  in (,   ) =  out ( (,   ), −  ); the function  (,   ) returns the point on the closest surface from  in direction   .Integrating over the sphere of directions  2 then gives the total radiance scattered towards ; this integral can be estimated with Monte Carlo as where  (  ) is the probability density function (PDF) with respect to solid angle used to sample the incident directions   .
As in Kajiya's formulation, sometimes it is more convenient to reformulate Equation 1 over surfaces.To keep the discussion independent of the choice of formulation, we use ∫ Ω  () d to generically represent the integral we want to evaluate with Ω as its domain.This integral can likewise be estimated using where   are independent random samples drawn from any source PDF  that is non-zero on the support of  .In rendering, one often draws samples proportional to individual terms of the rendering equation to reduce variance (e.g., the BSDF ).To perform even better importance sampling, ReSTIR instead uses RIS to draw samples approximately proportional to the product of multiple terms in the integrand (e.g., We review RIS and generalized RIS next (Sections 2.1 and 2.2); Section 2.3 discusses a streaming RIS implementation via reservoir sampling.Section 2.4 then describes how correlations arise within ReSTIR due to resampling.Section 2.5 discusses the Metropolis-Hastings algorithm we use in Section 3 to resolve correlation artifacts.

Resampled Importance Sampling (RIS)
RIS [Talbot et al. 2005;Lin et al. 2022] enables unbiased estimation and sample generation from a non-negative target function p () with an unknown normalization factor ∫ Ω p () d.It does so by rewriting the standard Monte Carlo estimator from Equation 3 as The normalization factor is estimated by generating  ≥ 1 candidate samples y = { 1 , . . .,   } from a source PDF  that may be suboptimal but easy to sample from (e.g.,  ∝ ), yielding . (5) The samples x = { 1 , . . .,   } in turn are selected by randomly choosing an index  ∈ {1, . . .,  },  times, from the candidate pool y with discrete probabilities: Here the resampling weight  for each candidate   is given by where  := 1/(  ) is called the (unbiased) contribution weight for   .The selected samples   are likewise given contribution weights that assume the role of a reciprocal PDF, though they are only unbiased estimates for elements of the resampled set x.This is because the parenthesized term for the normalization factor of p is itself an estimator that has variance.Each   ∈ x is also distributed only approximately in proportion to p (i.e., p is sampled perfectly only in the limit as  → ∞).Since we resample with replacement, the set x can contain duplicate samples, which reflects that samples are selected in proportion to p.With this setup, Talbot [2005] shows that the RIS estimator is unbiased as long as p and  are non-zero on the support of  , i.e., Combining with Multiple Importance Sampling (MIS).There are often several reasonable sampling strategies available in rendering, e.g., BSDF or light sampling.MIS [Veach and Guibas 1995b] allows multiple strategies to be combined robustly within RIS [Talbot 2005].When each candidate   has its own source PDF   , then MIS weights generalize the parenthesized term in Equation 5with Here,   ≥ 0 is the MIS weight for the th sampling technique.These weights must form a partition of unity, i.e.,  =1   () = 1.A common choice is the balance heuristic   () =   ()/  =1   () [ Veach and Guibas 1995b].With MIS, the resampling weight in Equation 7 becomes: Notice we recover   = 1/ when source PDFs are the same for each sample   .MIS weights play an important role in ReSTIR.Beyond reducing noise in the resampling weights, they also remove bias when the supports of the source and target distributions do not match integrand  's support (see Section 4 in Bitterli et al. [2020] and Section 5 in Lin et al. [2022] for further details).
In practice, using RIS with the balance heuristic is costly, as all sampling strategies (i.e., the source PDFs) must be evaluated for each candidate sample   .Bitterli [2022, Chapter 9.1.3]provides a similarly robust but more efficient heuristic called Pairwise MIS, which only requires  () PDF evaluations over the entire candidate pool.We use pairwise MIS when the number of sampling strategies  is greater than 2 (e.g., during spatial resampling in ReSTIR; see Section 2.4); otherwise we use the balance heuristic.

Generalized Resampled Importance Sampling (GRIS)
So far we assumed the resampling inputs   ∼   share a common integration domain Ω with integrand  .This assumption may no longer hold when reusing spatially or temporally across an image (as in ReSTIR), and depends on the integral formulation used for the rendering equation.For instance, ReSTIR applied to global illumination [Ouyang et al. 2021;Lin et al. 2022] generates samples from PDFs with respect to solid angle.Reuse across pixels therefore requires a change of integration domain, necessitating a correction term in the resampling weights [Ouyang et al. 2021, Equation 11].ReSTIR for direct lighting [Bitterli et al. 2020] instead integrates over the surface of all lights, ensuring Ω is fixed across samples.
Recent work by Lin et al. [2022] generalizes RIS to use candidate samples   originating from different domains Ω  .It achieves this via shift mapping, i.e., a bijective transformation of samples from one pixel to corresponding samples on another pixel [Lehtinen et al. 2013].In particular, if Ω denotes the domain of integration for  , and   : Ω  → Ω are shifts that map   ∈ Ω  to the modified sample  ′  ∈ Ω, then the resampling weight for   becomes where the Jacobian determinant | ′  /  | accounts for the change of integration domain from Ω  to Ω. (Jacobians also appear in MIS weights   ; see Appendix B).The rest of the RIS procedure in Section 2.1 remains unchanged-substituting these resampling weights to Equation 8 provides the contribution weight for the selected  ′  .Various shift mappings have been proposed to maximize the similarity between  ′  and   such that | ′  /  | ≈ 1 [Hua et al. 2019, Section 3].We describe the shift mappings we use in Section 4.

Weighted Reservoir Sampling (WRS)
WRS [Chao 1982] facilitates efficient RIS implementations using a single pass over elements in a stream { 1 , . . .,   } to select a  ← 0 ⊲number of samples seen so far 5: ← 0 ⊲contribution weight (set in Algorithm 2) 6: function update(, ) 7: if rand() < (/ sum ) then 10: ←  random sample.As in Section 2.1, each stream element has an associated resampling weight .The basic idea is to process the stream one element at a time, and to select-from the  <  elements processed so far-a sample   with probability  (  )/  =1  (  ).The next stream element  +1 then replaces   with probability  ( +1 )/ +1 =1  (  ).The stream length  need not be known ahead of time, and WRS can be used to select  > 1 samples if needed [Wyman 2021, Chapter 22.6].
WRS reduces the storage needed for resampling to  ( ).A lightweight data structure called a reservoir is typically used to process the stream and store the selected samples, the stream length  and the weight sum  =1  (  ); see Algorithm 1.

Reservoir-based Spatiotemporal Resampling
ReSTIR applies RIS and WRS in a chained fashion within and across pixels of an image.The first key idea is to approximately importance sample multiple terms in the rendering equation's integrand through a per-pixel target function p.The second is to reuse samples from neighboring pixels to exploit the similarity between their target functions.The algorithm performs four steps every frame: (1) (Initial resampling) Select  samples from a candidate pool of  samples at each pixel.Equations 12 and 8 provide the resampling and contribution weights for the candidate and selected samples respectively.A reservoir stores the selected samples and their contribution weights.(2) (Temporal resampling) Use Algorithm 2 to reuse samples across two corresponding pixels in consecutive frames  and  − 1.The resampling weight for each sample is computed using the contribution weight already stored in its reservoir.
(3) (Spatial resampling) For each pixel, select  random reservoirs from a small spatial neighborhood and merge them into the pixel's reservoir.The is similar to Algorithm 2 and can be repeated multiple times; for reference see Bitterli et al. [2020, Algorithm 4] and Ouyang et al. [2021, Algorithm 2].(4) (Final shading) Use Equation 9to compute each pixel's color.
Spatiotemporal reuse gives each pixel access to a large population of samples from its local neighborhood.As a result, ReSTIR quickly finds samples that make large contributions to pixels.Nonetheless, gains from sharing samples are not indefinite, and correlation artifacts may arise from undersampling, imperfect shift mappings, and wrongly set parameters.For instance, performing multiple rounds of spatial resampling with too small a pixel radius can lead to blotchy 12: ⊲Combine weighted samples into a single reservoir return  VPL-like artifacts.This happens when RIS cannot effectively importance sample the integrand, e.g., due to difficult-to-sample lighting.Likewise, inadequately designed shift mappings may introduce geometric singularities into a sample's resampling weight via the Jacobian determinant, causing the sample to be widely reused.
During temporal resampling, one must cap the stream length or confidence weight  of a temporally reused sample (Algorithm 2, line 3) to guarantee convergence [Lin et al. 2022, Section 6.4]-not doing so results in convergence to a wrong result.Unfortunately, the ideal  cap cannot always be determined in a scene-agnostic waysmall caps inadequately utilize the temporal history and result in higher variance (Lin et al. [2022, Figure 9]), while large caps increase correlation.In particular, increasing  cap decreases the relative weight and selection probability of newly proposed samples.As a result, an outlier sample's contribution has to decay for it to match a pixel's average value.Unfortunately, the outlier may be spread between neighboring pixels before it is eventually replaced.This can lead to visible correlation artifacts and sample impoverishment over multiple frames (Figure 4).We use the Metropolis-Hastings algorithm, described next, to address these issues in ReSTIR.

Metropolis-Hastings (MH)
Like RIS, the MH [Metropolis et al. 1953;Hastings 1970] algorithm generates a set of samples distributed proportionally to a non-negative and possibly unnormalized target function p.While RIS uses resampling to achieve this goal, MH instead constructs a Markov chain that has a stationary distribution proportional to p's probability density function p/ ∫ Ω p.In more detail, given an initial sample  0 ∈ Ω, MH incrementally constructs a sequence of random samples  0 ,  1 ,  2 , ... as follows: For instance, a large  cap value introduces correlations (right), while a small  cap inadequately exploits prior samples, leading to noise (left).Our approach offers greater leeway in setting parameter values that trade noise for correlation (see Figure 10).
(1) Generate a candidate sample   by applying a random mutation to the current sample   in the chain, i.e., sample   from a proposal density  (  →   ).
The acceptance probability (  →   ) ensures that samples are distributed proportional to the target function p.The detailed balance condition guarantees the existence of the Markov chain's stationary distribution by requiring the transition density between any two sample values to be equal: (15) To generate the correct distribution from all inputs, Markov chains must be ergodic.This can be guaranteed easily with mutations that always propose candidate samples over the entire support of p, i.e.,  (  →   ) > 0 for all   and   where p (  ) > 0 and p (  ) > 0.Even with this constraint, there is still much freedom in choosing mutation strategies-Section 4 describes the strategies we use.
Unlike RIS, MH does not estimate the value of integrals.It does however produce valid samples from its target function which can be used by a secondary estimator such as RIS for estimation (Section 3).
Eliminating start-up bias.MH assumes the initial sample  0 is generated with probability density proportional to p; using a sample not from this distribution results in start-up bias.A typical solution runs the Markov chain for numerous iterations until the initial state is "forgotten", i.e., discarding several early samples generated by MH.Sadly, determining the length of this burn-in period is tricky, as it depends on the initial sample value and its actual distribution.Veach [1998, Chapter 11.3.1]instead proposed resampling  0 from  candidate samples y = { 1 , . . .,   } generated using an easy-to-sample source PDF (much like Section 2.1).Equations 6 and 7 then provide the discrete probabilities and resampling weights (resp.)needed to  select a candidate, i.e.,  0 =   for some  ∈ {1, . . .,  }.Finally, contributions of mutated samples initialized from  0 are weighted by Equation 8 to guarantee unbiased estimation.Our method leverages ReSTIR's built-in resampling to automatically avoid start-up bias when performing mutations.

METHOD
RIS improves sample selection from a target distribution when given a large population of candidate samples.ReSTIR amasses a sizable, per-pixel candidate pool for resampling through spatiotemporal reuse, helping it quickly identify high-contribution samples via RIS.However, at times ReSTIR extensively reuses a few samples over multiple frames due to imperfect importance sampling and suboptimal parameters, as it has no mechanism to diversify its existing sample population.
Inspired by Sequential and Population Monte Carlo techniques (Section 6), we interleave reservoir resampling with MCMC mutations to mitigate correlations and sample impoverishment caused by spatiotemporal reuse.Our key observation is that mutating reservoir samples with the same target function as RIS helps to quickly decorrelate the resampled population, especially when it contains outliers.In Algorithm 3, we use Metropolis-Hastings to locally perturb the per-pixel temporal samples selected by Algorithm 2; interleaving resampling with mutations diversifies the samples ReSTIR shares between pixels.We discuss key aspects of our work next, starting with how to modify mutated samples' contribution weights to guarantee unbiased results.

Modified contribution weights.
A contribution weight  (Equation 8) estimates the reciprocal value of the target PDF p/ ∫ Ω p that a sample is approximately distributed according to. is needed to compute resampling weights for combining reservoirs (Algorithm 2, lines 7 and 11) and to estimate per-pixel shading (Equation 9).
Contribution weights are sample dependent.Thus, a sample that undergoes mutation cannot reuse the weight associated with its original state, i.e., a mutated sample's contribution weight should provide an unbiased estimate for the sample's reciprocal target PDF.Our key contribution is to show that the unbiased contribution weight for any mutated sample   , from a Markov chain  0 , ...,   , ..., can be computed via the relation return   Equation 16 does not depend on samples between  0 and   in the Markov chain and imposes no constraints on computing  ( 0 ), which can arise from prior resampling, runs of MH, or a mix of the two.This provides flexibility in where and when to mutate samples during ReSTIR (as long as mutations are confined to a given pixel).
One can get an intuitive feel for Equation 16 by substituting in the expression for  ( 0 ) from Equation 8: Notice that the estimated normalization factor for p, i.e., the sum of weights , remains unchanged for both the initial and mutated samples  0 and   .This normalization factor arises via RIS (e.g., Algorithm 2, lines 14-15) prior to performing mutations.Meanwhile, MH treats the resampling weights as fixed, simply redistributing a reservoir's sample population proportionally to the per-pixel target function p. Equation 16 then encodes any required correction to a sample's contribution weight to account for the sample mutation.
Start-up bias.Algorithm 3 does not require a burn-in period for mutations, even though the samples used to initialize MH are not distributed exactly according to p.This is because we use the unbiased contribution weights of mutated samples for subsequent steps in ReSTIR, including when computing shading and resampling weights for further reuse.This approach eliminates start-up bias completely for any mutated sample   by ensuring for any function  .Appendix A provides a formal proof.Note that avoiding start-up bias does not imply samples generated using MH are well-distributed according to p.However, since we initialize MH using reservoir samples roughly proportional to the target function, our method does not rely on MH to find important samples (see Figure 11)-rather it decorrelates and diversifies outlier samples by mutating them locally in proportion to p.
When to perform mutations?Temporal reservoirs often contain stale samples, as ReSTIR assigns higher relative importance to existing samples.We therefore mutate samples output by Algorithm 2 within each pixel (Figure 5), using the same per-pixel target function as RIS for the current frame.Mutating samples randomly after temporal resampling diversifies the inputs to spatial resampling, protecting against possibly escalating amounts of sample impoverishment caused by repeated reuse.
Applying Algorithm 3 to mutate samples after the initial or spatial resampling steps in ReSTIR (Section 2.4) is possible but not required.Like mutations, initial resampling serves to rejuvenate the sample population every frame (by introducing new independent samples into the population).Samples from spatial resampling are stored for future reuse; mutating them proportional to the current target function would cause them to lag by one frame.
Finally, Algorithm 3 places no restrictions on MH iteration count.To improve runtime performance, one could adaptively specify mutation counts per pixel (including no mutations) using, for instance, local correlation estimates.We leave development of such heuristics to future work and use a fixed, user-specified number of iterations.

IMPLEMENTATION DETAILS
We perform mutations for both direct and indirect illumination in ReSTIR using Kelemen et al. [2002]'s primary sample space (PSS) parameterization.This conveniently allows applying mutations directly to random number sequences used to generate light-carrying paths, while constraining path vertices to remain on the scene manifold.Moreover, it simplifies use of certain shift mappings in ReSTIR PT, e.g., the random replay shift [Lin et al. 2022, Section 7.2].

Direct Lighting
Our ReSTIR DI mutations perturb the directions of reservoir samples via their random numbers.For direct lighting, path ȳ = [y 0 , y 1 , y 2 ] and its PDF ( ȳ) equals   ()|cos  |/|y 2 − y 1 | 2 , where   is the PDF for importance sampling the BSDF ,  is the unit vector from y 1 to y 2 , and  is the angle between  and the geometric surface normal at y 2 .The PDF ( ȳ′ ) is defined analogously for the mutated  ′ , pointing from y 1 to y ′ 2 .Random numbers for the starting MH sample y 2 are recovered by inverting the sampling procedure for direction  [Bitterli et al. 2017].Since this mutation is symmetric, the transition kernels in Equation 19 cancel.

Indirect Illumination
For ReSTIR PT, our mutation strategies build on shift maps.Unlike a mutation, a shift mapping deterministically perturbs a base path x through one pixel into an offset path ȳ through another pixel.
For instance, a random replay shift reuses the random numbers that generate x to trace ȳ.Since tracing a full path is expensive, a reconnection is often used to connect the offset path to the base path at a given index , i.e., y  = x  for  ≥ .Connecting paths immediately with  = 2 is called the reconnection shift.Compared to random replay, reconnections are often better at producing paths with similar contributions for diffuse surfaces.But reconnecting y −2 , y −1 to x  on a glossy surface can introduce paths with nearzero throughput, or introduce geometric singularities when y −1 and x  are too close.We use Lin et al.'s [2022] hybrid shift strategy (see Figure 6) to evaluate mutations in ReSTIR PT.This shift mapping postpones reconnection using random replay until certain connectability conditions are met (e.g., surface roughness and distance between vertices).
Mutation strategies.As with direct lighting, one way to mutate a path is to perturb the random numbers used to generate it.Like a random replay shift, this approach expensively requires tracing a full path for each proposed mutation (which may be rejected).
A more computationally efficient approach mutates the offset path with random replay up to the reconnection vertex y  = x  , and then connects to the base path starting at x +1 instead.We observe this mutation strategy is not only faster, but also has higher acceptance (70% vs. 40% on the scene from Figure 1) as it minimizes changes to path geometry.Moreover, its paths have similar contributions to the offset paths it mutates.Note that mutating path vertices with random replay until the reconnection to x +1 can cause connectability conditions for the hybrid shift to fail.We reject such mutated samples by defining their transition PDF to be 0.
Taking a step further, our final strategy mutates only the reconnection vertex y  (Figure 6) while keeping the rest of the offset path unchanged, i.e., [y ′ 0 , . . ., y  strategy only slightly less effective at reducing correlations.It is, however, significantly faster when performing multiple mutations, as only rays from y −1 to y ′  and y ′  to x +1 need to be traced.We use this mutation strategy to generate results in Section 5, unless otherwise noted.
Finally, note that the transition kernels  ( ū′ → ū) and  ( ū → ū′ ) are no longer symmetric when offset paths contain a reconnection vertex.In Appendix C, we show that their ratio equals: where  −1 ,   and  +1 are unit vectors from y −1 to y  , y  to y +1 (= x +1 ) and y +1 to y +2 (= x +2 ) respectively,  is the angle between   and the surface normal at y +1 , and  is the solid angle PDF used to sample an outgoing direction.Primed quantities are defined similarly.Any mutations applied to random numbers for the subpath [y 0 , y 1 , . . ., y −1 ] do not factor in the ratio as they are symmetric.
Reservoir storage.Lin et al. [2022, Section 8.2] note ReSTIR PT stores additional data in the reservoir from Algorithm 1, specifically a seed for random replay and the resampled path's reconnection vertex.For the first two mutation strategies above, we need the path's entire random number sequence since PSS mutations transform this sequence-as a result, it cannot be regenerated from its original seed.This increases the reservoir size as path length grows.Luckily, our final mutation strategy avoids this overhead, only mutating random numbers that sample y ′  from fixed offset vertex y −1 .As in ReSTIR DI, we recover random numbers for y  by inverting the sampling of direction y  − y −1 .The only additional information we store is the offset vertex y +1 (which connects to mutated vertex y ′  ).

RESULTS AND DISCUSSION
We prototyped our method in the open-source Falcor rendering framework [Kallweit et al. 2022].All results use a GeForce RTX 3090 GPU at 1920 × 1080 resolution.Our implementation uses the settings (e.g., spatial neighborhood size and reuse radius) proposed in Bitterli et al. [2020] and Lin et al. [2022] for direct and indirect  illumination, with the exception of  cap = 50 in our ReSTIR PT tests.
Our supplementary videos show 1 spp results for all our scenes; Table 1 gives single frame timings.As shown in Figure 1 and Figures 7-10, short-range correlation artifacts are noticeably reduced in scenes with glossy materials and difficult lighting with just 1-5 mutations; further mutations have diminishing returns in improving image quality (Figure 9).Mutation cost overhead is generally less than simply increasing sample count (Figure 7), and recent denoisers [NVIDIA 2017] provide considerably better results with our decorrelated samples (see Figure 1).Figure 8 shows mutations greatly reduce sample impoverishment, with fewer reservoirs sharing the exact same sample realizations.
Compared to standard path tracing, ReSTIR is much faster at achieving equal-error via correlated sampling for real-time direct [Bitterli et al. 2020, Figure 8] and global illumination [Lin et al. 2022, Figure 13].Mutations however provide only marginal improvements in mean squared error in ReSTIR samplers (see Figures 9 and 11).Akin to blue-noise dithering [Georgiev and Fajardo 2016;Heitz and Belcour 2019], our image quality improves despite errors having similar magnitudes.The reason is mutating within a pixel leaves the sum of resampling weight unchanged in Equation 17, and these weights ultimately control RIS estimator variance (Equation 9).Mutations do slightly reduce variance, as they indirectly alter resampling weights of future samples thanks to spatiotemporal reuse of the new, more diverse sample population; the supplementary document has more details.In Figure 10 we also ablate  cap values to show the greater leeway our approach offers for this parameter, allowing use of larger values to trade noise for correlation.
Since ReSTIR often suffers from correlation artifacts, we quantify improvements in correlation by computing sample covariance  between pixels, which naturally generalizes sample variance.This metric measures the joint variability of two random variables (e.g., whether error in two pixels varies similarly).For pixels  and  in image  , the sample covariance    between  and  is given by where  is the number of images used to estimate covariance (we use  = 100), and Ī is the average of  images.To capture the joint variability of a pixel with its local neighborhood, in our experiments we average covariance estimates over boxes of a given radius centered at each pixel.We then further average over the entire image to get a single number.Figure 9 (bottom left) shows average radial covariance decreases with increasing spatial radius.This is expected as ReSTIR only reuses samples in local neighborhoods (so smallscale correlation artifacts are more pronounced); mutations reduce covariance in these short ranges.Table 1 lists the reduction in average covariance observed on our scenes, with pixel radius equal to 8. As correlations are typically localized, the reduction is even larger for the image insets in our figures compared to the results in Table 1.Ineffective shift mappings in ReSTIR often result in increased correlations; mutations compensate for this shortcoming.For instance, mutations typically have fewer correlation artifacts to resolve with a hybrid shift in ReSTIR PT compared to, e.g., random replay (Figures 11 and 8 respectively), which highlights the benefit of using good shift mappings.In contrast, mutations provide greater covariance reduction in the Forest scene rendered with ReSTIR DI (Figure 7, top), where higher covariance stems from vertex reconnections failing to preserve path contributions for low roughness surfaces.
Why mutations help?The supplemental document details why mutations reduce covariance, simplifying down to the following, arXiv, Vol. 1, No. 1, Article .Publication date: November 2022.somewhat unintuitive, phenomenon: without mutations, covariance between pixels  and  stems from mismatches between input sample distributions and the target functions at  and  (Equation 16 in the supplemental).However, in the limit of infinite mutations, covariance is determined by samples' mismatch with their own pixel's target function (Equation 10 in the supplemental) due to the ratio p ( 0 )/ p (  ) in the mutated contribution weight (Equation 16); this mismatch tends to be smaller.Though our analysis predicts that covariance does not vanish completely even with infinite mutations, our results show covariance is often reduced with just one mutation.

RELATED WORK
Our method builds directly on the recent ReSTIR family of algorithms for real-time direct [Bitterli et al. 2020] and global illumination [Ouyang et al. 2021;Lin et al. 2021Lin et al. , 2022]].We augment spatiotemporal reservoir resampling in ReSTIR with sample mutations, and demonstrate the complementary strengths of resampling and mutations in this framework.In graphics, our approach is most closely related to Metropolis Light Transport (MLT) [Veach and Guibas 1997] and associated techniques [Kelemen et al. 2002;Jakob and Marschner 2012;Lehtinen et al. 2013;Hachisuka et al. 2014;Otsu et al. 2018;Cline et al. 2005;Lai et al. 2007Lai et al. , 2009;;Bashford-Rogers et al. 2021].In the broader Monte Carlo landscape, our approach belongs to the class of algorithms that jointly use resampling and mutations for sampling problems, such as Sequential Monte Carlo (SMC) [Doucet et al. 2001] and Population Monte Carlo (PMC) [Cappé et al. 2004].We discuss the relation to MLT, SMC and PMC in more detail next; Table 2 provides a summary.We refer the reader to Bitterli et al. [2020, Section 7] and Lin et al. [2022, Section 9.3] for comparisons between ReSTIR and other rendering algorithms that exploit path reuse and spatial correlations.
Metropolis Light Transport.MLT uses statistically correlated samples generated by Metropolis-Hastings to solve the rendering equation.Unlike algorithms using independent samples, MLT is effective at finding difficult light paths by locally exploring the path space.It reuses samples by mutating high-contribution paths over the image.Algorithmically, our method resembles MLT in various ways.Both techniques require secondary estimators, respectively RIS and bidirectional path tracing (BDPT) [Lafortune and Willems 1993;Veach and Guibas 1995a], to normalize the MH target function.Samples used by these estimators are resampled into a smaller set to initialize MH (our Section 3 and Veach [1998, Chapter 11.3.1]),and contributions of mutated samples are effectively weighted by the same weights (Equation 17) to remain unbiased (our Appendix A and Veach [1998, Appendix 11.A]).
The crucial difference between our work and MLT lies in how samples are reused across pixels.MLT latches onto high-contribution paths and mutates them over the entire image while just resampling to eliminate start-up bias.Thus, MLT results often contain correlation artifacts caused by mutations, applying MH to both find important samples and redistribute them between pixels.In contrast, ReSTIR derives spatiotemporal reuse from resampling; in this paper, we mutate samples within each pixel to mitigate correlations and sample impoverishment from spatiotemporal resampling.As a result, our method does not require numerous MH iterations, as the primary purpose of mutations is not finding important paths (Figure 11).Further, our approach suits real-time rendering as it integrates seamlessly into ReSTIR.MLT can be adapted to mutate temporally, but unlike our work, the entire animated sequence must be available in advance [Van de Woestijne et al. 2017].Several features have recently been added to MLT, including sample stratification [Cline et al. 2005], MIS [Hachisuka et al. 2014] and enhanced mutation strategies [Jakob and Marschner 2012;Bitterli et al. 2017;Otsu et al. 2018;Kaplanyan et al. 2014].Though we mostly employ simple PSS-style mutations [Kelemen et al. 2002], many of these improvements can also be incorporated into our approach.[Doucet et al. 2001].As shown in the inset, the goal is maintaining a population of weighted samples distributed roughly proportional to an evolving target distribution (with unknown normalization factor).Sample weights are adjusted every iteration to reflect each sample's importance to the most recent distribution.Resampling discards samples with low weights and duplicates those with high weights.Mutations ensure the population does not contain identical samples.Unlike ReSTIR, which uses RIS, SMC methods use weighted importance sampling (WIS) to estimate correlated integrals in a chained fashion, i.e., the current step's sample weights and normalization factors are defined incrementally based on corresponding quantities from earlier steps [Del Moral et al. 2006 illumination of dynamic environment maps.Unlike ReSTIR, which derives its samples from spatiotemporally neighboring pixels, et al. [2006] instead maintain a fixed sample population per pixel that is resampled and mutated to be updated for each frame.Large populations are needed for effective importance sampling, as highcontribution samples are not shared between pixels; in contrast, ReSTIR often stores just a single sample per reservoir.SMC methods likely require MIS weights and shift mappings (like ReSTIR) to resolve bias and correctly derive effective spatiotemporal reuse from neighbors.Similar to our work, mutations mitigate sample impoverishment but do not provide reuse.
Population Monte Carlo.PMC methods also couple resampling and mutations to distribute weighted samples in proportion to a sequence of target functions [Cappé et al. 2004].The main added feature is they sample using parametric mixture models with simple source PDFs.Mixture probabilities are tuned for each target function using previously generated samples and their importance.
In rendering, the PMC framework has been used for direct lighting [Fan et al. 2007;Lai et al. 2015], global illumination [Lai and Dyer 2007;Lai et al. 2007] and animation [Lai et al. 2009].Lai et al.'s [2009] work is most relevant to ours: they derive sample reuse by mutating samples spatially and temporally across the image plane using Energy Redistribution Path Tracing (ERPT) [Cline et al. 2005].Resampling serves to select high-contribution samples while discarding those with small weights; it is also used to refresh the sample population (much like initial resampling in ReSTIR) and eliminate start-up bias from mutations.Unlike our method, they require knowing animated sequences in advance, precluding most real-time applications.

LIMITATIONS AND FUTURE WORK
In this paper, we provide an unbiased mechanism leveraging MCMC mutations to diversify ReSTIR's sample population.Often, just a single mutation per pixel effectively mitigates correlation artifacts in glossy scenes with complex lighting.However, as in most MCMC schemes, we cannot accurately predict the number of Metropolis-Hastings iterations needed to reduce correlations below a given threshold.Beyond the analysis in the supplemental document, further investigation is also needed to understand how mutations address sample impoverishment in ReSTIR-not just in terms of the number of duplicate samples (Figure 8), but also the discrepancy characteristics of the resulting sample population.
Mutating inside ReSTIR has a non-negligible run-time overhead.Though we demonstrate improvements on an equal-time covariance metric with simple mutation strategies in both ReSTIR DI and Re-STIR PT (Figure 7 and Table 1), more sophisticated mutations [Jakob and Marschner 2012;Bitterli et al. 2017;Otsu et al. 2018;Kaplanyan et al. 2014] could provide further gains.Our decision to mutate only after temporal (but not spatial) resampling is also informed in part by run-time considerations.As mentioned in Section 3, applying mutations selectively (i.e., not at each pixel every frame) based on e.g., local correlation heuristics could improve performance.As both mutations and ReSTIR's initial path candidates serve to rejuvenate the sample population, it may be interesting to carefully balance the costs of per-pixel mutations versus new path candidates.
Our proposed sample mutations reduce correlation between nearby pixels, leading to an error distribution (likely) closer to white noise.But blue noise error distributions are often superior with respect to human perception [Mitchell 1987]; perhaps our mutations could change to more directly optimize for blue noise characteristics.For example, when deciding mutation acceptance, we might consider both the target function and the neighboring pixel samples, preferring mutations that introduce differing sample values.A further improvement might apply the insights of Heitz and Belcour [2019] to optimize the image-space distribution of error rather than solely considering the sample values.
Like Metropolis Light Transport, mutating samples across pixels potentially unlocks further amortization by sharing samples over the entire image (e.g., using the expected values technique [Veach 1998, Section 11.5]).We leave such "cross-pixel" mutations to future work as they require adjusting a mutated sample's contribution weight (Equation 16) to account for varying integration domains.
More generally, by augmenting ReSTIR with mutations, our work establishes a closer correspondence between the RIS-based resampling techniques developed in graphics, and those in the broader statistics literature such as SMC and PMC.In particular, our approach stands to benefit from techniques such as annealed importance sampling [Neal 2001] used in SMC to reduce variance in the resampling weights [Ghosh et al. 2006, Section 4], as well as from adaptation strategies for mutation kernels developed in PMC to increase acceptance rates [Lai et al. 2007, Section 4.2].Moreover, as in these fields, mutations in ReSTIR open the door not just to artifact-free integration (of the rendering equation), but also to tracking and filtering problems-for instance using well-distributed sample populations generated by our approach as training data for path guiding.

Fig. 3 .
Fig. 3. Reservoir resampling suffers from sample impoverishment as it becomes more difficult to sample light-carrying paths.Top row, left to right: The Veach Ajar scene rendered using ReSTIR PT (random replay shift) at 1 spp with the door's angle decreasing.Bottom row: Heat maps visualize duplicate samples in 20 × 20 pixel neighborhoods.Black represents no duplicates, while white indicates the number of identical samples in a neighborhood.

Fig. 4 .
Fig.4.Parameters for ReSTIR sample reuse can be difficult to set in a scene agnostic way.For instance, a large  cap value introduces correlations (right), while a small  cap inadequately exploits prior samples, leading to noise (left).Our approach offers greater leeway in setting parameter values that trade noise for correlation (see Figure10).

Fig. 5 .
Fig. 5. Our approach introduces Metropolis-Hastings mutations as an additional block into the larger ReSTIR algorithm for spatiotemporal sample reuse.Samples are mutated within each pixel after temporal resampling (Algorithm 2) to mitigate correlation artifacts and sample impoverishment.

Fig. 6 .
Fig.6.The hybrid shift in ReSTIR PT reconnects the offset path to the base path when it encounters two consecutive diffuse vertices x 3 , x 4 ; prior to that it reuses random numbers from the base path to trace rays.Our mutation strategy perturbs the reconnection vertex y 4 in the offset path.

Fig. 7 .
Fig. 7. Correlation artifacts often do not disappear simply by using more samples, justifying the overhead of performing mutations.

Fig. 8 .
Fig. 8. Mutations mitigate sample impoverishment in ReSTIR by diversifying the sample population.The bottom row visualizes duplicate samples in 20×20 pixel neighborhoods on the scene from Figure 3.

ForestFig. 9 .
Fig. 9. Increasing sample mutations reduces short range correlation artifacts produced by ReSTIR, with even 1-5 mutations providing noticeable improvements in image quality (measured in the bottom left using average radial covariance).Mutations typically have little impact on mean squared error (shown in the bottom right at equal spp with  cap = 50) as we perturb samples only within each pixel.

Fig. 10 .
Fig. 10.By reducing correlation artifacts, mutations allow use of larger  cap values in ReSTIR to trade noise for correlation, yielding lower error in scenes with difficult to sample light-carrying paths (see Figure 4 for results with smaller  cap values).

Fig. 11 .
Fig. 11.Mutations do not reduce mean squared error in the Veach Ajar scene rendered using ReSTIR PT with the hybrid shift and  cap = 50.This suggests that in contrast to Metropolis Light Transport, resampling (and not mutations) finds important light-carrying paths in ReSTIR.Compared to the random replay shift in Figure 8, resampling with the superior hybrid shift does not introduce large correlation artifacts in this scene.
ALGORITHM 3: Mutate sample via Metropolis-Hastings Input: Pixel i, reservoir   from Algorithm 2, and iteration count Output: Reservoir   with its sample mutated in proportion to p 1: function mutateSample(,   , iters)

Table 2 .
Ghosh et al. [2006]rally reusing samples for estimation, instead of generating new samples every frame.SMC methods have found limited use in rendering;Ghosh et al. [2006]sample a sequence of per-pixel target functions for direct arXiv, Vol. 1, No. 1, Article .Publication date: November 2022.Overview of the role of resampling and mutations in MLT, PMC, SMC and ReSTIR.