Differentiable Systematic Resampling for Variational Sequential Monte Carlo

• Systematic resampling is widely adopted for its lower variance, but its discrete ancestor assignments block gradient flow and can yield biased or high-variance gradient estimators.Source: §1 · p. 1
• Soft multinomial resampling inherits higher variance, while EROT uses iterative Sinkhorn updates and can produce diffuse plans unlike systematic resampling’s sparse CDF order.Source: §1 · p. 2
• DSR uses sigmoid functions to soften CDF-level ancestor assignments, preserving systematic resampling’s CDF-ordered, banded structure while enabling full gradient flow.Source: §1 · p. 2
• As τ →0+, DSR converges to exact systematic resampling at a pointwise exponential rate under a positive-margin condition.Source: §1 · p. 2
• DSR achieves competitive or superior state estimation and dynamics learning in VSMC while training more than 2× faster than EROT.Source: §1 · p. 2

2.1 State space models and particle filters

• Particle filters approximate the filtering posterior by propagating and reweighting particles, then resample to replace the weighted set with equal-weight particles and avoid weight degeneracy.Source: §2.1 · p. 2

2.2 Variational Sequential Monte Carlo

• VSMC jointly learns transition and proposal models by maximizing a surrogate ELBO; its gradient includes a high-variance score-function term for discrete ancestors, which is dropped in practice, yielding a biased gradient.Continuous propagation and weighting use reparameterization; dropping the score-function term means gradients do not pass through ancestor selection. At finite temperature, a differentiable relaxation instead gives a surrogate objective that is not a lower bound on log p(y1:T).Source: §2.2 · p. 3

3.1 DSR formulation

• DSR replaces the hard systematic ancestor threshold with a smooth sigmoid surrogate, while leaving the random offset draw unchanged.The sigmoid approaches the step indicator as τ →0+, and the relaxed assignment is used to form convex combinations of ancestor particles.Source: §3.1 · p. 4
• At finite τ, DSR does not preserve expected replication counts; smaller τ sharpens the transport plan but produces larger gradients, while larger τ smooths gradients.The full transport matrix costs O(N^2) per resampling step, though batched GPU processing keeps the overhead small in practice.Source: §3.1 · p. 4

3.2 Particle ordering

• Morton sorting places spatially nearby particles next to one another in the CDF, so soft assignments mix nearby ancestors rather than distant particles.Source: §3.2 · p. 4
• At N = 200, Morton sorting reduces the L2 distance between DSR and hard systematic resampled particles by 35% versus random ordering.It also offers a modest filtering and transition-model advantage at low SNR.Source: §3.2 · p. 5

4 Theoretical Analysis

• As τ → 0+, DSR recovers systematic resampling for each resampled particle.The result establishes asymptotic exactness of the relaxation.Source: §4 · p. 5
• For positive margin ∆j, the pointwise error bound depends on particle diameter D, particle count N, and ∆j, with dominant term exp(−∆j/τ).Here ∆j is the distance between uj and the nearest CDF boundary.Source: §4 · p. 5
• The exponential rate is not uniform over offsets because ∆j can be arbitrarily small under some weight and offset configurations.The positive-margin condition holds almost surely, but does not yield a uniform rate.Source: §4 · p. 5
• Empirically, the per-particle L2 distance from hard systematic resampling decreases monotonically and accelerates as τ → 0+.This Appendix E observation is consistent with the exponential rate.Source: §4 · p. 5

5.1 Experimental Setup

• The evaluation spans a 4D linear-Gaussian model, Lorenz-63 with nonlinear camera observations, and CharacterTrajectories; most methods use 25 training particles and results are averaged over 10 runs.The common setup uses matched Adam settings, network architectures, and gradient clipping; the linear-Gaussian results use three runs.Source: §5.1 · pp. 5–6, §5.2 · p. 6

5.2 Linear Gaussian State Space Model

