Learning Closure of Dynamical Systems with Kernel Ridge Regression

Evan Habbershaw; John Harlim; Senwei Liang

math.NAOct 1, 2026

Learning Closure of Dynamical Systems with Kernel Ridge Regression

Authors: Evan Habbershaw, John Harlim, Senwei Liang

Organizations: Department of Mathematics, The Pennsylvania State University · Department of Mathematics and Institute for Computational and Data Sciences, The Pennsylvania State University · Department of Mathematics and Statistics, Texas Tech University

Abstract

We develop a closure modeling framework for identifying missing components of dynamical systems using Kernel Ridge Regression (KRR). The framework addresses two classes of closure problems: difference-equation closures arising in ODE and PDE settings, and algebraic closures arising from moment closure in kinetic equations. For the first class, we derive an error bound in an ODE setting that quantifies contributions from time integration, approximation of unresolved scales, and interpolation required to couple unresolved-scale effects to the resolved solver. Numerical experiments on the Lorenz-63 system and the Kuramoto-Sivashinsky equation demonstrate accurate long-horizon predictions and substantial improvements over an LSTM-based closure model. For the second class, we consider moment closure for a one-dimensional kinetic equation by modeling discrepancies between kinetic and macroscopic fluxes as a function of the resolved macroscopic variables. We compare global KRR models based on PCA coordinates with spatially local models. While the global model performs well for unimodal initial conditions, its accuracy deteriorates for bimodal initial conditions. Spatially local models with appropriate modeling inputs improve robustness and achieve higher predictive accuracy.

Figures & tables

Stage treatment	DM-KRR	RBF-KRR	LSTM
Heun / endpoint RK2	$4.1969\pm 0.0799$	$4.1783\pm 0.1165$	$4.4603\pm 1.3410$
Quadratic RK4	$5.9874\pm 0.0868$	$5.9242\pm 0.0980$	$4.4924\pm 1.2797$
Cubic RK4	$8.3815\pm 0.3321$	$7.9683\pm 0.4844$	$4.6008\pm 1.3046$
Quartic RK4	$8.4832\pm 0.4910$	$8.3477\pm 0.6612$	$4.6115\pm 1.3114$
Quintic RK4	$\mathbf{8.6505\pm 0.5239}$	$8.4533\pm 0.7001$	$4.6439\pm 1.2873$

Table 1: Mean test valid prediction time for the Lorenz–63 hybrid closure models with

N_{\textnormal{train}}=4096

. The KRR hyperparameters are selected using three rollout validations per candidate pair. Results are reported as the mean and standard deviation across independently trained-model test means.

$N_{\textnormal{train}}$	DM-KRR	RBF-KRR	LSTM
512	$8.1719\pm 0.4849$	$7.8863\pm 0.5426$	$4.2784\pm 1.3443$
1024	$7.4267\pm 1.3808$	$7.4753\pm 0.9587$	$4.2253\pm 1.3061$
2048	$8.2535\pm 0.4977$	$8.0014\pm 0.4216$	$4.2652\pm 1.2222$
4096	$\mathbf{8.3815\pm 0.3321}$	$7.9683\pm 0.4844$	$4.6008\pm 1.3046$

Table 2: Mean test valid prediction time for the Lorenz–63 hybrid models using cubic reconstruction of the closure at RK4 stage times. Results are reported as the mean and standard deviation across independently trained-model test means.

Panel of Figure 1 — Figure 1: Mean validation VPT (left) and a representative validation (right) over the diffusion-kernel hyperparameter grid for $N_{\textnormal{train}}=16384$ and $m=0$ , over ten training windows. The figure shows the mean landscape over ten independent training windows. The outlined square marker identifies the selected grid point.

Panel of Figure 2 — Figure 2: Left: Representative test rollout for the $N_{\textnormal{train}}=16384$ , $m=0$ , diffusion-kernel model. The signed physical-space difference between the full-order resolved field and the hybrid closure rollout is recorded. The dashed line marks the verified VPT, $2.608$ . Right: The time-averaged energy in the six retained Fourier modes is shown.

