Identifying ODEs from Unstructured Data with Causal Representation Learning

Alessandro Trenta; Riccardo Massidda; Davide Bacciu; Sara Magliacane

cs.LGSep 29, 2026

Identifying ODEs from Unstructured Data with Causal Representation Learning

Authors: Alessandro Trenta, Riccardo Massidda, Davide Bacciu, Sara Magliacane

Organizations: Department of Computer Science University of Pisa · Saarland Informatics Campus Saarland University & University of Amsterdam

Abstract

We study the problem of recovering the governing ODE of a dynamical system from unstructured, high-dimensional observations such as images. Existing methods for ODE discovery typically assume direct measurements of the variables, or do not provide theoretical guarantees on the learned variables and equations. While Causal Representation Learning (CRL) methods provide guarantees on identifying variables from high-dimensional observations up to component-wise diffeomorphisms, we show that in general these variables cannot be used directly as input to equation discovery methods, which typically assume that the variables will lead to sparse equations. So we introduce SParse Equivalent Equation Discovery AutoEncoder (SPEED-AE), a framework that combines a pretrained CRL method with a component-wise autoencoder that learns transformations of variables that are amenable to sparse ODE discovery. We show that for polynomial ODEs, this additional step allows us to restrict the identifiability of each variable from polynomial to monomial diffeomorphisms. Experiments on Lotka-Volterra, Lorenz, and a two-pendulum system show that SPEED-AE improves on the disentanglement of the CRL methods and that it recovers ODEs that are closest to the ground truth, while achieving state-of-the-art forecasting performance.

Figures & tables

Figure 1: SPEED-AE pipeline: we observe high-dimensional data from an invertible mapping ${\mathbf{x}}(t)=g({\mathbf{z}}(t))$ of the true latents ${\mathbf{z}}(t)$ . First, SPEED-AE employs a pretrained CRL model to obtain low-dimensional latents $\hat{{\mathbf{z}}}(t)$ ( red block). We assume that the recovered latents correspond to the ground-truth ${\mathbf{z}}$ up to a permutation and component-wise diffeomorphism $z_{i}=\hat{h}_{i}(\hat{z}_{\pi(i)})$ . For simplicity, in the figure we assume that the permutation is the identity function. SPEED-AE learns a set of component-wise autoencoders $(\phi^{\text{enc}}_{i},\phi^{\text{dec}}_{i})$ and an ODE by combining the functions in a library $\Theta(\tilde{{\mathbf{z}}})$ , imposing sparsity of this representation ( blue block). The output latents $\tilde{{\mathbf{z}}}(t)$ identify the ground-truth up to permutation and monomial and, in most cases, linear transformation.

Figure 2: Comparison between SINDyAE, SINDy applied to CRL latents (using CITRIS), and SPEED-AE (with CITRIS as CRL method and SINDy as equation discovery method) on the Lotka-Volterra dataset ( $\dot{z_{1}}=0.66z_{1}-1.33z_{1}z_{2}$ , $\dot{z_{2}}=-z_{2}+z_{1}z_{2}$ ). For each model, the left box reports scatter plots of the true versus learned latents. We then compare a test trajectory with the one predicted from only the initial condition and the discovered equation, after converting both to the true latent space, following the evaluation protocol in Section 6 . SPEED-AE is the only method that recovers a sparse equation and achieves strong long-term prediction performance.

Figure 3: Learned coefficients for Lotka-Volterra for each polynomial of the state variables (each row) and different methods (each column, where the last four are versions of SPEED-AE). The first column shows the ground truth coefficients and the black boxes highlight the non-zero coefficients.

Figure 4: Results for the Lotka-Volterra and Lorenz experiments with mean and standard deviation over

10

seeds. For all metrics, lower is better. We limit the height of the boxes for better readability, reporting the actual scales of the outlier results inside their bars.

