Learning the Infinitesimal Generator of Stochastic Diffusion Processes
Authors: Vladimir R Kostic
Hélène Halconruy
Timothée Devergne
Karim Lounici
Massimiliano Pontil
Venue: NeurIPS 2024 poster (2024)
Keywords: Stochastic Diffusion Processes
Infinitesimal Generator
RKHS
non-asymptotic learning bounds
Primary area: learning_theory
TL;DR: An energy-based formulation for learning  infinitesimal generator of stochastic diffusions with theoretical learning bounds.

Abstract: We address data-driven learning of the infinitesimal generator of stochastic diffusion processes, essential for understanding numerical simulations of natural and physical systems. The unbounded nature of the generator poses significant challenges, rendering conventional analysis techniques for Hilbert-Schmidt operators ineffective. To overcome this, we introduce a novel framework based on the energy functional for these stochastic processes. Our approach integrates physical priors through an energy-based risk metric in both full and partial knowledge settings. We evaluate the statistical performance of a reduced-rank estimator in reproducing kernel Hilbert spaces (RKHS) in the partial knowledge setting. Notably, our approach provides learning bounds independent of the state space dimension and ensures non-spurious spectral estimation. Additionally, we elucidate how the distortion between the intrinsic energy-induced metric of the stochastic diffusion and the RKHS metric used for generator estimation impacts the spectral learning bounds.

Learning the Infinitesimal Generator of Stochastic Diffusion Processes
Vladimir R. Kostic
Karim Lounici CMAP-Ecoue Polytechnique
CSML, Istituto Italiano di Tecnologia University of Novi Sad vladimir. kosticdiit.it
karin wwt..ocoooooollllecciiiuou.... edu
Hélène Halconruy
Timothée Devergne ATSIM, Istituto Italiano Tecnologia
SAMOVAR, Télécom Sud-Paris MODAL'X, Université Paris Nanterre
CSML
timothee. wewpoctee.deeeree@iii..oo
helene. tallcorruutteeeecc---uuppprrrr
Massimiliano Pontil
CSML, Istituto Italiano di Tecnologia
AI Centre, University College London massimiliano pon.aliiian.ponn1iiitt......
Abstract
We address data-driven learning of the infinitesimal generator of stochastic diffu- sion processes, essential for understanding numerical simulations of natural and physical systems. The unbounded nature of the generator poses significant chal- lenges, rendering conventional analysis techniques for Hilbert-Schnidt operators ineffective. To overcome this, we introduce a novel framework based on the energy functional for these stochastic processes. Our approach integrates physical priors through an energy-based risk metric in both full and partial knowledge settings. We evaluate the statistical performance of reduced-rank estimator in reproducing ker- nel Hilbert spaces (RKHS) in the partial knowledge setting. Notably, our approach provides learning bounds independent of the state space dimension and ensures non-spurious spectral estimation. Additionally, we elucidate how the distortion between the intrinsic energy-induced metric of the stochastic diffusion and the RKHS metric used for generator estimation impacts the spectral learning bounds.
1 Introduction
Continuous-time processes are often modeled using ordinary differential equations (ODEs), assuming deterministic dynamics. However, real systems in science and engineering that are modeled by ODEs subject to unfeasible-to-modde influences, necessitating the extension of deterministic models through stochastic differential equations (SDEs), see [40 42 and references therein. SDEs are advantageous for modeling inherently random phenomena. For instance, in finance, they specify the stochastic process governing asset behavior, a crucial step in constructing pricing models for financial derivatives [41]. Another compelling application arises in atomistic simulations, where SDEs are used to model the evolution of atomic systems subjected to thermal fluctuations through the Boltzmann distribution B7].
diverse range of SDEs can be represented as d.X, = a( ((X,)dt + b( b(X))WW, where X, = æ. Here, W denotes a (possibly multi-dimensional) Brownian motion, and the functions and b are commonly
38th Conference on Neural Information Processing Systems (NcurIPS 2024).

---

Aspect
1
19
43
Our work
Covers many SDEs
X (only Laplacian)
x Conly Langevin)
Risk metric
(20) metric
(2)) metrie
(2)) metric
energy
Physics-informed informed method
X
V (full info. needed)
K
V (partial info. needed)
Avoids spurious cigenvalues
X
x
X
IG crror bound
On x
of Var Ol
O(n i)
On R ) 2
Onn Spi ) <T
Spectral rates
X
x
x ((n)d))
Time complexity
Ola² + /d/
(vx
((rxAd))
Table 1: Comparison to previous kernel-based works on generator learning. Sample size is state-space dimension d, y is the regularization parameter of KKR and RRR and RRR rank parameter. Our learning bounds are derived in Theorem 2 where the parameters a, 3, quantify the intrinsic difficulty of the problem and impact of kernel choice on learning IG.
known as the drift and diffusion coefficients, respectively. Determining these coefficients from one or more trajectories, whether discretized or continuous, has been a key pursuit in "diffusion statisties" since the 1980s, as seen in works like [22 31 30 However, uncovering the drift and diffusion coefficients alone does not reveal all the intrinsic properties of a complex system, such as the metastable states of Langevin dynamics in atomistic simulations [see e.g. 49 and references therein]. Consequently, there has been shift and growing interest in the Infinitesimal Generator (IG) of an SDE, as its spectral decomposition offers a more in-depth understanding of system dynamics and behavior, thus providing a comprehensive picture beyond the mere identification of coefficients.
Recovering the spectral decomposition of the IG can theoretically be achieved by exploiting the well- studied Transfer Operators (TO) [see[2]] and references therein]. TOs represent the average evolution of state functions (observables) over time and, being linear and amenable spectral decomposition under certain conditions, they offer a valuable means of interpreting and analyze nonlinear systems. However, they require evenly sampled data high rate, which may be impractical. Additionally, TO approaches are purely data-driven, complicating the incorporation of partial or full knowledge of an SDE into the learning process. Thus, there growing interest in learning the IG directly from data, as it can handle uneven sampling and integrate SDE knowledge. The challenge lies in the IG being an unbounded operator, unlike TO which are often well-approximated by Seiheett-Shhiitt operators with comprehensive statistical theory (28). Unfortunately, the existing statistical theory collapses when applied to unbounded operators, prompting the need completely rethink the problem.
Related work Extensive research has explored learning dynamical systems from data [see the monographs 7 32] and references therein]. Analytical models for dynamics are often unavailable, motivating the need for data-driven methods. Two prominent approaches have emerged: deep learning-base methods (5) 12 (36), effective for learning complex data representations but lacking statistical analysis, and kernel methods [2. 6 10, 26 21 24| 27 28 50], offering solid statistical guarantees for the estimation of the TO but requiring careful selection of the kernel function. A related question, tackling the challenging problem of learning invariant subspaces of the TO, has recently led the development of several methodologies 21
52 51), some of which are based on deep canonical correlation analysis B 29]. In comparison, there have been significantly fewer works on learning the IG and only in very specific settings. In (55) deep learning approach is developed for Langevin diffusion, while (3 presents an extended dynamical mode decomposition method for learning the generator and clarifies its connection to Galerkin's approximation. However, neither of these two works provides any learning guarantees. In this respect the only previous work we are aware of are [ revising the Galerkin method for Laplacian diffusion, [19]. presenting a kernel approach for general diffusion SDE with full knowledge of drift and diffusion coefficients, and 43 addressing Langevin diffusion. As highlighted in Table I the bounds and analysis in these works are either restricted to specific SDES or are incomplete. Notably, none of these works proposed an adequate framework to handle the unboundedness of the IG, resulting in an incomplete theoretical analysis and suboptimal rates, sometimes explicitly depending on state space dimension. Moreover, the estimators proposed in these works are prone to spurious eigenvalues and do not offer guarantees on the accurate estimation of eigenvalues and eigenfunctions.
Contributions In summary, our main contributions are: 1) Proposing a fundamentally new idea to estimate the spectrum of adif-adjoint generators of stable Itō SDE from single trajectory. In contrast

---