Figure 3: Mean test VPT as a function of training-set size for the $m=0$ KS closure models. KRR error bars denote variability across independently trained models. The LSTM values are included as a qualitative comparison because the LSTM training-data convention and noise augmentation are not identical to those used for KRR. The dashed reference line marks the reported long-memory LSTM mean, $2.494688$ , whose reported standard deviation is $0.95904$ .

Data type	Phase shift	Varying parameters
I	$\mu=\eta\Delta x$	$A,\ \sigma,\ \eta$
II	$\mu_{1,2}=\phi\pm\frac{d}{2}$	$A_{1},A_{2},\sigma_{1},\sigma_{2},$ $d$ , $\phi$

Table 3: The two physical initial condition data types used in the kinetic experiments.

Model	Training data	KRR centers	$\mathcal{E}_{1}$	$\mathcal{E}_{2}$
Global	I	$8208$	$0.001285\pm 0.000200$	$0.028925\pm 0.005158$
Global	I+II	$8208$	$0.001005\pm 0.000234$	$0.008638\pm 0.007220$
zeroth-order local	I	$4032$	$0.001455\pm 0.000252$	$0.006088\pm 0.001457$
zeroth-order local	I+II	$4032$	$0.001891\pm 0.000341$	$0.003178\pm 0.000675$
first-order local	$(\textnormal{I},\,\textnormal{I})$	$(4032,8064)$	$0.000665\pm 0.000096$	$0.008098\pm 0.002001$
first-order local	$(\textnormal{I},\,\textnormal{I+II})$	$(4032,8064)$	$0.000793\pm 0.000105$	$0.003768\pm 0.001161$

Table 4: Comparison of the kinetic closure models under different closure models and training-data compositions. The columns

\mathcal{E}_{1}

and

\mathcal{E}_{2}

report the average prediction errors for the unimodal and bimodal testing data, respectively (lower values are better). For the first-order local model, the ordered pair in the training-data column gives the trajectory types used to train

(\widetilde{\mathcal{R}}_{0},\widetilde{\mathcal{R}}_{1})

, and the number of training centers is reported as

(N_{0},N_{1})

. The zeroth-order model

\widetilde{\mathcal{R}}_{0}

is held fixed between the two first-order comparisons.

Figure 4: Final-time density profiles for a Type II validation trajectory. Each panel compares the kinetic reference with one macroscopic model: the zero-correction baseline, the global PCA–KRR closure, the zeroth-order local closure, and the first-order local residual closure.

Figure 5: Space–time evolution of temperature for the same Type II validation trajectory as in Figure 4 . Panel (a) shows the kinetic reference temperature

\theta

, while panels (b)–(d) show the absolute temperature errors, for the global model, the zeroth-order local, and the first-order local model, respectively. The three error panels use a common linear color scale so that their magnitudes are directly comparable. The global model can reproduce the final macroscopic profile reasonably well while accumulating substantially larger structured transient errors during recursive rollout; both local models substantially reduce the structured rollout error and track the kinetic dynamics more faithfully throughout the trajectory.

Appendix figures & tables6 assets

Supplementary material from the paper’s appendix.

Appendix

$N_{\mathrm{train}}$	$m=0$		$m=1$
	RBF	Diffusion	RBF	Diffusion
2048	$1.0171\pm 0.1177$	$0.9121\pm 0.1183$	$0.9885\pm 0.1562$	$0.9161\pm 0.1567$
4096	$1.1241\pm 0.1884$	$1.1101\pm 0.1780$	$1.2156\pm 0.1276$	$1.1494\pm 0.1948$
8192	$1.6903\pm 0.3016$	$1.8164\pm 0.3150$	$1.6484\pm 0.3067$	$1.8659\pm 0.2142$
16384	$2.3650\pm 0.4693$	$\mathbf{2.6516\pm 0.3819}$	$2.1400\pm 0.4798$	$2.6490\pm 0.4311$