Model	CorrD	$\mathbf{x}$ Err	$\mathbf{z}$ Err	Recovered ODE
SINDyAE	$0.431$	$18.888$	$1.621$	$\begin{gathered}\ddot{z}_{1}-0.302\ddot{z}_{1}z_{1}+1.20\ddot{z}_{1}z_{1}^{2}+2.39\ddot{z}_{1}\dot{z}_{1}z_{1}=\\ =-0.32+0.75z_{2}-0.43z_{1}+0.17z_{1}^{3}+0.06z_{1}^{2}+\ldots\\ \ddot{z}_{2}+0.01\ddot{z}_{2}z_{2}+0.03\ddot{z}_{2}z_{2}^{2}+0.05\ddot{z}_{2}\dot{z}_{2}z_{2}=\\ =0.03+0.18z_{1}-0.07z_{1}^{3}-0.03z_{1}^{2}+\ldots\end{gathered}$
CITRIS+SINDy	$0.020$	$6.824$	$0.239$	$\begin{gathered}\ddot{z}_{1}-0.15\dot{z}_{1}^{2}+0.15\ddot{z}_{1}z_{1}+0.02\dot{z}_{1}^{2}z_{1}^{2}-0.012\dot{z}_{2}^{2}z_{2}^{2}+\ldots=\\ =0.64+11.31z_{1}-131z_{1}^{2}+\cos(0.14\dot{z}_{1}z_{1}-0.93z_{1})+\ldots\\ \ddot{z}_{2}-0.16\ddot{z}_{2}z_{2}-0.16\dot{z}_{2}^{2}+0.32\dot{z}_{2}z_{2}+\ldots=\\ =1.48+4.03z_{2}-0.04z_{2}^{2}+\cos(0.14z_{2}\dot{z}_{2}+0.84z_{2}-0.43)\ldots\end{gathered}$
SPEED-AE (C+S)	$0.004$	$6.155$	$0.177$	$\begin{gathered}\ddot{z}_{1}=-0.27z_{1}-0.61\sin(1.05z_{1})\\ \ddot{z}_{2}=-0.27z_{2}-0.61\sin(1.07z_{2})\end{gathered}$

Table 1: Results for the Pendulum experiment for the model.

Appendix figures & tables18 assets

Supplementary material from the paper’s appendix.

Appendix

Parameter	Lotka-Volterra	Lorenz	Pendulum
$n_{\text{train}}$	$20$	$1024$	$100$
$n_{\text{val}}$	$5$	$256$	$25$
$n_{\text{test}}$	$5$	$256$	$25$
$T$	$50$	$5$	$10$
$\Delta t$	$0.01$	$0.02$	$0.02$
num. steps	$5000$	$250$	$500$

Appendix

Table 2: Data generation parameters.

Model	Training Time (mins)
LSTM	$\approx 45$
Transformer	$\approx 62$
SINDyAE	$\approx 100$
MNNAE	$\approx 191$
Autoencoder (CITRIS-NF)	$\approx 22$
CITRIS-NF	$\approx 52$

Appendix

Table 5: Training time for each model on the Lotka-Volterra experiment.

Figure 5: Evaluation pipeline of our experiments.

Figure 6: Experimental results for the ablation on the independence from CRL method. Lotka-Volterra dataset with an additional version of SPEED-AE with an oracle CRL model that achieves perfect disentanglement and SINDy-based loss.

Panel of Figure 7 — Figure 7: Results on the synthetic experiments with perfect disentanglement, polynomial identifiability in the CRL latents, and SPEED-AE encoder and decoders implemented as polynomials. (a) Some of the runs show strong instabilities, while others (b) remain stable but with flat dynamics. In general, the polynomial functions are hard to optimize simultaneously to the dynamics.

Figure 8: Experimental results for the Lotka-Volterra and the Lorenz experiments.

Model	Lotka-Volterra	Lotka-Volterra	Lorenz (Stable)	Lorenz (Chaotic)
	$\alpha=\frac{2}{3},\beta=\frac{4}{3}$	$\alpha=\frac{6}{5},\beta=\frac{3}{2}$	$\rho=14$	$\rho=28$
LSTM	13.860 ±2.609	29.313 ±14.719	0.157 ±0.056	1.624 ±0.210
Transformer	16.001 ±2.076	25.999 ±6.919	0.422 ±0.076	1.902 ±0.210
SINDyAE	23.075 ±11.611	17.610 ±2.898	1.988 ±0.323	3.137 ±0.120
MNNAE	1.257 ±0.240	0.930 ±0.290	0.448 ±0.102	$\approx 10^{5}$ (failed)
SPEED-AE (C+S)	1.279 ±0.483	3.999 ±2.390	0.851 ±0.384	2.411 ±0.379

Appendix

Table 6:

{\mathbf{x}}

Err: mean and standard deviation across

10

seeds.

Panel of Figure 9 — Figure 9: Pearson cross correlation matrices between $\tilde{{\mathbf{z}}}$ and ${\mathbf{z}}$ for all experiments and models. The rightmost matrix in each subfigure corresponds to the true correlation matrix (i.e., between ${\mathbf{z}}$ and ${\mathbf{z}}$ ).