to all existing works, we exploit the geometry of the process via a novel energy (risk) functional; 2) In a certain sense we "fight fire (resolvent) with fire (generator)" to derive a new efficient learning method that is able to infer the best approximation of the resolvent of IG on the RKHS, independently of the time-sampling: 3) We prove the first IG spectral estimation finite sample bounds using (imperfect) partial knowledge, which notably, overcome the curse of dimensionality present in classical numerical methods; 4) Each important aspect of our learning method, especially in relation to the most relevant existing works, is empirically demonstrated to complement our theoretical analysis.
2 Background and the problem
Our drive to estimate the eigenvalues of the infinitesimal generator for an SDE like ① stems from its crucial role in characterizing dynamics in physical systems. This operator's closed form B relies on the drift and diffusion coefficient b, where we have partial knowledge: is assumed to be known, but a is not. To compensate for the lack of prior knowledge about a, we introduce the system's 's energy as an additional known quantity. Below, we detail the mathematical concepts framing the problem (generator, spectrum, energy) exemplified through the Langevin and Cox- Cox-Ingersoll -oosss processes. For detailed list of notation used throughout the paper we refer Table in the appendix.
Transfer operator and infinitesimal generator A variety of physical, biological, and financial systems evolve through stochastic processes X = (X+) where X₁ E X C Rd denotes the system's 's state at time t. A commonly employed model for such dynamics is captured by stochastic differential equations (SDEs) of the form
d.X,= = a(X,)dt + b(X,)dW, and Xo - xx
where I E X, W = (W?,. (..... ..... WPREER, is RR-dimensional (p IN) standard Brownian motion, the drift a X' Rdd and the diffussion b X Rdxp are assumed be globally Lipschitz and sub-linear, so that the SDE (T) admits an unique solution X = (X2)20 with values in (X, (())). Processes akin equations like ① are diverse, spanning models like Langevin and Cox-Ingersoll- Ross processes (see examples below), with broad applications in science and engineering. The process X is continuous-tin dynamics is described by a family of probability densities (pe )EER, and transfer operators (At)tER) such that for all / R4, E E B(2), x X and measurable function f : X R.
Markov process with almost surely continuous sample paths whose
P(X, E E|xO= x) = p1
f y)dy and Aef = S f(y)p(( v) = [[x(X) =. (2)
Evaluating A at x yields the expectation of f starting from x and evolving until time f, making the transfer operator crucial for understanding X dynamics. The family (Ar)EER satisfies the fundamental semigroup equation Atts = A, A., for s, I E R4. Here, we focus on the transfer operator's effect on the set L (X), choice driven by the existence of an invariant measure for A, on (X, (())) which satisfies T = т for all t E R_, Then, the process X is characterized by the infinitesimal generator L defined for every f E L (X) such that the limit Lf lim, -A,f-f)/+ exists in 2 The operator I is closed on its domain dom T.) which is equal to the Sobolev space
(L) = (s E 22 [I = || + |V < <x).
For SDE dynamics of the form ①. we can prove (see A.2 for details) that L is the second-order differential operator given, for any f (2), x E X, by
Lf(x) = V(()) (xx) 2 (M()) (V))))))())
1
where V² = (02, (eed denotes the Hessian matrix of f.
Spectral decomposition Knowing only the drift a and diffiusion is not enough to compute 2 or to understand quantitative aspects of dynamical phase transitions, such time scales and metastable states. The eigenvalues and eigenvectors of the generator are crucial for capturing these effects. To address the possible unbounded nature of L, one turn auxiliary operator, the resolvent, which, under certain conditions, shares the same eigenfunctions as L and becomes compact. When it exists and is continuous for some μ. E C, the operator L,, (nl L)- is the resolvent of L and the corresponding resolvent set is defined by
p(L) {μ E C| [2I - L) is bijective and L. is continuous).
3

---

We assume that L has a compact resolvent, meaning P(L) # 0 and there exists po E p(L) such that (HO -)- is compact. Under this assumption, and given that (L, (o)(L)) is self-adjoint, we can prove the spectral decomposition of the generator (see A.1 for details) as follows:
L L= -=> E
where (Ae)IEN are the cigenvalues of L, and the corresponding cigenfunctions f. E (X)), forming an orthonormal basis (fi)iEN, are also eigenfunctions of the transfer operator A,,
Dirichlet forms and energy To handle the initial lack of knowledge about the drift a, assume to
have access to another quantity, called the energy. defined as E(f) = lim,-so (x -e f,f))//td for all functions f C2 (2)) for which this limit exists, defining in the way the domain dom(E). The
associated Dirichlet form is the bilinear form defined by polarization for any f.g E dom(e) by
E(f.g) = 10:0--a --- f(Lgdd = 1 (-Lf)gdT.
(5)
For every f dom (), we have E(S) = f), As for every i E N, I E(f: fi) .xx fi(Adi)dd = ->;, we check that the eigenvalues of are negative. To relate to Dirichlet
form, we assume there exists a Dirichlet operator B = s s where = |s.|...| |spl E X' + ((r)=(s1(<) (..). (..x))/((((_(]]e Rdxp Rdxxp is a smooth function s.t. Lf = s(s ff) = sBB) and so that
f (-I)gdd We get that for any f E dom(E)
(s(x)s(1) (((c))sx)TV{(())) VV()))))) / (Bf(=))((((())))))))))))
e(f) / ||s TVI
-G-> - [|
(6)
x
which is reminiscent of the expected value of the kinetic energy in quantum mechanics [7].
In the following, while we discuss in detail the examples of Overdamped Langevin and the Cox- Ingersoll-Roo Ross processes, we briefly mention other ones that have Dirichlet form: the Wright-Fishi Fisher diffusion (in dimension one), which can be defined in the context of population genetics and can be adapted to model interest rates, see [16], the geometric Brownian motion which models the price process of a financial asset, the heat equation, the transport processes associated with advection-diffusio equation, see [33, and the process associated with Poisson's 's equation in electrostatics, see (25)
Brownian motion (a=0) that corresponds to the
Example (Langevin). Let ki, T E R7. The overdamped Langevin equation driven by a potential V Rdd + IR is given by dX,, -VV (V(x))d dtt + 22(kTT)dW; and Xo- I, where ko and T respectively represent the coefficient of friction and the temperature of the system. Its infinitesimal generator L is defined by Lf ===- -V Vf + (k,T) for f E w22(X), Since f(-Lf)g da
-S [
1 g()) = (R:,T) J Vg dr f[(((Lg) dn, generator ( (1
7
L is self- and associated to a gradient Dirichlet form with s(x) = (k,T) (2,)))))))))) Example (Cox-Ingersoll-Roos process). Let d 1, a, b, R, Ri. The Cox- Cox-Ingersoll--oo process is solution of the SDE dX,, = (a + bX, )dl a. X,dW, and Xo = x. Its infinitesimal f. By integration by parts, we can check that the generator L satisfies ((-Lf) = 0-x (x)gg(()) (x))gg(()(((a)) ((xx) = =s f(((Lg) dr, and it is associated to a gradient Dirichlet form with s(x) = /2/V2
generator L is defined f E (x) by Lf = (a bx) + Af.
Learning in reproducing kernel Hilbert spaces (RKHSs) Throughout the paper we let H be an RKHS and let k. : X' x X' + R be the associated kernel function. We let o X - H. be feature map (45 such that k(2 x²) (Φ(z), ((z')) for all x, x' E X. We consider RKHSs satisfying H. C. (2)) 45 Chapter 4.3], so that one can approximate L: (X') + (X)) with an operator G:H. + H. Notice that despite H C L (X'), the two spaces have different metric structures, that
is for all f.g E H, one in general has (f.9) R. F (J.9) In order to handle this ambiguity, we introduce the injection operator S + L2 (2)) such that for all f E H. the object S is the element of C3 (X)) which is pointwise equal to f E H. but endowed with the appropriate C3 norm.

---

and Then, the infinitesimal generator restricted to H. is simply LSx which can be estimated by S for some G E HS (H). This approach is based on the embedding tø: X + H. of the generator in the RKHS that can be defined for kernels & E (X' x X)) whenever one knows drift and diffusion coefficients, see App. Based on this observation, risk I|LS, - S, G, Hs(H..c3 = [=>~<|<((xx) G (t*(())) In scenarios where drift and diffusion coefficients are not known, then to becomes non-computable. However the process has the Dirichlet form (6). one can still empirically estimate the Galerkin projection " L.S, onto H. as considered in ([]. which in fact minimizes the same risk. Yet this approach is problematic due the unbounded nature of the generator and the associated estimators typically suffer from large number of spurious eigenvalues around zero, making the estimation of physically most relevant eigenfunctions unreliable even in the self-adjoint case. Conversely, classical numerical methods can compute the leading part of spectrum without spuriousness issues, suggesting that data-driven approaches should achieve similar reliability. In this paper, we address this problem by designing a novel notion of risk, leading to principled estimators designed to surmount these challenges.
B so that the reproducing property (fp((z),
(I.) holds true. [9] developed empirical estimators of LSS, that essentially minimize the
3 Novel statistical learning framework
In this section, we tackle the challenges in developing suitable generator estimators highlighted earlier. To this end, we introduce risk metric for resolvent estimation that can be efficiently minimized empirically, leading to good spectral estimation. Since L and (pl -)- share the same eigenfunctions, the main idea is to learn the (compact) resolvent, which can be effectively approximated by finite-rank operators, instead of learning the generator directly. However, this approach is challenging due to the lack of closed analytical forms for the action of the resolvent.
in order to approximate the action of the resolvent on the RKHS by some operator G: H + H. we introduce its embedding Xu= X' + H. via the reproducing property (x), l1)
First, given μ
formally given by (x) = So Elo )e Xo = x]dl, see App. B for details. Using this notation, we aim to estimate [(u] - S-)((r) = (Z),
[[al - S,,
E H. i.e. the objective is to estimate Xa(x)) x G* ((z) T-..e
An obvious metric for the risk would be the mean square error (MSE) w.r.t. distribution æ of the data. However, this becomes intractable since in general X₂. is not computable in closed form when either full or partial knowledge of the process is at hand. To mitigate this issue, we introduce a different ambient space in which we study the resolvent,
(X)) {f E dom(L) Ilee = C,1ff = (f. <0)),
where the norm now balances the energy of an observable f: X - R w.r.t. the invariant distribution Ilff and its energy w.r.t. the transient dynamics - (f, (f) Indeed
=
[=)(( = E (x] [|s(
(=))) [
(7)
where the last equality holds for Dirichlet gradient form 66). in which case (X) is simply weighted Sobolev space. Importantly, this energy functional can be empirically estimated from data sampled from ,, whenever full knowledge, that drift and diffusion coefficients of the SDE T.. or partial knowledge, i.e. the diffusion coefficient and Dirichlet operator B in (6). is at hand. With that mind, instead of the standard MSE, we introduce the energy- based risk functional as
min G: KHH
R(G) = E (.,()x(())) G* (())l)))-
(8)
Denoting by Zy: H WH(() the canonical injection, (8) can be equivalently written as
R(G) = [1 L) Zμ G (i [1( -)-1/25 (2I- S (9)
where we abbreviated W (x) and used z (1T-)) recalling that Hilbert- Schmidt and spectral norms for operators A: H + W. are (A
1/ ||||| --22
Therefore, 99) implies that the regularized energy norm, while dominating the classical L3 norm, exerts a balancing effect on the estimation of the resolvent. This leads us to the first general result regarding the well well-posedness of this framework.
5

---

Proposition 1. Given > O, let H C (x) be the RKHS associated to kernel k E (X x x) such that Z E HS (H, (X)), and let Ps. be the orthogonal projector onto the closure of Im(Z,) C W Then, for every ε > O, there exists a finite rank operator G: H - H such that R(G) S (1 - L) 112 W)) +.. Consequently, when k is universal, R(G) S E.
The previous proposition reveals that whenever the hypothetical domain is dense in the true domain and the injection operator is Hillbert-Schmice there is no irreducible risk and one can find arbitrarily good finite rank approximations of the generator's ' resolvent. Note that Zμ E HS (H. W (x)) is equivalent to Zz = S (μI (nI-L)S, being trace class operator, which is assured for our Examples 2 see the discussion in App.
Now, to address how minimization of the risk impacts the estimation of the spectral decomposition, let us define the operator norm error and the metric distortion functional, respectively, as
E(G) 11(1l - Z,, Ge HS (H), and 17(h) = ||hll// ||h lw, h E H. (10)
Proposition 2. Ler G = =E h4 @ 9. be the spectral decomposition of G: H + H.
where Re 2 X,+1 and let f: = Syhi /|S_ãhille. for forie [r]. Then every 0 and i E [r]
14g - 21
2E(G)n(hua)
S = and Ilfi - fill2 S [gap, - (ha)++ 1 Aill Xul where gap, is the difference between i-th and (i + 1)-th eigenvalue of (ul (
(11)
Note that the estimation of the eigenfunctions is first obtained in the norm with respect to the energy space, and then transformed to the C -orm,, as Il 2 //ll ffl f EW,y. Therefore, by controlling the operator norm error and metric distortion (see App. C), we can guarantee accurate spectral estimation. Consequently, this allows us to approximately solve the SDE (1 starting from an initial condition
E[h(x:) Xoo xo==]== = ~ Eierr dul (2., Fu)().
(12)
4 Empirical risk minimization
In this section we address empirical risk minimization (ERM), deriving two main estimators. The first one minimizes empirical risk with Tikhonov regularization, while the second introduces additional regularization in the form of rank constraints.
To present the estimator, we denote the covariance operator w.r.t. C by = S: = EE,=z0((2) & ((x)), the cross-covariance operator T S LSS,, capturing correlations between input and the outputs of the generator, and the covariance operator W,, = Z S (μT-L)), LL)S w.r.t. energy space W. All operators can be estimated from data, depending In this work we focus on the case when Dirichlet gradient form is known, i.c. we can define the embedding of the Dirichlet operator B s' V (x) + (X)) into RKHS do: X + H² via the reproducing property (do(x), h) H [z)
the available prior knowledge.
s(x) Dh(x) RD. h E H. More precisely, we have that dø(z) is p-dimensional vector with components dRo(I), where de0: X + H. is given via reproducing property (dző(z), h) = sx [ Dh(x), & E [p]. Hence, in this case we have that
-E (d((()) @ dø(x)] = EEw (x) ) @ dx((z)]. Moreover, defining wh: X' + H0+1 by wo(x) V/u((x), (//b(x),ddd6)) d10(x) dT E RP+11 we get
(13)
W = Fa (x) @ (xx) + dø(x)) @ do(x)) = Femz (x) @ wo(x))
(14)
Regularized risk Let us first introduce the regularized risk defined for some 7 > by - (ul--)'/2 Gliis(HH.))) + 7y "
GG) = ((uu-L)-1o
G E HS (H), (15)
which, after some algebra, can be written as
RG) = tr | / Gliss + (xIxa()) - |I 1/CC
variance term
bios term

---

where W,y = Wy + joyI is the regularized covariance w.r.t. W.
Hence, assuming the access the dataset D, (ar))e[jn) made of i.id. samples from the invariant
distribution т. replacing the regularized energy e₂ with its empirical estimate e, leads to the regularized empirical risk functional expressed as
(G) [
p. PPps + Ex(xx()]] 177
( | HS(H)
(16)
P.7
where Way and are estimated by their empirical counterparts Wny and C, respectively, via (T3)
Therefore, our regularized empirical risk minimization approach reduces to
min 1 Wh/TGllssc
(17)
and we analyze two different estimators, the first one G,,y is obtained by minimizing regularized empirical risk 17 over all HS (H), and, hence, the name Kernel Ridge Regression (KRR) of
the generators resolvent. The second one T minimizes (17) subject to the (hard) constraint that G most of (a priori fixed) rank r E N and, hence, is named Reduced Rank Regression (RRR) the generator's 's resolvent. Notice that when r = n, the two estimators coincide. After some algebra, one sees that both minimization problems have closed form solutions
G W p.7 Ō and
G W 1/2112
Ww W/W [I 1/20m M4
(18)
where [-I,, denotes the r-truncated SVD of a compact operator.
To conclude this section, we show how to compute the eigenvalue decomposition of (88) To this end, we define the sampling operators S: H + R" and Zy 7 R(L+p)) by
v h(xx), (x;)) Dh(a;)
h: = 0, [n],
(Sh); = vn h()s)
and (Z.l))aati=
k (h(()),
Further, let K = n- (k Ty JEJj(M) E R be kernel Gram matrix, and introduce the Gram matrices NE Rnxpp and M E [Ren xpn whose elements, for k [[+p], E [n]] are
dzz(); and M(a-T)owe((--)att =n (dd(()))
(19)
We note that although we have introduced the above matrices via inner products in H. they can be readily computed via the kernel and its gradients knowing the Dirichlet form, see
Theorem 1. Given > O and > O, let = Ja --NN +MM++yu)))1 - +yI. Let U))e(r) be the leading eigenpairs of the following generalized eigenvalue problem
f (J p.y - <e)Ke = 3, 0. Ku, = 88g, i,j E [r].
(20)
Denoting V, V,=[ui [...///JEEE Ruur and Σ, = diag(ō 1i.., G,),, (v,, wt, wier are eigentriplets of matrix V V2 E Brxx, then the eigenvalue decomposition the RRR estimator G 2, U,TS
is given by G = Exely 2,) 1 h4 @ g. where X μ 1/v,, g. =1; hi. = 22 for U, ( 1/2 NMM (KV, (KV,-NV,2, E (HHp)xx
-1/2 S.V V, of and
The main computational cost of our method, in view of (12). to solve SDE ① lies in the implicit inversion of Ju when solving 20). When computed with direct solvers this inversion is of the order Onn)p²), however leveraging on the fact that ربري is Schur's complement of the (1+p)m (1 p)m symmetric positive definite matrix and using classical iterative solvers, like Lanczos the generalized Davidson method, when r < n this cost can significantly be reduced to n("))" c.f. (18)
Spectral learning bounds
Recalling Prop. 2 in order to obtain the bounds on eigenvalues and eigenfunctions of the generator, it suffices to analyze the learning rates for the operator norm error & and metric distortion r7. For this purpose, we analyze the operator norm error of empirical estimator Gin using the decomposition
E(G) < (1 - -,, 2,, (G, GZa 11 22 ( (Gpo G G -ww
regularization bias
reduction
7

---

where Gp.7 = W is the minimizer of the full (i.e. without rank constraint) Tikhonov regularized risk and
1/2 AAr [W is the population version of the empirical estimator '4.7"
Note that, while the last two terms in the above decomposition depend on the estimator, the first term depends only on the choice of H. and the regularity of L w.r.t. H. In this work we focus on the classical kernel-based learning where one chooses a universal kernel [45] Chapter 4] and controls
the regularization bias with a regularity condition. For details see Rem. of App. E.2 Let ja > 0 be a prescribed parameter, we make following assumptions to quantify the difficulty of the learning problem:
(BK) Boundedness. There exists cw>0 such that ess supl| (x))]| S Sow, i.e. W6 (x, HPP): ---
(RC) Regularity. For some E (0, 2| there exists ca >> 0 such that C² 3 ca w
(SD) Spectral Decay. There exists 3 E (0, 1] and cg s.t X, S Sca j-3 for all j E J.
The above assumptions, discussed in more details in App E.1 are in the spirit of ttte-of-the-ar
analysis of statistical learning theory of classical regression RKHS spaces [14], recently extended
to regression of transfer operators [35] 28). In our setting, instead of relying on the injection into 2 as in the case of transfer operators, the relevant object is the injection Zp into the energy space W.
The first assumption 2ll =l]o(z)) + Exetol de dz)(r)) apart from needing the kernel to be bounded, we also need the Dirichlet form embedding to be bounded. Essentially, this means that the Dirichlet coefficients are not growing too fast w.r.t. the kernel's ' gradients decay. Having this, we assure that Zμ E HS (H. (X))
(BK)
is the main limiting factor of our approach. Indeed, since [()()))
which implies the operator norm error can be controlled.
Another key difference between generator and transfer operator regression is that the covariance w.r.t. the domain of the operator becomes W, μ Zn instead of C S, On the other hand, the "cross-covariance" that captures now RKHS cross-correlations of the resolvent w.r.t domain is simply - Z) = S (u)
L C. With this in mind, (RC) corresponds to the regularity condition in [28] and it quantifies the relationship between the hypothesis class (bounded operators in H-) and the object of interest (nd L)
1-11 Indeed. if has cigenfunctions that belong to a-interpolation space between H. and (X). (RC) holds truc. In particular, if f. E H for all i 2 2
(constant eigenfunction excluded), one has that a Z 1 (c.f. Proposition 7 in App. E.T.
Finally. (SD) quantifies the "size" of the hypothetical domain HL within the true domain (x)) via the effective dimension tr tr(W,,,,,, S (g((r)) 88 which, due to (BK) leads to another notion, known as the embedding property, quantifying the relationship between H' and essentially bounded functions in the domain of the operator
(KE) Kernel Embedding. There exists T E [8, 1] and such that
G = ess sup Oj [ale 22 (x)) 2, [[2, +00, where Z,. Ejer
(21)
グz 2h hy is the SVD of the injection operator Zμ. H + (x)...
Using the above assumptions, we prove generalization bound for the RRR estimator, notably address- ing the general case when the prior knowledge of the Dirichlet coefficient is inexact, i.e.
(DF)
Dirichlet form. For some eE(0, 1) there exists s.: Rd RP so that (1-e) Tsss ((+e)ss holds --.... where s: Rd is such that L= B* for B=s T V: 22(2') + (x)('))
Recalling the form of IG in Equation 3 the Dirichlet form in (FF) of a self-adjoint generator exists whenever the positive definite diffusion part satisfies uniform ellipticity conditions and the drift term allows integration by parts, leading <<x)=b(z:)/ /2. Thus, to obtain the partial knowledge
need, it is enough to estimate the diffusion function s, with the relative error bound in [DF] see App. E.
Theorem 2. Let (DF)
(RC)
(SD) and ( (KE) hold for some E 10. E (0.2) 3
0.
and T E 11. respectively, and let cl(Im[. (S,) = (2)) Denoting X = (s - - S.).
1 E IN. and given ö E 10,1) and r E [m), if RRR estimator is buili from the Dirichlet coefficients
Se: Rd + R and
7P2 a+8 and =;=== मलाश when ZT, Or YZn" the and E =n 213+7 when a <T (22)
then there exists a constant e>0, depending only on H and gap X 1 >.. 0, such that for large enough 2n with probability at least 1 - 8 in the i.i.d. draw of D' from it holds that
E(Gu,) (6r41 A (x++z) + to c 1 +c).
(23)
8

---

the rank reduction bias upper bounded by ($ -))) while in the exact knowledge case (<=0) the bounds on the variance terms critically rely on the well-known perturbation result for spectral projectors reported in Prop. 4 App. A The latter then chained to Pinelis- Sinelis-aakhnnenko s inequality and Minsker's s inequality for self-adjoint HS-operators, Props. 12 and [3] in App. [..31 respectively. When the knowledge is not exact, that is > O in (DF) this relative bound implies that e eW,, N < W
Proof sketch. The regularization bias is bounded by Co 70/2 by Prop. 9 of App. E.2
<e W,,, where the empirical covariance with the inexact Dirichlet coefficient 8, denoted by We. variance, paying the price of additive term c. Combining the bias and variance terms, we obtain the balancing equations for the regularization parameter and then the next result follows. O
μ.' This allows one to control the additional approximation error in the analysis of
First, note that the learning rate 23 implies the 2 -orm learning rate. Moreover, for a 2 T, it matches information theoretic lower bounds for transfer operator learning upon replacing parameters a, B and 5 related to the space W with their 23 analogues (28), see App. E.6 This motivates the development of the first mini-max optimality for the IG learning, for which our results are an important first step. Next, remark that Theorem @agrrantggerreeeee the reliability of fully data-driven methods when the diffusion coefficients are not known but estimated. Furthermore, when b is constant or linear (e.g. Overdamped Langevin and CIR), the classical estimation bounds coincide with the relative error bound of assumption (DF)
To conclude this section, we address the spectral learning bounds stemmiing from the Prop. The main task to achieve this is control the metric distortions, which we demonstrate in App. E.5 In this context, an additional assumption a: 2 1 is needed, since otherwise the metric distortions can blow-up due to eigenfunctions being out of the RKHS space. Importantly, our analysis reveals that
A S (3; 2//2441/;;) + ( In 11 + c), i [r].
(24)
where & = o. 7/ is the empirical spectral bias that informs how good the estimation of the particular eigenpair is (see Fig. a)). 0g being given in (20) and fT = =1 | / [W H- Importantly, (24) reveals that our data-driven method for spectral decomposition of differential operator L does not
suffer from the curse of dimensionality as present in the classical numerical methods, see App. E.6
6 Experiments
In this section, we showcase the key features of our method outlined in Table T We demonstrate that
our approach: (1) avoids the spurious effects noted in other IG methods [19 I. (2) is more effective
than transfer operator methods (27). and (3) validates our bounds a prediction task for a model with non-constant diffusion term. Further details are available in App.
One dimensional four well potential We first investigate the overdamped Langevin dynamics in a po- tential that presents four different wells, two principal wells and then in each of them two smaller wells, given by (zx) = 44x 00.e==00000000(((( exp((S0(22))))) (0((==)))=000ee((((((==))))) exp((800x=005))) - exp((40(x-00.5 This leads to three relevant eigenpairs: the slowest mode corresponds to the transition between the two principal wells, while the others two capture transitions between the smaller wells. Fig. I panel a), we show that the empirical bias $1 ã₁ fi allows us to choose the hyperparameters of the model, that is higher empirical bias coincides with unreliable estimation of the operator's eigenfunction. In panel b) we observe how it varies w.r.t. land-scale (y-axis) and regularization (x-axis) hyperparameters showcasing the robustness of the model, see also Fig. 3 of App. F for hyperparameter j.. Further, in panel e) we show the consistency of our model with the true Boltzmann distribution. Namely, we use our model to forecast the conditional probability density function (pdf) of the system. We perform the same procedure with prefect knowledge and imperfect diffusion coefficient estimated from data. We also report that if the same approach is used with the method described in 19 4 M. no dynamical quantity can be forecast due to the presence of numerous spurious eigenpairs, which prevent the system from relaxing towards the Boltzmann distribution. This issue is further illustrated in panel d) where we show how, contrary to KRR method of [19 we avoid spuriousness in the estimation of eigenvalues.
CIR model Next, with the CIR model we show that our method is not limited to Langevin process with constant diffusion. For this process, the conditional expectation of the state X, is analytically known. We can thus compare the prediction of our model with respect to this expectation using

---

a)
b)
0.00
C)
1.50
0.-
0.25
0.0
9a5
1:55
9
0.3
-000
-//5
1.00
51
-1.35
0.75
0.4
0.-
-1.50
-1.31
9.90
0.1
-2.00
0.23
-223
0.00
US Morsgage rates
d)
work
Iegersogg
10"
0.3
0.8
0.3
02-
011
-10
-2
2015
samplees
2.0
2.0
2.4
1.5
15
1.-
1.0
-15
3.2
0.5 0.
0.3 10
2.5
1.0
0.5
0.0
0.5
10
3.5 1.0
0.3
0.0
05
20
Figure 1: a) Empirical biases $1 = 01 77 and estimation of the first (nontrivial) eigenfunction of the IG of Langevin process under four-well potential. The ground truth is shown in black, our method RRR is red and blue for two different kernel lengthscales. d) Estimation by our method (black) of the eigenvalues for the
process (red) compared the methods in [19 1 for which eigenvalue histogram in blue shows spuriousness. e) Prediction RMSE for the CIR model w.r.t. number of samples. f) Performance of our data-driven method and fitted CIR model on the real data of US mortgage rates. g) The second eigenfunction of a Langevin process under Muller brown potential (white level lines) with its estimation by RRR h) and Transfer Operator (TO) in Observe that TO fails to recover the metastable state.
root mean squared error (RMSE) and compute it for different number of samples to validate our bounds. Conditional expectation were computed on 100 different simulations at In(2) /a which corresponds to the half life of the mean reversion. Results are shown in panel e) of Fig. |
US mortgage rates We have trained our method on real 30-year US mortgage rates dataset and contrasted it with the fitted CIR model using continuous ranked probability scores that are estimated from the forecasts obtained by of each of them, using data from Jamuary 2009 to December 2016. The initial condition was the last week December 2016 and the predictions were made for the years 2017 and 2018. Since the dataset real, we used the imperfect partial knowledge, that is, for our method, we estimated the diffusion coefficient only via a least squares calibration of CIR model over the training set. This allows more flexibility on the drift term in our model.
panel f) of Fig. D Each model has been trained
Muller-Brow potential We next study Langevin dynamics under more challenging conditions: the Muller-Brown potential. Panels gi- of Fig. depiic the second eigenfunction obtained by our method compared the ground truth one, as well the one found by the transfer operator approach, with the same number of samples. Notably, our physics informed approach outperforms transfer operator learning for this task. Note that with different lag times, were able to recover this second cigenfunction.
7 Conclusion
We developed a novel basee framework for learning the Infinitesimal Generator of stochastic diffusion SDEs using kernel methods. Our approach integrates physical priors, achieves fast error rates, and provides the first spectral learning guarantees for generator learning. A limitation is its computational complexity, scaling as n22 Future work will explore alternative methods to enhance computational efficiency and investigate broader suite of SDEs beyond stochastic diffusion.
10

---

Acknowledgements
This work was partially funded by EU Project ELIAS under grant agreement No. 101 101120237, and by the European Union NextGennratiinnUU initiative and the Italian National Recovery and Resilience Plan (PNRR) from the Ministry of University and Research (MUR), under Project PF0000013 CUP J53C22003010066 "Future Artificial Intelligence Research (FAIR)". We also extend our gratitude to the anonymous reviewers for their valuable feedback, which inspired us to enhance the quality and depth of our results.
References
[1] A Cabannes, V. and Bach, F. (2024). The Galerkin method beats graph-based approaches for spectral algorithms. In Proceedings of The 27th International Conference on Artificial Intelligence and Statistics, volume 238 of Proceedings of Machine Learning Research, pages 451-459. PMLR.
[2] Alexander, R. and Giannakis, D. (2020). Operator-theoretic framework for forecasting nonlinear time series with kernel analog techniques. Physica D: Nonlinear Phenomena, 109:132520.
[3] Andrew, G., Arora, R., Bilmes, J., and Livescu, K. (2013). Deep canonical correlation analysis. In International conference on machine learning, pages 1247-1255. PMLR.
[4] Bakry, D., Gentil, I,, and Ledoux, M. (2014). Analysis and Geomelry of Markov Diffusion Operators. Springer.
[5] Bevanda, P., Beier, M., Kerz, S., Lederer, A., Sosnowski, S., and Hirche, S. (2021). Koop- manizingFlows: Diffeomorphically Learning Stable Koopman Operators. arXiv preprint arXiv.2112.04085.
[6] Bouvrie, J. and Hamzi, B. (2017). Kernel Methods for the Approximation of Nonlinear Systems. SIAM Journal on Control and Optimization, 55(4):2460-2499)
[7] Brunton, S. L., Budišić, M., Kaiser, E., and Kutz, J. N. (2022). Modern Koopman Theory for Dynamical Systems. SIAM Review, 64(2):229-340.
[8] Caponnetto, and De Vito, E. (2007). Optimal rates for the regularized least-squares algorithm. Foundations of Computational Mathematics, 7(3):33--368.
[9] Comte, F. and Genon-Catalot, V. (2020). Nonparametric drift estimation for iid paths of stochastic differential equations. The Annals of Statistics, 4836):33363335..
[10] Das, S. and Giannakis, D. (2020). Koopman spectra in reproducing kernel Hilbert spaces. Applied and Computational Harmonic Analysis, 49(2)5573-607
[11] Davis, C. and Kahan, W. M. (1970). The rotation of eigenvectors by perturbation. iii. SIAM Journal Numerical Analysis, 7(1):1-46.
[12] Fan, F., Yi, B., Rye, D., Shi, G., and Manchester, R. (2021). Learning Stable Koopman Embeddings. arXiv preprint 1eXiv.210000633
[13] Fan, J. and Zhang, C.-H. (2003). A reexamination of diffusion estimators with applications to financial model validation. Journal of the American Statistical Association, 98(461) 1):1::11111 18-134.
[14] Fischer, S. and Steinwart, I. (2020). Sobolev norm learning rates for regularized least-squares algorithms. Journal of Machine Learning Research, 21(205):1-38.
[15] Florens-Zmirou, D. (1989). Approximate discrete-time schemes for statistics of diffusion processes. Statistics: A Journal of Theoretical and Applied Statistics, 20(4)554
[16] Grothaus, M. and Sauerbrey, M. (2023). Dirichlet form analysis of the jacobi process. Stochastic Processes and their Applications, 157:81-104.

