Koopman Kernel Regression
Authors: Petar Bevanda
Max Beier
Armin Lederer
Stefan Georg Sosnowski
Eyke Hüllermeier
Sandra Hirche
Venue: NeurIPS 2023 poster (2023)
Keywords: Kernel Methods
Regression
Statistical Learning Theory
Koopman Operator
Mode Decomposition
Dynamical Systems
Supervised Learning
TL;DR: We propose a novel Koopman-invariant and universal kernel on trajectories with generalization and convergence guarantees for LTI forecasting.

Abstract: Many machine learning approaches for decision making, such as reinforcement learning, rely on simulators or predictive models to forecast the time-evolution of quantities of interest, e.g., the state of an agent or the reward of a policy. Forecasts of such complex phenomena are commonly described by highly nonlinear dynamical systems, making their use in optimization-based decision-making challenging.
Koopman operator theory offers a beneficial paradigm for addressing this problem by characterizing forecasts via linear time-invariant (LTI) ODEs, turning multi-step forecasts into sparse matrix multiplication.
Though there exists a variety of learning approaches, they usually lack crucial learning-theoretic guarantees, making the behavior of the obtained models with increasing data and dimensionality unclear.
We address the aforementioned by deriving a universal Koopman-invariant reproducing kernel Hilbert space (RKHS) that solely spans transformations into LTI dynamical systems. The resulting Koopman Kernel Regression (KKR) framework enables the use of statistical learning tools from function approximation for novel convergence results and generalization error bounds under weaker assumptions than existing work. Our experiments demonstrate superior forecasting performance compared to Koopman operator and sequential data predictors in RKHS.

Koopman Kernel Regression
Petar Bevanda
Max Beier TU Munich max. wwwerereeeee de
Armin Lederer
TU Munich
TU Munich armin.lederee lemein.eddrrr@tummcco de
betar.bevanda@tum.co
Stefan Sosnowski
Eyke Hüllermeier
Sandra Hirche
TU Munich
LMU Munich eyke fifii.ii. fi. 1mu
TU Munich hi mirche@tunn de
sosnowski@tto de
Abstract
Many machine learning approaches for decision making, such reinforcement learning, rely on simulators or predictive models to forecast the time-evolution quantities of interest, e.g., the state of agent the reward of a policy. Forecasts of such complex phenomena are commonly described by highly nonlinear dynamical systems, making their in optimization-basel decision-making challenging. Koopman operator theory offers beneficial paradigm for addressing this problem by characterizing forecasts via linear time-invariant (LTI) ODEs, turning multi- step forecasts into sparse matrix multiplication. Though there exists a variety learning approaches, they usually lack crucial learning-theoretic guarantees, making the behavior of the obtained models with increasing data and dimensionality unclear. We address the aforementioned by deriving a universal Koopman invariant reproducing kernel Hilbert space (RKHS) that solely spans transformations into LTI dynamical systems. The resulting Koopman Kernel Regression (KKR) framework enables the use of statistical learning tools from function approximation for novel convergence results and generalization error bounds under weaker assumptions than existing work. Our experiments demonstrate superior forecasting performance compared to Koopman operator and sequential data predictors in RKHS.
-oopman-invariiat
Introduction
Dynamical systems theory is fundamental paradigm for understanding and modeling the time evolution of phenomenon governed by certain underlying laws. Such a perspective has been successful in describing countless real-world phenomena, ranging from engineering mechanics [1] and human movement modeling [2] to molecular and quantum systems [3, 4J. However, as the laws governing dynamical systems phenomena may have to rely on data rather than first principles. this regard, machine learning methods, which have shown immense potential in tackling complex tasks in domains such as language models [5] and computer vision [6], are coming to the fore. Though powerful, ofate-of-the-art neural vector fields [7] flows [8] commonly compose highly nonlinear maps for forecast, i.e. computing
often unknown, modeling and understanding the underlying
(f) x(0) J. f(x(t))
(1)
for, e.g. a scalar ODE = fx)) Hence, is often challenging to use such models in optimization- based decision making that relies simulators or predictive models, e.g., reinforcement learning [9-11]. particularly beneficial perspective for dealing with the aforementioned comes from Koopman operator theory [12-15]. Through point-spectral decomposition of Koopman operators, forecasts become superpositions of solution curves of a set of linear ODEs {& = Ay )
(0),
{I
lift
2
(2)
37th Conference on Neural Information Processing Systems (NcurIPS 2023).

---

where a vector valued function span( (g)) "lifts x onto a manifold Z := span (f2)) )
Throughout, we refer to these models linear invariint (LTI) predictors. The learning objective of such representations is twofold: spanning system trajectories by the learned manifold & and constraining the L.TI dynamics to The latter is long-standing challenge of Koopmanism [16-20],
as manifold dynamics of existing approaches "leak-out" [21] and limit predictive performance.
To tackle the aforesaid, we con- nect the representation theories of reproducing kernel Hilbert spaces (RKHS) and Koopman operators. As a first in the lit- erature, derive universal kernel whose RKHS exclusively spans manifolds invariant un- der the dynamics, as depicted in Figure 1. A key corollary of Figure 1: Illustration of on- m--maafifod dynamics of LTI predictors. unconstrained manifold dynamics is the lack of essential learning- guarantees, making the behavior of existing learned models unclear for increasing data and dimensionality. To address this, we utilize equivalences to function regression in RKHS to formalize a statistical learning framework for learning LTI predictors from sample trajectories of dynamical system. This, in turn, enables the use of statistical learning tools from function approximation for novel convergence results and generalization error bounds under weaker assumptions than before [22-24]. Thus, we believe that our Koopman Kernel Regression (KKR) framework takes the best of both RKHS and Koopmanism by leveraging modular kernel learning tools to build provably effective LTT predictors.
Our
ODES<)o()) x00)
eAit
200
Exisling
The remainder of this paper is structured as follows: We briefly introduce LTI predictors and discuss related work in Section The derivation of the KKR framework, including the novel Koopman RKHS, is presented in Section 3. In Section 4, we show the novel learning guarantees in terms of convergence and generalization error bounds. They are validated in comparison to the state-of-the-art through numerical experiments in Section 5.
Notation Lower/upper case bold symbols x x/X denote spatial vector/matrix valued quantities. A trajectory defines a curve By C traced out by the flow over time T= [0, IJ from any (r, x) E TxX. In discretizing T, collection of points 2H X from discrete time steps TH (to....... .n) is considered. IL' Y). with discrete-time analogues X (H, X) / YH C l (H, Y) with domain and co-domain separated by ".". The vector space of continuous functions on Xr endowed with the topology of uniform convergence on compact domain subsets is denoted C(Xr). The collection of bounded linear operators from Yr to YT is denoted as B(Tr). The adjoint of A E B(-) is A". Discrete-time eigenvalues read μ:= AAt X E C. random variable X defined on a probability space (12, A, p) has expectation E[X] = fg (())p(w).
The state/output trajectory spaces are denoted as [L X)
2
Problem Statement and Related Work
To begin, we formalize our problem statement and put our work into into context with existing work.
2.1 Problem Statement
Consider a forward-complete system comprising a nonlinear state-space model
a fax), xo x=()),
(3a)
= g(x),
(3b)
on a compact domain X C. Idd with a quantity of interest q C(X). The above system class includes
all systems with Lipschitz flow F' F'(xx)) 1/' (=(())) e.g., mechanical systems [27].
Inspired by the spectral decomposition of Koopman operators, we look to replace the nonlinear state-space model (3) by an LTI predictor
2 = Az, 20 = g(æ)),
(4a)
y = z,
(4b)
Although we outline the scalar output case for ease of exposition. expanding to a vector-valuee case is possible w.l.o.g. If required, forward completeness can be relaxed unboundedness observability [25, 26).

---