Panel of Figure 10 — Figure 10: Lotka-Volterra experiment with $\alpha=\frac{2}{3},\beta=\frac{4}{3}$ . Scatter plots of the found latents $\tilde{{\mathbf{z}}}$ against the true ones ${\mathbf{z}}$ . If their function is well defined, the plots well represent the maps $z_{i}=h_{i}(\tilde{z}_{i})$ . In the second row, we include the same plots for the CRL latents $\hat{{\mathbf{z}}}$ from CITRIS and DMSVAE. In most cases, SPEED-AE can go from a general diffeomorphism $h$ to a linear one.

Panel of Figure 11 — Figure 11: Lotka-Volterra experiment woth $\alpha=\frac{6}{5},\beta=\frac{3}{2}$ . Scatter plots of the found latents $\tilde{{\mathbf{z}}}$ against the true ones ${\mathbf{z}}$ . If their function is well defined, the plots well represent the maps $z_{i}=h_{i}(\tilde{z}_{i})$ . In the second row, we include the same plots for the CRL latents $\hat{{\mathbf{z}}}$ from CITRIS and DMSVAE. In most cases, SPEED-AE can go from a general diffeomorphism $h$ to a linear one.

Figure 12: Identified coefficients for the Lotka-Volterra experiment with transfer of CRL models and

\alpha=\frac{6}{5},\beta=\frac{3}{2}

. The true ones are in the leftmost column and highlighted by the black borders. For comparison, we include four baselines that consist of equation discovery directly applied to the CRL latents

\hat{{\mathbf{z}}}

Figure 13: Example test trajectory of the Lotka-Volterra experiment (first $10$ dimensions of ${\mathbf{x}}$ ). Ground truth (blue) against the predicted trajectory from the LSTM model (red). During rollout, the model collapses after around 500 steps.

Figure 14: Example test trajectory of the Lotka-Volterra experiment (first $10$ dimensions of ${\mathbf{x}}$ ). Ground truth (blue) against the predicted trajectory from the Transformer model (red). During rollout, the model collapses after a few iterations.

Panel of Figure 15 — Figure 15: Lorenz experiment with stable trajectories $\rho=14$ . Scatter plots of the found latents $\tilde{{\mathbf{z}}}$ against the true ones ${\mathbf{z}}$ . If their function is well defined, the plots well represent the maps $z_{i}=h_{i}(\tilde{z}_{i})$ . In the second row, we include the same plots for the CRL latents $\hat{{\mathbf{z}}}$ from CITRIS and DMSVAE. SPEED-AE models with DMSVAE as disentangler suffer more due to the CRL latents being less disentangled, although the SINDy-based loss is actually able to recover a linear map $h$ .

Panel of Figure 16 — Figure 16: Scatter plots of the found latents $\tilde{{\mathbf{z}}}$ against the true ones ${\mathbf{z}}$ . If their function is well-defined, the plots well represent the maps $z_{i}=h_{i}(\tilde{z}_{i})$ . In the second row, we include the same plots for the CRL latents $\hat{{\mathbf{z}}}$ from CITRIS and DMSVAE. SPEED-AE models with DMSVAE as disentangler suffer more due to the CRL latents being less disentangled, although the SINDy-based loss is actually able to recover a linear map $h$ .

Panel of Figure 17 — Figure 17: Identified coefficients for the Lorenz experiments in the stable ( $\rho=14$ , left) and chaotic case ( $\rho=28$ , right). The true ones are in the leftmost column and highlighted by the black borders. For comparison, we include four baselines that consist of equation discovery directly applied to the CRL latents $\hat{{\mathbf{z}}}$ .

Panel of Figure 18 — Figure 18: Example images from the experiment with two pendulums. Pendulums are colored differently to ensure invertibility of the mixing function.

Panel of Figure 19 — Figure 19: Pendulums experiment. Scatter plots of the found latents $\tilde{{\mathbf{z}}}$ against the true ones ${\mathbf{z}}$ . If their function is well defined, the plots well represent the maps $z_{i}=h_{i}(\tilde{z}_{i})$ . In the second plot, we include the same for the CRL latents $\hat{{\mathbf{z}}}$ from CITRIS. SPEED-AE is able to achieve linear identifiability from the general diffeomorphism of CITRIS. SINDyAE only identifies one variable, while the other likely represents a mapping between the two pendulums.

Explore similar work

Jun 25, 2026cs.LG

Recovering Governing Equations from Solution Data: Identifiability Bounds for Linear and Nonlinear ODEs