---

[17] Haag, R. (2012). Local quantum physies: Fields, particles, algebras. Springer Science & Business Media.
[18] Hogben, L., editor (2006). Handbook of Linear Algebra. CRC Press, Boca Raton, FL, USA.
[19] Hou, B., Sanjari, S., Dahlin, N., Bose, S., and Vaidya, U. (2023). Sparse learning of dynamical systems in RKHS: An operator-theoretic approach. In Proceedings of the 40th International Conference on Machine Learning. volume 202, pages 13352.
[20] Jacod, J. (2012). Discretization of Processes. Springer.
[21] Kawahara, Y. (2016). Dynamic mode decomposition with reproducing kernels for Koopman spectral analysis. In Advances in Neural Information Processing Systems, volume 29.
[22] Kessler, M., Lindner, A., and Sørensen, M. (2012). Statistical methods for stochastic differential equations. Monographs on Statistics and Applied Probability, 124:7-12.
[23] Klus, S., Nüske, F., Peitz, S., Niemann, J.-H., Clementi, C., and Schüttee C. (2020). Data-driven approximation of the koopman generator: Model reduction, system identification, and control. Physica D: Nonlinear Phenomena, 406: 132416.
[24] Klus, S., Schuster, I,, and Muandet, K. (2019). Eigendecompositions of transfer operators in reproducing kernel Hilbert spaces, Journal of Nonlinear Science, 30(1):283-315.
[25] Kolluru, S. H. (2008). Preliminary Investigations of a Stochastic Method to Solve Electrostatic and Electrodynamic Problems. PhD thesis, PhD Thesis.
[26] Kostic, V., Lounici, K., Inzerilli, P., Novelli, P., and Pontil, M. (2024). Consistent long-term forecasting of ergodic dynamical systems. In Proceedings of the 41st International Conference Machine Learning.
[27] Kostic, V., Novelli, P., Maurer, A., Ciliberto, C., Rosasco, L., and Pontil, M. (2022). Learning dynamical systems via Koopman operator regression in reproducing kernel Hilbert spaces. In Advances in Neural Information Processing Systems.
[28] Kostic, V. R., Lounici, K., Novelli, P., and Pontil, M. (2023a). Sharp spectral rates for Koopman operator learning. Thirty-seventh Conference on Neural Information Processing Systems.
[29] Kostic, V. R., Novelli, P., Grazzi, R., Lounici, K., and Pontil, M. (2023b). Learning invariant representations of homoenooogss stochastic dynamical systems. In The Twelfth International Conference on Learning Representations.
[30] Kutoyants, Y. A. (2013). Statistical Inference for Ergodic Diffussion Processes. Springer Science & Business Media.
[31] Kutoyants, Y. A. (2020). Parameter Estimation for stochastic processes. Automation and Remote Control, 8133):445-468
[32] Kutz, J. N., Brunton, S. L., Brunton, B. W., and Proctor, J. L. (2016). Dynamic Mode Decompo- sition. Society for Industrial and Applied Mathematics.
[33] Lejay, A., Lenôtre, L., and Pichot, G. (2019). Analytic expressions of the solutions of adyection- diffusion problems in one dimension with discontinuous coefficients. SIAM Journal on Applied Mathematics, 791511771
[34] Li, Q., Dietrich, F., Bollt, E. M., and Kevrekidis, G. (2017). Extended dynamic mode decomposition with dictionary learning: A data-driven adaptive spectral decomposition of the Koopman operator. Chaos: An Interdisciplinary Journal of Nonlinear Science, 27(10).
[35] Li, Z., Meunier, D., Mollenhauer, M., and Gretton, A. (2022). Optimal rates for regularized conditional mean embedding learning. In Advances in Neural Information Processing Systems.
[36] Lusch, B., Kutz, N., and Brunton, I.. (2018). Deep learning for universal linear embeddings of nonlinear dynamics. Nature Communications, 9(1).
12

---

[37] Mackey, M. C. (2011). Time's Arrow: The origins of thermodynamic behavior. Courier Corporation.
[38] Mardt, A., Pasquali, L., Wu, H., and Noé, (2018). VAMPnets for deep learning molecular kinetics. Nature Communications, 9(1..
[39] Minsker, S. (2017). On some extensions of Bernstein's inequality for self-adjoint operators. Statistics & Probability Letters, 127:111-119.
[40] Oksendal, B. (2013). Stochastic differential equations: an introduction with applications. Springer Science & Business Media.
[41] Pascucci, A. (2011). PDE and Martingale Methods in Option Pricing. Springer Milan.
[42] Pavliotis, G. A. (2014). Stochastic Processes and Applications. Springer New York.
[43] Pillaud-Vivien, L. and Bach, F. (2023). Kernelized diffusion maps. The Thirty Sixth Annual Conference Learning Theory, pages 5236-5259. PMLR.
[44] Reed, M. and Simon, B. (1980). I: Functional Analysis. Academic Press.
[45] Steinwart, I. and Christmann, (2008). Support Vector Machines. Springer New York.
[46] Tian, and Wu, (2021). Kernel embedding based variational approach for low-dimensional approximation of dynamical systems. Computational Methods in Applied Mathematics, 21(:)6635- 659.
[47] Trefethen, L. N. and Embree, M. (2020). Spectra and Pseudospectra: The Behavior of Nonnormal Matrices and Operators. Princeton University Press.
[48] Tropp, J. A. (2012). User-friendly tail bounds for sums of random matrices. Technical report.
[49] Tuckerman, M. E. (2023). Statistical mechanics: theory and molecular simulation. Oxford university press.
[50] Williams, M. O., Rowley, C. W., and Kevrekidis, G. (2015). A kernel-based method for data-driven Koopman spectral analysis. Journal of Computational Dynamics, 2(2):247-265.
[51] Wu, H. and Noé, F. (2019). Variational Approach for Learning Markov Processes from Time Series Data. Journal of Nonlinear Science, 30(1):23-66.
[52] Yeung, E., Kundu, S., and Hodas, N. (2019). Learning deep neural network representations for Koopman operators of nonlinear dynamical systems. In 2019 American Control Conference (ACC), pages 4832-4839. IEEE.
[53] Yosida, K. (2012). Functional analysis, volume 123. Springer Science Business Media.
[54] Zabezyk, (2020). Mathematical Control Theory: An Introduction. Systems Control: Foundations & Applications. Springer International Publishing.
[55] Zhang, W., Li, T., and Schütte, C. (2022). Solving eigenvalue pdes of metastable diffusion processes using artificial neural networks. Journal of Computational Physics, 465:1
[56] Zhou, D.-X. (2008). Derivative reproducing properties for kernel methods in learning theory. Journal of computational and Applied Mathematics, 220(1-2) 1-2))-56-446666
[57] Zwald, L. and Blanchard, G. (2005). On the Convergence of Eigenspaces in Kernel Principal Component Analysis. In Advances in Neural Information Processing Systems.
13

---