Appendix

Table 5: Mean test VPT for the KRR closure models with memory depths

m=0

and

m=1

. For each configuration, the reported statistics are the mean and standard deviation across ten independently trained-model test means, with

50

test rollouts per trained model.

$N_{\textnormal{train}}+2N_{V}$	Memory $m$	VPTs (ETDRK4)
$2.5\times 10^{5}$	19	$2.494688\pm 0.95904$
2048+4000	0	$1.41372390\pm 0.6669579$
4096+4000	0	$1.4588309\pm 0.71275799$
8192+4000	0	$1.5381616\pm 0.725299612$
16384+4000	0	$1.7413667\pm 0.8070518383$
2048+4000	1	$1.3993662\pm 0.6346926$

Appendix

Table 6: Test VPTs of true 6 modes and predicted 6 modes for the LSTM model using ETDRK4.

Representation	Type I	Type I+II (Mixed)
Global	$72\times 114=8208$	$72\times 57+72\times 57=8208$
Zeroth-order local	$72\times 56=4032$	$72\times 28+72\times 28=4032$
First-order local	$72\times 112=8064$	$72\times 56+72\times 56=8064$

Appendix

Table 7: KRR training-pair allocations used in the kinetic-model comparisons. In each table entry, the number

72

represents the number of Type I and Type II trajectories used for training. For the global model, the second factor denotes the number of temporal full-field snapshots selected per trajectory; for the local models, it denotes the number of local space–time interface samples selected per trajectory.

Figure 6: Final-time velocity profiles for the same bimodal rollout used in the main-text kinetic comparison in Figure 4 . Each panel compares the kinetic reference with one macroscopic model: the zero-correction baseline, the global closure, the zeroth-order local closure, and the first-order local residual closure.

Figure 7: Final-time temperature profiles for the same bimodal rollout used in the main-text kinetic comparison in Figure 4 . Each panel compares the kinetic reference with one macroscopic model, using the same model ordering as Figure 6 .

Figure 8: Final-time energy profiles for the same bimodal rollout used in the main-text kinetic comparison in Figure 4 . Each panel compares the kinetic reference with one macroscopic model, using the same model ordering as Figure 6 .

PDE Surrogate Modeling

PDE surrogate modeling builds fast approximations of partial differential equation solutions for simulation, prediction, and design.

See all 321 papers

Jun 30, 2026

Learning dynamical systems from noisy data with Weak-form Kernel Ridge Regression

Accurate prediction of complex dynamical systems from noisy measurements remains a significant challenge in scientific computing. Kernel ridge regression learning strategies are often effective when applied to clean data, but have limited success with noisy data. Recent work has observed that a weak formulation can act to filter noisy data, and different learning strategies have achieved increased noise robustness with a weak-form framework. In this manuscript, we give an overview of the filtering mechanism behind the weak formulation and provide a bias-variance error decomposition. Using these insights, we combine a weak formulation with a kernel learning strategy to propose Weak-form Kernel Ridge Regression (WKRR) for learning dynamical systems. The proposed framework is simple to implement, effective for both clean and noisy data, and outperforms several baseline methods. We demonstrate the performance of WKRR on chaotic benchmark systems in up to 64 dimensions, as well as 15,000-dimensional real-world fluid data.

Max Kreider, John Harlim, Daning Huang

Department of Mathematics The Pennsylvania State University, University Park, PA 16802, USA · Department of Mathematics, Institute for Computational and Data Sciences The Pennsylvania State University, University Park, PA 16802, USA · Department of Aerospace Engineering The Pennsylvania State University, University Park, PA 16802, USA

Kernel Ridge Regression Kernel Regression Nonlinear System Identification Chaotic Dynamical Systems Scientific ML

Ask PDF