• Hard achieved the better training-time posterior in the linear-Gaussian task (KLD 17.49 vs. 21.77 for DSR), but DSR recovered the true transition matrix more accurately (0.019 vs. 0.035 Frobenius error).The experiment uses the Kalman posterior as an oracle and isolates resampling’s effect on parameter learning.Source: §5.2 · p. 6, Table 1 · p. 6

5.3 Lorenz-63

• −29.94 dB NMSE at 10 dB SMNR was DSR’s best Lorenz-63 result; at 0 dB, DSR improved NMSE by 5.7 dB over EROT (−12.57 dB).Hard and Soft struggled more on longer horizons; DSR’s transition model was not best on every one-step metric at 0 dB, but its particles filtered more effectively at inference.Source: Table 2 · p. 7, §5.3 · p. 8

5.3.1 Gradient norm

• On Lorenz-63 at 10 dB, DSR maintained gradient norms around 4 × 10^3, roughly 3× lower than Hard and EROT, with similarly reduced within-epoch variance.Lower norms alone do not establish a better gradient signal, but the authors report this pattern alongside DSR’s improved performance.Source: §5.3.1 · p. 8, Figure 3 · p. 8, Table 3 · p. 8

5.3.2 Scaling with particle count

• Across N = 25, 50, 100, and 200, DSR was the only method to train stably; at N = 200, Hard and Soft NMSE degraded to −20.42 and −21.52 dB, respectively.DSR improved from −29.94 dB at N = 25 to −30.35 dB at N = 200, with standard deviation below 0.3 dB throughout; EROT did not complete at N = 200.Source: §5.3.2 · pp. 8–9

5.4 Character typing

• At 0 dB SMNR on CharacterTrajectories, DSR achieved −8.59 dB NMSE, a 0.82 dB improvement over EROT (−7.77 dB), corresponding to a 21% average MSE reduction.At 10 dB, all differentiable methods were comparable; DSR’s NMSE was −15.21 dB versus EROT’s −14.88 dB and Soft’s −14.69 dB.Source: Table 4 · p. 9, §5.4 · p. 9

5.5 Training time comparison

• On an NVIDIA T4, DSR training took 2.55 h versus 2.49 h for Hard, while EROT took 5.96 h.DSR’s single-pass sigmoid soft assignment avoids EROT’s 100 Sinkhorn iterations.Source: Table 5 · p. 9, §5.5 · pp. 9–10
• At finite temperature, DSR yields convex combinations rather than exact ancestor copies, which can place mass between modes in multimodal posteriors.Morton sorting limits this cross-mode mixing when modes are well-separated, but multimodal posteriors remain challenging for particle filters generally.Source: §5.5 · p. 10

6 Conclusion

• DSR preserves systematic resampling’s CDF-ordered, banded structure while enabling full gradient flow through the particle-filtering recursion.The conclusion presents it as a temperature-controlled soft relaxation that trades gradient flow against exactness.Source: §6 · p. 10

A Differentiable Systematic Sampling Algorithm

• MortonSort is typically optional when the filtering posterior is unimodal, and the algorithm can use batched matrix operations instead of explicit particle loops.The conclusion notes that particle ordering has less effect on soft assignments in this regime.Source: A Differentiable Systematic Sampling… · p. 13
• The algorithm computes a weight CDF and draws a uniform offset before assigning particles by systematic resampling or the DSR relaxation.Algorithm 3 labels its DSR and systematic-resampling assignment steps separately.Source: A Differentiable Systematic Sampling… · p. 13

Proofs of propositions

• The proof treats convergence of the soft transport plan and convergence of the resampled particles as separate parts.The proof of particle convergence justifies exchanging a limit and finite summation, with transport entries in [0, 1].Source: B Proofs of propositions · pp. 13–14
• For each resampled column, the proof defines its hard systematic ancestor and considers indices below and above that ancestor.The case analysis uses the CDF intervals and the sampled offset.Source: B Proofs of propositions · pp. 13–14
• The proof uses D, the diameter of the particle set, in bounding the distance between resampled particles.D is defined as the supremum of pairwise particle distances.Source: B Proofs of propositions · p. 15

C.1 Network Architectures

