Particle backward simulator

Fix a finite horizon $T\ge 0$ and a requested path count $S\ge 1$, used in either mode below. In the particle mode, the inputs include a particle count $N\ge 1$ and normalized filtering clouds $\{(x_t^i,w_t^i):1\le i\le N\},w_t^i\ge 0,{\sum }_{i=1}^N{w}_t^i=1,$ for $0\le t\le T$, produced, for example, by Sequential importance resampling for Feynman–Kac path targets. Each $x_t^i$ is the time-$t$ coordinate of its stored particle path. Write $M_{t+1}$ for the model transition kernel and require a product-measurable density $f_{t+1}(x^{\mathrm{\prime }}\mid x)$ relative to a declared sigma-finite measure ${\lambda }_{t+1}$, meaning $M_{t+1}(x,A)={\int }_A{f}_{t+1}(x^{\mathrm{\prime }}\mid x){\lambda }_{t+1}(d{x}^{\mathrm{\prime }})$ for each measurable $A$. For each $t<T$, define the predictive score with zero-weight terms omitted: $q_{t+1}^N(x^{\mathrm{\prime }})={\sum }_{i:w_t^i>0}{w}_t^i{f}_{t+1}(x^{\mathrm{\prime }}\mid x_t^i)\mathrm{.}$ For particle-mode path sampling, define the set of reachable indices recursively by $R_T=\{j:w_T^j>0\}$. Descending from $t=T-1$ to $0$, require $q_{t+1}^N(x_{t+1}^j)\in (0,\mathrm{\infty })$ for every $j\in R_{t+1}$, and put $R_t={\bigcup }_{j\in R_{t+1}}\{i:w_t^i>0,f_{t+1}(x_{t+1}^j\mid x_t^i)>0\}\mathrm{.}$ These are exactly the indices that can be drawn with positive probability at the corresponding backward step. For the all-marginals mode below, require instead that $q_{t+1}^N(x_{t+1}^j)\in (0,\mathrm{\infty })$ whenever the current smoothing weight ${\bar{w}}_{t+1\mathrm{\mid }T}^j$ is positive. Thus no denominator is required for an omitted zero-weight term. At $T=0$ there are no denominator conditions.

In particle mode, repeat independently for $s=1,\dots ,S$:

1. Draw $I_T^s$ from $\{1,\dots ,N\}$ with probabilities $w_T^i$ and set ${\tilde{X}}_T^s=x_T^{I_T^s}$.
2. For $t=T-1,T-2,\dots ,0$, set $b_t^i=0$ without evaluating the density whenever $w_t^i=0$. For $w_t^i>0$, set $b_t^i({\tilde{X}}_{t+1}^s)=\frac{w_t^i{f}_{t+1}({\tilde{X}}_{t+1}^s\mid x_t^i)}{q_{t+1}^N({\tilde{X}}_{t+1}^s)}\mathrm{.}$ Draw $I_t^s$ from this probability vector, and set ${\tilde{X}}_t^s=x_t^{I_t^s}$.
3. Return the tuple ${\tilde{X}}_{0:T}^s=(x_0^{I_0^s},\dots ,x_T^{I_T^s})\mathrm{.}$ This is a draw from the particle backward approximation to the smoothing path law. The coordinates are selected by separate backward indices and need not be the stored forward ancestor history of any one particle. The particle path-sampling output is the list ${\tilde{X}}_{0:T}^1,\dots ,{\tilde{X}}_{0:T}^S$.

If an exact filtering kernel ${\eta }_T$ and supplied backward conditional kernels $B_t$ from Bayes factorization of the smoothing backward kernel are available for general transition kernels, there is a separate exact-kernel mode: draw $S$ conditionally independent paths, each by drawing $X_T\sim {\eta }_T$ and then drawing $X_t\sim B_t(X_{t+1},\cdot )$ successively backward. This mode uses the actual kernels; it does not evaluate a transition density or a kernel on particle singletons. It returns exact smoothing paths under the supplied conditional versions wherever their factorization holds.

An optional all-marginals mode computes particle smoothing weights instead of drawing paths. Initialize ${\bar{w}}_{T\mathrm{\mid }T}^j=w_T^j$. For $t=T-1,\dots ,0$, let $J_{t+1}=\{j:{\bar{w}}_{t+1\mathrm{\mid }T}^j>0\}$. For each $j\in J_{t+1}$ and $i$ with $w_t^i>0$, set $b_t^{i,j}=\frac{w_t^i{f}_{t+1}(x_{t+1}^j\mid x_t^i)}{q_{t+1}^N(x_{t+1}^j)};$ set $b_t^{i,j}=0$ without evaluating the density if $j\mathrm{\notin }{J}_{t+1}$ or $w_t^i=0$. Update ${\bar{w}}_{t\mathrm{\mid }T}^i={\sum }_{j\in J_{t+1}}{\bar{w}}_{t+1\mathrm{\mid }T}^j{b}_t^{i,j}\mathrm{.}$ For every active $j$, $(b_t^{i,j}{)}_{i=1}^N$ is a probability vector, so the raw updated values are nonnegative and already sum to one; explicit normalization leaves them unchanged. Return ${\sum }_i{\bar{w}}_{t\mathrm{\mid }T}^i{\delta }_{x_t^i}$ as the approximate smoothed marginal at time $t$. The algorithm makes no finite-particle accuracy, consistency, or long-time stability claim.

When a required density denominator is zero or infinite, the corresponding mode is inadmissible; the algorithm does not replace it with uniform weights. Zero filtering or smoothing weights are omitted before density evaluation, so no term evaluates $0\cdot \mathrm{\infty }$ or $0\mathrm{/}0$. For a general transition kernel, use an exact supplied backward conditional kernel as in the second mode, or provide a valid dominating-density interface. A fixed-lag target remains a truncated-observation target; see Fixed-lag particle smoother.