Jun 3, 2026

Mamba-Assisted Non-Markovian Closure for Reduced-Order Modeling

Reduced-order modeling of high-dimensional dynamical systems is often hindered by closure effects arising from unresolved variables, which can introduce non-Markovian dependence into the resolved dynamics. Motivated by the history-dependent memory term arising in the Mori--Zwanzig formalism, we recast non-Markovian closure modeling as a sequence modeling problem and propose the Mamba-Assisted Closure (MAC) framework. MAC employs a Mamba-based sequence model to predict the closure from the resolved trajectory and couples the learned closure with the reduced-order governing equations through a numerical integrator to advance the resolved variables in time. During training, the selective scan mechanism in Mamba enables efficient parallel sequence processing with linear scaling in sequence length, while autoregressive inference proceeds through recurrent state updates at essentially constant per-step cost. We evaluate MAC on four benchmark systems with complementary characteristics: the viscous Burgers' equation, the chaotic two-scale Lorenz '96 system, the 3-bus DeMarco--Zheng power-grid system, and the dispersive Korteweg--de Vries equation. Across these benchmarks, MAC consistently improves predictive accuracy and long-time rollout stability relative to the comparison models, demonstrating an effective and computationally scalable approach to non-Markovian closure modeling.

Zhi-Feng Wei, Saad Qadeer, Panos Stinis

Advanced Computing, Mathematics, and Data Division, Pacific Northwest National Laboratory, Richland, WA 99354, USA · Department of Applied Mathematics, University of Washington, Seattle, WA 98195, USA · Division of Applied Mathematics, Brown University, Providence, RI 02912, USA

Reduced-Order Modeling

Ask PDF

May 29, 2026

The Dynamic-Probabilistic Consistency Gap in Chaotic Surrogate Modeling

Dynamical systems reconstruction (DSR) aims to learn surrogate models that capture the dynamics underlying time-series data. Reliably deploying these surrogates requires uncertainty estimates consistent with the learned dynamics. We expose a dynamic-probabilistic consistency (DPC) gap: the pursuit of finite-horizon probabilistic objectives can degrade dynamics or decouple predictive uncertainty from the local tangent dynamics it ought to reflect. We isolate three mechanisms behind this gap: core collapse, noise masking, and blind uncertainty. Specifically, we show that open-loop Gaussian rollout objectives can penalize Jacobian-generated covariance growth in chaotic systems, encouraging optimization shortcuts that weaken physical expansion or decouple uncertainty from it. To mitigate this gap, we propose KAFFEE (Kalman-Aware Framework For Ergodic Emulation), a differentiable extended Kalman filter-based training framework that evaluates likelihood on local predictive residuals (innovations) while transporting covariance through learned local Jacobians. On stochastic hyperchaotic Lorenz-96, KAFFEE reduces the identified failure modes, improves reconstruction of dynamical invariants relative to open-loop objectives, and maintains competitive predictive scores. We further show that the DPC gap appears when probabilistically adapting a DSR foundation model across 13 chaotic systems, where KAFFEE enables in-context Bayesian filtering while largely preserving zero-shot dynamics.

Andre Herz, Matthijs Pals, Daniel Durstewitz +1

1Interdisciplinary Center for Scientific Computing, Heidelberg University, Germany · Faculty of Mathematics and Computer Science, Heidelberg University, Germany · 3Dept. of Theoretical Neuroscience, Central Institute of Mental Health (CIMH), Mannheim, Germany +3

Dynamical Systems Surrogate Modeling Uncertainty Quantification Chaotic Dynamical Systems Extended Kalman Filtering

Ask PDF

Abstract

PDE Surrogate Modeling

Explore similar work

Learning dynamical systems from noisy data with Weak-form Kernel Ridge Regression

Mamba-Assisted Non-Markovian Closure for Reduced-Order Modeling

The Dynamic-Probabilistic Consistency Gap in Chaotic Surrogate Modeling