• Transition and proposal models use multilayer perceptrons with ReLU activations and predict Gaussian distributions with diagonal covariance matrices.The proposal model shares the transition model’s architecture but takes the previous state and observation as concatenated inputs.Source: §C.1 · p. 16

C.2 Hyperparameters

• The linear Gaussian experiment uses N = 25 particles and τ = 0.1, training on 1,000 trajectories of length 100 for 500 epochs.It uses a learnable transition matrix and diagonal covariance, with no separate proposal model; log-likelihood clipping is used to avoid extreme training instability.Source: §C.2 · p. 16
• Lorenz-63 training uses N = 25 particles for 500 epochs on 1,000 trajectories, while Character Typing uses 25 particles on 1,400 training sequences.The Lorenz-63 setting uses τ = 0.1 and 20-bit Morton curves; Character Typing uses an augmented state dimension of 9.Source: §C.2 · p. 16
• State-estimation evaluation uses hard systematic resampling and N = 1,000 particles.This evaluation setup otherwise matches the training setup.Source: §C.2 · p. 16

D.1 SMNR

• A low signal-to-measurement noise ratio (SMNR) corresponds to high measurement noise.The definition uses observation-noise variance per dimension, observation dimension, and test-trajectory count.Source: §D.1 · p. 16

D.2 State Estimation; D.3 Transition Model Quality

• State-estimation accuracy is measured by NMSE in dB, while transition-model quality includes one-step and five-step MSE, log-likelihood, and KL divergence.Lower NMSE, MSE, and KL divergence indicate better results; higher log-likelihood indicates better probabilistic modeling.Source: §D.2 · p. 17, §D.3 · p. 17

E.1 Transport Plan Visualization

• At τ = 0.01, DSR’s transport plan closely resembles hard systematic resampling; increasing τ diffuses mass across neighboring ancestors, while retaining CDF-ordered banding.Figure 5 compares the transport plans for hard systematic resampling, DSR at varying temperatures, and EROT.Source: §E.1 · pp. 17–18, Figure 5 · p. 18

E.2 Effective Number of Ancestors

• 1.70 effective ancestors per resampled particle at τ = 0.01 with N = 25; effective ancestor count grows with temperature and, at fixed temperature, with particle count.At τ = 0.1, DSR is sparser than EROT: 6.54 vs. 9.83 for N = 25 and 23.11 vs. 36.41 for N = 100.Source: §E.2 · p. 18

E.3 Correlation Between DSR and Systematic Resampling

• > 0.95 correlation with hard systematic resampling in every tested configuration; even at τ = 1.0, correlation stays above 0.956.The correlation indicates agreement in selected CDF regions, not numerical closeness of resampled particles.Source: §E.3 · pp. 18–19, Figure 6 · p. 19

E.4 L2 Distance and Effect of Morton Sorting

• 35% lower L2 distance at τ = 0.01 and N = 200 with Morton sorting than with random ordering (0.38 to 0.25).At τ ≥ 0.5, the transport plan is diffuse enough that ordering matters less.Source: §E.4 · p. 19

E.5 Bias and Variance of DSR Estimates

• Two to five orders of magnitude lower variance than hard systematic resampling for the linear test function, with lower MSE at every τ.For the quadratic test function, bias increases monotonically with τ; the favorable τ = 1.0 result for the linear function is not a general recommendation.Source: §E.5 · pp. 19–20, Table 6 · p. 20

E.6 Discussion

• At τ = 0.1, DSR is substantially sparser than EROT while remaining highly correlated (> 0.95) with hard systematic resampling; Morton sorting reduces approximation error, especially at low τ and high N.The authors report these patterns as supporting the theoretical predictions and the complementarity of locality-preserving ordering with the relaxation.Source: §E.6 · p. 20

F.1 Dependence on τ

