Source-linked AI summary
Symbolic Neural ODEs: Learning interpretable models from time-series data
Nibodh Boddupalli, Jeff Moehlis
TL;DR
Learning dynamical systems from time-series data requires models that remain coherent under repeated composition rather than merely fitting local transitions. The paper uses sparse neural ODEs trained with a finite-horizon multi-step loss, obtaining interpretable models with improved stability and agreement with both short-term dynamics and long-term statistics. For chaotic systems, it separates finite-horizon trajectory accuracy from long-term statistical fidelity and connects the two through theoretical bounds.
Problem
Data-driven dynamical modeling must address the limitation that one-step regression does not ensure correct behavior under repeated composition, especially for nonlinear and chaotic systems.
Method
The framework parameterizes the vector field with a sparse neural architecture and trains it using a finite-horizon multi-step prediction loss with progressively increasing horizon.
Results
Multi-step training produces sparse, interpretable models with improved stability and agreement with short-term dynamics and long-term statistics across fixed points, periodic orbits, and chaotic attractors.
Takeaways & Limitations
Aligning the training objective with repeated composition supports models that reproduce qualitative and statistical features of the underlying dynamics, not only local data fits.
Takeaways & Limitations
For chaotic systems, sensitivity to initial conditions makes long-term trajectory-level agreement inherently impossible, so evaluation must distinguish finite-horizon trajectories from long-term statistics.
Abstract
from arXiv · showhide
We present a machine learning framework for identifying sparse, interpretable models of dynamical systems directly from time-series data. Our approach parameterizes the underlying vector field using a neural architecture and trains it by minimizing a multi-step prediction loss over a finite horizon. To ensure numerical tractability, we optimize a mean absolute error objective averaged across prediction steps, and progressively increase the horizon during training. A key feature of this formulation is that it enforces consistency under repeated composition of the learned dynamics. As a result, the identified models exhibit significantly improved stability compared with approaches based on one-step regression of the vector field. When combined with sparsity-promoting regularization, this leads to parsimonious models that generalize beyond the training data. We demonstrate accurate recovery of systems exhibiting a wide range of behaviors, including stable and unstable fixed points, periodic orbits, and chaotic attractors. For chaotic systems, while long-term trajectory prediction is inherently limited by sensitivity to initial conditions, we show that multi-step training yields models with accurate short-term dynamics and strong agreement in long-time statistical properties, including mean, variance, and Lyapunov exponents. Moreover, we establish theoretical bounds linking trajectory error to statistical accuracy, providing a step toward a principled explanation for this behavior.
1. Introduction.
The paper addresses the challenge of learning interpretable dynamical models from data when local, one-step fitting may fail under repeated composition. It proposes multi-step training and demonstrates stable, sparse models across fixed points, periodic orbits, and chaotic attractors.
- Motivation: Existing data-driven approaches often minimize one-step prediction error, which does not guarantee correct behavior under repeated composition and can produce instability.Small local errors may accumulate over time, especially in nonlinear and chaotic systems.
- Contribution: Multi-step training optimizes trajectory agreement across a finite horizon, enforces compositional consistency, and progressively increases the horizon for tractability.The objective uses mean absolute error averaged across prediction steps.
- Contribution: Combining multi-step prediction with sparsity-promoting regularization yields compact, interpretable representations that generalize beyond the training data.The framework is designed to retain interpretability while learning nonlinear dynamics.
- Results: For chaotic systems, multi-step training supports accurate finite-horizon trajectories and agreement in statistical quantities such as means, variances, and Lyapunov exponents.The paper provides theoretical bounds linking trajectory error to statistical error.
- Results: Experiments recover systems with stable and unstable fixed points, periodic orbits, and chaotic attractors.Increasing the prediction horizon improves stability, robustness to noise, and agreement with short-term dynamics and long-term statistics.
2. Problem formulation.
The problem formulation represents sampled dynamical-system data through a stroboscopic map and neural ODE integration. A NODE horizon defines multi-step forecasting, allowing the learned vector field to be evaluated for local accuracy and repeated-composition stability.
- Data and dynamics: The formulation considers finite-dimensional autonomous dynamical systems observed at discrete time intervals, with states indexed by their components.The state space may include manifolds, while most examples use R^n.
- Data and dynamics: The dataset X0 contains initial states for which observations are available through H future time steps, and Xh collects the corresponding states at each horizon step.Sequential data from one trajectory can be treated as multiple initial-value problems.
- Data and dynamics: The stroboscopic map S advances the true system by one sampling interval, while the learned map ˜S is induced by numerical integration of the learned vector field.Forecasting multiple steps corresponds to repeated composition of the learned map.
- Stability under composition: Repeated composition can amplify small one-step discrepancies, so the formulation trains the model to account for error propagation beyond a single step.This issue is especially pronounced in nonlinear and chaotic systems.
- Optimization horizon: The NODE horizon H is the number of time steps over which the model is optimized for trajectory agreement and repeated-composition behavior.Increasing H makes optimization more sensitive to error growth, favoring stable long-term behavior.
- Neural ODE approach: SymNODE represents the learned vector field with a sparse, interpretable neural network and numerically integrates it using a fixed-step fourth-order Runge–Kutta scheme.The explicit numerical integration produces learned states as functions of the current state.
3. Symbolic Neural Network.
The SymNN uses stacked operational layers built from primitive functions to parameterize a symbolic neural ODE, which is recursively integrated and trained against multi-step trajectories. Increasing the optimization horizon penalizes compounding errors and favors models with accurate, stable behavior over longer intervals.
- Architecture: SymNN organizes K stacks of L operational layers, with each stack receiving outputs from preceding stacks and each layer generating four primitive-function outputs.Each operational layer takes the state, a constant scaling value, and previous-stack outputs as input.
- Neural ODE approach: SymNODE extends SymNN to continuous-time systems by numerically integrating the learned vector field with a fixed-timestep fourth-order Runge–Kutta scheme.The integrated state is initialized from x[0] and recursively produces predicted states ˜x[h].
- Optimization: Training minimizes mean absolute error averaged across the entire forecasting horizon rather than only the endpoint, encouraging trajectory tracking throughout prediction.Endpoint-only optimization can fit the final state while failing to follow the intervening trajectory.
- Effect of optimization horizon: Extending the NODE horizon makes the objective more sensitive to error growth, penalizing models with compounding trajectory errors and favoring stable long-term behavior.The theoretical results connect longer-horizon optimization with longer-term prediction optimality and control of average multi-step trajectory error.
- Regularization and limitations: The framework combines trajectory-based training with sparsity-promoting regularization to support parsimonious models, while its theoretical guarantees are tied to the corresponding training loss.Practical mAIC-based conclusions preserve the qualitative horizon dependence but do not retain the theoretical guarantees exactly.
4. Examples.
The examples test SymNODE on systems with fixed points, limit cycles, and nonlinear periodic dynamics, emphasizing whether multi-step training recovers globally consistent behavior from informative trajectories. Across these cases, identified models closely match reference dynamics, while data coverage determines whether structurally correct nonlinear terms can be recovered.
- Takens–Bogdanov normal form: The Takens–Bogdanov trajectories closely followed the true dynamics from the same initial condition, recovering the governing equation despite noisy data influenced by the saddle’s unstable manifold.The phase portrait and time-series comparison both showed close agreement.
- Hopf normal form: Sufficient transient behavior enabled exact recovery of the Hopf normal form, including exponential radial decay and linear phase growth toward the stable limit cycle.When data were restricted near the limit cycle, simpler linear models fit short-term behavior but missed the correct radial dynamics.
- Approximate periodic-orbit models: For stable periodic orbits, structurally different approximate models can still retain similar stability, dynamic range, and time-period to the exact model.The corresponding loss remains small over increasingly long horizons when the periodic orbit is only slightly perturbed.
- Simple pendulum: Multiple pendulum trajectories with different amplitudes supplied the variation needed to identify nonlinear amplitude dependence that a single trajectory could obscure.Single-trajectory data can support a locally accurate linear model with low short-term error but poor generalization across initial conditions.
4.5. Van der Pol oscillator.
The Van der Pol experiments examine whether SymNODE can identify stiff, multiscale relaxation oscillations. The method recovers the exact model across moderate and strongly separated time scales when data resolution and the NODE horizon are appropriate.
- Moderate stiffness: Approximately 600 data points from 1–3 trajectories with 1% NSR were sufficient to recover the exact Van der Pol model for moderate µ values between 3 and 6.These regimes were only mildly stiff.
- Strong scale separation: For µ = 15, 1200 noisy points from one trajectory and a 16-step horizon identified the exact model using 68 parameters after 800 epochs.The experiment used K = 2 stacks and L = 1 layer despite limited transient information.
- Strong scale separation: The µ = 15 model reproduced relaxation oscillations, nonlinear damping, transient behavior, and convergence to the stable periodic orbit.Agreement was evident in both the phase portrait and time-series comparison.
- Practical conditions: Accurate identification under strong nonlinearities and multiple time scales requires appropriately chosen data resolution and NODE horizon.In stiff regimes, numerical integration errors can accumulate rapidly under repeated composition.
- Practical conditions: The approach does not always recover the true model structure, although stable periodic orbits can preserve nearly matching dynamic range, stability characteristics, and time-period.This boundary is illustrated by structurally different models with similar asymptotic behavior.
4.6. Chaotic dynamics: sensitivity, horizon length, and statistical fidelity.
Chaotic systems expose the trade-off between finite-horizon trajectory accuracy and long-term statistical fidelity. Increasing the prediction horizon improves model stability and can recover accurate statistics even when trajectories eventually diverge.
- Sensitivity and horizon length: Chaotic dynamics cause rapid trajectory divergence through sensitivity to initial conditions and error growth under repeated composition.
- Sensitivity and horizon length: Longer horizons improve robustness to noise, with the Rössler model identified exactly from data containing up to 20% NSR.
- Sensitivity and horizon length: Increasing the horizon to H ≥32 leads to convergence to the exact model, whereas smaller horizons can become trapped in suboptimal minima.
- Chaotic trajectory fidelity: For Chua’s system, multi-step training reproduces attractor geometry and dynamics accurately despite structural differences from the true model.
- Chaotic trajectory fidelity: The Chua model tracks trajectories for about 10 Lyapunov time units before diverging, while its estimated means and standard deviations remain close to the true values.
- Statistical fidelity: For Lorenz dynamics, the increased-horizon model nearly matches the true means, standard deviations, and Lyapunov exponents.
- Limitations under noise: Higher noise can produce structurally different Lorenz models with small spurious constants, because averaged multi-step loss weakly penalizes small biases.
- Statistical fidelity: Across chaotic examples, finite-horizon tracking is followed by divergence, while statistical quantities remain accurate because learned models capture the attractor’s invariant measure.
5. Discussion and Summary.
The discussion argues that multi-step training aligns model fitting with repeated dynamical composition, improving stability, interpretability, and recovery of qualitative and statistical behavior. It also identifies computational cost, chaotic trajectory divergence, and noisy-data symmetry violations as important boundaries and future-work targets.
- Core perspective: Multi-step losses enforce accuracy under repeated composition rather than only local, one-step behavior.They penalize errors that accumulate across a finite prediction horizon.
- Core perspective: The LossH framework distinguishes local approximation from dynamical fidelity, especially when nonlinear or chaotic trajectories are iterated.Proposition 3.5 formalizes the connection between finite-horizon trajectory agreement and the training objective.
- Chaotic systems: For chaotic systems, the appropriate targets are finite-horizon trajectory accuracy and long-term statistical agreement rather than indefinite pointwise trajectory matching.The discussion names invariant statistics including means, standard deviations, and Lyapunov exponents.
- Model properties: Multi-step training acts as an implicit regularizer, improving robustness to noise and, with sparsity regularization, producing interpretable and dynamically consistent models.The regularization suppresses models that fit locally but behave unstably or non-physically under iteration.
- Limitations and future work: Multi-step training requires repeated integration and backpropagation through compositions, increasing computational graphs and training time.The trade-off can be managed through horizon length, architecture, and optimization strategy.
- Limitations and future work: Future work includes adaptive horizon selection, partial-observation formulations, and preserving symmetries in models estimated from noisy data.The proposed extensions include latent-variable and delay-embedding approaches.
A.1. Systems with Stable Fixed Points.
For systems with stable hyperbolic fixed points, the appendix analyzes bounded perturbations through linearization and Lyapunov arguments. It concludes that perturbed solutions asymptotically remain close to the unperturbed fixed point.
- Assumptions: A stable hyperbolic fixed point has a Jacobian whose eigenvalues all have negative real parts.The appendix assumes this stability condition before analyzing perturbations.
- Perturbation analysis: Linearizing the perturbed system near the fixed point yields dynamics driven by the stable matrix A and the perturbation term g(x∗).Higher-order terms are ignored for the leading-order analysis.
- Stability argument: Because A is Hurwitz, a quadratic Lyapunov function establishes asymptotic stability of the linearized unperturbed system.The derivative of the Lyapunov function is negative away from the origin.
- Conclusion: Solutions of the perturbed system asymptotically stay close to the fixed point of the unperturbed system.The initial-condition norm is absorbed into the bound's constant.
A.2. Systems with Stable Periodic Orbits.
For stable hyperbolic periodic orbits, the appendix separates phase motion from transverse isostable dynamics to analyze bounded perturbations. The resulting bounds keep solutions near the original orbit while allowing phase drift, consistent with persistence under small perturbations.
- Assumptions: A stable hyperbolic periodic orbit has non-trivial Floquet exponents with negative real parts.These exponents govern transverse stability around the orbit.
- Coordinate reduction: Phase and isostable coordinates describe perturbation effects along the orbit and transverse to it, respectively.The phase response curve and isostable response curves appear in the reduced equations.
- Perturbation analysis: The transverse dynamics reduce to a stable linear system with bounded forcing in the isostable coordinates.The vector form uses y and the diagonal matrix Λ, with forcing b(θ).
- Conclusion: The perturbed asymptotic solution remains in an O(∥g∥) tube around the unperturbed periodic orbit.Solutions may nevertheless drift along the neutrally stable phase direction at rate O(∥g∥).
- Numerical verification: Numerical examples for Hopf normal-form and Van der Pol models remain close to the true periodic orbit, except for short-lived fast-dynamics deviations in the Van der Pol case.The observation is consistent with a small perturbation g(x).
Appendix B. Heuristic analysis of trajectory error growth.
The heuristic analysis explains multi-step trajectory-error growth by separating direct modeling-error forcing from amplification through linearized dynamics. Repeated composition can therefore magnify small local errors in unstable or chaotic systems.
- Setup: The derivation is heuristic and uses shared initial conditions to study how modeling errors accumulate under repeated composition.It is intended to explain why multi-step errors may grow with forecast horizon.
- Linearization: The trajectory error is defined by comparing the true and predicted solutions, with the predicted state written as x + δ.The analysis formally linearizes the modeling error for small δ.
- Error mechanisms: The approximate error dynamics separate direct forcing from amplification of existing discrepancies by the linearized dynamics.This decomposition follows after neglecting higher-order terms and taking norms.
- Implication: In unstable or chaotic systems, repeated composition of the learned flow map can amplify small local modeling errors into longer-horizon prediction-error growth.The mechanism combines modeling-error forcing with dynamical amplification.
Appendix C. Bounds on Statistics.
Appendix C establishes bounds connecting trajectory RMSE with errors in time-averaged means, variances, and standard deviations. These results show that accurate trajectory reconstruction supports accurate first- and second-order statistics under bounded-trajectory assumptions.
- RMSE upper-bounds the ℓ2 error between the true and predicted trajectory means.The bound is ∥µ −˜µ∥2 ≤ RMSE.
- Under uniformly bounded trajectories, variance error is bounded by 2(M + ˜ M) RMSE.The result is stated in equation (C.9).
- The same variance bound appears after decomposing quadratic and mean terms and applying Cauchy-Schwarz.The proof uses a^2 −b^2 = (a −b)(a + b) and trajectory bounds.
- Standard-deviation control is weaker because the square-root mapping is not Lipschitz near zero.The resulting bound depends on trajectory bounds and dimension, and requires additional assumptions for stronger control.
- Trajectory RMSE therefore controls first- and second-order statistics without explicitly penalizing those statistics during training.The appendix relates accurate trajectory reconstruction to accurate time-averaged means and variances.
Appendix D. Computation details.
The computation details use five-fold training, with each fold allocating 80% of the data for training and 20% for testing. The authors report that these workloads can run on consumer-grade computing hardware, while larger workloads may benefit from stronger infrastructure.
- The reported network and dataset sizes allow the folds to run on consumer-grade CPUs and GPUs.Examples include laptop and desktop systems with Nvidia GPUs.
- Larger networks, longer horizons, larger datasets, or extensive model distillation can benefit from PyTorch and Python workloads on stronger infrastructure.