Learning governing equations from observed solution data is a fundamental challenge in scientific machine learning, yet the theoretical conditions under which a ground-truth ODE can be uniquely and stably identified from multiple solution observations remain largely undeveloped, and no quantitative analysis of the sample complexity of such learning tasks exists in the literature. To address this gap, we introduce the Hausdorff distance on solution sets as the natural metric for comparing differential equations, since it captures the worst-case separation between two equations over all admissible initial conditions and thus encodes the minimax structure of the identification problem. We establish identifiability bounds for governing ODEs across a wide class of structure equations--ranging from linear ODEs to nonlinear classes with Lipschitz (Hölder)-continuous vector fields--characterizing precisely when two distinct equations can be distinguished from solution data. Using this metric, we derive metric entropy estimates for the relevant ODE classes and analyze sample complexity bounds, quantifying how many solution observations are needed to reliably recover the governing equation.

Yang Pan, Helmut Bölcskei

Ordinary Differential Equations Identifiability

Ask PDF

Jun 23, 2026cs.LG

LLM-ACES: Closed-Loop Discovery of Dynamical Systems with LLM-Guided Adaptive Search

Recovering governing Ordinary Differential Equations (ODEs) from data is a central challenge in modeling dynamical systems across scientific domains. Existing approaches cast discovery as a static inference problem over fixed datasets, assuming that the observed trajectories are sufficiently informative. However, dynamical systems evolve over large state spaces, and limited data can make multiple equations observationally indistinguishable, leading to identifiability gaps and the recovery of incorrect governing equations. To address this, we introduce LLM-ACES, or LLM-guided Active Closed-loop Equation Search, a closed-loop framework that jointly optimizes symbolic hypothesis construction and adaptive data acquisition. In LLM-ACES, a large language model (LLM) proposes operator priors that partition the large search space into distinct regions, within which candidate equations are fit to the observed data. The disagreement among these candidates guides the acquisition of informative trajectories, creating a feedback loop that iteratively refines both the hypothesis space and the discovered dynamics. On 122 ODE systems spanning ODEBench and ODEBase, LLM-ACES achieves the lowest median NMSE, outperforming state-of-the-art baselines by several orders of magnitude while achieving a high symbolic accuracy of 46.2% and 52.4%, respectively. Our analysis further shows that LLM-ACES is sample-efficient, achieving better performance with one-tenth the data. Furthermore, LLM-ACES's feedback-driven data acquisition makes it robust to noise and recovers the correct symbolic structure, while baselines introduce spurious terms that fit the data locally but obscure the true governing relationships.

Nikhil Abhyankar, Sha Li, Sanchit Kabra +3

Virginia Tech

Dynamical Systems Closed Loop Discovery Identifiability

Ask PDF

Jul 15, 2026cs.AI

Automatic Ordinary Differential Equations Discovery For Biological Systems Using Large Language Model Powered Agentic System

Automatic scientific discovery has long been a goal of computational scholars - a machine that can discover nature's secrets on its own, moving computational systems beyond data-fitting tools toward the generation and refinement of mechanistic models of the universe. Recent advances in symbolic regression (SR) and large-language-model (LLM)-based agents suggest that such systems can recover equations from data, incorporate domain priors, and automate parts of the research workflow. However, most existing approaches either focus on narrow equation-discovery benchmarks or broad end-to-end automation pipelines, while biological systems remain comparatively underexplored. Here, we introduce the MEDA system, an LLM- and SR-powered agentic framework for discovering ordinary-differential-equation (ODE) models of biological and biologically inspired dynamical systems. MEDA retrieves background knowledge, defines admissible variables, generates mechanistic constraints, proposes candidate ODEs, and fits and evaluates them. We evaluate it across canonical model retrieval, reasoning-based extrapolation to unseen variants, and open-ended discovery, with and without experimental data. Across these settings, MEDA recovered the correct state variables, achieved strong structural recovery in retrieval and extrapolation tasks, and produced biologically plausible discovery-oriented models. Ablation and robustness analyses show that knowledge-guided formalization and mechanistic constraints are load-bearing components, whereas numerical fitting alone can preserve trajectory-compatible but biologically incorrect equations.

David Krongauz, Arad Zulti, Eran Segal +1

Department of Computer Science and Applied Mathematics, Weizmann Institute of Science, Rehovot, Israel · Department of Molecular Cell Biology, Weizmann Institute of Science, Rehovot, Israel · Mohamed bin Zayed University of Artificial Intelligence, Abu Dhabi, United Arab Emirates +2

Model Discovery Ordinary Differential Equations Symbolic Regression Drug Discovery Large Language Model Agents

Ask PDF