Switching stochastic differential equations (SSDEs) describe continuous-time dynamics whose parameters switch according to a latent regime process that follows a continuous-time Markov chain (CTMC). By allowing dynamics to change between regimes, SSDEs represent heterogeneous system behavior and have been applied across diverse fields. However, Bayesian inference for SSDEs remains difficult, and existing SSDE inference methods have limited applicability, with restrictions such as noise-free observations, univariate states, linear drift, or state-independent diffusion. In this study, we propose an approximate Markov chain Monte Carlo sampler for SSDEs using uniformization and factorized neural likelihood estimation (FNLE), a simulation-based inference method. Uniformization provides an exact representation of the CTMC but requires SDE transition densities over arbitrary time intervals. We approximate these densities by training a time-conditioned FNLE model. The resulting sampler is broadly applicable to SSDEs without requiring analytically tractable transition densities. In synthetic-data experiments, our method recovered regime paths and parameters for three SSDE models for which previous methods have limited applicability. We also applied our method to a real dataset and detected a regime transition.
Figures & tables
Figure 1: Illustration of a simulated Lotka–Volterra SSDE. The upper panel shows the two components of the continuous state and their noisy observations. Background shading indicates the active regime, and the lower panel shows the regime path.
Method
Observation error
Nonlinear drift
State-dependent diffusion
Intractable transition density
Multivariate states
Blackwell et al. [3]
△
△
△
×
✓
Hibbah et al. [13]
△
✓
✓
✓
△
Köhs et al. [19]
✓
△
×
△
✓
Köhs et al. [20]
✓
×
×
×
✓
Stumpf-Fétizon et al. [32]
△
✓
✓
✓
△
Proposed method
✓
✓
✓
✓
✓
Table 1: Applicability of SSDE Inference Methods
Figure 2: Comparison of inference results for the OU experiment using an exact MCMC sampler and two approximate MCMC samplers with FNLE trained on exact transitions and Euler–Maruyama simulations, respectively.
Model
Observation error
Nonlinear drift
State-dependent diffusion
Intractable transition density
Multivariate states
Ornstein–Uhlenbeck
✓
×
×
×
×
Lotka–Volterra
✓
✓
✓
✓
✓
Gene-expression CLE
✓
×
✓
✓
✓
Susceptible–infected–recovered epidemic model
✓
✓
✓
✓
✓
Table 2: Characteristics of the SSDE Models Used to Generate the Synthetic Datasets
Figure 3: Inference results for synthetic data from three SSDEs with regime-specific Lotka–Volterra, gene-expression chemical Langevin, and susceptible–infected–recovered epidemic models (left to right).
Figure 4: Real-data inference results for the Didinium – Paramecium time series using a Lotka–Volterra SSDE.
Appendix figures & tables13 assets
Supplementary material from the paper’s appendix.
Appendix
Experiment
nmax
nmaxh
Synthetic-data OU (approx. MCMC with Euler–Maruyama FNLE)
100
1.00
Synthetic-data LV
100
1.00
Synthetic-data gene-expression CLE
25
0.25
Synthetic-data susceptible–infected–recovered
20
0.20
Real-data LV
98
0.98
Appendix
Table 3: Experiment-Specific FNLE Training Settings.
Experiment
Component
Prior
Synthetic-data OU
scalar state
τ2∼InvGamma(10−3,10−3)
Synthetic-data LV
predator
τ12∼InvGamma(10−3,10−3)
Synthetic-data LV
prey
τ22∼InvGamma(10−3,10−3)
Synthetic-data gene-expression CLE
mRNA
τ12∼InvGamma(10−3,10−3)
Synthetic-data gene-expression CLE
protein fluorescence
τ22∼InvGamma(10−3,10−3)
Synthetic-data susceptible–infected–recovered
susceptible
τ12∼InvGamma(10−3,10−3)
Appendix
Table 4: Priors for the Observation Variances
Experiment
State path y1:L
Parameters Θ
Synthetic-data OU (exact MCMC)
2−8
2−5
Synthetic-data OU (approx. MCMC with exact FNLE)
2−8
2−5
Synthetic-data OU (approx. MCMC with Euler–Maruyama FNLE)
2−8
2−5
Synthetic-data LV
2−9
2−8
Synthetic-data gene-expression CLE
2−9
2−8
Synthetic-data susceptible–infected–recovered
2−9
2−8
Appendix
Table 5: Fixed NUTS Step Sizes by Experiment and Sampler
Experiment
Chains
Sweeps
Burn-in
Retained
Synthetic-data OU (exact MCMC)
4
10,000
5,000
5,000
Synthetic-data OU (approx. MCMC with exact FNLE)
4
10,000
5,000
5,000
Synthetic-data OU (approx. MCMC with Euler–Maruyama FNLE)
4
10,000
5,000
5,000
Synthetic-data LV
4
10,000
5,000
5,000
Synthetic-data gene-expression CLE
4
10,000
5,000
5,000
Synthetic-data susceptible–infected–recovered
4
10,000
5,000
5,000
Appendix
Table 6: MCMC Chain Counts and Lengths Used for Posterior Summaries. Sweep, burn-in, and retained-draw counts are per chain.
Model
Parameter
Regime 1
Regime 2
Structure
OU
κ
0.15
0.15
shared
μ
-1.2
1.2
regime-specific
σ
0.05
0.05
shared
LV
α
1.0
0.5
regime-specific
β
1.0
0.5
regime-specific
γ
1.0
0.2
regime-specific
Appendix
Table 7: Data-Generating SDE Parameters for the Synthetic Experiments. Regime-specific values are listed for regimes 1 and 2.
Model
Parameter
Structure
FNLE range
Prior
OU
κ
shared
[0.05,0.2]
LogNormal(0,1)
μ
regime-specific
[−2.0,2.0]
Normal(0,1)
σ
shared
[0.01,0.1]
LogNormal(0,1)
LV
α
regime-specific
[0.01,2]
LogNormal(0,1)
β
regime-specific
[0.1,3]
LogNormal(0,1)
γ
regime-specific
[0.01,2]
LogNormal(0,1)
Appendix
Table 8: SDE Parameter Inference Settings for the Synthetic Experiments. Parameter ranges and priors are on the physical scale.
Parameter
Structure
FNLE range
Prior
α
regime-specific
[0.1,4]
LogNormal(0,1)
β
regime-specific
[0.01,6]
LogNormal(0,1)
γ
regime-specific
[0.1,8]
LogNormal(0,1)
δ
regime-specific
[0.01,5]
LogNormal(0,1)
σ1
shared
[0.01,0.8]
LogNormal(−1,1)
σ2
shared
[0.01,2]
LogNormal(−1,1)
Appendix
Table 9: SDE Parameter Settings for the Real-Data LV Analysis. Parameter ranges and priors are on the physical scale.
Figure 5: Multi-chain diagnostics for the real-data Lotka–Volterra analysis. Left panels show chain-specific marginal posterior density estimates, and right panels show traces over all 10,000 sweeps. Colors distinguish the four chains; solid and dotted lines distinguish regimes 1 and 2 for regime-specific parameters. Vertical dashed lines mark the burn-in cutoff at 5,000 sweeps. The observation-noise panels show standard deviations τ1,τ2 , not variances.
OU sampler
Parameter
R^
Bulk ESS
Tail ESS
Exact MCMC
κ
1.015
198.4
1,004.2
μ1
1.004
1,187.8
6,199.0
μ2
1.005
1,478.3
5,050.6
σ
1.004
1,519.8
2,234.7
q12
1.000
16,311.0
19,035.3
q21
1.000
16,972.1
18,373.4
Appendix
Table 10: Per-parameter MCMC diagnostics for the three OU samplers on synthetic data.
Parameter
R^
Bulk ESS
Tail ESS
α1
1.006
421.2
1,090.5
α2
1.007
406.2
777.6
β1
1.006
395.9
973.0
β2
1.006
385.9
771.1
γ1
1.028
217.1
893.4
γ2
1.038
143.4
301.5
Appendix
Table 11: Per-parameter MCMC diagnostics for the synthetic Lotka–Volterra experiment.
Parameter
R^
Bulk ESS
Tail ESS
ρ1
1.003
805.9
1,595.1
ρ2
1.001
1,274.5
3,142.3
β
1.011
211.0
333.4
γ
1.007
334.2
691.6
δ
1.007
335.3
681.7
c
1.018
201.7
358.9
Appendix
Table 12: Per-parameter MCMC diagnostics for the gene-expression CLE experiment.
Parameter
R^
Bulk ESS
Tail ESS
β1
1.019
481.9
1,013.8
β2
1.054
61.6
338.2
γ1
1.065
66.2
485.9
γ2
1.003
645.6
1,278.9
q12
1.000
15,315.1
18,870.5
q21
1.001
10,209.2
14,825.5
Appendix
Table 13: Per-parameter MCMC diagnostics for the susceptible–infected–recovered experiment.
Parameter
R^
Bulk ESS
Tail ESS
α1
1.004
442.4
987.0
α2
1.005
929.8
909.5
β1
1.003
466.3
1,065.3
β2
1.005
871.0
888.9
γ1
1.011
298.2
732.9
γ2
1.007
1,160.2
2,167.1
Appendix
Table 14: Per-parameter MCMC diagnostics for the real-data Lotka–Volterra analysis.
Stochastic differential equations (SDEs) provide a flexible framework for modeling temporal dynamics in partially observed systems. A central task is to calibrate such models from data, which requires inferring latent trajectories and parameters from sparse, noisy observations. Classical smoothing methods for this problem are often limited by path degeneracy and poor scalability. In this work, we developed a novel method based on characterization of the posterior SDE in terms of conditional backward-in-time score defined as the gradient of a function solving a Kolmogorov backward equation with multiplicative updates at observation times. We learn this conditional score using neural networks trained to satisfy both the governing PDE and the observation-induced jump conditions, thereby integrating continuous-time dynamics with discrete Bayesian updates. The resulting score induces a posterior SDE with the same diffusion coefficient but a modified drift, enabling efficient posterior trajectory sampling. We further derive a likelihood-based objective for learning the SDE parameters, yielding an evidence lower bound (ELBO) for joint state smoothing and parameter estimation. This leads to a variational EM-style procedure, where the neural conditional score is optimized to approximate the smoothing distribution, followed by a maximization step over the SDE parameters using samples from the induced posterior. Experiments on nonlinear systems demonstrate accurate and stable inference with a very few observations demonstrating significant improved scalability compared to classical MCMC methods.
Yu Wang, Arnab Ganguly
Department of Mathematics Louisiana State University
Recovering dynamical systems from noisy observations is a recurring challenge across scientific domains, including neuroscience and physics. Latent stochastic differential equations (SDEs) address this by modeling the system as an unobserved state that evolves according to a learnable SDE and generates the observations. Variational inference (VI) provides a tractable objective for fitting latent SDEs. Traditional VI algorithms evaluate this objective by numerical simulation over a time discretization, trading fidelity for computational cost. A recent class of algorithms, simulation-free VI, sidesteps this tradeoff by parameterizing the posterior through its instantaneous marginals rather than its drift. In this work, we show that the efficiency of existing simulation-free VI algorithms comes at a price: their parameterizations restrict the approximate posterior to a subset of the SDEs available to simulation-based methods, degrading posterior inference and parameter learning. We propose Helmholtz-SDE, a simulation-free VI algorithm that closes this gap by optimizing over path laws compatible with a prescribed collection of marginals. Helmholtz-SDE recovers dynamics more faithfully than prior simulation-free methods, with the largest gains under high posterior uncertainty. It further matches the performance of simulation-based VI at a fraction of the runtime. Code is available at https://github.com/lindermanlab/helmholtz-sde .
Henry D. Smith, Brian L. Trippe, Scott W. Linderman
We introduce an amortized neural sampler that combines operator learning with flow methods for sampling. It maps SDE coefficient functions to pushforwards from a reference measure to the invariant measures, enabling efficient sampling across families of stochastic differential equations. Our framework shifts traditional sampling cost to an initial training phase, after which new SDE instances require only one encoder pass and a few ODE solver steps, independent of mixing time. To handle problems in high dimensions, we use Lagrangian trajectory sensors for the coefficient functions and cross attention in the architecture. We also theoretically establish the expressivity and resolution invariance of our framework. Experiments on 1D and 2D SDE families show competitive accuracy with substantial speedups over MCMC in regimes with slow mixing, transfer across sensor counts, and demonstration results on a 64D interacting particle SDE where traditional grid approaches are infeasible.
Ling Guo, Lei Li, Jingtong Zhang
Department of Mathematics, Shanghai Normal University, Shanghai, China · School of Mathematical Sciences, Shanghai Jiao Tong University, Shanghai, China · Institute of Natural Sciences, MOE-LSC, Shanghai Jiao Tong University, Shanghai, China