• −28.31 dB NMSE at τ = 0.01; performance plateaus for τ ∈ [0.05, 1.0], with NMSE varying by less than 0.7 dB across that range.The authors attribute unstable optimization at small τ to exploding gradient magnitudes and choose τ = 0.1 for their experiments.Source: §F.1 · pp. 20–21, Table 7 · p. 20
• As τ shrinks, per-particle L2 approximation error falls from 0.505 to 0.102 while median gradient norm and within-epoch standard deviation rise.For τ ≥ 0.05, the ratio of median within-epoch standard deviation to median gradient norm stays near 0.45; it rises to 0.55 at τ = 0.01.Source: §F.1 · p. 21, Table 8 · p. 21

F.2 Annealing τ during training

• Annealing τ from 1.0 to 0.1 gives small, consistent improvements across metrics and reduces NMSE standard deviation to 0.10 versus 0.28 with fixed τ.The schedule is tested on Lorenz-63 at SMNR = 10 dB under the camera measurement map.Source: §F.2 · p. 21, Table 9 · p. 21

F.3 Dependence on N

• DSR alone scales gracefully across tested particle counts, while Hard and Soft collapse at N = 200 and EROT exceeds available GPU memory.DSR's NMSE improves from −29.94 dB at N = 25 to around −30.4 dB from N = 50 onward, with run-to-run standard deviation ≤0.28 dB.Source: §F.3 · pp. 21–22, Table 10 · p. 22

F.4 Dependence on spatial sorting of CDF

• At 10 dB, Morton and random ordering are indistinguishable; at 0 dB, Morton improves MSE-5-step to 3.70 versus 6.55, while random ordering has NMSE −18.43 versus −18.27 dB.The NMSE difference is within one standard deviation; overall, Morton offers a modest transition-model advantage at low SNR.Source: §F.4 · pp. 22–23, Table 11 · p. 23

F.5 EROT regularization sensitivity at low SNR

• At 0 dB, increasing EROT’s ε from 0.5 to 1.0 stabilizes training and improves mean performance across metrics, but DSR still outperforms EROT on every measure.Training time is essentially unchanged between the two ε settings.Source: §F.5 · p. 23, Table 12 · p. 23

G Multimodal Posteriors and Inter-modal Mass

• The experiment deliberately separates modes along one coordinate; as dimension grows, informative bits from that coordinate are diluted among uninformative bits in the Morton code.The authors interpret this as loss of locality expected of space-filling curves in higher dimensions.Source: G Multimodal Posteriors and… · p. 23
• With Morton sorting, DSR places about 0.28 inter-modal mass, versus all resampled particles between modes without sorting; EROT places all mass there for d ≥3.The test uses a static two-component Gaussian with modes separated along coordinate 0 and measures resampled mass in |x0| < 3.Source: G Multimodal Posteriors and… · pp. 23–24, Table 13 · p. 24

H Additional Linear Gaussian Experiments

• The linear-Gaussian state-space model permits comparisons of evaluation-time marginal-likelihood decomposition and DSR gradients against an unbiased VSMC gradient.Its known transition model allows optimal filtering with the Kalman filter.Source: H Additional Linear Gaussian Experiments · p. 24

H.1 Marginal likelihood at evaluation time

• DSR’s Jensen’s gap is 2.196 against 3.352 for hard SR, while hard SR has smaller model error, 0.428 versus 0.575 for DSR.All four methods are within 0.7 of the true evidence; the Jensen’s-gap result indicates more accurate finite-particle approximation of the learned model’s evidence.Source: §H.1 · pp. 24–25, Table 14 · p. 24

H.2 Comparison against unbiased VSMC gradients

• The gradient comparison estimates an unbiased reference by averaging over bootstrap noise draws and exhaustively summing systematic-resampling ancestor vectors.The experiments use N = 4 and T = 4 on one- and two-dimensional linear-Gaussian models.Source: §H.2 · p. 25
• DSR has smaller gradient bias than hard SR at every tested temperature in both dimensions, and at τ = 1.0 bias falls from 2.324 to 1.256 for d = 2 and from 3.627 to 2.869 for d = 4.Further lowering τ reduces bias but sharply increases variance; at τ = 0.1, DSR is competitive on MSE, with the lowest MSE at d = 4 and second-lowest at d = 2.Source: §H.2 · p. 25, Table 15 · p. 25