Supplementary Material
This appendix includes additional background on SDE and RKHS, proofs of the results omitted in the main body and information about the numerical experiments.
notation
meaning
notation
meaning
minimum
V
maximum
r-truncated SVD of an operator
identity operator
HS
(H.
space of Hilbert- Scchmit operators H. + G
(H)
set of rank-r Hilbert-Schmide operators on H.
operator norm of an operator A
[A] HS
norm of operator A
(.(+)
i-th singular value of an operator
);(-)
i-th eigenvalue of an operator time hommeeooogeneuus Markov process invariant measure of the Markov process
X'
state space of the Markov process
(Xele20
p
transition kernel of the Markov process
π
0
drift of the Itō process
b
diffusion of the Itō process
(X)
1.2 space of functions on X' w.r.t. measure
A
transfer operator C (X)
(x*((x))
Sobolev space w.r.t. measure on X'
L
generator of the semigroup W2)2(())
8.
Dirichlet form diffusion
B
Dirichlet form operator on (xt)
1
shift parameter
(x)
regularized energy space canonical feature map
k(2,0) y)
kernel
H
reproducing kernel Hilbert space
S
canonical injection H + (X)
t6
generator embedding
Z
canonical injection P6 + (x) derivative
D
v-directional derivative
D
D"
embedding of the b-directional derivative
Do
derivative embedding
k-th component of Dirichlet operator embedding
ww
Dirichlet operator embedding
0j
j-th singular value of Z,,
J
countable index set of singular values of Zp
2g
j-th left singular function of Zp
hj
j-th right singular function of Z
1
function in (x) with the constant output
Y
regularization parameter
R
true risk
E
operator norm error irreducible risk
Ens
excess risk, i.e. HS norm error
Ro
Dn
dataset (w:) [E]n)
R
empirical risk
S
sampling operator w.r.t. (X))
Z,
sampling operator w.r.t. (x()
C
covariance operator w.r.t. (x))
Č
empirical covariance operator w.r.t. L (2)
C2
regularized covariance operator w.r.t. C (X²)
C
regularized empirical covariance operator w.E.t L (2
W
covariance operator w.r.t. W (X)
Wa
empirical covariance operator w.r.t. WH()''
W
regularized covariance operator w.r.t. (x)
W JLY
regularized empirical covariance operator w.r.t. W
T
derivative-covariange operator
T
empirical derivative-covariinn operator
K
kernel Gram matrix w.r.t. L (2')
K,
regularized kernel Gram matrix w.r.t. L (X)
F2
kernel Gram matrix w.r.t. (F(X)
F Pa
regularized kernel Gram matrix W.r.t WP(X)
M derivative-derivatiee kernel Gram matrix
N
crature-derivatii kernel Gram matrix
G
population estimator of (μl - -)-1 on 7
G
empirical estimator of (μI L) on H
Ga
population KRR estimator
Gp
empirical KRR estimator
G
population RRR estimator
G
empirical RRR estimator
p2
P
spectral projector
P
empirical spectral projector
7
metric distortion
11
empirical metric distortion
A
generator eigenvalue
X
eigenvalue of the empirical estimator
generator eigenfunction in (X)
F
empirical eigenfunction in (X)
h
right empirical cigenfunction
g
left empirical eigenfunction
Cw
boundness constant
PH
orthogonal projector W (X') onto Im (Z.)
regularity parameter
Ca
regularity constant
8
spectral decay parameter
C3
spectral decay constant
T
embedding parameter
Cr
embedding constant
Table 2: Summary of used notations.
14

---

A Background
A.1 Basics operator theory for Markov processes.
We provide here some basics on operator theory for Markov processes. Let X' C RD (d E IN) and (Xt)aeR, be v'-valued time-oomogeeeous Markov process defined filtered probability space (2, F, (Fi)eeR, PP) where F = aXx,, S t) is the natural filtration of (Xt) The dynamics of X can be described through of family of probability densities (pi)LER, such that for all E R4+ E' E (()),
P(X, E E - =) =)= = Pe(ss, yddy.
Let 9 be a set real valued and measurable functions on A. For any t E R+ the transfer operator (TO) A, 9 + G maps any measurable function f E to
(A.))a)= = J ((y)pe(,,
(25)
In theory of Markov processes, the family (At)cER, referred the Markov semigroup associated to the process X.
Remark A possible choice is G = (X). Here, we are interested in another choice of 9 related
to the existence of an invariant measure for Ag i.e., a o-finite measure on (X, B(x)) such that P = π. for any E 14,. In that case, we can choose G = (x) so that, for all La (X), Pif converges
in C (X) as t goes O. Note that Pof and Pao f fdx for suitable functiool f X R. integrable function
Within the existence of this invariant measure the process X is then characterized by the infinitesi- mal generator (IG) L C (x) + 2 (x) the family (At)rER, defined by
A -
L = lim
(26)
In other words, L characterizes the linear differential equation 0,A,f -AAf satisfied by the transfer operator. The domain of L denoted dom(L) coincides with the the Sobolev space (x))
(2) (f E (X)) Il = II + IVVlle <x)).
The spectrum of the IG can be difficult to capture due to the potential unboundedness of L. To circumvent this problem, one can focus on an auxiliary operator, the resolvent, which, under certain conditions, shares the same eigenfunctions as L and becomes compact. The following result can be found in Yosida's book (153), Chap. IX): For every μ > the operator (pl - L) admits an inverse Lp (uI -)-1 that is continuous operator on X' and
(Ml -) / A,dt.
The operator Ln the resolvent of L and the corresponding resolvent set of L is defined by P(L) = {H E C| (MI L) is bijective and L,, is continuous).
In fact, p(L) contains all real positive numbers and (μl - -) is bounded. Here, we assume that L has compact resolvent, i.c. p(L) * and there exists Ao p(L) such that Lno is compact. Note that, through the resolvent identity, this implies the compactness Let then μ p(L). As (L. (oo(L) is assumed to be aelf-anjoint, so does L,g, that L, is both compact and self-adjoint. Then, its spectrum Lp( = C \ p (L,,) is purely discrete and consists of isolated eigenvalues (aE)) such that | + associated with an orthonormal basis (ft)ie (see 553. chapter XI). In other words, the spectral decomposition of the resolvent writes
all resolvents L,, for ja E p(L).
L, Σ E
iEN
where the functions fa E (x) are also eigenfunctions of the operator we get
L
Σ E S
iEN
15

---

Example (Langevin). Let (k,T) E R. overdamped Langevin equation driven by a potential V: Rd - R given by dX, = =d
and Xo Its invariant measure (koT) dr where & denotes
the solution process is the Boltzman distribution π a normalizing constant. infinitesimal generator I. is defined by Lf = -VV (K,T) for f E W/.2(X),, and VV E (2) is positive and coercive, it has compact resolvent. Finally, since [(-Lf)g dm = J v (()))))(())= (1.TT--V( ( =)
g(x)dz (k,T)) J V f Vg dm
f(-Lg) dr, generator L is self-adjoint and associated gradient Dirichlet form with s(x) (k,, ((d)
Example 4 (Cox- Seoo--nngrrool--- -oss process). Let L.bER, cess is solution of the SDE dX, = (a bx,) + o //XX,d,, and Xo = x. Its invariant measure is Gamma distribution with shape parameter a/ 0² and rale parameter b/a². Its infinitesimal gen- erator L is defined for f E (X) by Lf v <²x Af. Note that by integration by parts, we can check that the generator L satisfies J(-Lf)g dv = 6²xx (()gg' x(dx)= f( du, form and it is associated to a gradient Dirichlet form with s(x) = o /2/12
The Cox- Cox-Ingersoll-Ron -ooss pro-
A.2 Infinitesimal generator for diffusion processes
After defining the infinitesimal generator for Markov processes (see A.1 for solution processes of equations like([]. Given smooth function (see for instance [4]. B, p.495) provides for E Rt,
we provide its explicit form C2 R), Itô's formula
/ Σ a.f(X + E'
- -
Σ (x,(xx))((xx
!. VS(X.) dx. + / Tr[x] TT[X](V2/)XX (V)
Recalling
① we get
+ f(X.) = f(Xo) + / 33x ) ds
S. ((x.))
(27)
Provided f and b are smooth enough, the expectation of the last stochastic integral vanishes so that we get
E[f(X;) = = f(x) + / ( (7X00 = x ds
a(x.)
Recalling that
lim (Aef - ft, we get for every x E X
+-++
[(f(X,))
| Xo = f(x)
Lf(x)
lim
1-:
= lim 1 // E[VF(X.) a(x.) + ( [
= which provides the closed formula for the IG associated with the solution process of (I.
-V(()) + + [T (V))))())
(28)
A.3 Spectral perturbation theory
Recalling that for bounded linear operator A on some Hilbert space H the resolvent set of the operator defined as ((A) {A E C A - AI is bijective and (A - Al) is continuous]. and its spectrum Sp(A) A) = C (p( ()) let S Sp(A) be isolated part of spectra, i.e. both and
Sp(A) X are closed in Sp(A) Than, the Riesz spectral projector PA H + H is defined by
P = 2T /
(IT - -1
dz,
(29)
16

---

where I' is any contour in the resolvent set Res(A) with À in its interior and separating À from μ. Indeed. we have that P = Px and H. Im(Px) Ker(Px) where lm(xx) and (er(P,) are both invariant under A and Sp(A (m)r)) = X, Sp( (KTC)))) = μ.. Moreover, Px + Pp 7 and
PaPp. PaPx- 0.
compact operator, then the Reszu-Schauder theorem, see e.g. is a discrete set having no limit points except possibly À = O. Moreover, for any nonzero À E Sp(T). then À is an eigenvauue (i.e. it belongs to the point spectrum) of finite multiplicity, and, hence, we can deduce the spectral decomposition in the form
Finally
441 assures that Sp(T)
A =
Σ Apx
(30)
AESp(A)
where geometric multiplicity of rx = raak is bounded by the algebraic multiplicity of X. If additionally A is normal operator, i.e. AA* A" A"A, then PA = Px is orthogonal projector for each A Sp) and Px = E'= 1/1 @ th,, where v/; are normalized eigenfunctions of A corresponding to and r'X both algebraic and geometric multiplicity of X.
We conclude this section with well-known perturbation bounds for eigenfunctions and spectral projectors of self-adjoint compact operators.
Proposition 3 (III). Let A be compact self-adjoint operator on a separable Hilbert space H. Given a pair () E C x H such that II = 1, let X be the eigenvalue of A that is closest to X and let be its normalized eigenfinction. I min{(A Vll Sp(A) \ {A} 0, then sin( sin(<((f. ()) & ||Af Xfll//.
Proposition 4 ((57). Let A and Å be two compact operators on separable Hilbert space. For nonempty index set J C N let
(A) = min - NXJ, (MA..
denote the spectral gap w.r.t J and let PJ and Po be the corresponding spectral projectors of A and A, respectively. If self-adjoint and for some ||A All (A), then
|| All
PP - Pyll
gap, (A)
B Representations in the RKHS
In section 2 we have defined the infinitesimal generator of a diffusion process and specified its form when associated with Dirichlet forms. These operators act on L3 (X) or specific subsets of it. To develop our learning procedure, we need to understand these operators actions when embedding into the RKHS, and define their versions for feature maps.
IG and Dirichlet operator in RKHS, As reminder, we consider H be an RKHS and let X'xX' + R be the associated kernel function. The canonical feature map is denoted by Ø(z) = k(x, -) for x E X k(r, x²) = (Ø(z), ) ((x')) for all 2r, x' E X. Assuming that i is square-integrable with respect to the measure π. we define the injection operator Sx H. + 2 (X') and its adjoint As a preliminary step, we need to define the
s; : (X) + H. by f .x ((I)(())((az) reproducing partial derivatives in RKHS, which we introduce via Mercer kernels.
Definition 1 (Mercer kernel). kernel function k X X - is called a Mercer kernel if is continuous and symmetric function such that for any finite set of points (I1,..., Nn) C X, the matrix (k((,, x, is positive semindefiitt.
Several standard kernels satisfy the Mercer property with & 2 1, including the Gaussian kernel which we will consider subsequently.
For s E NN and m E N, we define the index set 1Tdh = {a E N : |aa| S s} where |a| = Σ5= aj. for a (aI,..., 0m) E N For function : 1RT + IR and (I), 2m) E IRTH, we denote its partial derivative Da at point (if it exists) as
m
glal
D D²
II
D3 f(x)
Joiz Jammm f(x).
17

---

For function k E C22 x 2) with & C Rd and a E I. we define
ølal
D* ((0.0) k(,,) jmn dazz σ1 yam Jm
k(x, v)-
(31)
and
@la
((0.00((I,5) galyy) 0 ym
k(a, y).
Theorem 3 (Theorem in (56)) Let IN, X Rm and l. be Mercer kernel such that E C (X x X) with corresponding RKHS H. Then the following hold:
For any X, T
(D*h)a(y) = D*k(z, v) E H.
(32)
ii A partial derivative reproducing property holds for ITm (D" = ((D*)),, Yh E H.
(33)
Theorem 3 allows introduce the first and second order operators D and D that act on any feature map ((x) as
Dō(x) = ((D*)))))sse) and D²ő(z) ((D k)))
where the (D) and (D4+tsk)) can be defined via (32). Then, we define the operator that maps any feature map (bx:) to dø(x) = a(x) Do(a)) Denote [=1.... [al x E X + a(x)) [s1()) [{(xx)) Roxd We have
s(x) Dh(x)
Σ (c)))M()))
Σ (D)) 11 = (do(x),
so that we can define the embedding the Dirichlet operator B TV: (x) - (22 into RKHS dø: X' + HP via the reproducing property as (do(x), [BS s(x) Dh(x) e
RP. h H. In fact, dó is p-dimensional vector with components d.0: X - 76 given via (dv h) = sx(r) (z) E [p]. Then, we can define
T = - न (d)()) dø(=)) - Eze @ dx)(=))
which captures correlations between input and the outputs of the generator in the RKHS. Defining wo: X + H by wo(r)= = (//0((), d1o(z),... d)o(x)]T E Rd+ we can consider
Wμ = 2; = (z
= Eg-n (up(x) @ ((x) + dp(x) @ d = Eswaawp(p & wo(r)]
which corresponds the RKHS covariance operator w.r.t. energy space W.
Examples. One way to ensure that the essential assumption (KE) holds is to show that wo fulfils the boundedness condition (BK)). i.e. we E L (X, H Let's show that the property holds true for the Damped Langevin (see Example I and the CIR process (see Example (2) if we consider the
y where k > 0 is free parameter that sets the "spread" of the kernel. As reminder, for every x E X, we have
Radial Basis Function (RBF) kernel fx y) k2 (v)
||)((>)]1 = 1 /1)((=))1 Σ [d2 with ds = sx { Debaar.
Recalling that for the overdamped Langevin process s(x) = (KoT) (/2 )e(4))))))))) ss is diagonal get that for every
so that (s(x)* Dk. (k.(.),s(z)) DDx, = 0 for every E X. As |
x E X. [=) < f CH-
Consider now the CIR process. We have d p and s(x) = o / for any x E X. For the very same reasons, (s(x)) Dk: s(x) Dk., = for every x E X, so that [ |2 S pp== =: ez-- In both Langevin and CIR cases, have wo L (X, H Http when considering an RBF kerne..
18

---

C Statistical learning framework
C.1 Spectral perturbation bounds
In this section, we prove key perturbation result and discuss the properties of the metric distortion. We conclude this section with the approximation bound for arbitrary estimator G B, (H) that the basis of the statistical bounds that follow. This result is direct consequence of [27] and Davis-Khan spectral perturbation result for compact self-adjoint operators, [Τ]]
In the framework of Koopman operator learning [28). spectral bounds are expressed in terms of distortion metric between the RKHS H. and (X),, corresponding to the cost incurred from observing the operator's 's action on the H. rather than on its domain C (X)) Aligned with the risk
definition ( here we measure in certain way the distortion between the H. and W (X) as given in definitions in (TO)
Proposition 2. Let G G EiE th Ts be the spectral decomposition of G: H. where x; 2 R0t1 and let f; Snhi / /1S Sishlllle fori [r]. Then every μ and [r]
14 1 1 X 1p & ((G)m(h)))) and Il - fol2; S 1 Iggap, (ha))
(
(11)
where gap, is the difference between i-th and (i + 1)-th eigenvalue of ( L)
Proof. We first remark that
Then, from the first inequality, using that (ul - L)-1 1 is self-adjoint as operator (x) - (x()). obtain the first bound in (II).
So, observing that for every ( A) -l L) -()((---
(1=-1,)) (-- A)
-1 (1- -x) -1r- -
we conclude that (11- - X,) -1 (1- A) 11 (
E(6) (hx). and
min min(((p-).)) (p A) - (/- Sp((nl -) 11
gap, gap.-E(G)) ((a)).
So, applying Proposition 3 we obtain
|1 (u-
-,,0
sin(<((fi,)))) <
lgap, E n/11)]+
[gap, [(G) n(hi)]+
Igap: E(G))
Since, clearly IIft - fal' <5 2(1 cos( (2(( fo)) s 22inn( fi)), the proof of the second bound completed.
Next, we adapt the 288 Proposition 1] our setting follows.
Proposition S. Let & (,)... For all i E [r] the metric distortion of Thp w.r.l. energy space Wil WH(X)) can be tightly bounded as
1/ VWIW s (ho < ||/||/ o (2, ()
(34)
D Empirical estimation
The Reduced Rank Regression (RRR) estimator is the exact minimizer of 17 under fixed rank constraint. Specifically, RRR is the minimizer Gy of (G) within the set of bounded operators
19

---

(S) on H that have rank at most The regularization term y/l is added to ensure stability,
The closed form solution of the empirical RRR estimator is
G - "
(35)
while its population counterpart is given by G W 1/ [ In order to prove Theorem
1/2 pin Crr.
1
recall kernel matrices in (19) and define
pk
VTI and F P.
pKy VIN
VAN M + +p2d
(36)
V
M.
Now we provide the explicit forms of the matrices N and M in the case of Langevin (see Example and CIR (see Example 2) processes, considering an RBF kernel k(z, y) k2 (y) = xxp(-x|[z-0p yl ".
I
As a reminder (see (36). E RR"xpa and M RMaxp are Gram matrices whose elements, for lo E [1 + p). E [m] are given by (((( =n (0(x)) dao()))) and M (A-D))HH =n ( desd))))))
where dz0(z,) sx(2;) DD)xz,) and DØ(Ø;) Dk(zy,) is defined by (31). For [p], i.j [n],, we have
(d(x;), dx0(x;)) = (Kiz, 8h (x;)) Dkinsm
A(,,31 + thes) - k(,,)g))
-Kke, (a,,,,((3))) (n(g))
lim
1--0
h
(E(d)
lim ( (Ke. sk(T,)] T k((,113 + +hes) k((,,1z)
h-0
h
LE(d]
k(1, 2j + hel) k(zr); 2rj)
=sx(z;)
lim
1-+0
h
1E|d]
= sx,( (D(O,))))(,
= (x;)
((x! - =z )( x 23))
where we have used the continuity of the inner product to get the third line and the reproducing property obtain the following one. Similarly, for k, E lp) [n], we get
(dr((x;), deo(x,))) 1) - (x:) T (Dler, ((Dler
)) Vrag s<(x;)
lim
k(, + hei) ki,, 5r
h-50
LE/d
lim ( (2((x)) (Dlev. (e(e))(((3)))
k(,,0, + hee)) k((,
h
1E(4
k( k(-, + hhee) - -k(-,zz)
lim (Dler. 1rre sa (x()))((())))
h
LE(d)
Dlev k(1,, + her) Diev. (K
=x(()))s((3)) T lim h
h
M.tEld
(=)())((( T (Dr (Dr,k(ez 2,)r,
where we have used the partial derivative reproducing 32) property and where we define for I' + I,
k(v,))((I;, 23 + hea) D (l
Dr Dr,k(er,
lim
h0
(x - u/ )k(s, + hei) - (aft) - )
lim --2y
1-40
h
(x 2/) (()) -0,00)/((()00a0 (k(,,))))))
(37)
and for l [d].
Deer k( + hea) Die1.0)
Di,k(1s,
lim
h-10
h
(z) 2
(1)
(1) - -1)k(Ti, 2j + hez) (Ir x 100) )((),)))))))
= =-2 lim
h-30
h
= - - (x (1) -
() -
= 27(1 - 27(x!)) (9) 2000 (1),, (38)
20

---

For the overdamped Langevin (see Example 1 for k I d, sx(2,) = (KoT) and
(x(r))s((r)) = 7 so that
Ni,(k-1)n+5 = (o(x)), de(z))) (d(t)) 27(ko2)
(a) - xxx
(((1)))()))
and
((--))+ ((da)))))nn(((((( de = (N,T)) xk(,,,
where the elements of 2 are given by 37 and (38). For the CIR process (see Example 2 in dimension 1, we have s(x) = / 22 Then,
Ni,(k-1)nts ( (6(z;), (2
V20m Viyla; xz)) expe) xx,1²)
and
((--))=s(((-))na=== = (deo(( oamm / exp(-xlpz 2z1²).
Based on the previous formulas, using Theorem which we prove next, one can estimate generator's eigenpairs in practice.
Theorem 1. Given > 0 and 7. > O, let NNM leading eigenpairs of the following generalized eigenvalue problem
+yl. Let (07,vi)set be the
" (Jar - -1)Ke, = Un 0. Ku, 8 i,j E [r].
(20)
Denoting V - l.../vrl Rexr and Σ,, - dia Gr), (v,, w w MeEr are eigentriplets of matrix V V, EP E R"xx, then the eigenvalue decomposition the RRR estimator G = 2 U,V M.7
given by p.y = EvErr (1+ 7 hú @ Di where X, = 1- 1/v, Di = 11 -1/25' w and hi = U, for U, = ( 11-2121 BEEiII-N(M++MNI ( -KV,-=VVE, ( (Ht))xxx
Proof.
First, note that μ pJ > O is exactly Schurs's complement w.r.t. second diagonal block of Fug O, and that, due block inversion lemma 18] have that
JJJ JAY F F-11 = (A-1/2(M + +yyl)= N
1/2j--1
N(M + - M.7
(39)
A
where A is np X np matrix. The first step in computing the RRR estimator lies in computing truncating SVD W-1//O, that is CW n.7 problem formulation
Cq = u. E [r]. Now, using the low-rank eigenvalue 18] we have that qi Shu, U, and SW n.7 SW-1C$- v; = o, Now, recalling that S |1-1/2 2 we obtain and that (Z n 117 -1 + 1 10-1 2,,. we obtain
1/21 OF F 1/2 0] Kv; = 0770.
which after some algebra, using 39
- (I- Kit = =070%
ie. 11 -1- (1- KV, = V.
By normalizing right singular value functions qi of W/-l/2 C that is by asking that (g., 9: of Ke; 1, we obtain that [Wr1/²C], W 1/2 ŌQ, for Q [ai I...I gr]. In other words, we have
Gnn W 4. 22 1/2
SS v, V. = 2F n. Iu 1/2 0J KV, = 2,U,V.S. KV, 7 KV, ~VV ' -
where we used that J n.7 2] to obtain the result.
Once with this form, we apply Theorem
Finally, next result provides the reasoning for using empirical metric distortion of Theorem 23
21

---

Proposition 6. Under the assumptions of Theorem 1 for every i E [r]
7 |7H,ll
" F.U.U [[ U, wofl1
(40)
and
17 - nhis) - (n((t)) ^ n7:) (hr) i/s W Will
(41)
Proof. First, note that (40) follows directly from Theorem
Next, since for every i E [r].
(7:)- (n(hx))) 2 (Te, (W Ww)h)))
S S| Walb
h
we obtain
i (n(()))) -2
((((h)))-/
S (n(h)) A AR)) 1W - -Wall.
7
(n(h:)) ii²
E Learning bounds
E.1
Main assumptions
Recalling the form of IG in Equation 3 the Dirichlet form in (DF) for a self-adjoint L exists whenever the positive definite diffusion part satisfies uniform ellipticity conditions and the drift term allows integration by parts, leading to s(x) = Mx)/V2. Thus, to obtain the partial knowledge we need, is enough to estimate the inexact Dirichlet coefficient s, that satisfies the relative error bound in (DF), When the sample paths are observed continuously, such Se can be directly identified from these observations, making it a non-statistical problem. For discrete realizations of the process, the diffusion coefficient can be estimated non-paramn parametrically using various methods, such as pathwise estimation by computing the variance the increments over small intervals (20), kernel-based methods [I5), local polynomial regression
[3] Remark, however, that the estimation of the drift (needed for the full knowledge) is a much more demanding task. Different methods are reviewed in (30) and references therein. A more recent approach [9] drawing inspiration from particle systems, consists in constructing estimates from several i..dd paths the solution process.
Next, observe that (BK) implies Zμ, E HS (H, (x)),, which according to the spectral theorem for positive self-adjoint operators, has an SVD, i.e. there exists most countable positive sequence (o5))EJ, where J = {1, 2, 2..,,..., S N, and ortho-normal systems (2j)zEJ and (hy)seu of cl(Im (Z,))
and (er() respectively, such that Zyh = 022j and Z3 = o,h,, E J.
Now, given a 2 0, let us define scaled injection operator Zp : H + W (X) as
Zu.a Σoje3 2 hj.
(42)
JEJ
Clearly, we have that Zμ = Zp.11 while Im Z,. cl(Im Zn)). Next, we equip Im Im(Zu,a) az with norm 1 to build an interpolation space
M|a - { & (m(Zpa)) 11 = -o, 20-20(1 (f (f,zih < <00
JEJ
noting that the inner product in WH(X') given by bilinear energy functional
(f.9)w = (( - Lg)
We remark that for the space [M] just an RKHS seen as a subspace of (X)).. Moreover, we have the following injections
[M] + [H]1 [7] + [Ho (X)),
22

---

where 2120 02 2 0.
In addition, from (BK) we also have that RKHS H. can be embedded into
Wi (x) {ff W (2) II w ess suplff()]) ((()) <
that is, for some (0,1]
[H]] + [H]- - WH (X) - W
Now, according to 14. if Z 00, [7| + (X) denotes the injection operator, its boundedness
1/1, j E J, and the implies the polynomial decay of the singular values of Zy, i.e. (z)) S j-1 condition (KE) is assured.
Assumption (SD) allows one to quantify the effective dimension of H in ambient space (= allows one estimate the norms of whitened feature maps, in our generator setting vector-valued since they define rank(1+ +p) operators on H.
while the kernel embedding property (KE)allows
{(x) == p wozl H
(43)
This object plays a key role in deriving the learning rates for regression problems (see following result is bounding it.
29 and the
Lemma 1. Le (KE) hold for some E [3, 1] and er E (0). Then,
(17) , < 1, and |IEII3 ess supp(()))] <c((x))
E Eg [(t))
c+ (17)
, 1.
(44)
Proof. W.l.o.g. set == 1, observing that the only change in the proof is scaling 7 > 0. We first observe that for every j E J from definition of we and fact that h2 (z) [2,))](()) x-a.e., it holds that
m-..e.,
vE[++p]
implying that EiE(L+D) (wi((z), hzy)
(T() So, for every r 0.
[k(())
Σ Σ
(W-1/3 w() hT)
Σ
JEJ (EJiE(It)) of + , (wid(()))))) h,)²
Σ
1
jEJ [E(1+n)
2(1-))
(w;o(x), hy)
(09 (wid(a),h)) -- +1
Σ
oz
09 +
of
of
iEJ [E(1t)]
JEJ (E(L+p)
p (x))))))
LZ h,
LZ,,),((())))
S
27 oj =
-
of
JEJ
JEJ
and, due 21 we obtain || T On the other hand, we also have that tr(E,~>I((z) @ {(x))) = tr( W,
which is effective dimension of the RKHS H in W (2)) Therefore, following the proof of Fischer and Steinwart [14 Lemma 11] for classical covariances in L2, we show that the bound on the effective dimension is
tr( Ww) = Σ
<
1-87-8 C c-
,B<1 < 1,
(45)
of +7
, =1
JEN
For the case 3 = 1, it suffices see that
||2/|lw C_(2z] < |1asllwg
= ess ess 1 12)x))) 2, 21(x)(L2=2,1(00)
and, hence
tr(
Σ o3 D -
JEN
JEN
For 3 < 1 we can apply the same classical reasoning in the proof of Proposition 3 of g
23

---

Proposition 7. If the eigenfinctions of L belong to [H] then Condition (RC) is satisfied.
Proof. Note first that the resolvent (nI - -) admits the same cigenfunctions as the generator L, meaning that Im((ul - -) Z,.) C (I p,,)). But according to [54] Theorem 2.2], this last condition equivalent to Condition (RC)
Finally we prove here Proposition 8
Proposition 8. Given > O, let H (X)) be the RKHS associated to kernel k E C2 x x) such that Z HS (H, (x))), and let Pu be the orthogonal projector onto the closure of Im( ( S W (2). Then, for every > O, there exists a finite rank operator G: H - H such that R(G) E [1 - (x) L) Zp W) + Consequently, when k is universal, R(G) S &.
2
Proof. Recall first that since Zμ E HS (H. W (X)), according to the spectral theorem for positive self-adjoint operators, Zμ. admits an SVD. Its form is provided in (42 taking a = 1.
Now, recalling that [-r denotes the r-truncated SVD, i.e. [ ΣyEin jz2; @ hj, since 2p [ li = Ej>r o3. for every d 0 there exists E N such that 12,, 1ZH.DIIHS < pd/3. Consequently since all the eigenvalues of are non-positive, [1 - (Z, - [Z HS I 2| - [12]r/lHs// μ S 5/3. Next since (m(Pu( - L)- Zn,) S cl(lm((z)) for any j E [r]. there exists 95 E 7 s.t. - 2j - -Zug; || S 3r. and, denoting B, := ΣjEn o39; @ hy we conclude [|P - L) "[z,]- BrIl S 6/3. Finally we recall that the set of non-defective matrices is dense in the space of matrices [47]. implying that the set of non-defective rank-r linear operators is dense in the space of rank-r linear operators Hilbert space. Therefore, there exists non-defective G E (H) such that ||G B, llus < 6/(30, (z,)). So, we conclude
-
[1(d - ((cul-)-122 Zn GIHH
$ |1(I -p))(u) - -1)-1 Z lls + 11p1
1IpI (III-L)- L)- - Zplrllls
- - L()
Z,IH.S + IZ, -B,)|HS
$ II - -H] - L)) -Ll Zullls +
Example 5. These three conditions depend the process X (through its generator L) well the chosen RKHS. They are satisfied, for example, by choosing for k as a Gaussian kernel. Indeed, the sub-linearity conditions on A and B required to ensure the existence and uniqueness of the solution of the process of (1), also ensure that A and B are sufficiently "nice" to fulfil, notably, condition ( (BK).
E.2 Bounding the Bias
Recalling the error decomposition and passing to the and L -ormss, we have that
E(Ō) S L( 2,-1-11--- 1//
p.y ( (46)
regularization bias
rank reduction bias
estimator's 's variance
and continue prove the bound of the first term. Note that, while this proof technique is standard for operator learning (29 35), we present it here for the sake of completeness.
Proposition 9. Let G = WH.C for y > 0, and PH: (I) + WH(X)) be the orthogonal
projector onto cl(Im(Z,)) the assumptions (BK) (SD) and (RC)oold, then [|G S caoW
(a- 1)/2
for a E [1, ||G 5 ca (17)) 1)/2 for E (0, 11. and
- Zp Z.G S ca (mn) + II(I Pn) -
(47)
Proof. Recalling that PN = ΣJEJ Zj 3 2j, start by denoting the orthogonal projectors on the sub- space of l leading left singular functions of Zp, Ps EjE(N) 2j @ æj. respectively. Next, observe
24

---

that C-Z,(uI-L) -1 Z and ", + - Zp.. Therefore,
(WI- - Z -,, -,, 7n =(1
and, hence
(nI-L)-' 2p.)) -2
Σ jEJ
2j 82 22
of + pn
2j @ (Knts) (ul-L)- Zn
(0) + uソ)o,
jEJ
Σ JEJ
(o) + +n)oo 2j 6 hig C.
19
Therefore, for every & E J, (ul
- Zp - GnG 11 becomes
|(E
My
23 8 h,
C²
Σ (ož+u7)0,
h2 8 2j
JELN
SE(M]
HWW
which, due to (RC), implies that
170
Σ jE(K]
2j 2 2j
of + My
W-W
On the other hand,
Σ py JE(M] of +
Σ (o)(un)
2j 2 z 28
03 -1 -1 3j 2z 1 (17) $
2j 2 2j,
.E(A]
jE[M
where the inequality holds due to x² < r + 1 for all 2 and s E [0. 1]. Since the norm of the projector equals one, we get [[([u--L)= 2μ - G,, S Next, observe that
(Mr)
1227
Σ
(o7+ (2)2) @ (Z)33)
which is bounded by
22-2 (
1272 0z 2a
Σ
11 S Σ
(0) T0)+ S Σ of
jEj(
JEJ(M
]EJ((M]
Using triangular inequality, for every k E J, we have that bounded by
7,y G Glw
and, therefore,
|| PH((uT-L.)- -,,
S S
+ Σ (078)
JEJ((M)
25

---

and, hence, letting k + oo we obtain (P - - GyG S ca(x)) 3. Hence, 44 follows from triangular inequality.
Li and consider two cases.
To estimate |Guzyl note that (RC) implies |[ S Ca |
M4?
First, if (RC) holds for some [1, 2). then, clearly || 1Guu a E (0, 11. then
(a- 1)/2 c
On the other hand, if
o tt
(17)) (03(0))) 0f (2x) S (17)
140
(m7)= =
0f + Mn
and, thus, ||Gnyll S co (1-11//2
Remark 2. Inequality 47 says that the regularization bias is comprised of a term depending on the choice of 7/ and on term depending on the "aliinment" between H. and -m( Z The term [] Px) (uT 1 ll can be set to zero by Iwo different approaches. One option is to choose kernel that, in some way, minimizes | (n] L) Zull- Another is to choose universal kernel 45 Chapter 4], for which -)
S (lmmmS))). While we here develop theory for universal kernels, deep learning approaches that leverage on our approach can be
developed, which is the direction to pursue in future.
In order to proceed with bounding the bias due to rank reduction for both considered estimators, we first provide auxiliary result.
Proposition 10. Let B = WWT//2 C, let (RC) hold for some E (0, and for j E N denote
\j (1I- (x
(48)
Then for every E N,
a/2 5 N1 - c
S (B) s X7.
(49)
Proof. Start by observing that
B' B [nu () Zul" Z,W ",
Zz
-(uI-L)- (pI ( Z,, - Z,, (2,27 (pI -L)--D)-
implies that, using
and [ut-
(,,--))- -))
-)) -)-
Σ P4y (5 @ (8*2,) B' = 2,
jEJ of + pn
Next, similarly to the above, for every h E J, we have
2a
1n (8,2))
<a²
Σ
2082/
of + 4r
(17)
JE(N]
jE(N]
(03(4r)) (0002/))1)0/00 (22))012
=02
2z002z
//2 [ 1|/2
03 +1
SE(N]
and
7 (5 @ (5t=)) 0f + u jE][]]
Σ
S Σ
7 وي < 2 Σ
(0280013.
03 + + ur
JEJiK
So, before, letting 1s + get the result.
As consequence, we obtain the bound for the rank reduction bias of the RRR method.
26

---

Proposition 11 (RRR), Let (RC) hold for some E (0. Then the bias of G due to rank reduction, recalling
bounded as
a/1 x - co Sw (ur)a/4 2c, (1y) (1na)(1 < (G Gur)) & /itt (50)
Proof. Observe that
2| (G G/n))1 W1122 (G G/,))l
(15))1=110 = - [B| = or+1(B) I or+1((uI-L)) Zp)
while
|Z G/ 2 | 1 - Gprrl
2 ( ca w l (17)
2ca (1 (1Aa)/2
E.3
Bounding the Variance
E.3.1
Concentration Inequalities
All the statistical bounds present will relay two versions of Bernstein inequality. The first one is Pinelis and Sakhanenko inequality for random variables in a separable Hilbert space, see g Proposition 2].
Proposition 12. Let Ai, E [n] be i.dd copies of random variable A in a separable Hilbert space with norm ||-11-- If there exist constants >> and a > O such that for every m 2 2 F|[A||m < AmlAA m a², then with probability least 1 - 5,
1 Σ REIM
2 vn
2 σ
A2
A.
EA
log
02
(51)
On the other hand, we recall that in 39 dimension free version of the non-commutative Bernstein inequality for finte-diieensionl symmetric matrices proposed (see also Theorem 7.3.1 [48] for
an easier to read and slightly improved version) as well as extension to self-adjoint Hilbert-Schmidt operators on a separable Hilbert spaces.
Proposition 13. Let A,, i E [n]] be i.i.d copies of a Hilbert-Schmice operator on the separable Hilbert space. Let ||||| S almost surely, EA 0 and let E[A?] = V for some trace class operator V. Then with probability at least 8,
Σ sEin]
4
S
2c
(a + V
21[V|
Cx(
(52)
3n
where
Cx = log log tr(V)
||V]|
Proposition 14. Given 8 > 0, with probability in the i.i.d. draw of (ax) from n, it holds that
P(I,W Will S {n) 2 1-5.
where
2ew L(8) + Wall C()) and L(6) log
4 tr((W))
(6)
(53)
8 [W
Proof. Proof follows directly from Proposition applied to rank-(1+p) operators we((xr) 6o ww(x:) using the fact that W,, = E[wó(x:) @ wo(x.))
27

---

Proposition 15. Let KE)) hold T E [8, Given ó > 0, with probability in the i..dd draw of 1 n, it holds that
/// W W M.7 S (7) 21
1/21
(54)
p.7
where
(),
20, T 3n7
L' (r, 0) + V
2CTH- (1,,
ny -'
(55)
and
(7,,) In 4 In tr(W
[WW.WWyll
Moreover.
$ 1 - 14.7 1ir
2
(56)
Proof. The idea is apply Proposition 13 for operator E(x) 69 E(xr), where {(x) is defined in 443
Due to D. we have that ||||| & 1 S (47) ( cr. On the other hand, we have that
[g>=(E(x) @ E(=)]² = 15l1 E @ E(x)) = |1el W W W W-1/2 1/2 P.7
and, hence 54 follows. To complete the proof, observe that
1/ n.7 W n. W = ||
S (), 1/2 ( W W M47
and, hence for = (7, 8) smaller than one we obtain
/W1/2W-1 1 Mi7 = [I(W 1/2 1/2-1/2) 1
1 - || W
1/-1/22 W p., W 1/2
hold for some a E (0, 2), (0, 1] and E [8, 1]. Given 5 > O, with probability in the it.d. draw of (x:) (z:)²= from n, it holds
Proposition 16. Let (RC) (SD) and (KE)
P W 1/2 (W, Ww)W--CClla S 2 (7, 2 -,
P.7
where
2 (2.7) = 4 2 ca G7
cap-B
C+
(11)) (rr-a)/2 (a-1)/2 ow
, S
8 nyB +
(57)
n²gT
a 2
Proof. As before, w.l.o.g. set μ 1. First, recall that HS (7) equipped with I1-lls is sep- arable Hilbert space. Hence, we will apply Proposition 12 for A {(x) @ v(r), where
(xx) CW (wo) H2+P, To that end, observe that E is - (x)) 1m [
and that
2
E[I(()) S m! (1-1/2 //////))))))))))))) (V((Watw.))
Recalling Lemma we have that
the task to bound [[)(()) 12/1--- E
Using (RC)
[|+())) Sc.
Σ Σ
# (w hy)
1E[1+p]
JEJ iE|t
But, since
lta of + (wid(z),
( (W(1+a)/2 w w; (x), h2)
1:7
expanding as in proof of Lemma I and using (KE)
we have that
(2ta-1) 2
[ Sca Σ
(1) h2))
o7 <
7 7-(--a)
; S
0f +
of
cw
(a-T)
,a 2T.
JEJ [1(p)))
28

---

Therefore, we can set
A = cale y (2)
a S T,
{ (-(8Tr-)) , , T,
and o² = CacFCB
(a-t)
g , 2 T.
a 2
OW
in Proposition 12 to obtain (58)
Proposition 17. Let (KE) hold for E [8, 1]. Given o > O, with probability in the iid. draw of (ara) 2 from it holds
W 1.7
C) S (=)17,))} 2 1
where
4 422ow (7,
2 capr
+ CTMTT
In V ny n2a7
(58)
11
Proof. First, note that since [ W₁ imply analogous assumptions for C. Hence, we can apply Proposition 13 from 29
H. I 10- [ - and C I / W (KE)
and (SD) for
μ.-
which using
(uC,,-)) the observation that ||Wm,// =/= (uC, the proof.
/ |,,,,, ' |
completes
Next, we develop concentration bounds of some key quantities used to build RRR empirical estimator.
E.3.2 Variance and Norm of KRR Estimator
Proposition 18. Lef (RC) (SD) and KK)) hold for some (0, 2), (0, 1] and [8. Given 0 i < 1, then with probability at least 1 in the i.i.d. draw of (xz, R= from p
2. (7,88/3) + 2, 8/3)
P [//l/2(0 Mrr GHr)l
21-6
1- (7,
Proof. Note that W4y (G n.7 Gp,)
W/l/?(G
W (W W W and, hence,
W1l2 (Gpm Ga,z) W 10
W WHn W C C
WW/WW WW//W---1/ p.7
1/-1/2 C) - p.
1-1/2 -,,))= (59)
p.7 pi?
k.7
Thus, taking the norm and using ||,|'T|| I o (Sx) with the Propositions 16 and 15 we prove the first bound.
E.3.3 Variance of Singular Values
In this section we prove concentration of singular values computed in Theorem I necessary step to derive learining rates for RRR estimator.
Proposition 19. Let (RC)
(SD) and (KE) hold for some a E (0, 3 E (0, 1] and T E [8, 1). Let
B Wa and B Wi-t C. Given 8 0 of e (7, 1, then with probability at least
1 - in the i.i.d. draw of (I:)) π
|o² Se
B*Bl
(60)
where
(e. ( 2 6/33
(61)
1- (7,3
Proof. We start from the Weyl's inequalities for the square of singular values
[ - of < B B - B| ie [n].
29

---

But, since,
B*B-B'B=
m
denoting M = (-
MS W w,, and (G - Gp.)), ). we have
B - B = B BBMMMMMM+MM... MC B' = B" + (M B" (G = Ga,))
=3" + M BN (M' NW 147
=
(G (Wp=WW)Gmn + (GH,)"
Therefore, recalling that, due to 59
= R- W 11
M: P.y (M NGH,,),, we conclude
B--BBBBGGAA (C- C
(C- -G Gwn
+ (M +M WW -N
(62)
Next, observe that
is bounded by Proposi-
tion 17
JM-NG tions [6] and 17
M. (C- 1/2 (W
W W W Cl| is bounded by Proposi-
1Gp WW,
(W,, W M.
O'l is bounded by Proposition
Therefore, using additionally Proposition 15 result follows.
Remark that bound singular values can rely the fact
- - oB V (B) (BB) OB OB
< -
o
E.3.4 Variance of RRR Estimator
Recalling the notation B = W and B W W-1//2, let denote P. and P, denote the orthogonal projector onto the subspace of leading r right singular vectors of B and B, respectively. Then we have [B]. = BP, and [ = -P,, and, hence G = Gu. and T = Ga, Pr.
Proposition 20. Let (RC) (SD)
(EE) hold for some a E (0, 2), (0, 1] and E [B, Given then with probability at least 1 8 in the i.i.d. draw of (#:)
8 > 0 and > 0, if (7, 0) <
from n.
+ + + (7,
(B)
1 PM Gr,, $
e (1,11/88))
(63)
1-2,(,00/3
of( - o (B)
Proof. Start by observing that |2 Gr S [ W G(ani) and
(1/ HG -G - (WW1//WW1111 W- -
( // (W, Wa)G,,,PP pr
-(P. P.)).
Using that the norm of orthogonal projector P. is bounded by one and applying Propositions and 16 together with Propositions and 19 completes the proof.
30

---

E.4
Proof of Theorem 2
First, observe that (DF) implies that
W
where the empirical covariance with the inexact Dirichlet coefficient &, is denoted by W9 Since RRR algorithm now uses W when e > 0 the only change to our proof technique lies in the analysis of variance in section E..3 In particular,
only need adapt Propositions 15 and 16 which is
straightforward. Indeed, we have
|F(W - G = |F( + W, N.--W)) FGI
< FI, W,) II FIllGG
S ||(WW, +
for F = -W and G being either W C or I. Thus, since |I Fl - and ||G]] is either or bounded by /CH/ u., in conclusion the relative error of imperfect knowledge simply appears as the additive term the final guarantees. Thus, in what follows, we present that 0.
E.4.1 Operator Norm Error Bounds
Summarising previous sections, in order to prove Theorem 2 just need analyse the bounds e ex and ey Not that we fix the hyperparameter >> O, which affects the constants, but need to chose the decay rate of Tikhonov regularization parameter y to obtain balancing of bias and variance in the generalization bounds.
Let us first assume regime T. which covers "well-specified learning" when non-trivial eigenfunc- tions of L are inside H. Since a 2 t and B S т, we have that for large enough n one has
112 7r/2
(n21//2
(7,
1n
and ()
In i 2, 3,4.
(64)
7//2
7/122
But, since the bias term is S 20/2 and the slow term from Proposition
20 is 1/ nyy0 we can set
Tn n .ोग and obtain
n -1/2
n n-1/2
in
aa/2
8/2
= n 2000)
and
7/2
=
-+
2(0-8),
which, due to a 2 T, implies
lim (7,, 5/3) lim (1hn 5/3) 12+00
lim c2(155/3) n-+00
lim (7,,8/3) 0.
n-+30
Therefore, for this choice of regularization parameter Equation 44 with Propositions 9 and 20
assure that
E 0r+ l E(G)))
Aitt
(65)
Now, let us consider the more difficult learn case a < t when eignefunctions of the generator have only weaker norms than the RKHS one. Then, for large enough bounds in (64) hold for i = 3, but
n 1 1/2
n
(7(
V
(2r--)//
1n
(66)
Then, by balancing the slow terms with bias //2 2 we set Yn the and obtain
a/22
1/2
and
n- क[ナक,
2
(T+3 23/3
r/2
which, since 7 > 3, implies
lim o) (7n, 5/3) lim 2 (Y/, 8/3) R++00
lim
6/3) lim (7,888 O.
no
n-
31

---

and we obtain
E( M47)
xr+1(B) - ((G M,7) X,+1 S n गान.
(67)
Therefore, denoting
n 0 S T:
2(++))
0 S T, Q TT
and E
(68)
0+3
as a consequence the operator norm error bound in Theorem holds, and We have the following result on the estimation of the singular values of (μI- (.)- 7pp denoted by X using the singular values of B, denoted by Ts from Theorem 1
Proposition 21. Let (RC) (SD) and (KE) hold for some (0,2) and B E (0, 1] and [3. 1). and define (68). Then, there exists a constant > O such that for every given E (0, 1), large enough and with probability at least 1 - in the i.i.d. draw of (x=)==1 from a for all i [r].
o x E In
(69)
Proof. The proof is direct consequence of Propositions 10 and using (64)-66))
E,5 Spectral Learning Rates
Finally, we conclude the proof of Theorem 2 showing the concentration of eigenpairs. Recalling Proposition[] we need to combine the operator norm bound and metric distortion. Since, as indicated in Proposition the population KRR estimator grow in the operator norm whenever the regularity condition is violated a < 1, leading to possibly unbounded metric distortions w.r.t. increasing sample size, we restrict the case a 2
First, combining
665 67 and Proposition 21 we have that
E(G S (ōr+1 x2+1) + ce, In8
On the other hand, from Propositions 40 and 14 have that with failure probability 6,
(7:)² (n()))) S (o),
which, if we can prove that î and (h.) are bounded, concludes the proof for empirical spectral biases. To that end, recall that, e.f. Proposition 5 for RRR estimator we have
Grall o, (Z
||G
|G + +|G Guuill
m(aa) x
S or or((uT - - Zp) E(G
< |Guzll + (177 -1/2 11/2 [ (Gurr Gp.nil
Mi7
+ n 11/2--121
1 + n 213
X E(G).))
Vx; - n 11 1/2
Vx n1
where in the last inequality we have applied Propositions 9 and 18 proof of Theorem 2 is concluded.
Thus, using Proposition 40 the
E.6 Discussion of the learning rates
We first discuss the learning rates reported in Table Notice that although the papers we compare to employ different risk, our comparison remains meaningful because our energy- risk measure provides upper bound on their risk measure. This ensures that any upper bound derived with our risk also applies to theirs.
The learning bound for obtained is 1 covers only pure diffusion processes (Laplacian with constant weights). Their learning rate is non-parametric and depends the state space dimension a counter-intuitive way O(n xeी ) in [ Theorem 3]. highlighting a potential limitation of their approach. In comparison, when specify our RRR method with RBF kernel (i.e. B = 0), we achieve much faster parametric learning rate (n
32

---

The recent work of 43 covers Langevin processes via kernel approach, but they derive a sub-optimal learning bound for IG of order ( in [4] Theorem 4.4]. For Langevin diffusions, our RRR method with an RBF kernel achieves a faster parametric rate O(n -1/2). Moreover, the computational complexity of their method is On"'d"). which limits its application in realistic molecular dynamics scenarios.
As for 19 Theorem 7], although they considered general diffusions, they only derived suboptimal bound for the variance component of their risk, with an explicit dependence on the dimension of the state space.
Additionally, note that the mentioned works lack learning guarantees for eigenfunctions and eigen- values. Notably, their methods are prone to the spurious eigenvalue phenomenon, requiring expert manual review of each eigenpair to select plausible ones.
Next, we contrast our IG learning and well established TO learning. First, note that TO methods apply only to equally spaced data and the sampling frequency 1/At must be high enough to distinguish all relevant time-scales. Otherwise, since TO eigenvalues are edh.rt, small spectral gaps complicate learning (see [28] Thm. 3D.. Conversely, our IG method, which uses gradient information, is time- scale independent, handles irregularly spaced measurements, and does not rely time discretizations. Indeed, recalling the risk functional in Eq. (8), we see that the "label" of the model X (x( R G ((x) is the action of the resolvent. Since this "label" is not computable, we "fight fire (resolvent) with fire (generator)", i.e. we use the energy norm of Eq. (9) to rewrite regularized problem (8). Crucially, this allows to obtain estimators via energy covariance W in (14) that completely captures infinitesimal nature of the learning problem without needing time-lagged observations. This contrasts with TO methods, where choosing the time-lag AL is the major bottleneck in real applications. Indeed, to obtain data from an invariant distribution π one uses the trajectory data after some burn-in time needed ensure that ergodic mean approximates well æ. Then, the problem is reduced to studying only the dependence as is done e.g. in [27) using standard tools of 3-mixing and the method of blocks. This allows one to obtain non-parametric learning bounds for TO methods where the effective sample size suffers from the multiplicative effect of the time-lag to achieve approximate independence. So, TO methods "waste" a lot of data negatively impacting statistical accuracy. Contrary to this, our method can be applied to data with larger time-lags (even irregularly spaced) so that effective sample size is close to the true one. All this results in better generalization, as shown in Fig. I g)-h)-i) where our IG estimator captures ground truth significantly better than TO for the same sample size, this generalization across all time scales incurs quadratic computational complexity w.r.t. state dimension and not statistical accuracy. Lastly, with our additional discussion on imperfect knowledge, our IG method can also be applied a fully data-driven regime as TO methods.
To conclude this section, we contrast our method with classical numerical methods for the spectral decomposition of differential operators, such as finite-elements (FEM), For FEM, the approximation
error is A Ax1 S AAR. where p is the polynomial degree used to construct finite elements, h is the mesh size, and depends on the eigenfunctions' smoothness. As the number of mesh elements grows exponentially with d (i.e., ~ h hd),, reducing /t is the major bottleneck mitigated by computationally demanding adaptive higher order methods. On the other hand, our IG method that requires less or no knowledge has a quadratic impact of the d only on the computational complexity. Indeed, sample complexity depends domain that can be much lower than d. Therefore, RRR for IG's resolvent doesn" 't suffer from the curse of dimensionality as the FEM does.
the effective dimension of the equilibrium distribution on the
F Experiments
All the experiments were performed on a workstation with 125.6 GiB of memory and AMD@ Ryzen threadripper pro 3975wx 32-cores x 64 processor, graphics card was used. The ver- sion of python used is Python 3.9.18. The choice of the kernel in all experiments was Gaus- sian RBF with specified length-scales, and the hyperparameters were chosen via vooss-validtiion. The code to reproduce the results of the experiments can be found in the following repository: hawpr//boooooeeeeeeeeeeietaeeeee
_ernell
Four well potential For this experiment, we used an in-house code to simulate the system. The equations of motions were discretized using the Euler-Maruyama scheme with timestep of 10-2.
33

---

0.
01-
2.0
-0.7
-200
+2.2
200
022
535
141
-100
+0.0
000
-1.03
-0.73
-2.2
5000
2.35
Figure 2: Results of the RRR given by our method for two different length scales (blue and red) compared with ground truth (black) for the Langevin dynamics driven by a four well one dimensional potential.
2nd genfunction
3rd
Empirical
50
10
10
0.-
C6-
05-
3.00
0
-.25
1.00
4.75
02-
3.00
-.22
Figure 3: Panels a)-c): Test of the model's robustness with respect to the hyperparameter μs, tested for 30 different values between 10-3 and 5, compared to the ground truth result. Panel d): logarithm of the empirical bias function of the kernel length scale and the logarithm of regularization parameter 7.
RRR was fitted using 1000 points, = 5 and 10-5 The length scales used were 0.05 and 0.5. This experiment was reproduced 100 times leading to very small change in the estimation of the eigenfunctions. In Figure we report the result of one of them. The reduced rank regression was performed with a rank of 5. Further, in Figure we show the robustness of eigenfunctions w.r.t. choice of shift hyper parameter μl, as well as values of empirical bias for different values of length-scale and regularization hyperparameters. Concerning panel e) of Figure where we show the consistency of our model with the true Boltzmann distribution. Namely, we the conditional probability density function (pdf) of the system being in one bin, given it started at some point and after 100000 steps (so that it relaxes towards the equilibrium distribution). We bin the space into 50 bins and approximate the pdf via Eq. (12), where h is characteristic function of a accurate estimation,
model to forecast
bin. In order to achieve zero. Additionally, we perform the same procedure with imperfect diffusion coefficient estimated from data (k,T = 0.45 = 0.01 instead of true value 0.5 a.u.) by looking at the variance of increments over 10 steps. Note that the result (green dotted lines) unchanged and compares well with the analytical Boltzmann distribution (black lines). Finally, report that if the
the knowledge that the leading eigenvalue is
approach is used with the method described in [19 43 no dynamical quantity can be forecasted due to many spurious eigenpairs (see Figure la) which prevents the system from relaxing towards the Boltzmann distribution.
Muller Brown For this experiment, used in-house code simulate the system. The equations of motions were discretized using the Huler-Maruyama scheme with a timestep of 10 1 and a temperature of (arbitrary units). RRR fitted using 2000 points, μ 1 and y 10-5, The length scale used was 0.6. The reduced rank regression was performed with rank of 5
CIR model For this experiment, we reproduced 100 times the simulations in order to obtain statistical uncertainties. We used length scale 0.5, f = and y 10-6 The reduced rank regression was performed with rank of 2.
34

---

Ground truth EF
Generator EF (this work)
Transter operator
2.0
15
15
10
10
10
0.5
0.5
0.5
0.0
0.0
1.0
-0.5
100
2.0
15-
15
15-
10
10
10-
0.5
0.5
0.5
-0.5
-0.5
10
Figure 4: Results of the RRR given by our method (second column) compared to ground truth (first column) and transfer operator RRR (last column). Points are colored according to the value of the eigenfunction
US mortgage rates We have trained our method on real 30-year US mortgage rates dataset and contrasted with the fitted CIR model using continuous ranked probability scores that estimated from the forecasts obtained by of each of them. Each model has been trained using data from January 2009 December 2016. The initial condition was the last week of December 2016 and the predictions were made for the years 2017 and 2018. Since the dataset is real, we used the imperfect partial knowledge, that is, for our method, we estimated the diffusion coefficient only via least squares calibration of a CIR model over the training set. This allows more flexibility on the drift term in our model.
35

---

NeurIPS Paper Checklist
1. Claims
Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Question:
Answer: [Yes]
Justification:
we also provided summary table comparing our results to key related works.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No or
NA answer to this question will not perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results
be expected to generalize other settings.
is fine to include aspirational goals as motivation as long is clear that these goals are not attained by the paper.
2. Limitations
Question: Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes],
Justification: see the conclusion section.
Guidelines:
The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those not discussed in the paper.
The authors are encouraged to create separate "Limitations" section in their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification. asymptotic approximations only holding locally). The authors should reflect how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution low or images are taken in low lighting. Or speech-to-text system might be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and faimess.
While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't 't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory Assumptions and Proofs
Question:
For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
Answer: [Yes]
36

---

Justification:
assumptions are stated in Section 5 and proofs presented in the appendix.
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated referenced in the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but they appear in the supplemental material, the authors are encouraged to provide short proof sketch provide intuition.
Inversely, any informal proof provided in the core of the paper should complemented by formal proofs provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental Result Reproducibility
Question: Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that affects the main claims and/or conclusions of the paper (regardless whether the code and data are provided or not)?
Answer: [Yes]
Justification: the experimental setting is reported in Section 6 and expanded in the appendix.
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, a No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is novel architecture, describing the architecture fully might suffice, or if the contribution is a specific model and empirical evaluation, it may be necessary to either make possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g.. in the case of large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
While NeurIPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) the contribution is primarily a new algorithm, the paper should make clear how to reproduce that algorithm.
(b) the contribution is primarily a new model architecture, the paper should describe the architccture clearly and fully.
(c) If the contribution is new model (e.g., a large language model), then there should either way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome to describe the particular way they provide for reprodueibility. the case of closed-source models, may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers have some path to reproducing or verifying the results.
5. Open access to data and code
Question: Does the paper provide open access to the data and code, with sufficient instruc- tions to faithfully reproduce the main experimental results, as described in supplemental material?
37

---

Answer: [Yes]
Justification:
See App.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurIPS code and data submission guidelines (https://nips.cr
for more details.
While we encourage the release of code and data, we understand that this might not be possible, so "No" an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurIPS code and data submission guidelines (https
/ //ippe.cc/pullic/guiii wwwpc//ggue////eeeeceooc CoodeSubeiiiee ssi wsionnolli onnolicy)))) for more details.
The authors should provide instructions on data access and preparation, including how access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. If only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information possible in supplemental material (appended to the paper) recommended, but including URLs to data and code permitted.
6. Experimental Setting/Details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: please see App.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core of the paper to level of detail that is necessary to appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, or as supplemental material.
7. Experiment Statistical Significance
Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Question:
Answer: [Yes].
Justification: whenever appropriate we have reported these.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by error bars, confi- dence intervals, statistical significance tests, at least for the experiments that support main claims of the paper.
The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, overall run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call a library function, bootstrap, etc.)
The assumptions made should be given (e.g.. Normally distributed errors).
should be clear whether the error bar is the standard deviation the standard error of the mean.
38

---

It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of enrors is not verified.
For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out range (e.g. negative error rates).
If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text.
8. Experiments Compute Resources
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?
Question:
Answer: [Yes]
Justification: please App.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs well estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g.. preliminary or failed experiments that didn't make into the paper).
9. Code Of Ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https: /ate:://neuiissp http:p/ssccpppp.t.ttttttttttaeaa
Answer: [Yes]
Justification: We enjoy research and we respect other people work.
Guidelines:
The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.
If the authors answer No, they should explain the special cireumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., if there is a special consid- eration due to laws or regulations their jurisdiction).
10. Broader Impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [NA]
Justification: for other theoretical papers particular concerns are present.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
If the authors answer NA or No. they should explain why their work has no societal impact or why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that improvement in the quality of generative models could be used to
39

---

generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology.
there are negative socictal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification: please the previous answer.
Guidelines:
The answer NA means that the paper poses no such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards challenging, and many papers do require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Are the creators original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Question:
Answer: [NA]
Justification: we use synthetic data generated by ourselves.
Guidelines:
The answer NA means that the paper does not use existing assets.
The authors should cite the original paper that produced the code package dataset.
The authors should state which version of the asset is used and, if possible, include URL.
The name of the license (e.g., CC-BY 4.0) should be included for each asset.
For scraped data from a particular source (c.g., website), the copyright and terms of service of that source should be provided.
If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, wwppppswipteeeee www.reeswthhhe....c. com/ datasets has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that are re-packaged, both the original license and the license of the derived asset (if has changed) should be provided.
this information is not available online, the authors are encouraged to reach out to the asset's is creators.
13. New Assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
40

---

Answer: [NA]
Justification: see above.
Guidelines:
The answer NA means that the paper does release new assets.
Researchers should communicate the details the dataset/codeemoodll as part their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset used.
At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsourcing and Research with Human Subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well as details about compensation (if any)?
Answer: [NA]
Justification:
applicable.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Including this information the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then much detail as possible should be included the main paper.
According to the NeurIPS Code of Ethics, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional Review Board (IRB) Approvals or Equivalent for Research with Human Subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or equivalent approval/review based on the requirements of your country or institution) were obtained?
Answer: [NA]
Justification: applicable.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Depending on the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and expect authors to adhere the NeurIPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review.
41


Reviews:
[{"review_id": "5L6CJcbkrI", "reviewer": "Reviewer_nFzp", "summary": "This paper discusses learning the generator of stochastic diffusion processes in reproducing kernel Hilbert space. In particular,\nsection 2: background on the generator, Dirichlet form, energy, learning in RKHS, empirical risk in the Hilbert-Schmidt norm\nsection 3: introduce an energy-based risk functional for resolvent estimation which leads to spectral estimation\nsection 4: empirical risk minimization with Tikhonov regularization and rank constraints.", "questions": "1. In the paper, 'full knowledge' refers to given drift and diffusion coefficients, and 'partial knowledge' refers to the given diffusion coefficient and energy/Dirichlet operator. In practice, is it necessarily easier to obtain diffusion coefficient and energy/Dirichlet operator than to obtain drift and diffusion coefficients? \nIn a fully data-driven scenario, one can first estimate drift and diffusion coefficients, and then use physics-informed method with 'full knowledge' obtained. Does this make the 'full knowledge' methods cited in the paper more flexible than the 'partial knowledge' method proposed?\n\n2. As transfer operator A_t and the infinitesimal generator L satisfy the relation A_t = e^{Lt}, is there a consistency between the learned TO from cited TO methods with the learned IG from the proposed method?\n\n3. What are the kernels used for the 3 examples in experiments? How sensitive is the experimental result to the kernel selection, regularization, and rank parameter?\n\n4. Can the authors provide the learned dynamics of the 3 examples compared to the true dynamics, the invariant measure, and time scales?", "limitations": "Limitations have been discussed.", "rating": 5, "confidence": 3, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "The paper is well-written and has a thorough discussion on the background and the comparison of previous work. It proposes a novel energy risk functional and a reduced-rank estimator with dimension-free learning bounds. It also establishes spectral learning bounds for generator learning.", "weaknesses": "The weaknesses are mainly in the motivation (please see question 1 below) and the experiment (see questions 3 and 4 below). The authors did a great job in the theoretical comparison with previous work on TO and IG learning, however, the empirical comparison focuses on the spuriousness of eigenvalues and the recovery capability of the metastable states. To sufficiently demonstrate the practical advantage of the proposed method, other comparisons of the learned dynamics, phase transition, and invariant measure should also be included.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We appreciate the reviewer's insightful evaluation and valuable comments. In what follows, we aim to address the highlighted weaknesses and respond to the reviewer's questions.\n\n## Weaknesses:\n\nThank you for pointing this out. Indeed, clearer motivation beyond the Dirichlet form should help the reader understand our approach better. While some aspects of this are addressed in the general reply, in the following we try to address your concerns in detail.\n\n## Questions:\n\n 1. As we have discussed in the general reply, a large class of processes can be written in Dirichlet form and, importantly, for many of them statistical estimation of the Dirichlet coefficient (related to the diffusion coefficient) from data is much easier than recovering the drift term. This means that our method is possible to apply even when the challenging drift term is not estimated. That said, we are, to our knowledge, __the only existing method__ that is able to learn the infinitesimal generator model from __a single trajectory__ that is __theoretically consistent with the true dynamics__. Finally, as we show in the general reply, note that our method is robust under imperfectly estimated partial knowledge. This is not only proven, but also empirically demonstrated in Fig. 1 and 2 of the attached pdf. \n\n2. Indeed, there is the consistency of our IG method and the properly executed TO methods. This is theoretically backed by Thm. 2 in [22] for TO and our Thm. 2 for IG. However, note that TO methods apply only to equally spaced data and the sampling frequency  $1/\\Delta t$ must be high enough to distinguish all relevant time-scales. Otherwise, since TO eigenvalues are $e^{\\lambda_i \\Delta t}$, small spectral gaps complicate learning (see Thm. 3 [22]). Conversely, our IG method, which uses gradient information, is time-scale independent, handles irregularly spaced measurements, and does not rely on time discretizations. This results in better generalization, as shown in Fig. 1 c)-d)-e) where our IG estimator captures true phase transitions between meta-stable states significantly better than TO for the same sample size. Lastly, with our additional discussion on imperfect knowledge, our IG method can be safely applied in a fully data-driven regime as TO methods. \n\n3. Thank you for spotting this. Indeed, in App. F we didn’t specify that the RBF Gaussian kernel with specified length-scales was used in all cases. As reported in the generall reply, the robustness of kernel hyperparameters and Tikhonov regularization is essentially the same as for TO kernel estimators, while the IG estimator is very robust against the new shift hyperparameter, see Fig. 3 of the attached pdf. \n\n4. Thank you for this question.  Besides original comparisons of the time-scales (IG eigenvalues) and meta-stable states (IG eigenfunctions), providing the additional experiments, we demonstrate our method’ consistency with true dynamics from the perspective of equilibrium distribution recovery and  long term  forecasting of conditional cumulative density functions in Fig. 1 and 2 of the attached pdf, respectively. This might come as a surprise, but after a closer look, the cited works,  only showcase empirical evaluation of the leading eigenfunctions, which due to the reported effect of spuriousness in Fig, 1a requires picking proper eigenpairs from an abundance of wrongly estimated ones. Indeed, if one builds such models without expert knowledge needed to remove spurious estimation, the obtained dynamics is inconsistent with the true one, see caption of Fig. 1 of the attached pdf.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "Gz6SHNj8Xh", "reviewer": "Reviewer_6wyC", "summary": "In this paper, the authors consider the problem of learning the infinitesimal generator a Stochastic Diffusion Process (SDP). Compared to existing approaches such as [1] they tackle the unbounded nature of the generator by introducing a novel statistical framework which is based on the Dirichlet form associated with the SDP. In this framework, the authors estimate the resolvent of the generator which can be approximated with finite-rank operators. They consider a regularized and rank-truncated version of the regression loss. Theorem 1 provides a way to compute the eigenvalues and eigenvectors of the estimated operator. Spectral learning bounds are derived in Theorem 2. The theoretical results of the paper are completed with three examples: learning the overdamped Langevin generator for a one-dimensional four well potential, the overdamped Langevin generator for a Muller Brown potential and finally for a Cox-Ingersoll-Ross process (CIR process).\n\n[1] Cabannes et al. (2023) -- The Galerkin method beats Graph-Based Approaches for Spectral Algorithms", "questions": "See weaknesses", "limitations": "Limitations are addressed in Section 7 (\"Conclusion\")", "rating": 6, "confidence": 2, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "* Even though the paper is mathematically heavy and contains a lot of notation, I think it is well-written. I appreciate that (almost) all the assumptions needed before stating Theorem 2 are clearly laid out and explained. I also appreciate the rigour shown in the paper. \n\n* The obtained results regarding the spectral bounds are interesting and provide a strong theoretical grounding for the method.\n\n* From a methodological point of view the time complexity can be reduced compared to [1,2,3] if the rank is low. I think this is one of the strength of the method, although I would have appreciated more details on the choice of the hyperparameters and their interactions (see below).\n\n[1] Cabannes et al. (2023) -- The Galerkin method beats Graph-Based Approaches for Spectral Algorithms\n\n[2] Hou et al. (2023) -- Sparse learning of dynamical systems in RKHS: An operator-theoretic approach\n\n[3] Pillaud-Vivien et al. (2023) -- Kernelized Diffusion Maps", "weaknesses": "Before expressing my main concerns regarding this paper, I want to emphasize that I'm not an expert in this domain and therefore my understanding of the main competitors and methods used in the paper is lacking. \n\n* My main concern is regarding the experiments. I know this is a theoretical paper but a new statistical framework and a new learning procedure should be clearly validated against competing methods like [1,2,3]. Even though some qualitative conclusions are drawn I would have liked to see a more extensive study (for different choices of hyperparameters for both the target Langevin diffusion and the method introduced by the authors). \n\n* Indeed, there are little details in the paper about how to choose the crucial parameters of the algorithm like $\\mu$, $r$ and $\\gamma$. Establishing the robustness of the procedure regarding these hyperparameters seem crucial to validate the methodology.\n\n* As mentioned before, even though I appreciate the rigour of the paper it is notation heavy. Having to constantly rely on the (massive) table of page 13 to remember the notation hindered my understanding of the paper. Is there a way to either reduce the notational load or to incorporate a lightweight version of Table 2 in the main paper? This would greatly help the reader.\n\n* All assumptions are commented except (KE). What are the main limitations of this assumption? I am not familiar with the domain so a light explanation in the main text would also be appreciated.\n\n* Is the assumption of l.188 that there exists a Dirichlet operator realistic? It is true for Langevin operator and CIR (which already covers quite a lot of ground) but apart from these models? Does this impose anything on the diffusion?\n\nSome minor remarks:\n\n* l.178 \"Importantly, this energy functional can be empirically estimated from data sampled from π, whenever full knowledge, that is drift and diffusion coefficients of the SDE (1), or partial knowledge,\" -- This sentence is not really clear to me. How can one leverage the drift and diffusion coefficient here?\n\n* l.215 typo (spurious new line) \n\n[1] Cabannes et al. (2023) -- The Galerkin method beats Graph-Based Approaches for Spectral Algorithms\n\n[2] Hou et al. (2023) -- Sparse learning of dynamical systems in RKHS: An operator-theoretic approach\n\n[3] Pillaud-Vivien et al. (2023) -- Kernelized Diffusion Maps", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We appreciate the reviewer's insightful evaluation and valuable comments. In what follows, we aim to address the highlighted weaknesses and respond to the reviewer's questions.\n\n## Weaknesses:\n\n__Main contributions and broader impact.__ In the general reply we tried our best to provide more details on positioning of our main contributions in the larger context of solving SDEs, and help the reviewer evaluate the possible impact of the presented results. \n\n__Experiments.__ While we would like to stress out that we originally provided quantitative comparison to the cited methods ([1] is in fact just a special case of [2], while [3] is not applicable to Langevin dynamics) in Fig. 1 a) which implies the failure of this method to capture the true dynamics,  we follow the reviewer’s suggestion and provide additional experiments that make this point more clear for the reader. \n\n__Choice of hyper-parameters.__ As suggested, we provide the  empirical evaluation of the robustness of hyper-parameter choice, as indicated in the general answer. In summary, the robustness of the choice of kernel’s hyper-parameters and Tikhonov regularization is similar to general kernel methods for (vector-valued) regression. On the other hand when estimating IG via its resolvent the choice of additional shift hyper-parameter $\\mu$, as expected from the theory, has a much smaller impact on the performance, see Fig. 3 of the attached pdf.\nFinally, we wish to stress out that we do go beyond classical validation via prediction by introducing novel quantities, the empirical spectral biases $\\widehat s_i$ of eq. (24), which are (theoretically and empirically) good metrics to use for fine-tuning the model. \n\n__Notations.__ Thank you for your suggestion. We will prepare a smaller summary table that only introduces necessary notations for presenting the IG learning framework and the method, while omitting extensive notation related to proving the generalization bounds. We believe this is feasible, since, if the manuscript is accepted, the additional page can be used to present this table in the main body of the paper.\n\n__(KE) assumption.__ The kernel embedding (KE) assumption originates from the study of mini-max optimal excess risk bounds for classical kernel regression [12]. Whenever the kernel is bounded, (SD) and (KE) do _not impose any additional constraints_, they are just used to _describe the interplay between the data distribution and the chosen RKHS_. This is a formal way to quantify the impact of the kernel choice to the learning rate. We will additionally clarify this in the revised manuscript. \n\n__Dirichlet form.__ To complement the discussion on this issue from the general reply, we would like to stress out that many diffusion processes have an infinitesimal generator that can be expressed in the form (IG). In addition to the Langevin and the CIR processes already studied, we can mention: \n-  the Wright-Fisher diffusion (in dimension one), which can be defined in the context of population genetics and can be adapted to model interest rates, see [G], \n-  the geometric Brownian motion which models the price process of a financial asset, \n- the multi-dimensional Brownian motion (a=0) that corresponds to the heat equation, \n- the transport processes associated with advection-diffusion equation, see [H]\n- the process associated with Poisson's equation in electrostatics, see [I]. \n\n__Minor remarks.__ Thank you for pointing out the typo, we will fix it. Concerning the claim in line 178, it relates to the first equality in Equation (7). There, one can see that knowing drift and diffusion terms, we can evaluate $[L f] (x) $ and approximate the energy as empirical mean of $ \\mu |f (x) |^2 + f(x)[L f] (x) $ via samples $x$ from the invariant distribution.\n\n\n## Additional references:\n[A] Fukushima, Masatoshi, Yoichi Oshima, and Masayoshi Takeda. Dirichlet forms and symmetric Markov processes. Vol. 19. Walter de Gruyter, 2011\n\n[B] Jean Jacod. Discretization of processes. Springer, 2012 \n\n[C] Danielle Florens-Zmirou. Approximate discrete-time schemes for statistics of diffusion processes. Statistics: A Journal of Theoretical and Applied Statistics 20.4, 1989\n\n[D] Fan, J., & Zhang, C. (2003). A reexamination of diffusion estimators with applications to financial model validation. Journal of the American Statistical Association, 98(461), 118-134.\n\n[E] Yuri Kutoyants. Statistical inference for ergodic diffusion processes. Springer Science & Business Media, 2013.\nYuri Kutoyants. Parameter estimation for stochastic processes, 1984.\n\n[F] Fabienne Comte,Valentine Genon-Catalot Nonparametric drift estimation for iid paths of stochastic differential equations. The Annals of Statistics, 48(6), 2020\n\n[G] Martin Grothaus, Max Sauerbrey. Dirichlet form analysis of the Jacobi process. Stochastic Processes and their Applications, 157, 2023.\n\n[H] Antoine Lejay,  Lionel Lenôtre, Géraldine Pichot. Analytic expressions of the solutions of advection-diffusion problems in one dimension with discontinuous coefficients. SIAM Journal on Applied Mathematics, 79(5), 2019.\n\n[I] Sethu Hareesh Kolluru. Preliminary Investigations of a Stochastic Method to solve Electrostatic and Electrodynamic Problems. Phd Thesis, 2008.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "qSIAQq4tQS", "reviewer": "Reviewer_PUzn", "summary": "The paper considers a time-homogeneous Stochastic Differential Equation (SDE) with known diffusion part and known or unknown drift. The problem is to find properties of this equation from known data, in particular, to find a (low-rank) representation of its infinitesimal generator (IG). For this purpose, the resolvent of this operator is considered and the reduced-rank estimator in reproducing kernel Hilbert spaces (RKHS) together with the energy-based risk metric are used. Accuracy estimates for the found approximation to IG and time-complexity are given. Several model numerical experiments are conducted for proof-of-concept.", "questions": "- Can your method be extended to the case when neither drift nor diffusion coefficients are known?", "limitations": "Limitation are described in the paper. The paper is mostly theoretical.", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "- Detailed appendix with necessary reference material\n\n- Development of SDE theory with exact estimates on the result obtained", "weaknesses": "- The paper uses quite complex concepts and may be difficult to understand in a 9 page format.\n\n- There are not enough practical examples of the application of these results to machine learning\n\n\nMinor \n\nThe x-axis in Figure 1 (a) is not signed\n\nL120 comma is not needed\n\nL124 space after dot is missing\n\nL206: colon missing at the end of the line\n\nL215 A new line  is redundant (`this work we focus on the case when`)\n\nL275 space after dot is missing\n\nL915 \"anan\" -> \"and\"", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We appreciate the reviewer's insightful evaluation and valuable comments. In what follows, we aim to address the highlighted weaknesses and respond to the reviewer's questions.\n\n## Weaknesses:\n\n### Major:\n\n- Indeed, the topic of learning IG of a stochastic process with kernel-based methods, and, in particular, development of nonparametric learning rates contains several complex concepts. More so, since the paper is the first one to formalize a consistent way to learn IG from data. We tried our best to present the sophisticated theory and mitigate the complexities, believing that this work paves the path to efficient and reliable algorithms for important applications in ML (see also general reply).\n- In addition to our general reply, we would like to emphasise here that, from molecular dynamics to weather models, the ability to reliably learn eigenfunctions of differential operators that generate the dynamics is of capital importance for building trustworthy and physics informed AI. Hence, we strongly believe this work is of __high relevance__, if not for the most general ML community, definitely __for the AI-for-Science community__.  With that respect, we would like to further stress that our method is the first of its kind in machine learning to have error bounds that theoretically demonstrate the superiority of an ML approach for spectral estimation of differential operators over classical numerical methods like FEM in high-dimensional settings. While empirical evidence supporting this has been well documented in practice, particularly in fields such as molecular dynamics [42], to the best of our knowledge, this had never been theoretically proven until now. Following reviewer's suggestions, in the revised manuscript we will better emphasize possible applications and further improve readability and accessibility of our results.\n\n\n\n### Minor:\n\n- Thank you for spotting the typos, we will correct them in the revision. \n\n\n## Questions:\n\nAs we elaborate in the general reply, indeed we can show, both theoretically and empirically (Figures 1 and 2 in the attached pdf), that our method works in a fully data-driven regime. While empirically this might not come as a big surprise, we believe that __adding a theoretical proof__ motivated by the comments of reviewers is a very nice __additional contribution__. If the reviewer wishes to have more details on this aspect, we are happy to elaborate.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "Xe4yyRZMgr", "reviewer": "Reviewer_D2P3", "summary": "This paper proposes a relevant and sound approach for learning self-adjoint SDE generators via operator learning techniques.\n\nThe paper includes a compactification, a novel prior knowledge inclusion, and first-of-its-kind statistical learning guarantees that extend the known ones from discrete Markov processes.", "questions": "1. *Are Dirichlet forms commonly deduced from the diffusion part of the SDE?* It would help assess the ease of knowing these forms beforehand for readers unfamiliar with the notion.\n\n2. *Is it reasonable to think you could deal with a bias of the Dirichlet form using imperfect diffusion knowledge?*\n\n3. *Could you highlight what practical gains/losses are introduced w.r.t. discrete-time setting?* I would guess you need less data due to diffusion priors and not requiring as small of a discretization. On the other hand, IG regression complexity seems higher for high-dimensional systems.\n\n4. *Can a sample complexity gain be recognized using prior knowledge in this form for your theoretical results?* Could the performance gain w.r.t. TO learning be recognized using the derived bounds?\n\n5. *How arbitrary is the choice of $\\mu$?*\n\n6. *How does the time-sampled data come into play, and does the sampling rate of the data influence any aspect of the approach?* Does sampling the invariant distribution + Dirichlet form allow you to avoid time-derivative or \"1-step\" dynamics observations (e.g., compared to [A])? Based on the current writing, it is not immediately clear to me (but I might have missed something).\n\n7. *How do you compare (using full knowledge) to classical numerical methods using known drift and diffusion? Would your sample efficiency be better than classical FEM for a given precision?* You mention how, for known operators (via drift and diffusion knowledge), there would be no spuriousness in estimating eigenpairs. This is perhaps of separate interest.\n\n\n[A] - Meng, Y., Zhou, R., Ornik, M., & Liu, J. (2024). Koopman-Based Learning of Infinitesimal Generators without Operator Logarithm. http://arxiv.org/abs/2403.15688", "limitations": "1. *Possibly strong prior knowledge requirements.*\n2. *Limited to self-adjoint generators.* \n3. *Experiments are on 1D toy examples.*", "rating": 8, "confidence": 4, "soundness": 4, "presentation": 4, "contribution": 3, "strengths": "1. The paper is easy to read, with suitable mediation throughout.\n\n2. Inserting known diffusion effects into the Dirichlet forms and using the resolvent, provides an elegant way of including prior knowledge and well-posed estimation.\n\n3. These are the most complete statistical learning guarantees for spectra of self-adjoint generators (albeit reliant on partial knowledge).", "weaknesses": "1. Unclear extensions to entirely data-driven regimes (no partial knowledge) and requiring sampling from an invariant distribution.\n\n2. Experiments only consider toy examples that validate theory but do not showcase the practicality of the approach for more general cases.\n\n3. I am not sure that contribution **4)** is well supported by the current presentation, primarily due to hard-to-parse figures. A revised comparison and legends with sharper plots could make the results much easier to understand.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We appreciate the reviewer's insightful evaluation and valuable comments. In what follows, we aim to address the highlighted weaknesses and respond to the reviewer's questions.\n\n## Weaknesses:\n1. Thank you for emphasizing the realistic setting of imperfect partial knowledge. This motivated us to __show theoretically and empirically__ that indeed our IG learning method with estimated Dirichlet form __is guaranteed to work__. Concerning the sampling, our assumption to have i.i.d. data from the invariant distribution $\\pi$ can be relaxed as we reply in Q6. \n2. Concerning the experiments, please see our general reply.  \n3. We agree that more discussion of Fig. 1 is needed to enhance the understanding of the experiment section. If accepted, an extra page will allow us to improve it, as well as to  include additional content from the attached PDF and this rebuttal in the revision.\n\n## Questions:\n- Answers to both __Q1__ and __Q2__ are positive. Please see our general reply for the discussion. \n- __Q3:__ TO methods apply only to equally spaced data and the sampling frequency $1/\\Delta t$ must be high enough to distinguish all relevant time-scales. Otherwise, since TO eigenvalues are $e^{\\lambda_i \\Delta t}$, small spectral gaps complicate learning (see Thm. 3 [22]). Conversely, our IG method, which uses gradient information, is time-scale independent, handles irregularly spaced measurements, and does not rely on time discretizations (see reply to Q6). While it results in better generalization, as shown in Fig. 1 c)-d)-e) where our IG estimator captures ground truth significantly better than TO for the same sample size, this generalization across all time scales incurs quadratic computational complexity w.r.t. state dimension and not statistical accuracy (see also reply to Q7). Lastly, with our additional discussion on imperfect knowledge, our IG method can be safely applied in a fully data-driven regime.\n- __Q4:__ For i.i.d. data from $\\pi$, at this point, it is hard to answer this question, while for trajectory data, important improvement can be expected, see reply to Q6. While the standard assumptions (SD) and (KE) are the same for TO and IG error bounds, regularity assumptions in these two settings are not easily comparable. This motivates the development of non-standard regularity assumptions for TO and IG learning that could expose the benefit of the prior knowledge. This is an interesting problem in its own right, and it will be the subject of our future work. \n- __Q5:__  As illustrated in Fig. 3 in the attached pdf, estimators are quite robust w.r.t. value of the shift parameter $\\mu$. Theoretically, it is important that $\\mu$ is not too small.\n- __Q6:__ Thank you for raising this question, we will emphasize these important aspects of our method in the revision. \n  - Recalling the risk functional in eq. (8), we see that the “label” of the model $\\chi_{\\mu}(x) \\approx G^*\\phi(x)$ is the action of the resolvent. Since this “label” is not computable, we “__fight fire__ (resolvent) __with fire__ (generator)”, i.e. we use the energy norm of eq. (9) to rewrite regularized problem (8) as in line 222. Crucially, this allows us to obtain estimators via energy covariance $W_{\\mu}$ in (14) that __completely captures infinitesimal nature__ of the learning problem without needing time-lagged observations. This contrasts with TO methods, where choosing the time-lag $\\Delta t$ is the major bottleneck in real applications, [18,42]. \n  - Thank you for this reference (Meng et al. 2024)  that we will cite in the revision. Note that this method is exactly the Galerkin projection of the Yoshida approximation of IG on a finite-dimensional RKHS (dictionary of functions). It essentially corresponds to solving (8) where $\\phi$ is finite-dimensional, expectation is used instead of energy and $\\chi_{\\mu}$ is replaced with $\\mu^2\\chi_{\\mu}(x)-\\mu\\phi(x)$. To do this, authors estimate resolvent via equation in line 167 which requires approximating an operator integral via __equally spaced data with high sampling frequency (small $\\Delta t$)__, essentially suffering from the same bottleneck as TO methods. Since in (Meng et al. 2024) the provided learning theory and the empirical evaluation is limited, guarantees and limitations of this method are not fully clear to us.  \n  - In practical applications to obtain data from an invariant distribution $\\pi$ one uses the trajectory data after some burn-in time needed to ensure that ergodic mean approximates well $\\pi$. Then, the problem is reduced to studying only the dependence as is done e.g. in [21] using standard tools of $\\beta$-mixing and the method of blocks. This allows one to obtain non-parametric learning bounds for TO methods, c.f.  [22], where the effective sample size suffers from the multiplicative effect of the time-lag to achieve approximate independence. Recalling Q3, TO methods “waste” a lot of data negatively impacting statistical accuracy (we expect similar limitations for the method of  (Meng et al. 2024)). Contrary to this, our method can be applied to data with larger time-lags (even irregularly spaced) so that effective sample size is close to the true one.\n- __Q7:__ Compared to FEM, the key difference lies in the approximation error $ |\\lambda_k - \\widehat{\\lambda}_k| \\leq c h^{2p} |\\lambda_k|$, where $p$ is the polynomial degree used to construct finite elements, $h$ is the mesh size, and $c$ depends on the eigenfunctions' smoothness. As the number of mesh elements grows exponentially with $d$,  i.e., $\\sim h^{-d}$, reducing $h$ is the major bottleneck mitigated by computationally demanding adaptive higher order methods. On the other hand, our IG method that requires less or no knowledge has a quadratic impact of the $d$ only on the computational complexity. Indeed, sample complexity depends on the effective dimension of the equilibrium distribution on the domain that can be much lower than $d$.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:
We wish to thank all reviewers for their insightful evaluation of our paper. We appreciate all their comments and remarks, which we will incorporate in our revision. Before addressing each review in detail, we would like to point out some general remarks that apply to all of them.  

## __Assumptions on the process and availability of partial knowledge.__
The existence of Dirichlet forms for diffusion processes satisfying a stochastic differential equation like (1) is not too restrictive an assumption, as it encompasses many processes beyond Langevin and CIR such as advection-diffusion, Wright-Fisher diffusion, multidimensional Brownian motion, etc. (see reply to reviewer __6wyC__ for details and the list of additional references, there ref. [A] gives in-depth review of this topic). For the reviewers convenience, we briefly elaborate. Recalling the form of IG in eq. (3), the Dirichlet form in eq. (5) for a self-adjoint $L$ exists whenever the positive definite diffusion part $b b^\top$ satisfies uniform ellipticity conditions and the drift term $a$ allows integration by parts, leading to $s(x)=b(x)/\sqrt{2}$.  Thus, to obtain the partial knowledge we need, it is enough to estimate the diffusion function. When the sample paths are observed continuously, the diffusion coefficient $b$ can be directly identified from these observations, making  it a non-statistical problem. For discrete realizations of the process, __the diffusion coefficient can be estimated non-parametrically using various methods__. These methods include pathwise estimation by computing the variance of the increments over small intervals [B], kernel-based methods [C], local polynomial regression [D]. On the other hand, __the estimation of the drift is a much more demanding task__. Different methods are reviewed in [E] and references therein. A more recent approach [F] drawing inspiration from particle systems, consists in constructing estimates from $m$ i.i.d. paths of the solution process.

 ## __Imperfect partial knowledge and fully data-driven method.__ 
Inspired by the reviewers' questions, we show how to make our method fully data-driven. This theoretical discussion will be presented as a remark in the main body of the revision with detailed proofs in the appendix. In practice, whenever classical modeling from first principles is not feasible, we can first estimate the Dirichlet form and then apply our method. Hence, for a fully data-driven method, we need to analyze the impact of the imperfect knowledge. To this end, assume that we apply our model with some $\tilde s$ that incurs relative estimation error $\tilde \varepsilon>0$ of the form 
$$
\mathbf{(RE)} \quad \tilde \varepsilon s(x)^\top s(x) \preceq  \tilde s(x)^\top \tilde s(x) - s(x)^\top s(x) \preceq  \tilde \varepsilon s(x)^\top s(x),  
$$
where $\preceq$ is the standard Loewner ordering. Note that whenever $s$ is constant (e.g. Langevin) or linear (e.g. CIR), the estimation can be done from a single trajectory by estimating the variance of increments, and then __(RE)__ reduces to the standard relative error. So, our algorithm now uses empirical covariance $\widetilde W$ built from $\tilde s$, and, hence, the only change to our proof technique lies in the analysis of variance in App. E.3. In particular, we only need to adapt Prop. 15 and 16, which is straightforward since (RE) implies $\tilde\varepsilon \widehat W \preceq \widetilde W-\widehat W \preceq  \tilde\varepsilon \widehat W$. Indeed, we have 
$$
\Vert F(W-\widetilde W)FG \Vert = \Vert F(W\pm\widehat W-\widetilde W)F G \Vert \leq \Vert F(W-\widehat W)FG \Vert  + \tilde\varepsilon \Vert F(\widehat W \pm W)F\Vert \Vert G \Vert \leq \Vert F(W-\widehat W)FG \Vert + \tilde\varepsilon (\Vert F(W-\widehat W)F\Vert +1)\Vert G \Vert,
$$
for $F= W_{\mu,\gamma}^{-1/2}$ and $G$ being either $W_{\mu,\gamma}^{-1/2}C$ or $I$. So, in conclusion, __the relative error $\widetilde\varepsilon$ of imperfect knowledge simply appears as the additive term in the final guarantees.__ Namely, (23) and (24) hold upon replacing $\varepsilon_n^\star \ln \delta^{-1}$ with $\varepsilon_n^\star \ln \delta^{-1} + \tilde \varepsilon$. 

## __Experiments.__ 
Given the theoretically challenging nature of the problems addressed in this paper, we have intentionally limited the role of experiments to effectively highlight the main theoretical contributions. While we include two additional experiments in the attached pdf and test robustness on hyper-parameters choice, significantly extending empirical studies provided in the most related works on IG  [1, 15, 36], our focus remains consistent. We believe that developing realistic applications requires extensive work and contextual analysis, which  would constitute contribution in its own right.

## __Contributions.__ 
To conclude, let us clarify our contributions:
 1. We propose a __fundamentally new idea to estimate the spectrum of self-adjoint generators of stable Ito SDE from a single trajectory__. In contrast to all existing works, we exploit the geometry of the process via a novel energy (risk) functional. In a certain sense we “fight fire (resolvent ) with fire (generator) ” to 
2.  derive __a new efficient learning method__ that is able to learn the best approximation of the resolvent of IG on the RKHS __independently of the time-sampling__, (see reply Q6 to __D2P3__); 
3. We derive __the first IG spectral estimation finite sample bounds__ using (imperfect) partial knowledge, which __notably, overcome the curse of dimensionality present in classical numerical methods__ (see reply to __D2P3__, item Q7 for a quantitative discussion); 
4. __each important aspect of our learning method__, especially in relation to the most relevant existing works, __is empirically demonstrated__ to complement our theoretical analysis.  

We once again thank all the reviewers whose comments inspired the discussion above that significantly strengthens our paper and better demonstrates its broader impact.

Meta-review: The reviewers agree that this is an interesting paper that addresses a different and important problem from a typical NeurIPS submissions. Counterbalancing its strengths as a theoretically grounded estimation framework, the reviewers list some concerns about motivation, hyperparameter tuning, experimental evaluation and presentation. I think the authors have addressed the first two well enough (including additional results). I think there is some truth to the latter two points, and while I personally though the  presentation in the main body was fine, the over submission crams a lot of theory material into the supplementary material, and I do not think the reviewers have or can be expected to go carefully over it. As such, I think this paper sits just at the borderline.
Decision: Accept (poster)