with D E N and D-dimensional function approximator dense in C(X). Then, from initial
conditions X, S X that form non--eeurrent domain Xr. (4) admits universal approximation of
the flow of (3) such that Ve > 0, 3D so that = lur(() - e e g(x) < e, Vi E [0, 7 [28]². In this work, we aim to find a solution to the following constrained, functional optimization problem (OR) Output reconstruction:
min ll - c. grllvr c.g.A
(5a)
(KT) Koopman-invariance
such that g(x(1)) = g(x(0)), VI E [0,
(5b)
Although the sought-out model (4) is simple, the above problem is non-trivial and much of the existing body of work utilizes different simplifications that often lead to undesirable properties. In the following, we elaborate these properties and motivate our novel sample-based solution to (5), which remains relatively simple but nonetheless ensures a well-defined solution with strong learning guarantees.
2.2 Related Work
Koopman operator regression in RKHS Equipped with a rich set of estimators, operator regression in RKHS seeks a sampled-data solution
min [g(()(()) A (z)(((
(6)
with A a Hilbert-Schmid operator [29] - commonly known as KRR, and EDMD (PCR) or RRR when under different fixed-rank constraints [23, 24]. The choice of RKHS g is commonly one that is dense in suitable L space, e.g. that of the RBF kernel. By an additional projection, a quantity of interest can be predicted via mode decomposition of the estimated operator, leading to a model akin to (4). In the light of (5), the feature map g is predetermined while violating (KT) is merely minimized for single time-instant l. As a consequence, such approaches offering limited predictive power over the time interval [o, 7' of trajectory as displayed in Figure 1. The extent to which (KI) is violated due to spectral properties [30] or estimator bias [24] is known as spectral pollution [31]. The strong implications of this phenomena, motivate regularization [32] and spectral bias measures [24] to reduce its effects. Due to the above challenges, guarantees for Koopman operator regression (KOR) have only recently gained increased attention. Often, however, existing theoretical results [29, 33] are generally not applicable to nonlinear dynamics [34] due to the commonly unavoidable misspecification [35] of the problem (6) incurred by neglecting (KT). The first more general statistical learning results [23, 24] are derived in stochastic setting under the assumption that the underlying operator is compact and self-adjoint. In stark contrast, the same set of assumptions is restrictive for the deterministic setting [35]: compactness only holds for affine deterministic dynamics [36, 37] while self-adjointness is known to generally not hold for Koopman operators [13, 38, 39]. Regardless of the setting, however, the state-of-the-art exhibits alarming properties: forecasting error not necessarily vanishing with LTI predictor (4) rank [23, Theorem 1] and risk based on a single time-instant.
oblivious
to the time-series structure
Learning via Koopman eigenspaces Geared towards LTI predictors and closer to our own problem setting (2), another distinct family of approaches aims to directly learn the operator's invariant subspaces [28, 40-42). The goal is to fit g()) based on approximate Koopman operator eigenfunctions that still fit the output of interest (OR). However, existing data-driven approaches in this line of work rely on ad-hoc choices and lack essential carning-theoreet properties such as feasibility and uniqueness of solutions - prohibiting provably accurate and automated LTI predictor learning.
Kernels for sequential data Motivated by the lack of priors that naturally incorporate streaming and sequential data, there is an increasing interest in signatures [43]. They draw from the rich theory of controlled differential equations (CDEs) [44, 45] and build models that depend on a vareing observation history. An RKHS suitable for sequence modeling is induced by signature transformation of a base/static RKHS. Generally. if the latter is universal, so are the signature kernels [46]. While arguably more general and well-versed for discriminative and generative tasks [47], forecasting using signature kernels [48] comes at price, as their nonlinear dependence on observation streams leads to significant complexity increase compared to LTI predictors.
Background on prerequisite Koopman operator theory can be found in the supplementary material.
3

---

Motivated by the restrictions of existing Koopman-basee predictors, we propose a function approxi- mation approach that exploits exploits time-series data and Koopman operator theory to provably learn LTTI predictors. Through novel invariance transform we can satisfy (KT) by construction and directly minimize the forecasting risk over an entire time-interval (OR). In simple terms: Koopman operator regression fixes g(.) and regresses A and in (5), whereas our KKR approach selects A to jointly regress с and g(.). Similar in spirit to generalized Laplace analysis [21, 49], our approach allows the construction of eigenmodes from data without inferring the operator itself. Crucially, we demonstrate that selecting A requires no prior knowledge as confirmed by our theoretical results and experiments. To facilitate learning LTI predictors, we derive universal RKHSs that are guaranteed to satisfy (KI) over trajectories - a first in the literature. The resulting equivalences to function regression in RKHS allow for more general and complete learning guarantees in terms of consistency and risk bounds that are free of restrictive operator-theoretic assumptions.
3 Koopman Kernel Regression
With the optimization (5) being prohibitively hard due to nonlinear and possibly high dimensional constraints, we eliminate the constraints (5b) by enforcing the feature map g()) to have the dynamics of intrinsic LTI coordinates associated with Koopman operators, i... their (open) eigenfunctions [21].
Definition 1. A Koopman eigenfunction (xc:EC(X) satisfies xx(x)) e ex F (x)), 1L E [0, T).
is proven that Koopman eigenfunctions from Definition 1 are universal approximators of continuous functions [28] making them a viable replacement for the feature map g()) in (4). However, follow- ing their definition, it is evident that Koopman eigenfunctions are by means arbitrary due their inherent dependence on the dynamics' flow. Using the well-established fact that Koopman operators compose function with the flow, i.e., K. (-) = g( (FF"))) it becomes evident the eigenfunctions from Definition I are (semi)group invariants, they remain unchanged after applying {e At K') Thus, inspired by the seminal work of Hurwitz on constructing invariants [50], we can equivalently reformulate (5) as unconstrained problem and jointly optimize over eigenfunctions
Lemma 1 (Invariance transform). Consider a function E ((x)) over a set of initial conditions Xo C X that form non-recurrent domain KT. The invariance transform IJ transforms g into an Koopman eigenfunction Φx E C(KT). (3a) with LTI dynamics described by À E C
Φx = g(0o) = f c g( (2o))
7=0
The above Lemma is key stepping stone towards deriving a representer theorem for LTI predictors. However, it is also interesting in its own right as it provides an explicit expression for the flow of an eigenfunction from any point in the state space. Thus, it provides a recipe to obtain a function space that fulfills (KI) by construction. As we show in the following. a sufficiently rich set of eigenvalues (28) and Lemma will allow for a reformulation (5) into an unconstrained problem
min llyr ((FT)))s
(8)
where the operator M())=-1 [ox (-) Φx T-)T is universal and consisting of Koopman-invarinn functions.
3.1 Functional Regression Problem
Notice that the problem reformulation (8) is still intractable, as a closed-form expression for the flow map is generally unavailable even for known ODEs. This requires integration schemes that can introduce inaccuracies over time interval [0. T). Thus, to make the above optimization problem tractable, data samples are used ubiquitous in learning dynamical systems.
Assumption 1. A collection of N pairs of trajectories Dy=(=²) Ur Jim1 E ( N is available.
(s) Whh
By aggregating different invariance transformations (7) into the mode decomposition operator
M()=E)='>,>> Xy + Wr,
(9)
Proof's for all theoretical results can be found in the supplementary material.

---

we can formulate a supervised learning approach in the following.
Learning Problem With Assumption and Lemma 1, the sample-based approximation of problem (8) reduces to solving
min
VT 1 M(x)"))
(10)
while preserving the mode decomposition structure (9). To realize the above learning problem, we resort to the theory of reproducing kernels [51, 52] and look for an operator M E H. where H
is an RKHS. A well-established approach using RKHS theory is to select M as a solution to the regularized least squares problem
M arg min MEH
Σ lor M(a + 1l
(11)
with y R and ||-11 a corresponding RKHS norm. As our target a function-values mapping (-) an operator - 1l-l is induced by operator-valuee kernel K Xx B( space of bounded operators over the output space [53]. The salient feature of the above formulation (11) is its well-posedness: its solution exists and is unique for any H. expressed as
mapping to the
M(-= - K(,, 24)3, E
(12)
through representer theorem [54]. Still, due the Koopman-invarian structure (9) from Lemma 1, the choice of the RKHS for M arbitrary. Thus, the question is how to craft H. the solution M is decomposable into Koopman operator eigenfunctions (9), forming an LTI predictor.
Firstly, it is obvious that (9) consists of summands that may lie in different RKHS, denoted as {HA* 1/21- Then, H. constructed from the following direct sum of Hilbert spaces [55]:
H = A @ ......... ... HDD so that =
( +b SIEH> ....., ,............................... JDHH (13)
with S: F + H. (/1 ... S.)...)) + f1 + ..... + fD the summation operator [56]. Thus, to construct H. a specification of the RKHS collection { is required, so that it represents Koopman eigenfunctions from (9).
GAN from Assumption 1, Theorem 1 (Koopman eigenfunction kernel,, Consider trajectory data {
a X E C and a universal (base) kernel k: XxX + R.. Then, the kernel K XyxXx + B(Yr)
[T [T - MT-- ko (xT (T) (T), ((' drdr',
T
(xT,, x'T)
(14)
7=0 T'=0
(i) defines an RKHS H²,
(ii) universal for every eigenfunction of Definition corresponding X, (iii) induces a data- detteddenndent function space span { K (K (. ).... K K^ x())) that is
Koopman-invarian over trajectory-data {xT 1,WN
Ji=1-
In Theorem 1, we derive an eigenfunction RKHS by defining its corresponding kernel that embeds the invariance transformation (7) over data samples. Also, would like to highlight that the above result addresses a standinn open challenge in the Koopman operator community [19-21], ], i.e., defining universal function spaces that guaranteed to be Koopman invaaiann Now, are ready to introduce the Koopman kernel as the kernel obtained by combining "eigen-RKHS" as described in (13).
Proposition 1 (Koopman kernel). Consider trajectory data DN of Assumption I and a set of kernels {K^s] I= from Theorem 1. Then, the kernel K Kx XXT + B(Er) given by
K( x'7) = Ep- K (xr,
(15)
(i) defines an RKHS H == S(H $ ... ① HHAD),
(i) is universal for any outpul (3b), provided a sufficient amount of eigenspaces D.
Above, have derived the "Koopman-RKHS" H for solving the problem (11) with uni- versal RKHS spanning Koopman eigenfunctions. Thus, the sample-data solution for an eigenfunction flow follows from the functional regression problem (11) and takes the form Φx, (.) = EN
3 E Yr - providing a basis for the LTI predictor.
SSffficiett amount is rich enough set of eigenvalues (0.T10 from Bi (0) in C [57, Theorem 3.0.2).
5

---

3.2 Practicable LTI Predictor Regression
As functional approximation problem, the solution of (11) is not parameterized by vector-valued coefficients, but rather functions of time. Although there are a few options to deal with function-valued solutions [53], we consider a vector-valued solution. A common drawback of such a discretization involves the loss of the inter-sample relations along the continuous signal. Crucially, this problem does not apply in our case, as the inter-sample relationships remain modeled for the discrete-time "Koopman kernel" due to its causal structure. Importantly, the vector-valued solution allows us to preserve all of the desirable properties derived in the continuous case.
Consider sampling [0, 7 at H H-T/At regular intervals to yield a discrete-time dataset from Assump- tion discretized at points H {fo- (fo..-t }.
As discretization of a function over time, with a slight solving the
abuse of notation, we denote the target vectors yи [y((t) -.. (tn)) T Thus, time- and data-discretized version of the problem (5) that takes the form of a linear coregionalization model [58, 59].
{ Corollary 1 (Time-discrete Koopman kernel), Consider trajectory data { i==1 and let pe:: e 11 [u]. -. Then, the scalar-induced matrix kernel KP:; XH x XH + B(HH)
K ("I 2H') = 13p33 (1+1)² Em=
E E" In=op 1 (tm) ah (tn)) p3""
(16)
satisfies the properties (i) (i)- Theorem I over H., so that it defines an RKHS HDs is universal
per Definition I over H. with span ("( (", (1) ( KM m (2))) (KI) over {= (8hN
Ji=1- Given
collection of kernels {KP) (K)))== the matrix Koopman kernel Kd;; Xи*XH + B(Yn)
K(EH, <(oH, = =n')=
Σド= (ZH,
(17)
satisfies the properties (i)-(ii from Proposition I over HI, defining RKHS = B.....
Now, we are fully equipped to obtain the time-discrete solution to our initial problem (5) provided dataset of trajectories. Before presenting the solution to Koopman Kernel Regression, we introduce
some helpful shorthand notation. We use the following kernel matrix abbreviations: kxx
[l)( ( Ma, k(x, X)= = [k x 1 Kxx- [K(a, (ax, x: 1
K [((x()) x [ and K(x, X): =
Proposition 2 (KKR). Consider discrete-time dataset of Assumption 1, D {a (N)T with @ the Kronecker product. Then, B=(Kxnxi ( YH (18)
H VH Ji=1 and
(t) OhN
let y== - [V/ ... ..H
a; = kxxxxxx Xaxa (IN @ pz; B,
defines unique time-sampled solution to (11) in terms of eigenfunctions (b(zr) = [k , ax ,1
determining an LTI predictor with A diag([u ... .Dl),
z = Az, 20 = ②0xx))
(19a)
9 = 1 1z.
(19b)
Notice how in (18), we re-scale the trajectory domain to that of the state-space. This enables us to write the forecast of (19), with slight abuse of notation, using extended observability matrix [60]
tu To)(z)),
T [ A 1 1TA
(20)
The confinement to a non-recurrent domain plays a crucial role in making the base kernel RKHSs isometric to "eigen-RKHSs" Hd s H via invariance transforms, guaranteeing feasible return from the time-series domain Xy to the state-space domain Xg C for evaluating the model over initial conditions.
The salient feature of our proposed KKR framework compared to existing methods is the fact that Koopman- invariance (KI) over data samples is independent from the outcome of an optimization algorithm, e.g. minimizing the forecusting risk to compute 3 (18). Thus, we are able to directly optimize for a downstream task (forecasting) (OR) given suitably rich set of eigenvalues.
Remark 1.
s For discrete-time predictors, we omit the sime-step specification and denote the next state with "(-)tn.

---

3.3 Selecting Eigenvalues
Until now, we have used the sufficient cardinality D E of an eigenvalue set that encloses [61] or the true spectrum. However, we have provided no insight regarding the selection of D spectral components or how they can be estimated. Here, we go beyond the learning-independent and non- constructive existence result of [28] and provide consistency guarantee and relate it to sampling cigenvalues without the knowledge of the true spectrum.
Proposition 3. Consider the oracle Koopman kernel K ( oh') and a dense set {1z} in (1).). Then, | (#H, "H') Σ}=1 K (IH, FH (PH')l((Z))) 0, V HH Vy' BH as
0o.
As shown in Proposition 3, even we do not know the oracle kernel, we can arbitrarily approximate it by sampling from a dense set supported on the closed complex unit disk IP (0) [57, Theorem 3.0.2] with the emror vanishing the limit Dc There is loss of generality when considering the unit disk as any finite radius disk can be scaled in the interval [0, T. Furthermore, approximation of the oracle kernel by sampling a distribution over (0) leads to an almost sure (1/vD) convergence rate. It is conceivable that faster rates can be obtained in practice by including prior knowledge to subspace
shape the spectral distribution, e.g. using well-known concepts such leverage-scores orthogonality [62, 63]. Based on spectral priors one can include a more biased sampling technique by precomputing components of the operator spectrum, e.g. computing Fourier averages [64], to determine the phases Wj of complex-conjugeen
pairs μ = [4 | eااي and sample the modulus from physics-informed informed distribution. However, rigorous considerations of optimized and efficient sampling are beyond the scope of this paper and rather a topic of future work.
another physics
3.4 Numerical Algorithm and Time-Complexity
Algorithm 1 Regression and LTI Forecasts using KKR
For a better overview, the pseudocode for regression and forecasting of our method are shown in Algorithm We also put the time-compp teee-explletty of our algorithm into per- spective w.r.t Koopman operator regres- sion of PCR/RRR [23] and ridge regres- sion using state-of-the-art signature ker- nels [48] (RR-Sig-PDE) in Table 1. The training complexity of our KKR is com- parable to that of RR- PDE regression and generally better than that of PCR/RRR. Given that accurate LTI forecasts require higher-rank predictors, the seemingly mild quadratic dependence makes D2 > NH and leads to a more costly matrix inversion. Furthermore, our LTI predictor also has slightly better forecast complexity due to not depending trajectory length. Obviously, due mere matrix multiplication after initial nonlinea map, LTI predictors have a significantly lower evaluation complexity than the nonlinear predictor of Sig-PDE's. Due to requiring updated observation sequences as inputs, Sig-PDE kernels introduce a raw evaluation complexity that is also quadratic in sequence length.
(8) L Data (- Vit 1i=1. Eigenvalues {
function REGRESS(ID, {1g})=))
form Gramians xx (A Xn 1D K XnnX
fit mode operator (-): XH + Yн (18, right)
recover eigenfunctions ((.): ; Xo + Zn (18, left)
construct r: Zo ~ W (20, right)
return LTI predictor T (-) Xo + Yи
end function
function FORECAST(2))
"lift" 2g = Φ (x))
rollout Uи = T T2
return trajectory ÚH
end function
Method
Training O(N DD H' O'D'N"H*+N'IP' O(N O
H-step forecast
N
trajectories
Table 1:
KKR (ours)
O(DH DNd)
H
trajectory length
Time com- plexities.
PCR/RRR
HD
D
predictor rank
I
time-delays
RR-Sig-PDE
O(NHPP)
d
dim(input data)
4 Learning Guarantees
With a completely defined KKR estimator, we assess its essential the behavior of the learned functions w.r.t. to the ground truth with increasing dataset size.
properties, i.e.,
7

---

4.1 Consistency
Although estabbitall in most function approximation settings [65-67], the setting of Koopman- based LTI predictor learning for nonlinear systems is void of consistency guarantees. Here we use definition of universal consistency from [68] that describes the uniform convergence of the learned function to the target function as the sample size goes to infinity for any compact input space X and every target function q (X). The existing convergence results for Koopman-basee L.TI predictors [69] are in the sense of strong operator topology - allowing the existence of empirical eigenvalues that are not guaranteed to be close to true ones even with increasing data [70]. This lack of spectral convergence has cascaded effect in Koopman operator regression as, in turn, the convergence of eigenfunctions and mode coefficients is not guaranteed. Here, the convergence of modes is replaced by the convergence of eigenfunctions, and convergence of spectra is replaced by the convergence of (20) to the mode decomposition operator M Tá + M =T with the estimate denoted by (.). Theorem 2 (Universal consistency). Consider universal kernel K (17) and data distribution supported Xμ Y Then, as N + co, ||M-M - 0 and 11p - - lle - 0,
4.2
Generalization Gap: Uniform Bounds
Due to formulating the LTTI predictor learning problem as a function regression problem in an RKHS, we can utilize well-established concepts from statistical learning to provide bounds on the generalization capabilities of KKR. Given dataset of trajectories, the following empirical risk is minimized
Rx(M) = EiE(N) Ilvi
M (aff))
which is "in-sample" mean square error (MSE) w.r.t. a trajectory-data generating distribution pp of i.i.d initial conditions. The true risk/generalization error of estimator is the "out-of-sample MSE of the model on the entire domain and denoted as R(.). Those quantities are, in essence, the model's performance on test and training data, respectively. Allowing for statements on the test performance with an increasing amount of data by means of training performance is a desirable feature in data-driven learning. Hence, we analyze our model in terms of the generalization gap
(R( Rx (M) = [B wa) [llyn - M (#n) [2 1/ EN 11vl (0 M (m
(21)
To ensure well-specified problem, we require models in the hypothesis to admit a bounded norm. Assumption 2 (Bounded RKHS Norm). The unknown function M has bounded norm in the RKHS H attached to the Koopman kernel K (.,,). )), ie., || S. B for some B E R..
The above smoothness assumption is mild, e.g., satisfied by band-limited continuous trajectories [71] and computable from data [72, 73). In stark contrast, well-specified Koopman operator regression [23] requires the operator map the RKHS onto itself, which is very strong assumption [34, 35).
To derive the main result of this section, we utilize the framework of Rademacher random variables for measuring complexity of our model's hypothesis space, concept generally explored in [74] and more particularly for classes of operator- valued kernels in [75]. Conveniently, the derivation is, in terms of the RKHS HA, similar to standard methods IKHS-based complexity bounds [74]. We use well- well-known results based on concentration inequalities to provide high probability bounds on a model's 'ss generalization gap in terms of those complexities. Finally, we upper bound any constant with quantities specified in our assumptions and can state the following result.
( Theorem 3 (Generalization Gap of KKR), Let DA = [ar Ji=l be dataset as in Assumption I consistent with a Lipschitz system on a non-recurrent domain. Then the generalization gap (21) of model M from Proposition 2 under Assumption 2 is, with probability 1 1- ô, upper bounded by
KHH
8log 2
(
H VT
(R(M) - RN(( < 4RB
E COL
(22)
N
where an upper bound on the loss in the domain, and k. the supremum of the base kernel.
We observe an overall dependence of order O(//N) w.r.t. data points, resembling the regular Monte Carlo rate to be expected when working with Rademacher complexities. Remarkably, an increase in the order of the predictor D cannot widen the generalization gap but will eventually decrease the
8

---

Sample Trejectories
Time Horizon Dependence
Ours
Sig,PDE
PCR
Irue
((f))
Generalizztioon Gap
lest Za
10
training ænp
10
1
0.2 0.1
0.6
0.8 12
10 12
14
Figure 2: Forecasting performance (48 i.i.d. runs) for the bi-stable system for H=14 and N=50 for respectively optimal D0KR=100, DrcR=10 and 15 delays for Sig-PDEs. Left: Exemplary trajectories showing the advantage of learning with time-series kernels. Right: The generalization gap with an increasing forecast horizon, demonstrating generalization advantages of KKR.
empirical risk due the consistency of eigenspaces (Proposition 3). Combined, our findings are a substantial improvement, both quantitatively and in terms of interpretability, over existing risk bounds on forecasting error [23, Theorem 1J. Additionally, our intuitive non-recurrence requirement is easily verifiable from data. In contrast, the Koopman operator regression in RKHS comes with various strong assumptions [35] that require commonly unavailable expert knowledge. Also, the generalization of existing goopman-basee statistical learning approaches depends on rank while ours rank-independent. The significant implications of our results are demonstrated in the following.
5 Numerical Experiments
In our experiments", report the squared error of the forecast vector for the length of data trajectories averaged over multiple repetitions with corresponding min-max intervals. We validate our theoretical guarantees and compare to state-of-the-art operator and time-series approaches in RKHS. For fairness, the same kernel and hyperparameters are chosen for our KKR, PCR (EDMD), RRR [23] and regression with signature kernels (Sig-PDE) [48]. Note, PCR and RRR are provided with the same trajectory data split into one-step data pairs while the time and observation dimedsllas are fed as data to the Sig-PDE regressor due to its recurrent structure, Along with code for reproduction of our experiments, we provide a JAX [76] reliant Python module implementing a sklearn [77] compliant KKR estimator at https:/gitthub.com wttp:/hhhw..cooTM-uttu TOLIIM/UM--EEE IUo-oo//ooo
Bi-stable system Consider an ODE & = ar + bz that arises in modeling of nonlinear friction. The parameters are 4, b -4, making for bi-stable system fixed points +1. The numerical results are depicted in Figure 2. Sample trajectories both on training and testing data indicate the utility of the forecast risk minimization of KKR. While EDMD correctly captures the initial trend of most trajectories it fails to match the accuracy of Sig-PDE or our KKR predictors that utilize time-series structure. Furthermore, the behavior of KKR's generalization gap for an increasing time horizon T = H At, At = 1/14s closely matches our theoretical analysis.
Van der Pol oscillator Consider ODE 5 #(2 10x²) 0.8z describing dissipative system whose nonlinear damping induces a stable limit cycle - a phenomenon present in various dynamics. In Figure 3 two fundamental effects are validated: the generalization gap with increasing data and consistency with test risk that does not deterio- rate for increasing eigenspace cardi-
Table 2: Average risk (20 runs) [x10-2] for Van der Pol for various spectral sampling and lengthscales, 2-200, H=14.
p(u)
uniform
boundary-
physics inforseed
16
200
16
200
16
200
R4
13.7
5.38
11.2
5.38
5.60
5.58
nality. The performance of PCR/RRR strongly tied to predictor rank while Sig-PDE's less so w.r.t. delay length. R1-10-1
R1-100
6.46
0.78
4.10
0.78
0.97
0.92
7.12
1.74
4.33
1.74
1.83
1.80
Crucially, our KKR approach does not require careful choice of the eigenspace cardinality to perform for specific amount of data. Although the eigenvalues that determine the eigenspaces are randomly chosen from a uniform distribution in the unit ball, KKR consistently outperforms PCR/RRR. In Table 2 we show the spectral sampling and hyperparameter effects. We employ the following strategies: uniform- uniform distribution on the complex unit disk, boundary-biased a eigenvalues
distribution on the complex unit disk skewed towards the unit circle, of various vector field Jacobians. As expected, physics-informer performs well with lower rank
Additional details the numerical experiments can be found in the supplementary material.

---

Training Sample Cardinality Dependence
Delay Length
Generalization Gap
10"
10 100
Test Rlak
10-1
10
Curs Siy-PDE
RAR
10
10-2
10
10a
10 Eigenspace Cardinality
10a
Figure 3: Forecasting risks (20 i..d. runs) for the Van der Pol system over a time-horizon
14
(T = 1s). Left: Generalization gap for the best D / (ours 500, PCR 62, RRR 100, RR-Sig-PDE 10) is depicted with a growing number of data points. Right: Test risk behavior with an increasing amount of eigenspaces is shown for N 200. Shaded areas depict min-max risk intervals.
Σ Absolute Error
60
train
Ours test
10
40
PCR train
PCR test
10
Risk
20
10
10-2
0
20
40
60
10
10a
10a
10
timesteps
RBF kernel lengthscale
Figure 4: Cumulative eтror and forecast risks train-test splits) for flow past cylinder data and greater usable l-range and accuracy. Left: Cumulative absolute error for the best D/E (ours 200/70, PCR 200/35) is depicted over timesteps. Right: Forecast risk for 99 steps within range of RBF lengthscales. Shaded areas depict min-max intervals.
H 99. Our KKR with onders-of-magninee
compared to uninformed approaches. However, it is outperformed by unit-ball sampling approaches for higher rank due to a lack of coverage. Table 3 includes CPU timings for completeness.
#data = N x = = 200 x 14
KKR
PCR
RRR
Sig-PDE
Table 3: Computation
Training | Forecast ms)
8.0/ 54
90/ 84 88/150
8.6/5900
times for Van der Pol.
Flow past a cylinder We consider high- high-dimensional data velocity magnitudes a illustrated in Figure 5.
Kármán vortex street under varying initial cylinder placement,
The cylinder position is varied on a 7x7 grid in a 50x100- dimensional space and the flow is recorded over = The quantity of interest a velocity magnitude sensor placed in the wake of the cylinder. In forecasting from an initial ve- locity field, KKR outperforms PCR by orders-of-mant mfaggnuude shown in Figure 4. We omit Sig-PDE regression due persistent divergence after 2 20 steps. The latter is hardly surprising, given that Sig-PDE models iterate one step pre- dictions based only on the shapes of time-delays while LTI models directly output time-series from initial conditions.
Figure 5: Flow illustration. Area of initial cylinder positions shaded,
Conclusion
We presented a novel statistical learning framework for learning LTI predictors using trajectories of dynamical system. The method is rooted in the derivation of novel RKHS over trajectories. which solely consists of universal functions that have LTT dynamics. Equivalences with function regression in present in previous literature. Another key con-
RKHS allow provide consistency guarantees tribution is a novel rank-independent generalization bound for i.i.d. sampled trajectories that directly describes forecasting performance. The significant implications of the proposed approach are con- firmed in experiments, leading to superior performance compared to Koopman operator and sequential data predictors in RKHS. this work, we confined our forecasts cific length of trajectory data, where the choice of spectra is arbitrary. However, exploring more effica- cious spectral sampling schemes a natural next step for extending our results to asymptotic regimes that include, e.g., periodic and quasi-periodic behavior. It has to be noted that rector-valued kernel methods have limited scalability with growing number training data and output dimensionality. Therefore, exploring solutions that improve scalability is an important topic for future work. Further- more, to enable the use of LTT predictors in safety-critich domains, the quantification of the forecast- ing error essential. Hence, deriving uniform prediction error bounds for KKR is of great interest.
non-recurrent domain for spe-
10

---

Acknowledgements
The authors acknowledge the financial support of the EU Horizon 2020 research and innovation programme "SeaClear" (ID 87 $2125 and the ERC Consolidator grant "CO-MAN" (ID 864686). Petar Bevanda also thanks Jan Brüdigam for feedback, and Vladimir Kostié and Pietro Novelli for useful discussions operator regression.
References
[1] J. L. Meriam, L. G. Kraige, and N. Bolton, Engineering Mechanics: Dynamics. & Sons, 2020.
John Wiley
[2] A. Billard, S. Mirrazavi, and N. Figueroa, Learning for Adaptive and Reactive Robot Control:
A Dynamical Systems Approach. MIT Press, 2022.
[3] V. May and O. Kühn, Charge and Energy Transfer Dynamics Molecular Systems. 2011.
Wiley,
[4] J. Johansson, P. Nation, and F. Nori, "QuTiP: An open-source Python framework for the dynamics of open quantum systems," Computer Physics Communicutions, vol. 183, no. 8, pp. 1760-1772, 8 2012.
[5] T. Brown, B. Mann, N. Ryder, M. Subbiah, J. D. Kaplan, P. Dhariwal, A. Neelakantan, Shyam, G. Sastry, A. Askell, Agarwal, A. Herbert- Vosss, G. Krueger, Henighan, R. Child, A. Ramesh, D. Ziegler, J. Wu, C. Winter, C. Hesse, M. Chen, E. Sigler, M. Litwin, S. Gray, B. Chess, Clark, C. Berner, McCandlish, A. Radford, I. Sutskever, and D. Amodei, "Lan- guage Models are Few-Shot Learners," in Advances in Neural Information Processing Systems, 2020, pp. 1877-1901.
[6] A. Radford, W. Kim, C. Hallacy, A. Ramesh, G. Goh, Agarwal, G. Sastry, A. Askell, P. Mishkin, J. Clark, G. Krueger, and I. Sutskever, "Learning Transferable Visual Models From Natural Language Supervision," Proceedings of the 38th International Conference on Machine Learning, 2021, pp.
[7] R. T. Q. Chen, Y. Rubanova, J. Bettencourt, and D. K. Duvenaud, "Neural Ordinary Differential Equations," in Advances in Neural Information Processing Systems, 2018
[8] M. Biloš, J. Sommer, S. S. Rangapuram, Januschowski, and S. Günnemann, "Neural Flows: Efficient Alternative to Neural ODEs," Advances in Neural Information Processing Systems, 2021.
[9] M. Janner, Q. Li, and S. Levine, "Offline Reinforcement Learning One Big Sequence Modeling Problem," in Advances in Neural Information Processing Systems, 2021, pp. 1273- 1286.
[10] K. Chua, R. Calandra, R. McAllister, and S. Levine, "Deep reinforcement learning in a handful of trials using probabilistic dynamics models," in Advances in Neural Information Processing Systems, 2018
[11] E. A. Theodorou, Buchli, and S. Schaal, "A Generalized Path Integral Control Approach to Reinforcement Learning," Journal of Machine Learning Research, vol. 11, pp. 3137-3181, 2010.
[12] M. Budišić, R. Mohr, and I. Mezić, "Applied Koopmanism," Chaos, vol. 22, 4, 10 2012.
[13] A. Mauroy, I. Mezić, and Y. Susuki, The Koopman Operator in Systems and Control, ser. Lecture Notes in Control and Information Sciences. 2020, vol. 484.
Cham: Springer International Publishing,
[14] I. Mezić and A. Banaszuk, "Comparison of systems with complex behavior," Physica D: Nonlinear Phenomena, vol. 197, no. 1, pp. 101-133, 2004.
[15] Mezić, "Spectral Properties of Dynamical Systems, Model Reduction and Decompositions." Nonlinear Dynamics, vol. 41, no. 1-3, pp. 309-325, 8 2005.

---

[16] S. L. Brunton. B. W. Brunton, J. L. Proctor, and J. N. Kutz, "Koopman Invariant Subspaces and Finite Linear Representations of Nonlinear Dynamical Systems for Control," PLOS ONE, vol. 11, no. 2, p. e0150171, 2 2016.
[17] S. L. Brunton and N. Kutz, Iata-Driven Science and Engineering. Cambridge University Press, 2019.
[18] S. E. Outo and C. W. Rowley, "Koopman Operators for Estimation and Control of Dynamical Systems," Annual Review of Control, Robotics, and Autonomous Systems, vol. 4, p. 2021, 2021.
[19] P. Bevanda, S. Sosnowski, and S. Hirche, "Koopman operator dynamical models: Learning, analysis and control," Annual Reviews in Control, vol. 52, pp. 197-212, 2021.
[20] S. L. Brunton, M. Budišić, E. Kaiser, and J. N. Kutz, "Modern Koopman Theory for Dynamical Systems," SIAM Review, vol. 64, no. 2, pp. 229-340, 2022.
[21] I. Mezić, "Spectrum of the Koopman Operator, Spectral Expansions in Functional Spaces, and State-Space Geometry," " Journal of Nonlinear Science, vol. 30, no. 5, pp. 2091-2145, 2020.
[22] S. Klus, I. Schuster, and K. Muandet, "Eigendecomppsitions of Transfer Operators in Repro- ducing Kernel Hilbert Spaces," Journal of Nonlinear Science, vol. 30, no. 1, pp. 283-315, 2020.
[23] V. Kostic, P. Novelli, A. Maurer, C. Ciliberto, L. Rosasco, and M. Pontil, "Learning Dynamical Systems via Koopman Operator Regression in Reproducing Kernel Hilbert Spaces," in Advances Neural Information Processing Systems, 2022, pp. 4017-4031.
[24] V. Kostic, K. Lounici, P. Novelli, and M. Pontil, "Koopman Operator Learning: Sharp Spectral Rates and Spurious Eigenvalues." " 2 2023. [Online]. Available: http://arrivoorraaa//33220.00000
[25] D. Angeli and E. D. Sontag, "Forward completeness, unboundedness observability, and their Lyapunov characterizations." Systems & Control Letters, vol. 38, no. 4-5, pp. 209-217, 12 1999.
[26] V. Andrieu and L. Praly, "On the existence of a Kazantzis-Kravari// Kravaris/Luenbergee observer," SIAM Journal Control and Optimization, vol. 45, no. 2, pp. 422 456, 2006.
[27] M. Krstic, "Fomward-Complete Systems," in Delay Compensation for Nonlinear, Adaptive, and PDE Systems. Boston: Birkhäuser, 2009, pp. 171-190.
[28] M. Korda and I. Mezic, "Optimal Construction of Koopman Eigenfunctions for Prediction and Control," IEEE Transactions on Automatic Control, vol. 65, no. 12, pp. 5114-5129, 12 2020.
[29] Klus, I. Schuster, and K. Muandet, "Eigendecomposstioon of transfer operators in reproducing kernel hilbert spaces," Journal of Nonlinear Science, vol. 30, pp. 283-315, 2 2020.
[30] M. J. Colbrook, L. J. Ayton, and M. Szöke, "Residual Dynamic Mode Decomposition: Robust and verified Koopmanism,"" " Tech. Rep.
[31] M. Lewin and E. Sere, "Spectral pollution and how avoid it," Proceedings of the London Mathematical Society, vol. 100, no. 3, pp. 864-900, 2010.
[32] M. Khosravi, "Representer theorem for learning koopman operators," " IEEE Transactions on Automatic Control, vol. 68, no. 5, pp. 2995-3010, 2023.
[33] F. Philipp, M. Schaller, K. Worthmann, S. Peitz, and F. Nüske, "Error bounds for kernel-based approximations of the Koopman operator," 2023. [Online]. Available: http://arxiv.or//ab//23..... .08637
[34] S. Das and D. Giannakis, "Koopman spectra in reproducing kernel Hilbert spaces," Applied and Computational Harmonic Analysis, vol. 49, no. 2. pp. 573-607, 9 2020.
[35] C. Valva and Giannakis, "Consistent spectral approximation of Koopman operators using resolvent compactification." 9 2023. [Online]. Available:
12

---

[36] R. K. Singh and A. Kumar, "Compact Composition Operators," Journal of the Australian Mathematical Society, vol. 28, no. 3, pp. 309-314, 1979.
[37] M. Ikeda, Ishikawa, and Y. Sawano, "Boundedness of composition operators on reproducing kernel Hilbert spaces with analytic positive definite functions," Journal of Mathematical Analysis and Applications, vol. 511, no. 1, p. 126048. 7 2022.
[38] P. Cvitanović, R. Artuso, R. Mainieri, G. Tanner, and G. Vattay, Chaos: Classical and Quantum. Copenhagen: Niels Bohr Inst., 2016. [Online]. Available: http://ChaosBook.org
[39] I. Mezić, "On Numerical Approximations of the Koopman Operator," Mathematics, vol. 10, no. 7, 1180, 2022.
[40] P. Bevanda, J. Kirmayr, S. Sosnowski, and S. Hirche, "Learning the Koopman Eigendecompo- sition: A Diffeomorphic Approach," in 2022 American Control Conference (ACC), 2022, pp. 2736-2741.
[41] Bevanda, M. Beier, S. Kerz, Lederer. S. Sosnowski, and S. Hirche, "Diffeomorphically Learning Stable Koopman Operators." IEEE Control Systems Letters, vol. 6, pp. 3427-3432 2022.
[42] E. M. Bollt, "Geometric considerations of good dictionary for Koopman analysis of dynamical systems: Cardinality, "primary eigenfunction," and efficient representation." Communications in Nonlinear Science and Numerical Simulation, vol. 100, 9 2021.
[43] F. J. Király and H. Oberhauser, "Kernels for Sequentially Ordered Data," Journal of Machine Learning Research, vol. 20, pp. 1-45, 2019.
[44] J. Lyons, "Differential equations driven by rough signals." Revista Matemática Iberoameri- cana, vol. 14, no. 2. pp. 215-310, 1998.
[45] P. K. Friz and N. B. Victoir, Multidimensional Stochastic Processes as Rough Paths. University Press, 2 2010.
Cambridge
[46] D. Lee and H. Oberhauser, "The Signature Kernel," 5 2023. [Online]. Available: http://arxxiv.or/abbs//200 04665
[47] M. Lemercier, C. Salvi, T. Cass, E. V. Bonilla, T. Damoulas, and T. J. Lyons, "SigGPDE: Scaling Sparse Gaussian Processes Sequential Data," Proceedings of the 38th International Conference on Machine Learning, vol. 139. PMLR, 10 2021, pp. 6233-6242.
[48] C. Salvi, T. Cass, Foster, T. Lyons, and W. Yang. "The Signature Kernel Is the Solution of Goursat PDE," SIAM Journal on Mathematics of Data Science, vol. 3, no. 3, pp. 873-899, 2021.
[49] I. Mezić, "Analysis fluid flows via spectral properties the koopman operator," " Annual Review of Fluid Mechanics, vol. 45, no. pp. 357-378, 2013.
[50] A. Hurwitz, "über die Erzeugung der Invarianten durch Integration," Nachrichten von der Gesellschaft der Wissenschaften zu Göttingen, siahhematisch-Plyyiaaissht Klasse, vol. 1897, pp. 71-72, 1897.
[51] B. Schölkopf and A. J. Smola, Learning with Kernels.
The MIT Press, 2018.
[52] Ingo Steinwart and Andreas Christmann, Support Vector Machines, 1st ed., Information Science and Statistics. New York, NY: Springer, 2008.
[53] H. Kadri, E. Duflos, P. Preux, S. Canu, A. Rakotomamonjy, and J. Audiffren, "Operator-value valued Kernels for Learning from Functional Response Data," Journal of Machine Learning Research, vol. 17, no. 20, pp. 1-54, 2016.
[54] C. A. Micchelli and M. Pontil, "On Learning Vector- Functions," " Neural Computation, vol. 17, no. 1, pp. 177-204, 2005.
13

---

[55] N. Aronszajn, "Theory of Reproducing Kernels," " Transactions of the American Mathematical Society, vol. 68, no. 3, p. 337, 1950.
[56] T. Hotz and F. E. Telschow, "Representation by Integrating Reproducing Kernels," " 2 2012. [Online]. Available: http://arii..rr/aa///0024444
[57] K. Küster, "The Koopman Linearization of Dynamical Systems." 2015. [Online]. Available: https:/hhomepagss.ass.tiiiii hww.h/sssaaastt///eeioooomoooooo eokkarikkustteng
[58] M. A. Álvarez, L. Rosasco, and N. D. Lawrence, "Kernels for vector-valued functions: A review," pp. 195-266, 2011.
[59] A. Lederer, A. Capone, T. Beckers, Umlauft, and Hirche, "The Impact of Data on the Stability of Learning-Based Control," in Proceedings of the 3rd Conference on Learning for Dynamics and Control, 2021, pp. 623-635.
[60] J. Qin, "An overview of subspace identification," Computers and Chemical Engineering, vol. 30, no. 10-12, pp. 1502-1513, 9 2006.
[61] N. Dunford, "Spectral Theory. I Convergence to Projections," Transactions of the American Mathematical Society, vol. 54, no. 2. p. 185, 9 1943.
[62] Bach, "On the equivalence between kernel quadrature rules and random feature expansions." Journal of Machine Learning Research, vol. 18, no. 21, pp. 1-38, 2017.
[63] Z. Li, J.-F. Ton, D. Oglic, and D. Sejdinovic, "Towards a unified analysis of random fourier features," " Journal of Machine Learning Research, vol. 22, no. 1, 2021.
[64] A. Mauroy and I. Mezić, "On the use of Fourier averages to compute the global isochrons of (quasi)periodic dynamics," " Chaos, vol. 22, no. 3, 7 2012.
[65] A. Caponnetto and E. De Vito, "Optimal rates for the regularized least-squares algorithm." Foundations of Computational Mathematics, vol. 7, no. 3, pp. 331-368, 7 2007.
[66] C. Carmeli, E. De Vito, and A. Toigo, "Vector Valued Reproducing Kernel Hilbert Spaces of Integrable Functions and Mercer Theorem," Analysis and Applications, vol. 4, no. 4, pp. 377-408, 2006.
[67] C. Carmeli, E. De Vito, A. Toigo, and V. Umanità, "Vector Valued Reproducing Kernel Hilbert Spaces and Universality," Analysis and Applications, vol. 08, no. 01, pp. 19-61, 1 2010.
[68] A. Caponnetto, C. A. Micchelli, M. Pontil, and Y. Ying, "Universal Multi-Task Kernels," Journal of Machine Learning Research, vol. 9, no. 52, pp. 1615-1646, 2008.
[69] M. Korda and I. Mezić, "On Convergence of Extended Dynamic Mode Decomposition to the Koopman Operator," Journal of Nonlinear Science, vol. 28, no. 2, pp. 687-710, 2018.
[70] J. A. Rosenfeld, R. Kamalapurkar, L. F. Gruss, and T. T. Johnson, "Dynamic mode decompo- sition for continuous time systems with the liouville operator," Journal of Nonlinear Science, vol. 32, p. 5, 2 2022.
[71] M. Kanagawa, Hennig, D. Sejdinovic, and B. K. Sriperumbudur, "Gaussian Processes and Kernel Methods: A Review on Connections and Equivalences," 7 2018. [Online]. Available: http://arxxvvorr//abs/ https////88000000 1807.820555 007.02582
[72] Scharnhorst, T. Maddalena, Y. Jiang, and C. N. Jones, "Robust Uncertainty Bounds in Reproducing Kernel Hilbert Spaces: A Convex Optimization Approach," 4 2021. [Online]. Available:
[73] B. C. Csaji and B. Horvath, "Nonparametric, Nonasymptotic Confidence Bands With Paley- Wiener Kernels for Band-Limited Functions," IEEE Control Systems Letters, vol. 6, pp. 3355- 3360, 2022.
[74] I.. Bartlett and Mendelson, "Rademacher and Gaussian Complexities: Risk Bounds and Structural Results." Journal of Machine Learning Research, vol. 3, pp. 463-482, 3 2002.
14

---

[75] R. Huusari and H. Kadri, "Entangled Kernels Beyond Separability," Journal of Machine Learning Research, vol. 22, no. 24, pp. 1-40, 2021.
[76] Bradbury, R. Frostig, P. Hawkins, J. Johnson, C. Leary, D. Maclaurin, G. Necula, A. Paszke,
J. VanderPlas, S. Wanderman-Miine, and Q. Zhang, "JAX: composable transformations of Python+NummPy programs," 2018. [Online]. Available: http://gittbb.comggogl////
[77] Pedregosa, G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel, M. Blondel, P. Prettenhofer, R. Weiss, V. Dubourg, J. Vanderplas, A. Passos, D. Cournapeau, M. Brucher, M. Perrot, and E. Duchesnay, "Scikit-leam: Machine learning in Python," " Journal of Machine Learning Research, vol. 12, pp. 2825-2830, 2011.
15


Reviews:
[{"review_id": "InHuNbZhUV", "reviewer": "Reviewer_aayt", "summary": "Learning time-series obeying nonlinear dynamics is essential for forecasting. Often Koopman operator theory is used to replace the non-linear dynamics with a linear system where the state evolution is given by a high-dimensional matrix that has to be estimated/learned. In this paper, the authors bring the kernel methods together with Koopman operators to obtain a mathematical framework where learning guarantees can be studied. Eigenfunctions of the Koopman operator associated to a dynamical system can be approximated by the basis functions of a Reproducing Kernel Hilbert Space (RKHS) and learning these eigenfunctions reduces to functional kernel regression. The paper derives risk bounds using statistical learning theory and demonstrate their method with experiments.", "questions": "1. What is the computational complexity/cost of this method, and how does it relate to the existing methods?", "limitations": "Limitations section is missing.", "rating": 7, "confidence": 3, "soundness": 4, "presentation": 4, "contribution": 3, "strengths": "The paper is well-written with a clear presentation. The problem of Koopman linearization is well-defined and techniques introduced are novel in terms of applying kernel methods for function approximation in this framework. Also, supplementary material is well-written and very helpful!\n\nTheorem 1 clearly defines the operator-valued Koopman Kernel for a single eigenfunction in terms of a chosen base kernel (e.g. Gaussian RBF) and the rank of the full Koopman kernel is obtained by direct-summation of $\\bar D$ such RKHS's. The main problem in this formalism is the choice of eigenvalues for the Koopman Kernel. The authors propose a method for sampling eigenvalues, while it is not clear how to me (See weaknesses).", "weaknesses": "1. In my opinion, previous methods are not explained well, and it is not easy to understand how this method differs from the previous works.\n\n2. It could be helpful to comment on the spectral biases of the Koopman Kernels based on the base kernel and also the choice of Koopman eigenvalues. To me, it is still not clear how $\\mu_i$'s are chosen, and it is essential to understand that point because Koopman kernel requires both the choice of base kernel and the spectrum.\n\n3. Experiments might be extended to cover how the kernel hyperparameters affect the regression risk. For example, it is well-known in kernel methods community that certain spectral decay is essential for generalization. Since there are two spectra involved in this case (spectrum of the base kernel and spectrum of Koopman operator), it would be interesting to see how both affect generalization.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their recognition of our work.\n### **On weaknesses**\n1. Guiding the reader through the limitations of existing work to motivate our own approach seems to have been a less fortunate decision. Admittedly, we see how the explicit difference to the operator regression approach of EDMD might get lost in that. We will elaborate more on the fundamental difference between direct eigenfunction learning or KKR compared to operator learning (e.g. EDMD) in the updated paper.\nSpecifically: Our KKR approach is not a DMD method, as there is no operator regression involved. KKR is a function regression method that uses an RKHS comprised of Koopman eigenspaces to fit an output of interest. Conversely,  EDMD-like methods regress a matrix surrogate of an infinite-dimensional operator and only subsequently project it on an output of interest. We will formally support the above elaboration with a corresponding optimization objective.\n2. Tracing back to operator theoretic results we rely on, it is sufficient to sample a distribution with support on the closed unit disk. On a non-recurrent domain, any and all eigenvalues are legitimate while the eigenfunctions are so rich that they are universal approximators. Additionally, our consistency results prescribe a finite amount of eigenfunctions that leads to arbitrarily accurate forecasting while sampling from any distribution with support on the closed unit disk is the recipe to achieve that. In addressing the reviewer's concerns, we now\nalso show the experimental performance of different sampling distributions that are\nbiased towards specific regions of the unit disc, demonstrating the effect different\nsampling distributions have on the learning outcome. The figures supporting these\nexperiments can be found in the the global response PDF.\n3. We appreciate this interesting remark. Based on the reviewer's request we extend our experiments with a study of the effects different hyperparameters and spectral sampling schemes have on KKR's generalization (attached PDF in the global response, Figure 2). \n\n### **On limitations**\nAn explicit limitations section in the main paper turned out to be infeasible due to spatial constraints. We stress, however, that the conclusion section very much includes the limitations of our approach (also acknowledged by reviewer \\#GSvu). \n\n### **Answers to questions**\n1. We appreciate the reviewer's question. Therefore, we now provide training and forecasting complexities of our KKR method.\nThe complete comparison to Koopman operator time-series regression state-of-the-art is provided in the global response and the attached PDF therein.\n   - TL;DR: The methods are not directly comparable, due to their different solution approaches: in different regimes different methods would prevail.\nHowever, when constructing higher-rank predictors - required for forecast accuracy - our KKR method has a more favorable time-complexity than EDMD/RRR (linear vs quadratic in rank) and generally slightly better forecast complexity."}, {"review_id": "jwNF40S9W4", "reviewer": "Reviewer_xMDr", "summary": "The presented manuscript suggests a novel algorithm for learning dynamical systems trajectories called Koopman Kernel Regression. This method based on the Koopman operator theory proposes a unified way to construct observables such that the transformed dynamical system becomes LTI. Moreover, this study presents convergence results and the generalization risk bounds for the KKR method. The experiments confirm the given theory and demonstrate the performance of the KKR over existing kernel-based methods to learn the Koopman operator based on the observed snapshots.", "questions": "1) what is a particular form of the universal kernel $k$ used in experiments? Do EDMD and KKR use the same kernel inside?\n2) what about the computational complexity of the proposed algorithm? How does it scale with the increasing dimension of the ODE that generates training trajectories? Experiments use only 1d trajectories.\n3) please provide a complete step-by-step description of the proposed method in the algorithm LaTeX environment. It is hard to extract the final algorithm steps from the presented list of theorems and remarks.\n4) the experiment section should be extended to more challenging dynamical systems like the Lorenz system, and systems with high-dimensional states. The runtime comparison with EDMD is very interesting, too.\n5) please provide a more detailed description of the experimental setup. In particular, what spectrum sampling strategy is used and why? What happened if one will use another strategy and how to select the best one?", "limitations": "No limitations are provided in the manuscript.", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 2, "contribution": 3, "strengths": "1) convergence results and generalization risk bounds are developed and supported by experimental results\n2) also experiments illustrate the more accurate and robust w.r.t. prediction horizon compared to EDMD\n3) the introductory part with a description of the Koopman operator approach and the corresponding challenges is well-written and clearly states the objective of this study \n4) kernel for Koopman eigenfunction is derived and theoretically analyzed", "weaknesses": "1) the key difference with the EDMD approach is not explicitly stated in the text, only the issues of the latter\n2) it is unclear how to compose the complete algorithm with data samples as input and prediction as output from the presented theoretical results. The readability and accents in sections 3.2-3.3 should be fixed (from theory to practical usage)\n3) how to select the kernel $k$ is not discussed in sufficient detail, only required properties are mentioned.\n4) the discussion of the sampling strategy's importance for the obtained good experimental results is missing. In addition, how this step affects the performance of the KRR for arbitrary input data is interesting.\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their constructive comments.\n### **On weaknesses**\n1. It seems that the decision to guide the reader through the limitations of existing work, in the hope of better motivating our own approach, was less fortunate. Admittedly, we realize how the explicit difference to the operator regression approach of EDMD might thus get lost. We want to highlight that our KKR is a function regression method that uses an RKHS comprised of Koopman eigenspaces to fit an output of interest. Conversely, EDMD-like methods regress a matrix surrogate of an infinite-dimensional operator and subsequently project it on an output of interest. In the updated paper, we will formally support the above elaboration with a corresponding optimization objective.\n2. We understand that our focus on learning-theoretic results can be disadvantageous in understanding how to use our method in practice. Based on the reviewer's comment, we will provide the pseudocode for our algorithm. Additionally, we will also provide ``scikit-learn`` compatible open-source code with the updated paper.\n3. It seems that we miss a reference to Appendix C of the supplemental in the main text. In the updated version, we will point to it in a footnote. We believe the additional details and experiments therein can alleviate the reviewer's concerns. Based on reviewer's feedback, we now also provide an extended study of kernel hyperparameters and sampling distributions in the global response comment and PDF.\n4. Given the evident importance, we will address these aspects more clearly in the updated version of the main paper. We now also show the experimental performance of different sampling distributions that are biased towards specific regions of the unit disc, demonstrating their effect on the learning outcome. The figures supporting these experiments can be found in the global response and accompanying PDF.\n### **On limitations**\nAn explicit limitations section in the main paper turned out to be infeasible due to spatial constraints. We stress, however, that the conclusion section very much includes the limitations of our approach (also acknowledged by reviewer #GSvu).\n\n### **Answers to questions**\n\n1. We sided with the popular choice of the squared exponential kernel. Indeed, the same kernel and kernel hyperparamters are used for both methods. More details are found in the supplemental, Appendix C. \n2. We appreciate the reviewer's relevant question on complexity. Therefore, we now provide training and forecasting complexities of our KKR method together with a comparison to Koopman operator as well as time-series modeling state-of-the-art is provided in the global response and the accompanying PDF. \\\nWe would like to point out that in the Van der Pol experiment we consider a vector-valued observable which amount to solving 2 regressions. Also, the content of Appendix C demonstrates that our method performs well in a sparse data regime for high-dimensional fluid flow data.\n3. We welcome this suggestion. As mentioned, we will try to fit an algorithm in the form of pseudocode for the updated version. \n4. Based on the reviewer's request, we conducted a run-time comparison to EDMD in the global rebuttal response (attached PDF, Table 2). \\\nRegarding challenging dynamics, the Appendix C of the supplemental already provides additional experiments. In particular, we demonstrate that our KKR approach performs well in a sparse data regime on a high-dimensional dataset in flow-over-cylinder setting.\n5. In our statistical experiments we sampled from a uniform distribution over the closed unit disk to not introduce any prior knowledge or preconceptions.  As any eigenvalues from the closed unit-disc are legitimate (based on Koopman operator theory), we sampled from a uniform distribution over the closed unit disk in our statistical experiments. The idea was that such a \"flat\" distribution would not introduce any prior knowledge or preconceptions. Also, the experiments were randomized such that no single run had the same set of eigenvalues, confirming the theoretical results, which prescribe that the more of them we take, the better it gets. \\\nBased on the reviewer's comment, in the attached single-page PDF of the global response, we now assess different sampling distributions - demonstrating their effect on the learning outcome. \\\nWe appreciate raising the interesting question on optimal spectral sampling. As part of future work, we are currently investigating the problem of optimal spectral sampling."}, {"review_id": "hOfQQkTEMi", "reviewer": "Reviewer_GSvu", "summary": "This papers introduces a, to my knoweldge, new kernel on trajectories for dynamical system learning derived from Koopman operator theory. Consistency guarantees and generalisation bounds are derived using classical RKHS arguments. Numerical experiments on Bi-stable system and Van der Pol oscillator are performed using the proposed technique.", "questions": "- For readers unfamiliar with the theory of dynamical systems it is worth reminding the definition of forward-complete system in line 74, perhaps as a footnote.\n- Equation (1b) is slightly ambiguous as the function $q$ acts pointwise on the path $x: t \\mapsto x(t)$; it would be clearer to write $y(t) = q(x(t))$ for any $t \\in \\mathbb{T}$.\n- It is worth specifying what are the minimal regularity conditions required from $f$ and $q$ in (1a) and (1b) so that the flow in line 76 exists and is unique.\n- The map $g$ in equation (2a) is vaguely introduced as a \"feature map\"; can the authors provide more details regarding this map?\n- I think combining equations (3a) and (3b) into a single equation would make the optimisation objective clearer.\n- In lines 87-88 it is claimed that (3) *is generally not analytically solvable or tractable*; this is not true in general, as (3) is nothing but an ODE and numerous methods for fitting parameters of ODEs are available, including the recent Neural ODE framework [1] and many others.\n-  In lines 137-138 it is mentioned that \"the Koopman eigen-functions in Definition 1 are universal approximators of continuous functions\"; can the authors give more details about this statement? What's the hypothesis space of continuous functions mentioned in this sentence? \n- Lemma 1 is vague and not rigorously stated: what is the space $\\mathbb X_0$? Is $C(\\mathbb X_0)$ the space of continuous functions from $\\mathbb X_0$ to $\\mathbb R$? What's the regularity of $g$? To talk about continuity one needs a topology on input and output spaces; what are the choice of topologies in this Lemma? What does the wording \"... are assigned...\" in line 146 mean precisely? The function $M$ implicitly depends on $c,g,A$ so I don't understand the claim of the Lemma. The norm $|| \\cdot ||_{\\mathbb Y_T}$ has not been been defined.\n- Point (ii) in Theorem 1: what do you mean by universal here?\n- It seems to me that the kernel defined in equation (10) does not take its values on the space $\\mathcal{B}(\\mathbb Y_T)$ of bounded linear operators on $\\mathbb Y_T$. \n- Is the kernel defined in Proposition 1 positive (semi)definite?\n- Theorem 2 needs to be written in much higher level of details, recalling the definitions of $M, \\mathbf \\Gamma, \\hat \\phi$.\n- The link to the github repository in footnote 6 breaks the anonymity of submission rule as it reveals the identity of the some of the authors.\n\n**References**\n\n[1] Chen, Ricky TQ, et al. \"Neural ordinary differential equations.\" Advances in neural information processing systems 31 (2018).", "limitations": "The authors adquately discuss limitations and future work.", "rating": 5, "confidence": 4, "soundness": 2, "presentation": 2, "contribution": 3, "strengths": "The idea of making links between Koopman operator theory and kernel methods is interesting. The proposed kernel is to my knowledge new.", "weaknesses": "There is a general lack of mathematical rigour throughout the paper (see my remarks and questions in the next section).  Some statements are unjustified and at times untrue. For example in line 49-50 it is claimed that the RKHS associated to the derived kernel is the first in the literature that posseding universal features that have linear dynamics. *Signature kernels* [1] have such property and are universal for the much broader class of continuous functions on (rough) paths. Furthermore, the numerical examples are very simplistic and the experimental details are poorely described.\n\nThe paper has good potential. I am willing to increase my rating if the authors are able to address my points below reagarding the mathematical rigour and provide a fair theoretical and numerical comparison with [1].\n\n**References**\n\n[1] Salvi, Cristopher, et al. \"The Signature Kernel is the solution of a Goursat PDE.\" SIAM Journal on Mathematics of Data Science 3.3 (2021): 873-899.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer on the extensive constructive feedback that certainly improves our paper.\n### **On weaknesses**\nWe respectfully disagree with the reviewer on unjustified and untrue statements as we elaborate below. \n##### **Difference to signature kernels**\n- Although the provided reference [1] suggests a valuable kernel for sequence modeling, the Koopman and signature kernels solve different problems (Koopman eigenfunctions vs Goursat PDEs), have different properties and are applied in different contexts (ODEs vs controlled DEs [2]). Our KKR approach exploits the Koopman kernel as means to an LTI state-space model determined by its initial condition. In contrast, when used for forecasts, the signature kernel processes a sequence of past observations driving the predictor - leading to a more local yet nonlinear model. Which model type should be preferred is, of course, debatable and ultimately depends on the downstream task. For optimization-based decision, the computational simplicity of forecasting is paramount. In that case, our LTI state-space model is advantageous as it reduces forecasting to matrix multiplication or a linear constraint.\nFollowing the reviewer’s suggestion, **we support the above theoretical comparison by experiments (cf. global response and accompanying PDF).**\n##### **Linear dynamics & universality**\n- Signature kernels have linear dynamics in the restrictive case when sequences themselves are linear in time [1, Sec. 3.1]. For nonlinear dynamics, the suggested Sig-PDE kernel has piece-wise linear dynamics.\nIn contrast, the Koopman kernel is universal and its dynamics linear, regardless of whether the time-series are generated by a linear or nonlinear dynamical system. We will briefly clarify this in the updated version.\n#### **Simplicity and lack of experimental detail**\n- It seems that we miss a reference to Appendix C of the supplemental in the main text. In the updated version, we will point to it in a footnote. We believe the additional details and experiments therein (including a flow-over-cylinder dataset) can alleviate the reviewer's concerns. \n### **Answers to questions**\n- For the sake of brevity we hoped to cover it by delegating to a reference but try to build this into the main text.\n- We believe the suggested change is possibly misleading at that instant in the paper. The continuous-time model has no explicit time dependence and separates into an ODE and an measurement function - ubiquitous in dynamical systems modeling [3]. Thus, we would prefer to keep our notation.\n- We require Lipschitz continuous maps $\\boldsymbol{f}$ and ${q}$. We clarify this in the updated version.\n- We intended to state this vaguely and set up an abstraction layer to discuss various methods in the literature. This feature/``lifting'' map is an important aspect in learning LTI predictors, and the literature is filled with different choices of this mapping. Finding a Koopman-invariant (KI) choice for it that fits the system at hand is the open challenge we are aiming to address. \n- We can understand how this might be tempting but believe that it would compromise clarity. Satisfying (3b) is the very property that existing works on constructing LTI models are unable to guarantee. \n- The fact that (3) is a constrained functional regression problem does not lead to it being generally tractable, let alone analytically solvable. Approximations need to be put in place to make it tractable using discretization schemes, e.g., time-discretization or space-discretization. We do not introduce such discretizations immediately, which allows us to preserve the core structure of (3), such as Koopman-invariance. We will expand on this aspect in the updated version.\n-  It is considers $\\mathbb{C}$-valued continuous functions over a compact set equipped with the sup-norm. Due to spatial constraints, we refer to our supplemental (Appendix A) that includes auxiliary results.\n-  The set $\\mathbb{X}_0$ is introduced on line 81 and it represents initial conditions that form a non-recurrent domain. It suffices to take $g \\in C(\\mathbb{X}_0)$ of complex-valued continuous functions equipped with the sup-norm. This will be added in the notation paragraph.\nWe thank the reviewer for pointing out a lack of clarity, we will restate Lemma 1 as follows: \"Under the flow of (1a), a function $g(\\cdot)$ is transformed into a Koopman eigenfunction of Definition 1 via...\". \n    * TL;DR for the claim of the lemma: (1) given a sufficiently rich amount of eigenvalues, $g$'s can be replaced by eigenfunctions, (2) by extension $\\boldsymbol{A}$ is diagonal, (3)  $\\boldsymbol{c} =\\boldsymbol{1}$ w.l.o.g. We will explicitly require the sufficient richness in the updated version. The norm corresponds an $L^2$ norm, as implied by the definition of the output space.\n- We mean that any open eigenfunction of Definition 1 can be approximated to arbitrary precision using this kernel.\n- This is due to the regression goal being an operator $M$ which is a bounded operator over the output domain, well-established in operator-valued regression [64-66]. We will include the reasoning in the updated paper.\n- Indeed it is. Considering conjugate pairs of eigenvalues (possible w.l.o.g.), it becomes clear that the kernel is real-valued and positive definite. By polarizing their integrands, it follows that is p.s.d. as the integrals and sums of kernels preserve the base kernel's definiteness. It is probably interesting to note here that even without the usage of conjugate pairs, we did not observe any issues related to non-positive definiteness in our experiments.\n-  We will provide more details in the updated paper.\n-  Footnote 6 provides a publicly available link to the code accompanying [24] so **there is no breach of anonymity**.\n\n[2] P.K. Friz et al., Multidimensional stochastic processes as rough paths: theory and applications, Cambridge University Press, 2010.\\\n[3] H.K. Khalil, Nonlinear Systems, Pearson Education, 2002."}, {"review_id": "iB8fpN74QR", "reviewer": "Reviewer_oAtn", "summary": "This paper proposes a Koopman eigenfunction kernel that induces a Koopman-invariant RKHS and the corresponding analytical solution of the regression problem with this kernel. They further provide theoretical guarantees and numerical experiments to demonstrate the efficacy of the proposed method.", "questions": "[1] In Lemma 2, how to guarantee the Koopman eigenfunctions span the continuous functions space $C(\\mathbb{X})$, or you cannot obtain the universality of eigenfunctions? \n\n[2] The way of omitting spatial modes in Eq. (5) is not standard, and even it can simplify the formulation in your mathematical deduction, the spatial modes are still required in numerical experiments to approximate any measurement function $q$. How do you represent any measurement $y_T=q(x_T)$ without coefficients in Eq. (5)?\n\n[3] In Lemma 1, the authors define $g$ as $g:\\mathbb{X}_0\\to C(\\mathbb{X}_0)$, in this way $g$ is a functional not a function. However, $\\phi_\\lambda$ is the Koopman eigenfunction related to the original dynamics $f$ and $g$ is any functional, how do you guarantee the Eq. (4) stand? If the Lemma 1 was wrong, the Koopman-invariant property in Theorem 1 (iii) cannot stand either.\n\n[4] Why not compare your method with the kernel based Koopman operator learning methods such as [R1-3]? The EDMD and PCR are classic methods but their performance is now lagging behind.\n\n__References__\n\n[R1] Baddoo, P. J., Herrmann, B., McKeon, B. J., & Brunton, S. L. (2022). Kernel learning for robust dynamic mode decomposition: linear and nonlinear disambiguation optimization. Proceedings of the Royal Society A, 478(2260), 20210830.\n\n[R2] Kostic, V., Novelli, P., Maurer, A., Ciliberto, C., Rosasco, L., & Pontil, M. (2022). Learning dynamical systems via koopman operator regression in reproducing kernel hilbert spaces. Advances in Neural Information Processing Systems, 35, 4017-4031.\n\n[R3] Meanti, G., Chatalic, A., Kostic, V. R., Novelli, P., Pontil, M., & Rosasco, L. (2023). Estimating Koopman operators with sketching to provably learn large scale dynamical systems. arXiv preprint arXiv:2306.04520.", "limitations": "[1] In reality, the Koopman eigenvalues and eigenfunctions are unknown and are exactly the learning objectives.  However, this paper requires the Koopman eigenvalues to construct the Koopman kernel, which is impractical in real-world scenarios.\n\n[2] The authors only test their method in toy models instead of real-world dataset such as PIV data for flow over a cylinder or molecular dynamics.", "rating": 4, "confidence": 3, "soundness": 2, "presentation": 2, "contribution": 2, "strengths": "[1] Designing the kernel that induces a Koopman-invariant RKHS is a challenging problem, and this paper contributes to this problem.\n\n[2] The theoretical results in Section 4 guarantee the performance of the proposed method.", "weaknesses": "[1] The presentation is poor, for example, there should be more clearly mathematical description about the Koopman approximation theory in RKHS and the claimed open challenge in the Koopman operator community; in Line 13, a noun is missing after 'aforementioned'.\n\n[2] The authors only test the proposed method in low-dimensional toy models.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "First of all, we wish to thank the reviewer for their insight and the catching a notational oversight.\n### **On weaknesses**\n1.  We thank the reviewer for their insight. To clarify, we do not address the problem of learning Koopman operators but LTI predictors and try to convey that in Section 2.  Although popular, EDMD-like approaches are by no means the only way to obtain LTI predictors. We take a function approximation approach and thus present non-DMD method (at least not in the classical sense) for learning LTI predictors. We will expand on this difference in the updated paper and support it by a corresponding optimization objective.\n2. We thank the reviewer for voicing their concerns but we believe that the flow-over-cylinder experiment from Appendix C is no toy example. There, our KKR-based\nLTI predictor performs well in a sparse data regime on a high-dimensional dataset.\n### **On limitations**\n1.  We definitely see how this could be true for much of the existing work. However,  we stress that, **our work has no requirement of knowing the spectra, which is confirmed in both theory and experiments**. Why this is true traces back to operator theoretic results we rely on. On a non-recurrent domain any and all eigenvalues are legitimate [1] while the eigenfunctions are so rich that they are universal approximators [2]. When learning, different outputs of interest bring out different eigenfunctions for the same eigenvalues - just like different modes correspond to different measurement functions in the mode analysis of an estimated Koopman operator. Finally, our consistency results prescribe a finite amount of eigenfunctions that leads to arbitrarily accurate forecasting while sampling from any distribution with support on the closed unit disk is the recipe to achieve the latter. Importantly, our experiments are randomized in a way that no single run had the same set of eigenvalues, confirming that knowing or learning the spectra is not required.\n2. We thank the reviewer for voicing their concern. However, **there is nothing to suggest that KKR only works on toy examples based both on theory and experiments**, cf. Appendix C of the supplemental. There, KKR performs quite well in a sparse data regime when constructing an LTI predictor in a flow-over-cylinder setting.\n### **Answers to questions**\n1. The Lemma 2 of the Appendix A is an auxiliary result that connects existing results on non-recurrence and operator theory that include universality. For our approach, non-recurrence and a universal base kernel provide the ingredients for obtaining universality, cf. Proposition 1, (ii). Also, in Section 3.3 of our work we provide a \"recipe\" for constructing a universal RKHS using Koopman eigenfunction RKHSs.\n2. Indeed, omitting spatial modes is not standard but possible w.l.o.g. in our setting. Contrary to Koopman operator learning approaches, our KKR directly learns the eigenfunctions for a measurement function. Therefore, the modes are absorbed by individual eigenfunctions so that a normal form, e.g., modes equal to one, can be subsumed. Although we remark on this in footnote of Lemma 2, we will include the reasoning in the updated main text. \n3. We are very thankful to the reviewer for catching this unfortunate typo in our manuscript. It was never our intent to introduce such an ambiguity. We want to clarify that we do not use or require functionals for the maps $q(\\cdot),\\boldsymbol{g}(\\cdot)$ or $\\boldsymbol{\\phi}(\\cdot)$ and related theoretical results. Accordingly, all occurrences of the expression \"$:\\mathbb{X}_0 \\rightarrow C\\left(\\mathbb{X}_0\\right)$\" are going to be replaced with \"$\\in C\\left(\\mathbb{X}_0\\right)$'' as functions are what we consider. **Thus, the correctness of our results is unaffected.**\n4. We appreciate bringing works [R1-3] to our attention. When it comes to problem setting and scope, the work in [R1] is unrelated to our approach. The reasoning is threefold: 1) no LTI predictor learning, 2) missing a statistical learning framework, 3) lack of any statistical learning guarantees. \\\nActually, we compare our own to the results of [R2] both theoretically and empirically. From the two estimators (PCR and RRR) that are shown in [R2] we implemented one, as the paper showed that they are quite similar in terms of performance. Due to the reviewer's request, we will add the RRR estimator to our comparison. The new comparison (cf. global response PDF) with RRR fully confirms our intuition as our KKR estimator clearly outperforms RRR by a similar margin observed for PCR (EDMD). \\\nFinally, there was no possibility of being aware of [R3] as the arXiv preprint appeared weeks after our submission. Although improving the scalability of PCR and RRR, the forecasting performance is bounded by that of PCR and RRR estimators [R2] for the same data. Therefore, we have decided to not include this method for comparison in our paper, but instead show the results for ridge regression using state-of-the-art *signature kernels* for time series modeling as suggested by the reviewer \\#GSvu. \n\n[1] Bollt, E. M. (2021). Geometric considerations of a good dictionary for Koopman analysis of dynamical systems: Cardinality, “primary eigenfunction,” and efficient representation. Communications in Nonlinear Science and Numerical Simulation, 100. \\\n[2] Korda, M., & Mezic, I. (2020). Optimal Construction of Koopman Eigenfunctions for Prediction and Control. IEEE TAC, 65(12), 5114–5129."}]

Author rebuttal:
We thank the reviewers for taking the time to review our manuscript and providing constructive comments that certainly improve our manuscript. We have addressed the reviewers’ comments and will update our paper accordingly. Please find below a joint response to some common concerns of the reviewers as well as the content description of the attached single-page PDF. Please let us know if there are any further suggestions. We look forward to hearing back from you!

We would like to point out that there is a lot of useful information already available in our supplemental that includes, but is not limited to:
   - background on (Koopman) operator theory - Appendix A
   - proofs of theoretical results - Appendix B
   - implementation details and additional experiments (e.g., including high-dimensional fluid flow data) - Appendix C \
 In the updated paper we plan to also add a statistical comparison to other methods (see additional experiments below) in the flow-over-cylinder experiment.
### **Difference to existing work**
We want to highlight that we do not directly address the problem of learning Koopman operators or time-series modeling. While both can be used for forecasting, they ought to be analyzed with rather different tools. Our work **exploits time-series structure and Koopman operator theory to provably learn *LTI predictors***. \
Although popular and supported to statistical learning concepts, Koopman operator learning (e.g. EDMD) is by no means the only way to address the problem of *statistical learning of LTI predictors*. One of our primary goals was to **provide a novel approach for statistical learning of LTI predictors** *rooted in function (not operator) approximation* theory that comes **with stronger theoretical guarantees** than existing work. The practical significance of these guarantees is demonstrated in our experimental evaluation.\
We stress that under the condition of non-recurrence, **knowing or learning Koopman operator eigenvalues is not required** for forecasting performance using LTI predictors - *which is confirmed both in theory and experiments*.

### **Additional results in the attached PDF**
- Following the requests of reviewers #GSvu and #oAtn, we include:
   -  comparison to RRR and Sig-PDE kernel regression in a forecasting task (attached PDF, Figure 1) 
- Based on the suggestions of reviewers #aayt and #xMDr, we provide:
   - complexity analysis in comparison to existing work (attached PDF, Table 1)
   - runtime comparison (attached PDF, Table 2)
   - an empirical study of the effects  different hyperparameters and spectral sampling schemes have on KKR's generalization (attached PDF: Figure 2, Table 3)

#### **Forecast risk comparison to RRR and RR-Sig-PDE (Figure 1)**
For a fair comparison using the suggested signature kernels, we ran ``scikit-learn`` Ridge Regression (RR) using Sig-PDE kernels
https://github.com/crispitagorico/sigkernel. All methods used the same base kernel and the same lengthscales. For fairness w.r.t. RR-Sig-PDE (due to requiring time-delays), $l$-step backward simulation is done from $N$ i.i.d. trajectories of length $H$ used for PCR, RRR and KKR. Shown in Figure 1 of the attachment.

#### **Time-complexity (Table 1)**
Based on the reviewers' suggestions we include: 
- *KKR train complexity:*$~\mathcal{O}(\text{inverse})+\mathcal{O}(\text{creation})=\mathcal{O}(\text{inverse})+D(\mathcal{O}(\text{Grammian})+\mathcal{O}(\text{Kronecker})+\mathcal{O}(\text{sum}))=\mathcal{O}(N^3H^3+D(N^2H^2d+N^2H^2+N^2H^2)) \in \mathcal{O}(N^3H^3+DN^2H^2d)$
- *KKR H-step forecast complexity:*$~\mathcal{O}(\text{matrix product})+\mathcal{O}(\text{kernel sections})=\mathcal{O}(DH)+\mathcal{O}(DNd)$ 

The complete comparison to PCR/RRR and RR-Sig-PDE is found in Table 1 of the attached PDF. As the comparison methods do not use multi-step/forecasting targets, we introduced $n=NH$ to make classical complexities for regression using PCR/RRR and RR + Sig-PDE kernels estimators recognizable.\
Our training complexity is comparable to RR-Sig-PDE and not directly comparable to PCR/RRR due to different solution approaches.
However, when constructing higher-rank predictors - often required for forecast accuracy - our KKR method scales more favorably than PCR/RRR (linear vs quadratic in $D$) with slightly better forecast complexity. The aforementioned LTI predictors scale favorably against RR-Sig-PDE's iterated nonlinear mapping evaluation due to a data-independent matrix multiplication after an initial nonlinear mapping. 

#### **Runtime comparison (Table 2)**
We include these results based on the suggestion of reviewer \#xMDr. However, the results are not entirely comparable due to our KKR model being implemented in JAX https://github.com/google/jax. We adapt PCR/RRR from https://github.com/csml-iit-ucl/kooplearn and Sig-PDE kernels from https://github.com/crispitagorico/sigkernel. Note that, in the short rebuttal time, we did not manage to have the other methods running on the GPU or with precompiled code.
#### **Spectral sampling and hyperparameters (Figure 2, Table 3)**
In Figure 2 and Table 3 we show the influence of different base kernel lengthscales as well as spectral sampling strategies on generalization. \
For lengthscale selection, it is clear that a very small lengthscale does not generalize well quickly due to very slow spectral decay, as expected for kernel-based methods. \
Also confirming our theoretical results, as more eigenvalues are sampled, there is no significant difference observed between spectral sampling strategies and the only effect left influencing generalization is the spectral decay of the base kernel. 

#### We will release ``scikit-learn``-compatible open-source code for the final submission, if the paper gets accepted.

#### We thank the reviewers and AC for their time and would happily provide additional clarification on our rebuttal as well as answer any further questions.

Meta-review: Authors propose a new kernel derived from Koopman operator theory. Compared to prior methods based on Koopman operators, e.g. EDMD, it focuses on making a linear predictor. The paper provides generalization risk bounds, which are supported by experiments. However, the reviewers pointed out that the writing, assumption, and details on the experiments are missing. I strongly encourage the authors to make the final manuscript complete.
Decision: Accept (poster)