# Learning the Infinitesimal Generator
of Stochastic Diffusion Processes

Vladimir R. Kostic
Istituto Italiano di Tecnologia
University of Novi Sad
<EMAIL_ADDRESS>
Karim Lounici
CMAP-Ecole Polytechnique
<EMAIL_ADDRESS>
Helene Halconruy
Telecom Sud-Paris
<EMAIL_ADDRESS>
Timothee Devergne
Istituto Italiano di Tecnologia
<EMAIL_ADDRESS>
Massimiliano Pontil
Istituto Italiano di Tecnologia
University College London
<EMAIL_ADDRESS>
(May 2024)

###### 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.

##  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 are subject to unfeasible-to-
model influences, necessitating the extension of deterministic models through
stochastic differential equations (SDEs), see [33, 35] 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 [34]. In
atomistic simulations, where SDEs are used to model the evolution of atomic
systems subjected to thermal fluctuations via the Boltzmann distribution [30].

A diverse range of SDEs can be represented as
$dX_{t}=a(X_{t})dt+b(X_{t})dW_{t}$, where $X_{0}=x$. Here, $W$ denotes a
(possibly multi-dimensional) Brownian motion, and the functions $a$ and $b$
are commonly 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 statistics"
since the 1980s, as seen in works like [18, 24, 25]. 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. 42, and references therein].
Consequently, there has been a 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 21, and references
therein]. TOs represent the average evolution of state functions (observables)
over time and, being linear and amenable to spectral decomposition under
certain conditions, they offer a valuable means to interpret and analyze
nonlinear systems. However, they require evenly sampled data at a 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 is 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 Hilbert-Schmidt operators with comprehensive
statistical theory [22]. Unfortunately, the existing statistical theory
collapses when applied to unbounded operators, prompting the need to
completely rethink the problem

Related work  Extensive research has explored learning dynamical systems from
data [see the monographs 7, 26, and references therein]. Analytical models for
dynamics are often unavailable, motivating the need for data-driven methods.
Two prominent approaches have emerged: deep learning-based methods [5, 11,
29], effective for learning complex data representations but lacking
statistical analysis, and kernel methods [2, 6, 9, 16, 17, 20, 21, 22, 43],
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 to the development of several methodologies [17, 27, 31, 39,
45, 44], some of which are based on deep canonical correlation analysis [3,
23]. In comparison, there have been significantly fewer works on learning the
IG and only in very specific settings. In [48] a deep learning approach is
developed for Langevin diffusion, while [19] 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 [1], revising the Galerkin method for Laplacian diffusion,
[15], presenting a kernel approach for general diffusion SDE with full
knowledge of drift and diffusion coefficients, and [36] addressing Langevin
diffusion. As highlighted in Appendix E.6, 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) Underpinning a
framework for learning the generator of stochastic diffusion processes in
reproducing kernel Hilbert space, using a novel energy risk functional; 2)
Proposing a reduced-rank estimator with dimension-free learning bounds; 3)
Establishing first-of-its-kind spectral learning bounds for generator
learning; 4) Demonstrating the practical advantages of our approach over
previous methods.

##  2 Background and the problem

Our drive to estimate the eigenvalues of the infinitesimal generator for an
SDE like (1) stems from its crucial role in characterizing dynamics in
physical systems. This operator’s closed form (3), relies on the drift $a$ and
diffusion coefficient $b$, where we have partial knowledge: $b$ is assumed to
be known, but $a$ is not. To compensate for the lack of prior knowledge about
$a$, we introduce the system’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-Ingersoll-Ross processes. For
detailed list of notation used throughout the paper we refer Table 1 in the
appendix.

Transfer operator and infinitesimal generator  A variety of physical,
biological, and financial systems evolve through stochastic processes
$X=(X_{t})_{t\in\mathbb{R}^{+}}$, where
$X_{t}\in\mathcal{X}\subset\mathbb{R}^{d}$ denotes the system’s state at time
$t$. A commonly employed model for such dynamics is captured by stochastic
differential equations (SDEs) of the form

$dX_{t}=a(X_{t})dt+b(X_{t})dW_{t}\quad\text{and}\quad X_{0}=x,$ (1)

where $x\in\mathcal{X}$, $W=(W_{t}^{1},\dots,W_{t}^{p})_{t\in\mathbb{R}_{+}}$
is a $\mathbb{R}^{p}$-dimensional ($p\in\mathbb{N}$) standard Brownian motion,
the drift $a:\mathcal{X}\to\mathbb{R}^{d}$, and the diffusion
$b:\mathcal{X}\to\mathbb{R}^{d\times p}$ are assumed to be globally Lipschitz
and sub-linear, so that the SDE (1) admits an unique solution
$X=(X_{t})_{\geqslant 0}$ with values in
$(\mathcal{X},\mathcal{B}(\mathcal{X}))$. Processes akin to equations like (1)
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 a continuous-time Markov process with almost surely continuous
sample paths whose dynamics is described by a family of probability densities
$(p_{t})_{t\in\mathbb{R}_{+}}$ and transfer operators
$(A_{t})_{t\in\mathbb{R}_{+}}$ such that for all $t\in\mathbb{R}_{+}$,
$E\in\mathcal{B}(\mathcal{X})$, $x\in\mathcal{X}$ and measurable function
$f:\mathcal{X}\to\mathbb{R}$,

$\mathbb{P}(X_{t}\in
E|X_{0}=x)=\int_{E}p_{t}(x,y)dy\;\,\text{and}\;\,A_{t}f=\int_{\mathcal{X}}f(y)p_{t}(\cdot,y)dy=\mathbb{E}\big{[}f(X_{t})\,|\,X_{0}=\cdot\big{]}.$
(2)

Evaluating $A_{t}f$ at $x$ yields the expectation of $f$ starting from $x$ and
evolving until time $t$, making the transfer operator crucial for
understanding $X$ dynamics. The family $(A_{t})_{t\in\mathbb{R}_{+}}$
satisfies the fundamental semigroup equation $A_{t+s}=A_{t}\circ A_{s}$, for
$s,t\in\mathbb{R}_{+}$. Here, we focus on the transfer operator’s effect on
the set $\mathcal{L}^{2}_{\pi}(\mathcal{X})$, a choice driven by the existence
of an invariant measure $\pi$ for $A_{t}$ on
$(\mathcal{X},\mathcal{B}(\mathcal{X}))$ which satisfies $A_{t}^{*}\pi=\pi$
for all $t\in\mathbb{R}_{+}$. Then, the process $X$ is characterized by the
infinitesimal generator $L$ defined for every
$f\in\mathcal{L}^{2}_{\pi}(\mathcal{X})$ such that the limit
$Lf=\lim_{t\rightarrow 0^{+}}(A_{t}f-f)/t$ exists in
$\mathcal{L}^{2}_{\pi}(\mathcal{X})$. The operator $L$ is closed on its domain
$\mathrm{dom}(L)$ which is equal to the Sobolev space

$\mathcal{W}^{1,2}_{\pi}(\mathcal{X})=\\{f\in\mathcal{L}^{2}_{\pi}(\mathcal{X})\;|\;\lVert
f\rVert_{\mathcal{W}}^{2}=\lVert f\rVert_{\mathcal{L}^{2}_{\pi}}+\lVert\nabla
f\rVert_{\mathcal{L}^{2}_{\pi}}<\infty\\}.$

For SDE dynamics of the form (1), we can prove (see A.2 for details) that $L$
is the second-order differential operator given, for any
$f\in\mathcal{L}^{2}_{\pi}(\mathcal{X})$, $x\in\mathcal{X}$, by

$Lf(x)=\nabla
f(x)^{\top}a(x)+\frac{1}{2}\mathrm{Tr}\big{[}b(x)^{\top}(\nabla^{2}f(x))b(x)\big{]},$
(3)

where $\nabla^{2}f=(\partial_{ij}^{2}f)_{i\in[d],j\in[p]}$ denotes the Hessian
matrix of $f$.

Spectral decomposition  Knowing only the drift $a$ and diffusion $b$ is not
enough to compute (2) or to understand quantitative aspects of dynamical phase
transitions, such as 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 can turn to an 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 $\mu\in\mathbb{C}$, the operator $L_{\mu}=(\mu I-L)^{-1}$ is the
resolvent of $L$ and the corresponding resolvent set is defined by

$\rho(L)=\big{\\{}\mu\in\mathbb{C}\,|\,(\mu I-L)\,\text{is bijective
and}\,L_{\mu}\,\text{is continuous}\big{\\}}$

We assume $L$ has a compact resolvent, meaning $\rho(L)\neq\emptyset$ and
there exists $\mu_{0}\in\rho(L)$ such that $(\mu_{0}I-L)^{-1}$ is compact.
Under this assumption, and given that $(L,\mathrm{dom}(L))$ is self-adjoint,
we can prove the spectral decomposition of the generator (see A.1 for details)
as follows:

$L=\sum_{i\in\mathbb{N}}\lambda_{i}\,f_{i}\otimes f_{i}$ (4)

where $(\lambda_{i}){i\in\mathbb{N}}$ are the eigenvalues of $L$, and the
corresponding eigenfunctions $f_{i}\in\mathcal{L}^{2}_{\pi}(\mathcal{X})$,
forming an orthonormal basis $(f_{i})_{i\in\mathbb{N}}$, are also
eigenfunctions of the transfer operator $A_{t}$.

Dirichlet forms and energy  To handle the initial lack of knowledge about the
drift $a$, we assume to have access to another quantity, called the energy,
defined as $\mathfrak{E}(f)=\lim_{t\to
0}\int_{\mathcal{X}}(f(f-A_{t}f))/td\pi$ for all functions
$f\in\mathcal{L}^{2}_{\pi}(\mathcal{X})$ for which this limit exists, defining
in the way the domain $\mathrm{dom}(\mathfrak{E})$. The associated Dirichlet
form is the bilinear form defined by polarization for any
$f,g\in\mathrm{dom}(\mathfrak{E})$ by

$\mathfrak{E}(f,g)=-\int_{\mathcal{X}}f(Lg)d\pi=\int_{\mathcal{X}}(-Lf)gd\pi.$
(5)

For every $f\in\mathrm{dom}(\mathfrak{E})$, we have
$\mathfrak{E}(f)=\mathfrak{E}(f,f)$. As for every $i\in\mathbb{N}$,
$0\leq\mathfrak{E}(f_{i},f_{i})=-\int_{\mathcal{X}}f_{i}(\lambda_{i}f_{i})d\pi=-\lambda_{i}$,
we check that the eigenvalues of $L$ are negative. To relate $L$ to Dirichlet
form, we assume there exists a Dirichlet operator $B=s^{\top}\nabla$ where
$s\\!=\\![s_{1}\,|\ldots|\,s_{p}]:x\\!\in\\!\mathcal{X}\mapsto
s(x)\\!=\\![s_{1}(x)\,|\ldots|\,s_{p}(x)]\\!\in\\!\mathbb{R}^{d\times p}$ is a
smooth function s.t. $Lf=s(s^{\top}\nabla f)=s(Bf)$ and so that

$\int_{\mathcal{X}}(-Lf)gd\pi=\int_{\mathcal{X}}(s(x)s(x)^{\top}\nabla
f(x))^{\top}\nabla
g(x)\pi(dx)=\int_{\mathcal{X}}(Bf(x))^{\top}(Bg(x))\pi(dx).$

We get that for any $f\in\mathrm{dom},(\mathfrak{E})$

$\mathfrak{E}(f)=\int_{\mathcal{X}}\lVert s^{\top}\nabla
f\rVert^{2}d\pi=\mathbb{E}_{x\sim\pi}[\lVert Bf(x)\rVert^{2}],$ (6)

which is reminiscent of the expected value of the kinetic energy in quantum
mechanics [13]. To illustrate this, we specify the generator and the
associated Dirichlet gradient form for two processes: Overdamped Langevin and
Cox-Ingersoll-Ross.

######  Example 1 (Langevin).

Let $k_{b},T\in\mathbb{R}_{+}^{*}$.The _overdamped Langevin equation_ driven
by a _potential_ $V:\mathbb{R}^{d}\rightarrow\mathbb{R}$ is given by
$dX_{t}=-\nabla V(X_{t})dt+\sqrt{2(k_{b}T)}dW_{t}\quad\text{and}\quad
X_{0}=x$, where $k_{b}$ and $T$ respectively represent the coefficient of
friction and the temperature of the system. Its infinitesimal generator $L$ is
defined by $Lf=-\nabla V^{\top}\nabla f+(k_{b}T)\Delta f$, for
$f\in\mathcal{W}^{1,2}_{\pi}(\mathcal{X})$. Since
$\int(-Lf)g\,d\pi=-\int\Big{[}\nabla\Big{(}(k_{b}T)\nabla
f(x)\frac{e^{-(k_{b}T)^{-1}V(x)}}{Z}\Big{)}\Big{]}g(x)dx=(k_{b}T)\int\nabla
f^{\top}\nabla g\,d\pi=\int f(-Lg)\,d\pi$, generator $L$ is self-adjoint and
associated to a gradient Dirichlet form with
$s(x)=(k_{b}T)^{1/2}(\delta_{ij})_{i\in[d],j\in[p]}$.

######  Example 2 (Cox-Ingersoll-Ross process).

Let $d=1$, $a,b\in\mathbb{R}$, $\sigma\in\mathbb{R}_{+}^{*}$. The _Cox-
Ingersoll-Ross_ process is solution of the SDE
$dX_{t}=(a+bX_{t})dt+\sigma\sqrt{X_{t}}dW_{t}\quad\text{and}\quad X_{0}=x$.
Its infinitesimal generator $L$ is defined for
$f\in\mathcal{L}^{2}_{\pi}(\mathcal{X})$ by $Lf=(a+bx)\nabla
f+\frac{\sigma^{2}x}{2}\Delta f$. By integration by parts, we can check that
the generator $L$ satisfies
$\int(-Lf)g\,d\pi=\int\frac{\sigma^{2}x}{2}f^{\prime}(x)g^{\prime}(x)\,\pi(dx)=\int
f(-Lg)\,d\pi$, and it is associated to a gradient Dirichlet form with
$s(x)=\sigma\sqrt{x}/\sqrt{2}$.

Learning in reproducing kernel Hilbert spaces (RKHSs)  Throughout the paper we
let $\mathcal{H}$ be an RKHS and let
$k:\mathcal{X}\times\mathcal{X}\to\mathbb{R}$ be the associated kernel
function. We let $\phi:\mathcal{X}\to\mathcal{H}$ be a feature map [38] such
that $k(x,x^{\prime})={\langle\phi(x),\phi(x^{\prime})\rangle}$ for all
$x,x^{\prime}\in\mathcal{X}$. We consider RKHSs satisfying
$\mathcal{H}\subset\mathcal{L}^{2}_{\pi}(\mathcal{X})$ [38, Chapter 4.3], so
that one can approximate
$L:\mathcal{L}^{2}_{\pi}(\mathcal{X})\to\mathcal{L}^{2}_{\pi}(\mathcal{X})$
with an operator $G:\mathcal{H}\to\mathcal{H}$. Notice that despite
$\mathcal{H}\subset\mathcal{L}^{2}_{\pi}(\mathcal{X})$, the two spaces have
different metric structures, that is for all $f,g\in\mathcal{H}$, one in
general has ${\langle f,g\rangle}_{\mathcal{H}}\neq{\langle
f,g\rangle}_{\mathcal{L}^{2}_{\pi}}$. In order to handle this ambiguity, we
introduce the injection operator
$S_{\pi}:\mathcal{H}\to\mathcal{L}^{2}_{\pi}(\mathcal{X})$ such that for all
$f\in\mathcal{H}$, the object $S_{\pi}f$ is the element of
$\mathcal{L}^{2}_{\pi}(\mathcal{X})$ which is pointwise equal to
$f\in\mathcal{H}$, but endowed with the appropriate $\mathcal{L}^{2}_{\pi}$
norm.

Then, the infinitesimal generator restricted to $\mathcal{H}$ is simply
$LS_{\pi}$ which can be estimated by $S_{\pi}G$ for some
$G\in{\rm{HS}}\left(\mathcal{H}\right)$. This approach is based on the
embedding $\ell\phi\colon\mathcal{X}\to\mathcal{H}$ of the generator in the
RKHS that can be defined for kernels
$k\in\mathcal{C}^{2}(\mathcal{X}\,{\times}\,\mathcal{X})$ whenever one knows
drift and diffusion coefficients, see B, so that the reproducing property
${\langle\ell\phi(x),h\rangle}_{\mathcal{H}}\,{=}\,[LS_{\pi}h](x)$ holds.
Based on this observation [15] developed empirical estimators of $LS_{\pi}$
that essentially minimize the risk $\lVert
LS_{\pi}-S_{\pi}G\rVert^{2}_{{\rm{HS}}\left(\mathcal{H},\mathcal{L}^{2}_{\pi}\right)}=\mathbb{E}_{x\sim\pi}\lVert\ell\phi(x)-G^{*}\phi(x)\rVert^{2}_{\mathcal{H}}$.
In scenarios where drift and diffusion coefficients are not known, then
$\ell\phi$ becomes non-computable. However if the process has the Dirichlet
form (6), one can still empirically estimate the Galerkin projection
$(S_{\pi}^{*}S_{\pi})^{\dagger}S_{\pi}^{*}LS_{\pi}$ onto $\mathcal{H}$, as
considered in [1], which in fact minimizes the same risk. Yet this approach is
problematic due to the unbounded nature of the generator and the associated
estimators typically suffer from a large number of spurious eigenvalues around
zero, making the estimation of physically most relevant eigenfunctions
unreliable even for self-adjoint generators. Conversely, classical numerical
methods can compute the leading part of a 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.

##  3 Novel statistical learning framework

In this section, we tackle the challenges in developing suitable generator
estimators highlighted earlier. To this end, we introduce a risk metric for
resolvent estimation that can be efficiently minimized empirically, leading to
good spectral estimation. Since $L$ and $(\mu I\,{-}\,L)^{-1}$ 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.

First, given $\mu>0$, in order to approximate the action of the resolvent on
the RKHS by some operator $G\colon\mathcal{H}\to\mathcal{H}$, we introduce its
embedding $\chi_{\mu}\colon\mathcal{X}\to\mathcal{H}$ via the reproducing
property ${\langle\chi_{\mu}(x),h\rangle}_{\mathcal{H}}=[(\mu
I\,{-}\,L)^{-1}S_{\pi}h](x)$, formally given by
$\chi_{\mu}(x)=\int_{0}^{\infty}\mathbb{E}[\phi(X_{t})e^{-\mu
t}dt\,|\,X_{0}=x]$, see B for details. Using this notation, we aim to estimate
$[(\mu I\,{-}\,L)^{-1}S_{\pi}h](x)\approx[Gh](x)$, $h\in\mathcal{H}$, i.e. the
objective is to estimate $\chi_{\mu}(x)\approx G^{*}\phi(x)$ $\pi$-a.e.

An obvious metric for the risk would be the mean square error (MSE) w.r.t.
distribution $\pi$ of the data. However, this becomes intractable since in
general $\chi_{\mu}$ 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

$\mathcal{W}^{\mu}_{\pi}(\mathcal{X})=\\{f\in\operatorname*{\text{\rm
dom}}(L)\;|\;\lVert
f\rVert_{\mathcal{W}^{\mu}_{\pi}(\mathcal{X})}^{2}\equiv\mathfrak{E}_{\mu}[f]={\langle
f,(\mu I{-}L)f\rangle}_{\mathcal{L}^{2}_{\pi}}<\infty\\},$

where the norm now balances the energy of an observable
$f\colon\mathcal{X}\to\mathbb{R}$ w.r.t. the invariant distribution $\lVert
f\rVert_{\mathcal{L}^{2}_{\pi}}^{2}$ and its energy w.r.t. the transient
dynamics $-{\langle f,Lf\rangle}_{\mathcal{L}^{2}_{\pi}}$. Indeed

$\mathfrak{E}_{\mu}[f]=\mathbb{E}_{x\sim\pi}[\mu\lvert
f(x)\rvert^{2}-f(x)[Lf](x)]=\mathbb{E}_{x\sim\pi}[\mu\,\lvert
f(x)\rvert^{2}+\lVert s(x)^{\top}\nabla f(x)\rVert^{2}],$ (7)

where the last equality holds for Dirichlet gradient form (6), in which case
$\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$ is simply a weighted Sobolev space.
Importantly, this energy functional can be empirically estimated from data
sampled from $\pi$, whenever full knowledge, that is drift and diffusion
coefficients of the SDE (1), or partial knowledge, i.e. the diffusion
coefficient and Dirichlet operator $B$ in (6), is at hand. With that in mind,
instead of the standard MSE, we introduce the energy-based risk functional as

$\min_{G\colon\mathcal{H}\to\mathcal{H}}\mathcal{R}(G)=\mathfrak{E}_{\mu}[\lVert\chi_{\mu}(\cdot)-G^{*}\phi(\cdot)\rVert_{\mathcal{H}}].$
(8)

Denoting by $Z_{\mu}\colon\mathcal{H}\to\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$
the canonical injection, (8) can be equivalently written as

$\mathcal{R}(G)=\lVert(\mu
I\,{-}\,L)^{-1}Z_{\mu}-Z_{\mu}G\rVert_{{\rm{HS}}\left(\mathcal{H},\mathcal{W}\right)}^{2}=\lVert(\mu
I\\!-\\!L)^{-1/2}S_{\pi}-(\mu
I\\!-\\!L)^{1/2}S_{\pi}G\rVert_{{\rm{HS}}\left(\mathcal{H},\mathcal{L}^{2}_{\pi}\right)}^{2},$
(9)

where we abbreviated $\mathcal{W}=\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$ and
used $Z_{\mu}^{*}=S_{\pi}^{*}(\mu I{-}L)$, recalling that Hilbert-Schmidt and
spectral norms for operators $A\colon\mathcal{H}\to\mathcal{W}$, are $\lVert
A\rVert_{{\rm{HS}}\left(\mathcal{H},\mathcal{W}\right)}{=}\sqrt{\operatorname{\text{\rm
tr}}(A^{*}(\mu I{-}L)A)}$ and $\lVert
A\rVert_{\mathcal{H}\to\mathcal{W}}{=}\sqrt{\lambda_{1}(A^{*}(\mu
I{-}L)A)}\geq\mu^{-1/2}\|A\|_{\mathcal{H}{\to}\mathcal{L}^{2}_{\pi}}$.

Therefore, (9) implies that the regularized energy norm, while dominating the
classical $\mathcal{L}^{2}_{\pi}$ one, exerts a balancing effect on the
estimation of the resolvent. This leads us to the first general result
regarding the well-posedness of this framework.

######  Proposition 1.

Given $\mu>0$, let
$\mathcal{H}\\!\subseteq\\!\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$ be the RKHS
associated to kernel $k\in\mathcal{C}^{2}(\mathcal{X}\\!\times\\!\mathcal{X})$
such that
$Z_{\mu}\in{\rm{HS}}\left(\mathcal{H},\mathcal{W}^{\mu}_{\pi}(\mathcal{X})\right)$,
and let $P_{\mathcal{H}}$ be the orthogonal projector onto the closure of
$\operatorname*{\text{\rm
Im}}(Z_{\mu})\subseteq\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$. Then for every
$\varepsilon>0$ there exists a finite rank operator
$G\colon\mathcal{H}\,{\to}\,\mathcal{H}$ such that
$\mathcal{R}(G)\,{\leq}\,\lVert(I\,{-}\,P_{\mathcal{H}})(\mu
I\,{-}\,L)^{-1}Z_{\mu}\rVert^{2}_{{\rm{HS}}\left(\mathcal{H},\mathcal{W}\right)}\,{+}\,\varepsilon$.
Consequently, when $k$ is universal, $\mathcal{R}(G)\leq\varepsilon$.

The previous proposition reveals that whenever the hypothetical domain is
dense in the true domain and the injection operator is Hilbert-Schmidt, there
is no irreducible risk and one can find arbitrarily good finite rank
approximations of the generator’s resolvent. Note that
$Z_{\mu}\in{\rm{HS}}\left(\mathcal{H},\mathcal{W}^{\mu}_{\pi}(\mathcal{X})\right)$
is equivalent to $Z_{\mu}^{*}Z_{\mu}\\!=\\!S_{\pi}^{*}(\mu I\\!-\\!L)S_{\pi}$
being a trace class operator, which is assured for our Examples 1 and 2, see
the discussion in E.

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

$\mathcal{E}(G)=\lVert(\mu
I\,{-}\,L)^{-1}Z_{\mu}-Z_{\mu}G\rVert_{\mathcal{H}\to\mathcal{W}},\,G\in{\rm{HS}}\left(\mathcal{H}\right),\;\text{
and }\;\eta({h})=\lVert
h\rVert_{\mathcal{H}}/\lVert{h}\rVert_{\mathcal{W}},\,{h}\in\mathcal{H}.$ (10)

######  Proposition 2.

Let
$\widehat{G}=\sum_{i\in[r]}(\mu\\!-\\!\widehat{\lambda}_{i})^{-1}\,\widehat{h}_{i}\otimes\widehat{g}_{i}$
be the spectral decomposition of $\widehat{G}\colon\mathcal{H}\to\mathcal{H}$,
where $\widehat{\lambda}_{i}\geq\widehat{\lambda}_{i+1}$ and let
$\widehat{f}_{i}=S_{\pi}\widehat{h}_{i}\,/\,\lVert
S_{\pi}\widehat{h}_{i}\rVert_{\mathcal{L}^{2}_{\pi}}$, for $i\in[r]$. Then for
every $\mu>0$ and $i\in[r]$

$\frac{\lvert\lambda_{i}-\widehat{\lambda}_{i}\rvert}{|\mu-\lambda_{i}||\mu-\widehat{\lambda}_{i}|}\leq\mathcal{\mathcal{E}}(\widehat{G})\eta(\widehat{h}_{i})\quad\text{
and
}\quad\lVert\widehat{f}_{i}-f_{i}\rVert^{2}_{\mathcal{L}^{2}_{\pi}}\leq\frac{2\,\mathcal{E}(\widehat{G})\eta(\widehat{h}_{i})}{\mu\,[\text{\rm
gap}_{i}-\mathcal{E}(\widehat{G})\eta(\widehat{h}_{i})]_{+}},$ (11)

where $\text{\rm gap}_{i}$ is the difference between $i$-th and $(i+1)$-th
eigenvalue of $(\mu I\,{-}\,L)^{-1}$.

Note that the estimation of the eigenfunctions is first obtained in the norm
with respect to the energy space, and then transformed to the
$\mathcal{L}^{2}_{\pi}$-norm, as
$\|f\|_{\mathcal{W}}\,{\geq}\,\sqrt{\mu}\|f\|_{\mathcal{L}^{2}_{\pi}}$,
$f\,{\in}\,W_{\mu,\gamma}$. Therefore, by controlling the operator norm error
and metric distortion (refer to C), we can guarantee accurate spectral
estimation. Consequently, this allows us to approximately solve the SDE (1)
starting from an initial condition.

$\mathbb{E}[h(X_{t})\,|\,X_{0}\\!=\\!x]\\!=\\![e^{Lt}S_{\pi}h](x)\approx\textstyle{\sum_{i\in[r]}}\,e^{\widehat{\lambda}_{i}t}\,{\langle\widehat{g}_{i},h\rangle}_{\mathcal{H}}\widehat{h}_{i}(x).$
(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.
$\mathcal{L}^{2}_{\pi}$ by
$C=S_{\pi}^{*}S_{\pi}={\mathbb{E}_{x\sim\pi}}[\phi(x)\otimes\phi(x)]$, the
cross-covariance operator $T=S_{\pi}^{*}LS_{\pi}$, capturing correlations
between input and the outputs of the generator, and the covariance operator
$W_{\mu}=Z_{\mu}^{*}Z_{\mu}=S_{\pi}^{*}(\mu I{-}L)S_{\pi}$ w.r.t. energy space
$\mathcal{W}$. All operators can be estimated from data, depending on the
available prior knowledge. In this work we focus on the case when

Dirichlet gradient form is known. Then we can define the embedding of the
Dirichlet operator
$B=s^{\top}\nabla\colon\mathcal{L}^{2}_{\pi}(\mathcal{X})\to[\mathcal{L}^{2}_{\pi}(\mathcal{X})]^{p}$
into RKHS $d\phi\colon\mathcal{X}\to\mathcal{H}^{p}$ via the reproducing
property as ${\langle
d\phi(x),h\rangle}_{\mathcal{H}}=[BS_{\pi}h](x)=s(x)^{\top}Dh(x)\in\mathbb{R}^{p}$,
$h\in\mathcal{H}$. More precisely, we have that $d\phi(x)$ is a
$p$-dimensional vector with components $d_{k}\phi(x)$, where
$d_{k}\phi\colon\mathcal{X}\to\mathcal{H}$ is given via reproducing property
${\langle d_{k}\phi(x),h\rangle}_{\mathcal{H}}=s_{k}(x)^{\top}Dh(x)$,
$k\in[p]$. Hence, in this case we have that

$T=-\mathbb{E}_{x\sim\pi}[d\phi(x)\otimes
d\phi(x)]=-\textstyle{\sum_{k\in[p]}}\mathbb{E}_{x\sim\pi}[d_{k}\phi(x)\otimes
d_{k}\phi(x)].$ (13)

Moreover, defining $w\phi\colon\mathcal{X}\to\mathcal{H}^{p+1}$ by
$w\phi(x)=[\sqrt{\mu}\phi(x),d_{1}\phi(x),\ldots
d_{p}\phi(x)]^{\top}\in\mathbb{R}^{p+1}$, we get

$W_{\mu}=\mathbb{E}_{x\sim\pi}[\mu\phi(x)\otimes\phi(x)+d\phi(x)\otimes
d\phi(x)]=\mathbb{E}_{x\sim\pi}[w\phi(x)\otimes w\phi(x)].$ (14)

Regularized risk  Let us first introduce the regularized risk defined for some
$\gamma>0$ by

$\mathcal{R}_{\gamma}(G)=\lVert(\mu I\\!-\\!L)^{-1/2}S_{\pi}-(\mu
I\\!-\\!L)^{1/2}S_{\pi}G\rVert_{{\rm{HS}}\left(\mathcal{H},\mathcal{L}^{2}_{\pi}\right)}^{2}+\gamma\mu\lVert
G\rVert_{\mathrm{HS(\mathcal{H})}}^{2},\quad
G\in{\rm{HS}}\left(\mathcal{H}\right),$ (15)

which, after some algebra, can be written as

$\displaystyle\mathcal{R}_{\gamma}(G)$ $\displaystyle\,=\,$
$\displaystyle\operatorname{\text{\rm tr}}\big{[}G^{*}(\mu C{-}T{+}\mu\gamma
I)G{-}2CG{+}S_{\pi}^{*}(\mu I\,{-}\,L)^{-1}S_{\pi}\big{]}$
$\displaystyle\,=\,$ $\displaystyle\underbrace{\lVert
W_{\mu,\gamma}^{-1/2}C-W_{\mu,\gamma}^{1/2}G\rVert_{{\rm{HS}}\left(\mathcal{H}\right)}^{2}}_{\text{variance
term}}+\underbrace{\mathfrak{E}_{\mu}[\chi_{\mu}(\cdot)]-\lVert
W_{\mu,\gamma}^{-1/2}C\rVert_{{\rm{HS}}\left(\mathcal{H}\right)}^{2}}_{\text{bias
term}}$

where $W_{\mu,\gamma}=W_{\mu}+\mu\gamma I$ is the regularized covariance
w.r.t. $\mathcal{W}$.

Hence, assuming the access to the dataset $\mathcal{D}_{n}=(x_{i})_{i\in[n]}$
made of i.i.d. samples from the invariant distribution $\pi$, replacing the
regularized energy $\mathfrak{E}_{\mu}$ with its empirical estimate
$\widehat{\mathfrak{E}}_{\mu}$ leads to the regularized empirical risk
functional expressed as

$\widehat{\mathcal{R}}_{\gamma}(G)=\lVert\widehat{W}_{\mu,\gamma}^{-1/2}\widehat{C}-\widehat{W}_{\mu,\gamma}^{1/2}G\rVert_{{\rm{HS}}\left(\mathcal{H}\right)}^{2}+\widehat{\mathfrak{E}}_{\mu}[\chi_{\mu}(\cdot)]-\lVert\widehat{W}_{\mu,\gamma}^{-1/2}\widehat{C}\rVert_{{\rm{HS}}\left(\mathcal{H}\right)}^{2},$
(16)

where $W_{\mu,\gamma}$ and $C$ are estimated by their empirical counterparts
$W_{\mu,\gamma}$ and $C$, respectively, via (13).

Therefore, our regularized empirical risk minimization approach reduces to

$\min_{G}\lVert\widehat{W}_{\mu,\gamma}^{-1/2}\widehat{C}-\widehat{W}_{\mu,\gamma}^{1/2}G\rVert_{{\rm{HS}}\left(\mathcal{H}\right)}^{2},$
(17)

and we analyze two different estimators, the first one
$\widehat{G}_{\mu,\gamma}$ is obtained by minimizing regularized empirical
risk (17) over all $G\in{\rm{HS}}\left(\mathcal{H}\right)$, and, hence, the
name Kernel Ridge Regression (KRR) of the generators resolvent. The second one
$\widehat{G}^{r}_{\mu,\gamma}$ minimizes (17) subject to the (hard) constraint
that $G$ is at most of (a priori fixed) rank $r\in\mathbb{N}$ and, hence, is
named Reduced Rank Regression (RRR) of the generator’s resolvent. Notice that
when $r=n$, the two estimators coincide. After some algebra, one sees that
both minimization problems have closed form solutions

$\widehat{G}_{\mu,\gamma}=\widehat{W}_{\mu,\gamma}^{-1}\widehat{C}\quad\quad\text{
and
}\quad\quad\widehat{G}^{r}_{\mu,\gamma}=\widehat{W}_{\mu,\gamma}^{-1/2}[\\![\widehat{W}_{\mu,\gamma}^{-1/2}\widehat{C}]\\!]_{r},$
(18)

where $[\\![\cdot]\\!]_{r}$ denotes the $r$-truncated SVD of a compact
operator.

To conclude this section, we show how to compute the eigenvalue decomposition
of (18). To this end, we define the sampling operators
$\widehat{S}\colon\mathcal{H}\to\mathbb{R}^{n}$ and
$\widehat{Z}_{\mu}\colon\mathcal{H}\to\mathbb{R}^{(1+p)n}$ by

$\displaystyle(\widehat{S}h)_{i}=\tfrac{1}{\sqrt{n}}h(x_{i}),i\\!\in\\![n],$
and
$\displaystyle(\widehat{Z}_{\mu}h)_{kn+i}=\begin{cases}\tfrac{\sqrt{\mu}}{\sqrt{n}}\,h(x_{i}),&k=0,i\in[n],\\\
\tfrac{1}{\sqrt{n}}s_{k}(x_{i})^{\top}Dh(x_{i}),&k\\!\in\\![p],i\\!\in\\![n].\end{cases}$

Further, let
$\textsc{K}\\!=n^{-1}[k(x_{i},x_{j})]_{i,j\in[n]}\\!\in\\!\mathbb{R}^{n\times
n}$ be kernel Gram matrix, and introduce the Gram matrices
$\textsc{N}\\!\in\\!\mathbb{R}^{n\times pn}$ and
$\textsc{M}\\!\in\\!\mathbb{R}^{pn\times pn}$ whose elements, for
$k\in[1\\!+\\!p],i,j\in[n]$ are

$\textsc{N}_{i,(k-1)n+j}\\!=\\!n^{-1}{\langle\phi(x_{i}),d_{k}\phi(x_{j})\rangle}_{\mathcal{H}}\;\text{
and
}\;\textsc{M}_{(k\\!-\\!1)n\\!+\\!i,(\ell\\!-\\!1)n\\!+\\!j}\\!=\\!n^{-1}{\langle
d_{k}\phi(x_{i}),d_{\ell}\phi(x_{j})\rangle}_{\mathcal{H}}.$ (19)

We note that although we have introduced the above matrices via inner products
in $\mathcal{H}$, they can be readily computed via the kernel and its
gradients knowing the Dirichlet form, see D.

######  Theorem 1.

Given $\mu>0$ and $\gamma>0$, let
$\textsc{J}_{\mu,\gamma}\\!=\\!\textsc{K}\\!-\\!\textsc{N}(\textsc{M}\\!+\\!\gamma\mu
I)^{-1}\textsc{N}^{\top}\\!+\\!\gamma I$. Let
$(\widehat{\sigma}_{i}^{2},v_{i})_{i\in[r]}$ be the leading eigenpairs of the
following generalized eigenvalue problem

$\mu^{-1}(\textsc{J}_{\mu,\gamma}-\gamma
I)\textsc{K}v_{i}=\widehat{\sigma}_{i}^{2}\textsc{J}_{\mu,\gamma}v_{i},\quad
v_{i}^{\top}\textsc{K}v_{j}\\!=\\!\delta_{ij},\;i,j\in[r].$ (20)

Denoting
$\textsc{V}_{r}\\!=\\![v_{1}\,|\\!\ldots\\!|\,v_{r}]\\!\in\\!\mathbb{R}^{n\times
r}$ and $\Sigma_{r}\\!=\\!\operatorname{\text{\rm
diag}}(\widehat{\sigma}_{1},\ldots,\widehat{\sigma}_{r})$, if
$(\nu_{i},w^{\ell}_{i},w^{r}_{i})_{i\in[r]}$ are eigentriplets of matrix
$\textsc{V}_{r}^{\top}\textsc{V}_{r}\Sigma_{r}^{2}\in\mathbb{R}^{r\times r}$,
then the eigenvalue decomposition the RRR estimator
$\widehat{G}^{r}_{\mu,\gamma}=\widehat{Z}_{\mu}\textsc{U}_{r}\textsc{V}_{t}^{\top}\widehat{S}$
is given by
$\widehat{G}^{r}_{\mu,\gamma}=\sum_{i\in[r]}(\mu\\!-\\!\widehat{\lambda}_{i})^{-1}\,\widehat{h}_{i}\otimes\widehat{g}_{i}$,
where $\widehat{\lambda}_{i}=\mu\\!-\\!1/\nu_{i}$,
$\widehat{g}_{i}\\!=\\!\nu_{i}^{-1/2}\widehat{S}^{*}\textsc{V}_{r}w^{\ell}_{i}$
and $\widehat{h}_{i}\\!=\\!\widehat{Z}_{\mu}^{*}\textsc{U}_{r}w^{r}_{i}$ for
$\textsc{U}_{r}\\!=\\!(\mu\gamma)^{\\!-\\!1}[{\mu}^{\\!-\\!1/2}I\,|\,\\!-\\!\textsc{N}(\textsc{M}\\!+\\!\gamma\mu
I)^{\\!-\\!1}]^{\top}(\textsc{K}\textsc{V}_{r}\\!-\\!\mu\textsc{V}_{r}\Sigma_{r}^{2})\\!\in\\!\mathbb{R}^{(1\\!+\\!p)n\times
r}$.

The main computational cost of our method, in view of (12), to solve SDE (1)
lies in the implicit inversion of $\textsc{J}_{\mu,\gamma}$ when solving (20).
When computed with direct solvers this inversion is of the order
$\mathcal{O}(n^{3}p^{3})$, however leveraging on the fact that
$\mu\textsc{J}_{\mu,\gamma}$ is Schur’s complement of the
$(1\\!+\\!p)n{\times}(1\\!+\\!p)n$ symmetric positive definite matrix and
using classical iterative solvers, like Lanczos or the generalized Davidson
method, when $r\ll n$ this cost can significantly be reduced to
$\mathcal{O}(r\,n^{2}p^{2})$, c.f. [14].

##  5 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 $\mathcal{E}$ and metric distortion $\eta$. For this
purpose, we analyze the operator norm error of empirical estimator
$\widehat{G}^{r}_{\mu,\gamma}$ using the decomposition

$\mathcal{E}(\widehat{G})\leq\underbrace{\lVert(\mu
I\,{-}\,L)^{-1}Z_{\mu}\\!-\\!Z_{\mu}G_{\mu,\gamma}\rVert_{\mathcal{H}\to\mathcal{W}}}_{\text{regularization
bias}}\\!+\\!\underbrace{\lVert
Z_{\mu}(G_{\mu,\gamma}\\!\\!-\\!G^{r}_{\mu,\gamma})\rVert_{\mathcal{H}\to\mathcal{W}}}_{\text{rank
reduction bias}}\\!+\\!\underbrace{\lVert
Z_{\mu}(G^{r}_{\mu,\gamma}\\!\\!-\\!\widehat{G}^{r}_{\mu,\gamma})\rVert_{\mathcal{H}\to\mathcal{W}}}_{\text{estimator's
variance}},$

where $G_{\mu,\gamma}\,{=}\,W_{\mu,\gamma}^{-1}C$ is the minimizer of the full
(i.e. without rank constraint) Tikhonov regularized risk and
$G^{r}_{\mu,\gamma}\,{=}\,W_{\mu,\gamma}^{-1/2}[\\![W_{\mu,\gamma}^{-1/2}C]\\!]_{r}$
is the population version of the empirical estimator
$\widehat{G}^{r}_{\mu,\gamma}$.

Note that, while the last two terms in the above decomposition depend on the
estimator, the first term depends only on the choice of $\mathcal{H}$ and the
regularity of $L$ w.r.t. $\mathcal{H}$. In this work we focus on the classical
kernel-based learning where one chooses a universal kernel [38, Chapter 4] and
controls the regularization bias with a regularity condition. For details see
Rem. 2 of App. E.2. Let $\mu>0$ be a prescribed parameter, we make following
assumptions to quantify the difficulty of the learning problem:

  1. (BK)

_Boundedness._ There exists $c_{\mathcal{W}}{>}0$ such that
$\displaystyle{\operatorname*{\text{\rm ess\,sup}}_{x\sim\pi}}\lVert
w\phi(x)\rVert^{2}\\!\leq\\!c_{\mathcal{W}}$, i.e.
$w\phi\\!\in\\!\mathcal{L}^{\infty}_{\pi}(\mathcal{X},\mathcal{H}^{1\\!+\\!p})$;

  1. (RC)

_Regularity._ For some $\alpha\in(0,2]$ there exists $c_{\alpha}>0$ such that
$C^{2}\preceq c_{\alpha}^{2}W_{\mu}^{1+\alpha}$;

  1. (SD)

_Spectral Decay._ There exists $\beta\,{\in}\,(0,1]$ and $c_{\beta}\,{>}\,0$
s.t. $\lambda_{j}(W_{\mu})\,{\leq}\,c_{\beta}\,j^{-1/\beta}$, for all $j\in
J$.

The above assumptions, discussed in more details in App E.1, are in the spirit
of state-of-the-art analysis of statistical learning theory of classical
regression in RKHS spaces [12], recently extended to regression of transfer
operators [28, 22]. In our setting, instead of relying on the injection into
$\mathcal{L}^{2}_{\pi}$ as in the case of transfer operators, the relevant
object is the injection $Z_{\mu}$ into the energy space $\mathcal{W}$.

The first assumption (BK) is the main limiting factor of our approach. Indeed,
since $\lVert
w\phi(x)\rVert^{2}=\mu\lVert\phi(x)\rVert^{2}+\sum_{k\in[p]}\lVert
d_{k}\phi(x)\rVert^{2}$, 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_{\mu}\in{\rm{HS}}\left(\mathcal{H},\mathcal{W}^{\mu}_{\pi}(\mathcal{X})\right)$
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_{\mu}=Z_{\mu}^{*}Z_{\mu}$ instead of $C=S_{\pi}^{*}S_{\pi}$. On the other
hand, the "cross-covariance" that captures now RKHS cross-correlations of the
resolvent w.r.t domain is simply $Z_{\mu}^{*}(\mu
I\,{-}\,L)^{-1}Z_{\mu}=S_{\pi}^{*}(\mu I{-}L)(\mu I\,{-}\,L)^{-1}S_{\pi}=C$.
With this in mind, (RC) corresponds to the the regularity condition in [22]
and it quantifies the relationship between the hypothesis class (bounded
operators in $RKHS$) and the object of interest $(\mu I\,{-}\,L)^{-1}$.
Indeed, if $L$ has eigenfunctions that belong to $\alpha$-interpolation space
between $\mathcal{H}$ and $\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$, (RC) holds
true. In particular, if $f_{i}\in\mathcal{H}$ for all $i\geq 2$ (constant
eigenfunction excluded), one has that $\alpha\geq 1$ (c.f. Proposition 7 in
App. E.1.

Finally, (SD) quantifies the "size" of the hypothetical domain $\mathcal{H}$
within the true domain $\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$ via the
effective dimension $\operatorname{\text{\rm
tr}}(W_{\mu,\gamma}^{-1}W_{\mu})\leq c_{\beta}(\mu\gamma)^{-\beta}$, which,
due to ((BK)), leads to another notion, known as the embedding property,
quantifying the relationship between $\mathcal{H}$ and essentially bounded
functions in the domain of the operator

  1. (KE)

_Kernel Embedding._ There exists $\tau\in[\beta,1]$ and such that

$c_{\tau}=\operatorname*{\text{\rm
ess\,sup}}_{x\sim\pi}\textstyle{\sum_{j\in\mathbb{N}}}\sigma^{2\tau}_{j}[\mu|z_{j}(x)|^{2}-z_{j}(x)[Lz_{j}](x)]<+\infty,$
(21)

where $Z_{\mu}=\sum_{j\in J}\sigma_{j}z_{j}\otimes h_{j}$ is the SVD of the
injection operator
$Z_{\mu}\colon\mathcal{H}\to\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$.

Using the above assumptions, we prove a generalization bound for the RRR
estimator of the resolvent.

######  Theorem 2.

Let (RC), (SD), and (KE) hold for some $\alpha\in(0,2]$, $\beta\in(0,1]$ and
$\tau\in[\beta,1]$, respectively, and let $\operatorname*{\text{\rm
cl}}(\operatorname*{\text{\rm
Im}}(S_{\pi}))=\mathcal{L}^{2}_{\pi}(\mathcal{X})$. Denoting
$\lambda_{k}^{\star}=\lambda_{k}(S_{\pi}^{*}(\mu I\,{-}\,L)^{-1}S_{\pi})$,
$k\in\mathbb{N}$, and given $\delta\in(0,1)$ and $r\in[n]$, if RRR estimator
is build from the Dirichlet form and

$\text{if }\,\alpha\,{\geq}\,\tau,\;\gamma\asymp
n^{-\frac{1}{\alpha+\beta}}\text{ and
}\varepsilon^{\star}_{n}=n^{-\frac{\alpha}{2(\alpha+\beta)}},\,\text{ or if
}\,\alpha\,{<}\,\tau,\;\gamma\asymp n^{-\frac{1}{\beta+\tau}}\text{ and
}\varepsilon^{\star}_{n}=n^{-\frac{\alpha}{2(\beta+\tau)}}$ (22)

then there exists a constant $c\,{>}\,0$, depending only on $\mathcal{H}$ and
gap $\lambda_{r}^{\star}{-}\lambda_{r+1}^{\star}\,{>}\,0$, such that for large
enough $n\geq r$ with probability at least $1\,{-}\,\delta$ in the i.i.d. draw
of $\mathcal{D}_{n}$ from $\pi$ it holds that

$\mathcal{E}(\widehat{G}^{r}_{\mu,\gamma})\leq(\widehat{\sigma}_{r+1}\,\wedge\,\textstyle{\sqrt{\lambda_{r+1}^{\star}})}+c\,\varepsilon^{\star}_{n}\,\ln\delta^{-1}.$
(23)

_Proof sketch._ The regularization bias is bounded by
$c_{\alpha}\,\gamma^{{\alpha/2}}$ by Prop. 8 of App. E.2, the rank reduction
bias is upper bounded by $\lambda_{r+1}(S_{\pi}^{*}(\mu
I\,{-}\,L)^{-1}S_{\pi})$, while 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 is then chained to Pinelis-Sakhanenko’s inequality
and Minsker’s inequality for self-adjoint HS-operators, Props. 11 and 12 in
App. E.3.1, respectively. Combining the bias and variance terms, we obtain the
balancing equations for the regularization parameter and then the next result
follows. $\Box$

Note that the learning rate (23) implies the $\mathcal{L}^{2}_{\pi}$-norm
learning rate. Moreover, for $\alpha\geq\tau$, it matches information
theoretic lower bounds for transfer operator learning upon replacing
parameters $\alpha$, $\beta$ and $\tau$ related to the space $\mathcal{W}$
with their $\mathcal{L}^{2}_{\pi}$ analogues [22], see App. E.6. This
motivates the development of the first mini-max optimality for the IG
learning, for which our results an important first step.

To conclude this section, we address the spectral learning bounds steaming
from the Prop. 2. The main task to do so is to control the metric distortions,
which we show how in App. E.5. In this context, an additional assumption
$\alpha\geq 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

$\tfrac{\lvert\lambda_{i}-\widehat{\lambda}_{i}\rvert}{|\mu-\lambda_{i}||\mu-\widehat{\lambda}_{i}|}\leq(\widehat{s}_{i}\,\wedge\,\textstyle{2\sqrt{\lambda_{r+1}^{\star}/\lambda_{r}^{\star}})}+c\,\varepsilon^{\star}_{n}\,\ln\delta^{-1},\;i\in[r],$
(24)

where $\widehat{s}_{i}=\widehat{\sigma}_{i}\,\widehat{\eta}_{i}$ is the
empirical spectral bias that informs how good is the estimation of the
particular eigenpair is (see Fig. 1 a)), $\widehat{\sigma}_{i}$ being given in
(20) and
$\widehat{\eta}_{i}\\!=\\!\lVert\widehat{h}_{i}\rVert_{\mathcal{H}}/\|\widehat{W}_{\mu}^{1/2}\widehat{h}_{i}\|_{\mathcal{H}}$.

##  6 Experiments

In this section, we showcase the key features of our method outlined in Table
LABEL:tab:summary. We demonstrate that our approach: (1) avoids the spurious
effects noted in other IG methods [15, 1], (2) is more effective than transfer
operator methods [21], and (3) validates our bounds in a prediction task for a
model with a non-constant diffusion term. Further details are available in
Appendix F.

Figure 1: a) Empirical biases
$\hat{s}_{1}=\widehat{\sigma}_{1}\,\widehat{\eta}_{1}$ and estimation of the
first (nontrivial) eigenfunction of the IG of a Langevin process under a four
well potential. Ground truth is black, our method RRR is red and blue for two
different kernel lengthscales. b) Eigenvalue estimation for the same process
compared to the methods in [15, 1], for which eigenvalue histogram in blue
shows spuriousness. c) Estimation of the second eigenfunction of a Langevin
process under Muller brown potential (white level lines) by RRR, Transfer
Operator (TO) in d) and ground truth in e). Observe that TO fails to recover
the metastable state. f) Prediction RMSE for the CIR model.

One dimensional four well potential  We first investigate the overdamped
Langevin dynamics in a potential that presents four different wells, two
principal wells and then in each of them two smaller wells, given by
$V(x)\,{=}\,4(x^{8}{+}0.8\exp(-80(x^{2})){+}0.2\exp(-80(x{-}0.5)^{2}){+}0.5\exp(-40(x{+}0.5)^{2}))$.
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. In Figure 1, we show that the empirical
bias $\widehat{s}_{1}=\widehat{\sigma}_{1}\,\widehat{\eta}_{1}$ allows us to
choose the hyperparameters of the model: higher empirical bias leads to wrong
operator. Compared to KRR [15, 1], there is no spuriousness in the estimation
of our eigenvalues, as shown in panel b).

Muller Brown potential  We next study Langevin dynamics under more challenging
conditions: the Muller brown potential. Figure 1(c-e) depicts the second
eigenfunction obtained by our method compared to the ground truth one and that
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, we were able to recover
this second eigenfunction.

CIR model  Finally, we show with the CIR model that our method is not limited
to Langevin process with constant diffusion. For this process, the conditional
expectation of the state $X_{t}$ is analytically known. We can thus compare
the prediction of our model with respect to this expectation using 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 $t=\ln(2)/a$ which corresponds to the half life of the mean
reversion. Results are shown in panel d) Figure 1.

##  7 Conclusion

We developed a novel energy-based 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 $n^{2}d^{2}$. Future work will explore
alternative methods to enhance computational efficiency and investigate a
broader suite of SDEs beyond stochastic diffusion.

## References

  * A Cabannes and Bach, [2024] 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.
  * Alexander and Giannakis, [2020] Alexander, R. and Giannakis, D. (2020).  Operator-theoretic framework for forecasting nonlinear time series with kernel analog techniques.  Physica D: Nonlinear Phenomena, 409:132520.
  * Andrew et al., [2013] Andrew, G., Arora, R., Bilmes, J., and Livescu, K. (2013).  Deep canonical correlation analysis.  In International conference on machine learning, pages 1247–1255. PMLR.
  * Bakry et al., [2014] Bakry, D., Gentil, I., and Ledoux, M. (2014).  Analysis and Geometry of Markov Diffusion Operators.  Springer.
  * Bevanda et al., [2021] Bevanda, P., Beier, M., Kerz, S., Lederer, A., Sosnowski, S., and Hirche, S. (2021).  KoopmanizingFlows: Diffeomorphically Learning Stable Koopman Operators.  arXiv preprint arXiv.2112.04085.
  * Bouvrie and Hamzi, [2017] Bouvrie, J. and Hamzi, B. (2017).  Kernel Methods for the Approximation of Nonlinear Systems.  SIAM Journal on Control and Optimization, 55(4):2460–2492.
  * Brunton et al., [2022] 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.
  * Caponnetto and De Vito, [2007] Caponnetto, A. and De Vito, E. (2007).  Optimal rates for the regularized least-squares algorithm.  Foundations of Computational Mathematics, 7(3):331–368.
  * Das and Giannakis, [2020] Das, S. and Giannakis, D. (2020).  Koopman spectra in reproducing kernel Hilbert spaces.  Applied and Computational Harmonic Analysis, 49(2):573–607.
  * Davis and Kahan, [1970] Davis, C. and Kahan, W. M. (1970).  The rotation of eigenvectors by a perturbation. iii.  SIAM Journal on Numerical Analysis, 7(1):1–46.
  * Fan et al., [2021] Fan, F., Yi, B., Rye, D., Shi, G., and Manchester, I. R. (2021).  Learning Stable Koopman Embeddings.  arXiv preprint arXiv.2110.06509.
  * Fischer and Steinwart, [2020] Fischer, S. and Steinwart, I. (2020).  Sobolev norm learning rates for regularized least-squares algorithms.  Journal of Machine Learning Research, 21(205):1–38.
  * Haag, [2012] Haag, R. (2012).  Local quantum physics: Fields, particles, algebras.  Springer Science & Business Media.
  * Hogben, [2006] Hogben, L., editor (2006).  Handbook of Linear Algebra.  CRC Press, Boca Raton, FL, USA.
  * Hou et al., [2023] Hou, B., Sanjari, S., Dahlin, N., Bose, S., and Vaidya, U. (2023).  Sparse learning of dynamical systems in RKHS: An operator-theoretic approach.  In Krause, A., Brunskill, E., Cho, K., Engelhardt, B., Sabato, S., and Scarlett, J., editors, Proceedings of the 40th International Conference on Machine Learning, volume 202 of Proceedings of Machine Learning Research, pages 13325–13352. PMLR.
  * Inzerilli et al., [2023] Inzerilli, P., Kostic, V., Lounici, K., Novelli, P., and Pontil, M. (2023).  Consistent long-term forecasting of ergodic dynamical systems.
  * Kawahara, [2016] Kawahara, Y. (2016).  Dynamic Mode Decomposition with Reproducing Kernels for Koopman Spectral Analysis.  In Lee, D., Sugiyama, M., Luxburg, U., Guyon, I., and Garnett, R., editors, Advances in Neural Information Processing Systems, volume 29. Curran Associates, Inc.
  * Kessler et al., [2012] Kessler, M., Lindner, A., and Sørensen, M. (2012).  Statistical methods for stochastic differential equations.  Monographs on Statistics and Applied Probability, 124:7–12.
  * Klus et al., [2020] Klus, S., Nüske, F., Peitz, S., Niemann, J.-H., Clementi, C., and Schütte, C. (2020).  Data-driven approximation of the koopman generator: Model reduction, system identification, and control.  Physica D: Nonlinear Phenomena, 406:132416.
  * Klus et al., [2019] 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.
  * Kostic et al., [2022] 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.
  * [22] Kostic, V. R., Lounici, K., Novelli, P., and massimiliano pontil (2023a).  Sharp spectral rates for koopman operator learning.  In Thirty-seventh Conference on Neural Information Processing Systems.
  * [23] Kostic, V. R., Novelli, P., Grazzi, R., Lounici, K., et al. (2023b).  Learning invariant representations of time-homogeneous stochastic dynamical systems.  In The Twelfth International Conference on Learning Representations.
  * Kutoyants, [1984] Kutoyants, Y. A. (1984).  Parameter estimation for stochastic processes.  Berlin: Heldermann.
  * Kutoyants, [2013] Kutoyants, Y. A. (2013).  Statistical inference for ergodic diffusion processes.  Springer Science & Business Media.
  * Kutz et al., [2016] Kutz, J. N., Brunton, S. L., Brunton, B. W., and Proctor, J. L. (2016).  Dynamic Mode Decomposition.  Society for Industrial and Applied Mathematics.
  * Li et al., [2017] Li, Q., Dietrich, F., Bollt, E. M., and Kevrekidis, I. 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).
  * Li et al., [2022] 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.
  * Lusch et al., [2018] Lusch, B., Kutz, J. N., and Brunton, S. L. (2018).  Deep learning for universal linear embeddings of nonlinear dynamics.  Nature Communications, 9(1).
  * Mackey, [2011] Mackey, M. C. (2011).  Time’s Arrow: The origins of thermodynamic behavior.  Courier Corporation.
  * Mardt et al., [2018] Mardt, A., Pasquali, L., Wu, H., and Noé, F. (2018).  VAMPnets for deep learning of molecular kinetics.  Nature Communications, 9(1).
  * Minsker, [2017] Minsker, S. (2017).  On some extensions of Bernstein’s inequality for self-adjoint operators.  Statistics & Probability Letters, 127:111–119.
  * Oksendal, [2013] Oksendal, B. (2013).  Stochastic differential equations: an introduction with applications.  Springer Science & Business Media.
  * Pascucci, [2011] Pascucci, A. (2011).  PDE and Martingale Methods in Option Pricing.  Springer Milan.
  * Pavliotis, [2014] Pavliotis, G. A. (2014).  Stochastic Processes and Applications.  Springer New York.
  * Pillaud-Vivien and Bach, [2023] Pillaud-Vivien, L. and Bach, F. (2023).  Kernelized diffusion maps.  In The Thirty Sixth Annual Conference on Learning Theory, pages 5236–5259. PMLR.
  * Reed and Simon, [1980] Reed, M. and Simon, B. (1980).  I: Functional Analysis.  Academic Press.
  * Steinwart and Christmann, [2008] Steinwart, I. and Christmann, A. (2008).  Support Vector Machines.  Springer New York.
  * Tian and Wu, [2021] Tian, W. and Wu, H. (2021).  Kernel embedding based variational approach for low-dimensional approximation of dynamical systems.  Computational Methods in Applied Mathematics, 21(3):635–659.
  * Trefethen and Embree, [2020] Trefethen, L. N. and Embree, M. (2020).  Spectra and Pseudospectra: The Behavior of Nonnormal Matrices and Operators.  Princeton University Press.
  * Tropp, [2012] Tropp, J. A. (2012).  User-friendly tail bounds for sums of random matrices.  Technical report.
  * Tuckerman, [2023] Tuckerman, M. E. (2023).  Statistical mechanics: theory and molecular simulation.  Oxford university press.
  * Williams et al., [2015] Williams, M. O., , Rowley, C. W., and Kevrekidis, I. G. (2015).  A kernel-based method for data-driven Koopman spectral analysis.  Journal of Computational Dynamics, 2(2):247–265.
  * Wu and Noé, [2019] Wu, H. and Noé, F. (2019).  Variational Approach for Learning Markov Processes from Time Series Data.  Journal of Nonlinear Science, 30(1):23–66.
  * Yeung et al., [2019] 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.
  * Yosida, [2012] Yosida, K. (2012).  Functional analysis, volume 123.  Springer Science & Business Media.
  * Zabczyk, [2020] Zabczyk, J. (2020).  Mathematical Control Theory: An Introduction.  Systems & Control: Foundations & Applications. Springer International Publishing.
  * Zhang et al., [2022] 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:111377.
  * Zhou, [2008] Zhou, D.-X. (2008).  Derivative reproducing properties for kernel methods in learning theory.  Journal of computational and Applied Mathematics, 220(1-2):456–463.
  * Zwald and Blanchard, [2005] Zwald, L. and Blanchard, G. (2005).  On the Convergence of Eigenspaces in Kernel Principal Component Analysis.  In Advances in Neural Information Processing Systems.

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 $\wedge$ minimum $\vee$ maximum
$[\\![\,\cdot\,]\\!]_{r}$ $r$-truncated SVD of an operator $I$ identity
operator ${\rm{HS}}\left(\mathcal{H},\mathcal{G}\right)$ space of Hilbert-
Schmidt operators $\mathcal{H}\to\mathcal{G}$ ${\rm{B}}_{r}({\mathcal{H}})$
set of rank-$r$ Hilbert-Schmidt operators on $\mathcal{H}$ $\lVert A\rVert$
operator norm of an operator $A$ $\lVert A\rVert_{\rm{HS}}$ Hilbert-Schmidt
norm of operator $A$ $\sigma_{i}(\cdot)$ $i$-th singular value of an operator
$\lambda_{i}(\cdot)$ $i$-th eigenvalue of an operator $\mathcal{X}$ state
space of the Markov process $(X_{t})_{t\geq 0}$ time-homogeneous Markov
process $p$ transition kernel of the Markov process $\pi$ invariant measure of
the Markov process $a$ drift of the Itô process $b$ diffusion of the Itô
process $\mathcal{L}^{2}_{\pi}(\mathcal{X})$ L2 space of functions on
$\mathcal{X}$ w.r.t. measure $\pi$ $A$ Transfer operator on
$\mathcal{L}^{2}_{\pi}(\mathcal{X})$ $\mathcal{W}^{1,2}_{\pi}(\mathcal{X})$
Sobolev space w.r.t. measure $\pi$ on $\mathcal{X}$ $L$ Generator of the
semigroup on $\mathcal{W}^{1,2}_{\pi}(\mathcal{X})$ $s$ Dirichlet form
diffusion $B$ Dirichlet form operator on
$\mathcal{W}^{1,2}_{\pi}(\mathcal{X})$ $\mu$ shift parameter
$\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$ regularized energy space $k(x,y)$
kernel $\phi$ canonical feature map $\mathcal{H}$ reproducing kernel Hilbert
space $S_{\pi}$ canonical injection
$\mathcal{H}\hookrightarrow\mathcal{L}^{2}_{\pi}(\mathcal{X})$ $\ell\phi$
generator embedding $Z_{\mu}$ canonical injection
$\mathcal{H}\hookrightarrow\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$ $D^{v}$
$v$-directional derivative $D$ derivative $D^{v}\phi$ embedding of the
$v$-directional derivative $D\phi{}$ derivative embedding $w_{k}\phi$ $k$-th
component of Dirichlet operator embedding $w\phi$ Dirichlet operator embedding
$\sigma_{j}$ $j$-th singular value of $Z_{\mu}$ $J$ countable index set of
singular values of $Z_{\mu}$ $z_{j}$ $j$-th left singular function of
$Z_{\mu}$ $h_{j}$ $j$-th right singular function of $Z_{\mu}$ $\mathbb{1}$
function in $\mathcal{L}^{2}_{\pi}(\mathcal{X})$ with the constant output 1
$\gamma$ regularization parameter $\mathcal{R}$ true risk $\mathcal{E}$
operator norm error $\mathcal{E}_{\rm HS}$ excess risk, i.e. HS norm error
$\mathcal{R}_{0}$ irreducible risk $\mathcal{D}_{n}$ dataset
$(x_{i})_{i\in[n]}$ $\widehat{\mathcal{R}}$ empirical risk $\widehat{S}$
sampling operator w.r.t. $\mathcal{L}^{2}_{\pi}(\mathcal{X})$
$\widehat{Z}_{\mu}$ sampling operator w.r.t.
$\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$ $C$ covariance operator w.r.t.
$\mathcal{L}^{2}_{\pi}(\mathcal{X})$ $\widehat{C}$ empirical covariance
operator w.r.t. $\mathcal{L}^{2}_{\pi}(\mathcal{X})$ $C_{\gamma}$ regularized
covariance operator w.r.t. $\mathcal{L}^{2}_{\pi}(\mathcal{X})$
$\widehat{C}_{\gamma}$ regularized empirical covariance operator w.r.t.
$\mathcal{L}^{2}_{\pi}(\mathcal{X})$ $W_{\mu}$ covariance operator w.r.t.
$\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$ $\widehat{W}_{\mu}$ empirical
covariance operator w.r.t. $\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$
$W_{\mu,\gamma}$ regularized covariance operator w.r.t.
$\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$ $\widehat{W}_{\mu,\gamma}$ regularized
empirical covariance operator w.r.t. $\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$
$T$ derivative-covariance operator $\widehat{T}$ empirical derivative-
covariance operator K kernel Gram matrix w.r.t.
$\mathcal{L}^{2}_{\pi}(\mathcal{X})$ $\textsc{K}_{\gamma}$ regularized kernel
Gram matrix w.r.t. $\mathcal{L}^{2}_{\pi}(\mathcal{X})$ $\textsc{F}_{\mu}$
kernel Gram matrix w.r.t. $\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$
$\textsc{F}_{\mu,\gamma}$ regularized kernel Gram matrix w.r.t.
$\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$ M derivative-derivative kernel Gram
matrix N feature-derivative kernel Gram matrix $G$ population estimator of
$(\mu I\,{-}\,L)^{-1}$ on $\mathcal{H}$ $\widehat{G}$ empirical estimator of
$(\mu I\,{-}\,L)^{-1}$ on $\mathcal{H}$ $G_{\mu,\gamma}$ population KRR
estimator $\widehat{G}_{\mu,\gamma}$ empirical KRR estimator
$G^{r}_{\mu,\gamma}$ population RRR estimator $\widehat{G}^{r}_{\mu,\gamma}$
empirical RRR estimator $P$ spectral projector $\widehat{P}$ empirical
spectral projector $\eta$ metric distortion $\widehat{\eta}$ empirical metric
distortion $\lambda$ generator eigenvalue $\widehat{\lambda}$ eigenvalue of
the empirical estimator $f$ generator eigenfunction in
$\mathcal{L}^{2}_{\pi}(\mathcal{X})$ $\widehat{f}$ empirical eigenfunction in
$\mathcal{L}^{2}_{\pi}(\mathcal{X})$ $\widehat{h}$ right empirical
eigenfunction $\widehat{g}$ left empirical eigenfunction $c_{\mathcal{W}}$
boundness constant $P_{\mathcal{H}}$ orthogonal projector in
$\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$ onto $\operatorname*{\text{\rm
Im}}(Z_{\mu})$ $\alpha$ regularity parameter $c_{\alpha}$ regularity constant
$\beta$ spectral decay parameter $c_{\beta}$ spectral decay constant $\tau$
embedding parameter $c_{\tau}$ embedding constant

Table 1: Summary of used notations.

##  Appendix A Background

###  A.1 Basics on operator theory for Markov processes.

We provide here some basics on operator theory for Markov processes. Let
$\mathcal{X}\subset\mathbb{R}^{d}$ ($d\in\mathbb{N}$) and
$(X_{t})_{t\in\mathbb{R}_{+}}$ be a $\mathcal{X}$-valued time-homogeneous
Markov process defined on a filtered probability space
$(\Omega,\mathcal{F},(\mathcal{F}_{t})_{t\in\mathbb{R}_{+}},\mathbb{P})$ where
$\mathcal{F}_{t}=\sigma(X_{s},\,s\leq t)$ is the natural filtration of
$(X_{t})_{t\in\mathbb{R}_{+}}$. The dynamics of $X$ can be described through
of a family of probability densities $(p_{t})_{t\in\mathbb{R}_{+}}$ such that
for all $t\in\mathbb{R}_{+}$, $E\in\mathcal{B}(\mathcal{X})$,

$\mathbb{P}(X_{t}\in E|X_{0}=x)=\int_{E}p_{t}(x,y)dy.$

Let $\mathcal{G}$ be a set of real valued and measurable functions on
$\mathcal{X}$. For any $t\in\mathbb{R}_{+}$ the _transfer operator_ (TO)
$A_{t}:\mathcal{G}\rightarrow\mathcal{G}$ maps any measurable function
$f\in\mathcal{G}$ to

$(A_{t}f)(x)=\int_{\mathcal{X}}f(y)p_{t}(x,y)dy.$ (25)

In theory of Markov processes, the family $(A_{t})_{t\in\mathbb{R}_{+}}$ is
referred to as the _Markov semigroup_ associated to the process $X$.

######  Remark 1.

A possible choice is $\mathcal{G}=\mathcal{L}^{\infty}(\mathcal{X})$. Here, we
are interested in another choice of $\mathcal{G}$ related to the existence of
an invariant measure for $A_{t}$ i.e., a $\sigma$-finite measure $\pi$ on
$(\mathcal{X},\mathcal{B}(\mathcal{X}))$ such that $P_{t}^{*}\pi=\pi$ for any
$t\in\mathbb{R}_{+}$. In that case, we can choose
$\mathcal{G}=\mathcal{L}^{2}_{\pi}(\mathcal{X})$ so that, for all
$f\in\mathcal{L}^{2}_{\pi}(\mathcal{X})$, $P_{t}f$ converges to $f$ in
$\mathcal{L}^{2}_{\pi}(\mathcal{X})$ as $t$ goes to $0$. Note that $P_{0}f=f$
and $P_{\infty}f=\int fd\mathcal{X}$ for a suitable integrable function
$f:\mathcal{X}\rightarrow\mathbb{R}$.

Within the existence of this invariant measure $\pi$, the process $X$ is then
characterized by the _infinitesimal generator_ (IG)
$L:\mathcal{L}^{2}_{\pi}(\mathcal{X})\rightarrow\mathcal{L}^{2}_{\pi}(\mathcal{X})$
of the family $(A_{t})_{t\in\mathbb{R}_{+}}$ defined by

$L=\underset{t\rightarrow 0^{+}}{\lim}\frac{A_{t}-I}{t}.$ (26)

In other words, $L$ characterizes the linear differential equation
$\partial/\partial_{t}A_{t}f=LA_{t}f$ satisfied by the transfer operator. The
domain of $L$ denoted $\mathrm{dom}(L)$ coincides with the the Sobolev space
$\mathcal{W}^{1,2}_{\pi}(\mathcal{X})$

$\mathcal{W}^{1,2}_{\pi}(\mathcal{X})=\\{f\in\mathcal{L}^{2}_{\pi}(\mathcal{X})\;|\;\lVert
f\rVert_{\mathcal{W}}^{2}=\lVert f\rVert_{\mathcal{L}^{2}_{\pi}}+\lVert\nabla
f\rVert_{\mathcal{L}^{2}_{\pi}}<\infty\\}.$

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 ([46], Chap. IX) : For every $\mu>0$, the operator
$(\mu I-L)$ admits an inverse $L_{\mu}=(\mu I-L)^{-1}$ that is a continuous
operator on $\mathcal{X}$ and

$(\mu I-L)^{-1}=\int_{0}^{\infty}e^{-\mu t}A_{t}dt.$

The operator $L_{\mu}$ is the resolvent of $L$ and the corresponding resolvent
set of $L$ is defined by

$\rho(L)=\big{\\{}\mu\in\mathbb{C}\,|\,(\mu I-L)\,\text{is bijective
and}\,L_{\mu}\,\text{is continuous}\big{\\}}.$

In fact, $\rho(L)$ contains all real positive numbers and $(\mu I-L)^{-1}$ is
bounded. Here, we assume that $L$ has compact resolvent, i.e.
$\rho(L)\neq\emptyset$ and there exists $\mu_{0}\in\rho(L)$ such that
$L_{\mu_{0}}$ is compact. Note that, through the resolvent identity, this
implies the compactness of all resolvents $L_{\mu}$ for $\mu\in\rho(L)$. Let
then $\mu\in\rho(L)$. As $(L,\mathrm{dom}(L))$ is assumed to be self-adjoint,
so does $L_{\mu}$, so that $L_{\mu}$ is both compact and self-adjoint. Then,
its spectrum $\operatorname*{\text{\rm
Sp}}(L_{\mu})=\mathbb{C}\setminus\rho(L_{\mu})$ is purely discrete and
consists of isolated eigenvalues $(\lambda_{i}^{\mu})_{i\in\mathbb{N}}$ such
that $|\lambda_{i}^{\mu}|\to 0$ associated with an orthonormal basis
$(f_{i})_{i\in\mathbb{N}}$ (see [46], chapter XI). In other words, the
spectral decomposition of the resolvent writes

$L_{\mu}=\sum_{i\in\mathbb{N}}\lambda_{i}^{\mu}\,f_{i}\otimes f_{i}$

where the functions $f_{i}\in\mathcal{L}^{2}_{\pi}(\mathcal{X})$ are also
eigenfunctions of the operator $L$ we get

$L=\sum_{i\in\mathbb{N}}\lambda_{i}\,f_{i}\otimes f_{i}.$

######  Example 3 (Langevin).

Let $(k_{b}T)\in\mathbb{R}$.The _overdamped Langevin equation_ driven by a
_potential_ $V:\mathbb{R}^{d}\rightarrow\mathbb{R}$ is given by
$dX_{t}=-\nabla V(X_{t})dt+\sqrt{2(k_{b}T)}dW_{t}\quad\text{and}\quad
X_{0}=x$. Its invariant measure of the solution process $X$ is the _Boltzman
distribution_ $\pi(dx)=z^{-1}e^{-(k_{b}T)^{-1}V(x)}dx$ where $z$ denotes a
normalizing constant. Its infinitesimal generator $L$ is defined by
$Lf=-\nabla V^{\top}\nabla f+(k_{b}T)\Delta f$, for
$f\in\mathcal{W}^{1,2}_{\pi}(\mathcal{X})$, and if $\nabla
V\in\mathcal{L}^{2}_{\pi}(\mathcal{X})$ is positive and coercive, it has
compact resolvent. Finally, since
$\int(-Lf)g\,d\pi=-\int\Big{[}\nabla\Big{(}(k_{b}T)\nabla
f(x)\frac{e^{-(k_{b}T)^{-1}V(x)}}{Z}\Big{)}\Big{]}g(x)dx=(k_{b}T)\int\nabla
f^{\top}\nabla g\,d\pi=\int f(-Lg)\,d\pi$, generator $L$ is self-adjoint and
associated to a gradient Dirichlet form with
$s(x)=(k_{b}T)^{1/2}(\delta_{ij})_{i\in[d],j\in[p]}$.

######  Example 4 (Cox-Ingersoll-Ross process).

Let $d=1$, $a,b\in\mathbb{R}$, $\sigma\in\mathbb{R}_{+}^{*}$. The _Cox-
Ingersoll-Ross_ process is solution of the SDE
$dX_{t}=(a+bX_{t})dt+\sigma\sqrt{X_{t}}dW_{t}\quad\text{and}\quad X_{0}=x$.
Its invariant measure $\pi$ is a Gamma distribution with shape parameter
$a/\sigma^{2}$ and rate parameter $b/\sigma^{2}$. Its infinitesimal generator
$L$ is defined for $f\in\mathcal{L}^{2}_{\pi}(\mathcal{X})$ by
$Lf=(a+bx)\nabla f+\frac{\sigma^{2}x}{2}\Delta f$. Note that by integration by
parts, we can check that the generator $L$ satisfies
$\int(-Lf)g\,d\pi=\int\frac{\sigma^{2}x}{2}f^{\prime}(x)g^{\prime}(x)\,\pi(dx)=\int
f(-Lg)\,d\pi$, and it is associated to a gradient Dirichlet form with
$s(x)=\sigma\sqrt{x}/\sqrt{2}$.

###  A.2 Infinitesimal generator for diffusion processes

After defining the infinitesimal generator for Markov processes (see A.1), we
provide its explicit form for solution processes of equations like(1). Given a
smooth function $f\in\mathcal{C}^{2}(\mathcal{X},\mathbb{R})$, Itô’s formula
(see for instance [4], B, p.495) provides for $t\in\mathbb{R}_{+}$,

$\displaystyle f(X_{t})-f(X_{0})$
$\displaystyle=\int_{0}^{t}\sum_{i=1}^{d}\partial_{i}f(X_{s})dX_{s}^{i}+\tfrac{1}{2}\int_{0}^{t}\sum_{i,j=1}^{d}\partial_{ij}^{2}f(X_{s})d\langle
X^{i},X^{j}\rangle_{s}$ $\displaystyle=\int_{0}^{t}\nabla
f(X_{s})^{\top}dX_{t}+\tfrac{1}{2}\int_{0}^{t}\mathrm{Tr}\big{[}X_{s}^{\top}(\nabla^{2}f)(X_{s})X_{s}\big{]}ds.$

Recalling (1), we get

$\displaystyle f(X_{t})$ $\displaystyle=f(X_{0})+\int_{0}^{t}\bigg{[}\nabla
f(X_{s})^{\top}a(X_{s})+\tfrac{1}{2}\mathrm{Tr}\big{[}b(X_{s})^{\top}(\nabla^{2}f(X_{s}))b(X_{s})\big{]}\bigg{]}ds$
$\displaystyle\quad+\int_{0}^{t}\nabla f(X_{s})^{\top}b(X_{s})dW_{s}.$ (27)

Provided $f$ and $b$ are smooth enough, the expectation of the last stochastic
integral vanishes so that we get

$\mathbb{E}[f(X_{t})|X_{0}=x]=f(x)+\int_{0}^{t}\mathbb{E}\Big{[}\nabla
f(X_{s})^{\top}a(X_{s})+\tfrac{1}{2}\mathrm{Tr}\big{[}b(X_{s})^{\top}(\nabla^{2}f(X_{s}))b(X_{s})\big{]}\Big{|}X_{0}=x\Big{]}ds$

Recalling that $L=\underset{t\rightarrow 0^{+}}{\lim}(A_{t}f-f)/t$, we get for
every $x\in\mathcal{X}$,

$\displaystyle Lf(x)$ $\displaystyle=\lim_{t\rightarrow
0}\frac{\mathbb{E}[f(X_{t})|X_{0}=x]-f(x)}{t}$
$\displaystyle=\lim_{t\rightarrow
0}\frac{1}{t}\bigg{[}\int_{0}^{t}\mathbb{E}\Big{[}\nabla
f(X_{s})^{\top}a(X_{s})+\tfrac{1}{2}\mathrm{Tr}\big{[}(X_{s})^{\top}(\nabla^{2}f(X_{s}))b(X_{s})\big{]}\Big{]}ds\Big{|}X_{0}=x\bigg{]}$
$\displaystyle=\nabla
f(x)^{\top}a(x)+\tfrac{1}{2}\mathrm{Tr}\big{[}b(x)^{\top}(\nabla^{2}f(x))b(x)\big{]},$
(28)

which provides the closed formula for the IG associated with the solution
process of (1).

###  A.3 Spectral perturbation theory

Recalling that for a bounded linear operator $A$ on some Hilbert space
$\mathcal{H}$ the resolvent set of the operator $A$ is defined as
$\rho(A)=\left\\{\lambda\in\mathbb{C}\,|\,A-\lambda I\text{ is bijective and
}(A-\lambda I)^{-1}\text{ is continuous}\right\\}$, and its spectrum
$\operatorname*{\text{\rm Sp}}(A)=\mathbb{C}\setminus\\{\rho(A)\\}$, let
$\lambda\subseteq\operatorname*{\text{\rm Sp}}(A)$ be isolated part of
spectra, i.e. both $\lambda$ and $\mu=\operatorname*{\text{\rm
Sp}}(A)\setminus\lambda$ are closed in $\operatorname*{\text{\rm Sp}}(A)$.
Than, the Riesz spectral projector
$P_{\lambda}\colon\mathcal{H}\to\mathcal{H}$ is defined by

$P_{\lambda}=\frac{1}{2\pi}\int_{\Gamma}(zI-A)^{-1}dz,$ (29)

where $\Gamma$ is any contour in the resolvent set $\operatorname*{\text{\rm
Res}}(A)$ with $\lambda$ in its interior and separating $\lambda$ from $\mu$.
Indeed, we have that $P_{\lambda}^{2}=P_{\lambda}$ and
$\mathcal{H}=\operatorname*{\text{\rm
Im}}(P_{\lambda})\oplus\operatorname{\text{\rm Ker}}(P_{\lambda})$ where
$\operatorname*{\text{\rm Im}}(P_{\lambda})$ and $\operatorname{\text{\rm
Ker}}(P_{\lambda})$ are both invariant under $A$ and $\operatorname*{\text{\rm
Sp}}(A_{|_{\operatorname*{\text{\rm Im}}(P_{\lambda})}})=\lambda$,
$\operatorname*{\text{\rm Sp}}(A_{|_{\operatorname{\text{\rm
Ker}}(P_{\lambda})}})=\mu$. Moreover, $P_{\lambda}+P_{\mu}=I$ and
$P_{\lambda}P_{\mu}=P_{\mu}P_{\lambda}=0$.

Finally if $A$ is compact operator, then the Riesz-Schauder theorem, see e.g.
[37], assures that $\operatorname*{\text{\rm Sp}}(T)$ is a discrete set having
no limit points except possibly $\lambda=0$. Moreover, for any nonzero
$\lambda\in\operatorname*{\text{\rm Sp}}(T)$, then $\lambda$ is an eigenvalue
(i.e. it belongs to the point spectrum) of finite multiplicity, and, hence, we
can deduce the spectral decomposition in the form

$A=\sum_{\lambda\in\operatorname*{\text{\rm Sp}}(A)}\lambda\,P_{\lambda},$
(30)

where geometric multiplicity of $\lambda$,
$r_{\lambda}=\operatorname*{\text{\rm rank}}(P_{\lambda})$, is bounded by the
algebraic multiplicity of $\lambda$. If additionally $A$ is normal operator,
i.e. $AA^{*}=A^{*}A$, then $P_{\lambda}=P_{\lambda}^{*}$ is orthogonal
projector for each $\lambda\in\operatorname*{\text{\rm Sp}}(A)$ and
$P_{\lambda}=\sum_{i=1}^{r_{\lambda}}\psi_{i}\otimes\psi_{i}$, where
$\psi_{i}$ are normalized eigenfunctions of $A$ corresponding to $\lambda$ and
$r_{\lambda}$ is both algebraic and geometric multiplicity of $\lambda$.

We conclude this section with well-known perturbation bounds for
eigenfunctions and spectral projectors of self-adjoint compact operators.

######  Proposition 3 ([10]).

Let $A$ be compact self-adjoint operator on a separable Hilbert space
$\mathcal{H}$. Given a pair
$(\widehat{\lambda},\widehat{f})\in\mathbb{C}\times\mathcal{H}$ such that
$\lVert\widehat{f}\rVert=1$, let $\lambda$ be the eigenvalue of $A$ that is
closest to $\widehat{\lambda}$ and let $f$ be its normalized eigenfunction. If
$\widehat{g}=\min\\{\lvert\widehat{\lambda}-\lambda\rvert\,|\,\lambda\in\operatorname*{\text{\rm
Sp}}(A)\setminus\\{\lambda\\}\\}>0$, then
$\sin(\sphericalangle(\widehat{f},f))\leq\lVert
A\widehat{f}-\widehat{\lambda}\widehat{f}\rVert/\widehat{g}$.

######  Proposition 4 ([50]).

Let $A$ and $\widehat{A}$ be two compact operators on a separable Hilbert
space. For nonempty index set $J\subset\mathbb{N}$ let

$\text{\rm
gap}_{J}(A)=\min\left\\{\lvert\lambda_{i}(A)-\lambda_{j}(A)\rvert\,|\,i\in\mathbb{N}\setminus
J,\,j\in J\right\\}$

denote the spectral gap w.r.t $J$ and let $P_{J}$ and $\widehat{P}_{J}$ be the
corresponding spectral projectors of $A$ and $\widehat{A}$, respectively. If
$A$ is self-adjoint and for some $\lVert A-\widehat{A}\rVert<\text{\rm
gap}_{J}(A)$, then

$\lVert P_{J}-\widehat{P}_{J}\rVert\leq\frac{\lVert
A-\widehat{A}\rVert}{\text{\rm gap}_{J}(A)}.$

##  Appendix 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 $\mathcal{L}^{2}_{\pi}(\mathcal{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 a reminder, we consider $\mathcal{H}$ be
an RKHS and let $k:\mathcal{X}\times\mathcal{X}\to\mathbb{R}$ be the
associated kernel function. The canonical feature map is denoted by
$\phi(x)=k(x,\cdot)$ for $x\in\mathcal{X}$
$k(x,x^{\prime})={\langle\phi(x),\phi(x^{\prime})\rangle}$ for all
$x,x^{\prime}\in\mathcal{X}$. Assuming that $k$ is square-integrable with
respect to the measure $\pi$, we define the injection operator
$S_{\pi}:\mathcal{H}\hookrightarrow\mathcal{L}^{2}_{\pi}(\mathcal{X})$ and its
adjoint $S_{\pi}^{*}:\mathcal{L}^{2}_{\pi}(\mathcal{X})\rightarrow\mathcal{H}$
by $S_{\pi}^{*}f=\int_{\mathcal{X}}f(x)\phi(x)\pi(dx).$ As a preliminary step,
we need to define the reproducing partial derivatives in RKHS, which we
introduce via Mercer kernels.

######  Definition 1 (Mercer kernel).

A kernel function $k\,:\,\mathcal{X}\times\mathcal{X}\rightarrow\mathbb{R}$ is
called a Mercer kernel if it is a continuous and symmetric function such that
for any finite set of points $\\{x_{1},\ldots,x_{n}\\}\subset\mathcal{X}$, the
matrix $(k(x_{i},x_{j}))_{i,j=1}^{n}$ is positive semi-definite.

Several standard kernels satisfy the Mercer property with $s\geq 1$, including
the Gaussian kernel which we will consider subsequently.

For $s\in\mathbb{N}$ and $m\in\mathbb{N}$, we define the index set
$I_{s}^{m}=\\{\alpha\in\mathbb{N}^{m}\,:\,|\alpha|\leq s\\}$ where
$|\alpha|=\sum_{j=1}^{s}\alpha_{j}$, for
$\alpha=(\alpha_{1},\ldots,\alpha_{m})\in\mathbb{N}^{m}$. For a function
$f:\mathbb{R}^{m}\to\mathbb{R}$ and $x=(x_{1},\ldots,x_{m})\in\mathbb{R}^{m}$,
we denote its partial derivative $d^{\alpha}f$ at point $x$ (if it exists) as

$D^{\alpha}f(x)=\prod_{j=1}^{m}D_{j}^{\alpha_{j}}f(x)=\frac{\partial^{|\alpha|}}{\partial^{\alpha_{1}}x_{1}\cdots\partial^{\alpha_{m}}x_{m}}f(x).$

For a function $k\in\mathcal{C}^{2s}(\mathcal{X}\times\mathcal{X})$ with
$\mathcal{X}\subset\mathbb{R}^{d}$ and $\alpha\in I_{s}^{d}$, we define

$D^{\alpha}k(x,y)=D^{(\alpha,0)}k(x,y)=\frac{\partial^{|\alpha|}}{\partial^{\alpha_{1}}x_{1}\cdots\partial^{\alpha_{m}}x_{m}}k(x,y).$
(31)

and

$D^{(0,\alpha)}k(x,y)=\frac{\partial^{|\alpha|}}{\partial^{\alpha_{1}}y_{1}\cdots\partial^{\alpha_{m}}y_{m}}k(x,y).$

######  Theorem 3 (Theorem 1 in [49]).

Let $s\in\mathbb{N}$, $\mathcal{X}\subseteq\mathbb{R}^{m}$ and $k$ be a Mercer
kernel such that $k\in\mathcal{C}^{2s}(\mathcal{X}\times\mathcal{X})$ with
corresponding RKHS $\mathcal{H}$. Then the following hold:

  * i.

For any $x\in\mathcal{X}$, $\alpha\in I^{m}_{s}$,

$(D^{\alpha}k)_{x}(y)=D^{\alpha}k(x,y)\in\mathcal{H}.$ (32)

  * ii.

A partial derivative reproducing property holds for $\alpha\in I_{s}^{m}$

$(D^{\alpha}h)(x)=\langle(D^{\alpha}k)_{x},h\rangle_{\mathcal{H}},\quad\forall
h\in\mathcal{H}.$ (33)

Theorem 3 allows us to introduce the first and second order operators $D$ and
$D^{2}$ that act on any feature map $\phi(x)$ as

$D\phi(x)=((D^{e_{i}}k)_{x})_{i\in[d]}\quad\text{and}\quad
D^{2}\phi(x)=((D^{e_{i}+e_{j}}k)_{x})_{i,j\in[d]}$

where the $(D^{e_{i}}k)_{x}$ and $(D^{e_{i}+e_{j}}k)_{x}$ can be defined via
(33). Then, we define the operator $d$ that maps any feature map $\phi(x)$ to
$d\phi(x)=s(x)^{\top}D\phi(x)$. Denote
$s^{\top}=[\bar{s}_{1}|\dots|\bar{s}_{d}]:x\in\mathcal{X}\mapsto
s^{\top}(x)=[\bar{s}_{1}(x)|\dots|\bar{s}_{d}(x)]\in\mathbb{R}^{p\times d}$.
We have

$s(x)^{\top}Dh(x)=\sum_{i=1}^{d}\bar{s}_{i}(x)\partial_{i}h(x)=\sum_{i=1}^{d}\bar{s}_{i}(x)\langle
D^{e_{i}}\phi(x),h\rangle_{\mathcal{H}}=\langle
d\phi(x),h\rangle_{\mathcal{H}},$

so that we can define the embedding of the Dirichlet operator
$B=s^{\top}\nabla\colon\mathcal{L}^{2}_{\pi}(\mathcal{X})\to[\mathcal{L}^{2}_{\pi}(\mathcal{X})]^{p}$
into RKHS $d\phi\colon\mathcal{X}\to\mathcal{H}^{p}$ via the reproducing
property as ${\langle
d\phi(x),h\rangle}_{\mathcal{H}}=[BS_{\pi}h](x)=s(x)^{\top}Dh(x)\in\mathbb{R}^{p}$,
$h\in\mathcal{H}$. In fact, $d\phi$ is a $p$-dimensional vector with
components $d_{k}\phi\colon\mathcal{X}\to\mathcal{H}$ given via ${\langle
d_{k}\phi(x),h\rangle}_{\mathcal{H}}=s_{k}(x)^{\top}Dh(x)$, $k\in[p]$. Then,
we can define

$T=-\mathbb{E}_{x\sim\pi}[d\phi(x)\otimes
d\phi(x)]=-\textstyle{\sum_{k\in[p]}}\mathbb{E}_{x\sim\pi}[d_{k}\phi(x)\otimes
d_{k}\phi(x)].$

which captures correlations between input and the outputs of the generator in
the RKHS. Defining $w\phi\colon\mathcal{X}\to\mathcal{H}^{p+1}$ by
$w\phi(x)=[\sqrt{\mu}\phi(x),d_{1}\phi(x),\ldots
d_{p}\phi(x)]^{\top}\in\mathbb{R}^{d+1}$, we can consider

$\displaystyle W_{\mu}$ $\displaystyle=Z_{\mu}^{*}Z_{\mu}=S_{\pi}^{*}(\mu
I{-}L)S_{\pi}$
$\displaystyle=\mathbb{E}_{x\sim\pi}[\mu\phi(x)\otimes\phi(x)+d\phi(x)\otimes
d\phi(x)]=\mathbb{E}_{x\sim\pi}[w\phi(x)\otimes w\phi(x)]$

which corresponds to the RKHS covariance operator w.r.t. energy space
$\mathcal{W}$.
Examples. One way to ensure that the essential assumption ((KE)) holds is to
show that $w\phi$ fulfils the boundedness condition ((BK)), i.e.
$w\phi\\!\in\\!\mathcal{L}^{\infty}_{\pi}(\mathcal{X},\mathcal{H}^{1\\!+\\!p})$.
Let’s show that the property holds true for the Damped Langevin (see Example
1) and the CIR process (see Example 2) if we consider the Radial Basis
Function (RBF) kernel $k(x,y)=k_{x}(y)=\exp(-\kappa\lVert x-y\rVert^{2})$
where $\kappa>0$ is a free parameter that sets the “spread” of the kernel. As
a reminder, for every $x\in\mathcal{X}$, we have

$\lVert w\phi(x)\rVert^{2}=\mu\lVert\phi(x)\rVert^{2}+\sum_{k=1}^{p}\lVert
d_{k}\phi(x)\rVert^{2}\quad\text{with}\quad d_{k}=s_{k}(x)^{\top}D\phi(x).$

Recalling that for the overdamped Langevin process
$s(x)=(k_{b}T)^{1/2}(\delta_{ij})_{i\in[d],j\in[p]}$, $ss^{\top}$ is diagonal
so that $\langle s(x)^{\top}Dk_{x}(\cdot),s(x)^{\top}Dk_{x}(\cdot)\rangle=0$
for every $x\in\mathcal{X}$. As $\lVert\phi(x)\rVert^{2}=1$, we get that for
every $x\in\mathcal{X}$, $\lVert w\phi(x)\rVert^{2}\leq\mu=:c_{\mathcal{H}}$.
Consider now the CIR process. We have $d=p=1$ and
$s(x)=\sigma\sqrt{x}/\sqrt{2}$ for any $x\in\mathcal{X}$. For the very same
reasons, $\langle s(x)^{\top}Dk_{x}(\cdot),s(x)^{\top}Dk_{x}(\cdot)\rangle=0$
for every $x\in\mathcal{X}$, so that $\lVert
w\phi(x)\rVert^{2}\leq\mu=:c_{\mathcal{H}}$.
In both Langevin and CIR cases, we have
$w\phi\\!\in\\!\mathcal{L}^{\infty}_{\pi}(\mathcal{X},\mathcal{H}^{1\\!+\\!p})$
when considering an RBF kernel.

##  Appendix 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\in{\rm{B}}_{r}({\mathcal{H}})$ that is the
basis of the statistical bounds that follow. This result is a direct
consequence of [21] and Davis-Khan spectral perturbation result for compact
self-adjoint operators, [10].

In the framework of Koopman operator learning [22], spectral bounds are
expressed in terms of a distortion metric between the RKHS $\mathcal{H}$ and
$\mathcal{L}^{2}_{\pi}(\mathcal{X})$, corresponding to the cost incurred from
observing the operator’s action on the $\mathcal{H}$ rather than on its domain
$\mathcal{L}^{2}_{\pi}(\mathcal{X})$. Aligned with the risk definition (9),
here we measure in a certain way the distortion between the $\mathcal{H}$ and
$\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$ as given in definitions in (10).

See 2

###### Proof.

We first remark that

$\|((\mu
I\,{-}\,L)^{-1}-\widehat{\lambda}_{i}\,I_{\mathcal{L}^{2}_{\pi}(\mathcal{X})})^{-1}\|_{\mathcal{H}{\to}\mathcal{W}}^{-1}\leq\|((\mu
I\,{-}\,L)^{-1}Z_{\mu}-Z_{\mu}\widehat{G})\widehat{h}_{i}\|_{\mathcal{H}{\to}\mathcal{W}}/\lVert
Z_{\mu}\widehat{h}_{i}\rVert\leq\mathcal{E}(\widehat{G})\eta(\widehat{h}_{i}).$

Then, from the first inequality, using that $(\mu I\,{-}\,L)^{-1}$ is self-
adjoint as operator
$\mathcal{W}^{\mu}_{\pi}(\mathcal{X})\to\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$,
we obtain the first bound in (11).

So, observing that for every
$(\mu{-}\lambda)^{-1}{\in}\operatorname*{\text{\rm Sp}}((\mu
I\,{-}\,L)^{-1}){\setminus}\\{(\mu{-}\lambda_{i})^{-1}\\}$,

$\lvert(\mu{-}\widehat{\lambda}_{i})^{-1}{-}(\mu{-}\lambda)^{-1}\rvert{\geq}\lvert(\mu{-}\lambda_{i})^{-1}{-}(\mu{-}\lambda)^{-1}\rvert{-}\lvert(\mu{-}\widehat{\lambda}_{i})^{-1}{-}(\mu{-}\lambda_{i})^{-1}\rvert$

we conclude that
$\lvert(\mu{-}\widehat{\lambda}_{i})^{-1}{-}(\mu{-}\lambda)^{-1}\rvert{\geq}\lvert(\mu{-}\lambda_{i})^{-1}{-}(\mu{-}\lambda)^{-1}\rvert{-}\mathcal{E}(\widehat{G})\,\eta(\widehat{h}_{i})$,
and

$\min\\{\lvert(\mu{-}\widehat{\lambda}_{i})^{-1}-(\mu{-}\lambda)^{-1}\rvert\,|\,(\mu{-}\lambda)^{-1}\in\operatorname*{\text{\rm
Sp}}((\mu
I\,{-}\,L)^{-1})\setminus\\{(\mu{-}\lambda_{i})^{-1}\\}\\}\geq\text{\rm
gap}_{i}{-}\mathcal{E}(\widehat{G})\,\eta(\widehat{h}_{i}).$

So, applying Proposition 3, we obtain

$\displaystyle\sin(\sphericalangle(\widehat{f}_{i},f_{i}))$
$\displaystyle\leq\frac{\lVert(\mu
I\,{-}\,L)^{-1}\widehat{f}_{i}-(\mu{-}\widehat{\lambda}_{i})^{-1}\,\widehat{f}_{i}\rVert}{[\text{\rm
gap}_{i}-\mathcal{E}(\widehat{G})\,\eta(\widehat{h}_{i})]_{+}}\leq\frac{\lVert((\mu
I\,{-}\,L)^{-1}Z_{\mu}-Z_{\mu}\widehat{G})\widehat{h}_{i}\rVert/\lVert
Z_{\mu}\widehat{h}_{i}\rVert}{[\text{\rm
gap}_{i}-\mathcal{E}(\widehat{G})\,\eta(\widehat{h}_{i})]_{+}}$
$\displaystyle\leq\frac{\mathcal{E}(\widehat{G})\,\eta(\widehat{h}_{i})}{[\text{\rm
gap}_{i}-\mathcal{E}(\widehat{G})\,\eta(\widehat{h}_{i})]_{+}}.$

Since, clearly $\lVert\widehat{f}_{i}-f_{i}\rVert^{2}\leq
2(1-\cos(\sphericalangle(\widehat{f}_{i},f_{i}))\leq
2\sin(\sphericalangle(\widehat{f}_{i},f_{i}))$, the proof of the second bound
is completed. ∎

Next, we adapt the [22, Proposition 1] to our setting as follows.

######  Proposition 5.

Let $\widehat{G}\,{\in}\,{\rm{B}}_{r}({\mathcal{H}})$. For all $i\in[r]$ the
metric distortion of $\widehat{h}_{i}$ w.r.t. energy space
$\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$ can be tightly bounded as

$1\,\,/\sqrt{\lVert
W_{\mu}\rVert}\,\leq\,\eta(\widehat{h}_{i})\,\leq\,\lVert\widehat{G}\rVert\,/\,\sigma_{\min}^{+}(Z_{\mu}\widehat{G}).$
(34)

##  Appendix D Empirical estimation

The _Reduced Rank Regression_ (RRR) estimator is the exact minimizer of (17)
under fixed rank constraint. Specifically, RRR is the minimizer
$\widehat{G}^{r}_{\mu,\gamma}$ of $\widehat{\mathcal{R}}_{\gamma}(G)$ within
the set of bounded operators $\mathrm{HS}_{r}({\mathcal{H}})$ on $\mathcal{H}$
that have rank at most $r$. The regularization term $\gamma\lVert
G\rVert_{\rm{HS}}^{2}$ is added to ensure stability. The closed form solution
of the empirical RRR estimator is

$\widehat{G}^{r}_{\mu,\gamma}=\widehat{H}_{\gamma}^{-1/2}[\\![\widehat{H}_{\gamma}^{1/2}C]\\!]_{r},$
(35)

while its population counterpart is given by
$G^{r}_{\mu,\gamma}=W_{\mu,\gamma}^{-1/2}[\\![W_{\mu,\gamma}^{1/2}C]]\\!]_{r}$.

In order to prove Theorem 1, recall kernel matrices in (19) and define

$\textsc{F}_{\mu}=\left[\begin{matrix}\mu\textsc{K}&\sqrt{\mu}\textsc{N}\\\
\sqrt{\mu}\textsc{N}^{\top}&\textsc{M}\end{matrix}\right]\quad\text{ and
}\textsc{F}_{\mu,\gamma}=\left[\begin{matrix}\mu\textsc{K}_{\gamma}&\sqrt{\mu}\textsc{N}\\\
\sqrt{\mu}\textsc{N}^{\top}&\textsc{M}+\gamma\mu I\end{matrix}\right].$ (36)

Now we provide the explicit forms of the matrices N and M in the case of
Langevin (see Example 1) and CIR (see Example 2) processes, considering an RBF
kernel $k(x,y)=k_{x}(y)=\exp(-\kappa\lVert x-y\rVert^{2})$. As a reminder (see
(36)), $\textsc{N}\\!\in\\!\mathbb{R}^{n\times pn}$ and
$\textsc{M}\\!\in\\!\mathbb{R}^{pn\times pn}$ are Gram matrices whose
elements, for $k\in[1+p],i,j\in[n]$ are given by

$\textsc{N}_{i,(k-1)n+j}\\!=\\!n^{-1}{\langle\phi(x_{i}),d_{k}\phi(x_{j})\rangle}_{\mathcal{H}}\text{
and
}\textsc{M}_{(k\\!-\\!1)n\\!+\\!i,(\ell\\!-\\!1)n\\!+\\!j}\\!=\\!n^{-1}{\langle
d_{k}\phi(x_{i}),d_{\ell}\phi(x_{j})\rangle}_{\mathcal{H}},$

where $d_{k}\phi(x_{j})=s_{k}(x_{j})^{\top}D\phi(x_{j})$ and
$D\phi(x_{j})=Dk(x_{j},\cdot)$ is defined by (31). For $k\in[p],i,j\in[n]$, we
have

$\displaystyle{\langle\phi(x_{i}),d_{k}\phi(x_{j})\rangle}_{\mathcal{H}}$
$\displaystyle={\langle
k_{x_{i}},s_{k}(x_{j})^{\top}D{k_{x_{j}}}\rangle}_{\mathcal{H}}$
$\displaystyle=\Big{\langle}k_{x_{i}},s_{k}(x_{j})^{\top}\Big{(}\underset{h\to
0}{\lim}\frac{k(\cdot,x_{j}+he_{l})-k(\cdot,x_{j})}{h}\Big{)}_{l\in[d]}\Big{\rangle}$
$\displaystyle=\underset{h\to
0}{\lim}\Big{\langle}k_{x_{i}},s_{k}(x_{j})^{\top}\Big{(}\frac{k(\cdot,x_{j}+he_{l})-k(\cdot,x_{j})}{h}\Big{)}_{l\in[d]}\Big{\rangle}$
$\displaystyle=s_{k}(x_{j})^{\top}\Big{(}\underset{h\to
0}{\lim}\frac{k(x_{i},x_{j}+he_{l})-k(x_{i},x_{j})}{h}\Big{)}_{l\in[d]}$
$\displaystyle=s_{k}(x_{j})^{\top}\big{(}D^{(0,e_{l})}k(x_{i},x_{j})\big{)}_{l\in[d]}$
$\displaystyle=2\gamma
s_{k}(x_{j})^{\top}\big{(}(x_{i}^{(l)}-x_{j}^{(l)})k(x_{i},x_{j})\big{)}_{l\in[d]}$

where we have used the continuity of the inner product to get the third line
and the reproducing property to obtain the following one. Similarly, for
$k,\ell\in[p],j\in[n]$, we get

$\displaystyle{\langle
d_{k}\phi(x_{i}),d_{\ell}\phi(x_{j})\rangle}_{\mathcal{H}}$
$\displaystyle{=}\Big{\langle}s_{k}(x_{i})^{\top}\big{(}D^{(e_{l^{\prime}},0)}k(x_{i},\cdot)\big{)}_{l^{\prime}\in[d]},s_{\ell}(x_{j})^{\top}\Big{(}\underset{h\to
0}{\lim}\frac{k(\cdot,x_{j}+he_{l})-k(\cdot,x_{j})}{h}\Big{)}_{l\in[d]}\Big{\rangle}$
$\displaystyle{=}\underset{h\to
0}{\lim}\Big{\langle}s_{k}(x_{i})^{\top}\big{(}D^{(e_{l^{\prime}},0)}k(x_{i},\cdot)\big{)}_{l^{\prime}\in[d]},s_{\ell}(x_{j})^{\top}\Big{(}\frac{k(\cdot,x_{j}+he_{l})-k(\cdot,x_{j})}{h}\Big{)}_{l\in[d]}\Big{\rangle}$
$\displaystyle{=}\underset{h\to
0}{\lim}\Big{\langle}\big{(}D^{(e_{l^{\prime}},0)}k(x_{i},\cdot)\big{)}_{l^{\prime}\in[d]},s_{k}(x_{j})s_{\ell}(x_{j})^{\top}\Big{(}\frac{k(\cdot,x_{j}+he_{l})-k(\cdot,x_{j})}{h}\Big{)}_{l\in[d]}\Big{\rangle}$
$\displaystyle{=}s_{k}(x_{i})s_{\ell}(x_{j})^{\top}\Big{(}\underset{h\to
0}{\lim}\frac{D^{(e_{l^{\prime}},0)}k(x_{i},x_{j}+he_{l})-D^{(e_{l^{\prime}},0)}k(x_{i},x_{j})}{h}\Big{)}_{l^{\prime},l\in[d]}$
$\displaystyle{=}s_{k}(x_{i})s_{\ell}(x_{j})^{\top}(\mathfrak{D}_{l^{\prime},l}k(x_{i},x_{j}))_{l^{\prime},l\in[d]},$

where we have used the partial derivative reproducing (33) property and where
we define for $l^{\prime}\neq l$,

$\displaystyle\mathfrak{D}_{l^{\prime},l}k(x_{i},x_{j})$
$\displaystyle=\underset{h\to
0}{\lim}\frac{D^{(e_{l^{\prime}},0)}k(x_{i},x_{j}+he_{l})-D^{(e_{l^{\prime}},0)}k(x_{i},x_{j})}{h}$
$\displaystyle=-2\gamma\underset{h\to
0}{\lim}\frac{(x_{i}^{(l^{\prime})}-x_{j}^{(l^{\prime})})k(x_{i},x_{j}+he_{l})-(x_{i}^{(l^{\prime})}-x_{j}^{(l^{\prime})})k(x_{i},x_{j})}{h}$
$\displaystyle=-4\gamma^{2}(x_{i}^{(l^{\prime})}-x_{j}^{(l^{\prime})})(x_{i}^{(l)}-x_{j}^{(l)})k(x_{i},x_{j}),$
(37)

and for $l\in[d]$,

$\displaystyle\mathfrak{D}_{l,l}k(x_{i},x_{j})$ $\displaystyle=\underset{h\to
0}{\lim}\frac{D^{(e_{l},0)}k(x_{i},x_{j}+he_{l})-D^{(e_{l},0)}k(x_{i},x_{j})}{h}$
$\displaystyle=-2\gamma\underset{h\to
0}{\lim}\frac{(x_{i}^{(l)}-x_{j}^{(l)}-h)k(x_{i},x_{j}+he_{l})-(x_{i}^{(l)}-x_{j}^{(l)})k(x_{i},x_{j})}{h}$
$\displaystyle=\big{[}2\gamma\
-4\gamma^{2}(x_{i}^{(l)}-x_{j}^{(l)})^{2}\big{]}k(x_{i},x_{j})=2\gamma[1-2\gamma(x_{i}^{(l)}-x_{j}^{(l)})^{2}]k(x_{i},x_{j}).$
(38)

For the overdamped Langevin (see Example 1), for $k\leq d$,
$s_{k}(x_{i})=(k_{b}T)^{1/2}e_{k}$ and
$s_{k}(x_{i})s_{\ell}(x_{j})^{\top}=(k_{b}T)I$ so that

$\textsc{N}_{i,(k-1)n+j}=n^{-1}{\langle\phi(x_{i}),d_{k}\phi(x_{j})\rangle}_{\mathcal{H}}=2\gamma(k_{b}T)^{1/2}n^{-1}(x_{i}^{(k)}-x_{j}^{(k)})k(x_{i},x_{j})$

and

$\textsc{M}_{(k\\!-\\!1)n\\!+\\!i,(\ell\\!-\\!1)n\\!+\\!j}\\!=\\!n^{-1}{\langle
d_{k}\phi(x_{j}),d_{\ell}\phi(x_{j})\rangle}_{\mathcal{H}}=(k_{b}T)n^{-1}(\mathfrak{D}_{l^{\prime},l}k(x_{j},x_{j}))_{l^{\prime},l\in[d]},$

where the elements of $\mathfrak{D}$ are given by (D) and (D). For the CIR
process (see Example 2) in dimension $d=1$, we have
$s(x)=\sigma\sqrt{x}/\sqrt{2}$. Then,

$\textsc{N}_{i,(k-1)n+j}=n^{-1}{\langle\phi(x_{i}),d_{k}\phi(x_{j})\rangle}_{\mathcal{H}}=\sqrt{2}\sigma\gamma
n^{-1}\sqrt{x_{j}}(x_{i}-x_{j})\exp(-\gamma|x_{i}-x_{j}|^{2})$

and

$\textsc{M}_{(k\\!-\\!1)n\\!+\\!i,(\ell\\!-\\!1)n\\!+\\!j}\\!=\\!n^{-1}{\langle
d_{k}\phi(x_{i}),d_{\ell}\phi(x_{j})\rangle}_{\mathcal{H}}=\sigma^{2}\gamma
n^{-1}\sqrt{x_{i}}\sqrt{x_{j}}\exp(-\gamma|x_{i}-x_{j}|^{2}).$

Based on the previous formulas, using Theorem 1, which we prove next, one can
estimate generator’s eigenpairs in practice.

See 1

###### Proof.

First, note that $\mu\,\textsc{J}_{\mu,\gamma}\succ 0$ is exactly Schurs’s
complement w.r.t. second diagonal block of $\textsc{F}_{\mu,\gamma}\succ 0$,
and that, due to block inversion lemma [14], we have that

$\textsc{F}_{\mu,\gamma}^{-1}=\left[\begin{matrix}\mu^{-1}\textsc{J}_{\mu,\gamma}&{\mu}^{-1/2}\textsc{J}_{\mu,\gamma}^{-1}\textsc{N}(\textsc{M}+\gamma\mu
I)^{-1}\\\ {\mu}^{-1/}(\textsc{M}+\gamma\mu
I)^{-1}\textsc{N}^{\top}\textsc{J}_{\mu,\gamma}^{-1}&\textsc{A}\end{matrix}\right].$
(39)

where A is some $np\times np$ matrix. The first step in computing the RRR
estimator lies in computing a truncating SVD of
$\widehat{W}_{\mu,\gamma}^{-1/2}\widehat{C}$, that is
$\widehat{C}\widehat{W}_{\mu,\gamma}^{-1}\widehat{C}q_{i}=\widehat{\sigma}_{i}^{2}q_{i}$,
$i\in[r]$. Now, using the low-rank eigenvalue problem formulation [14], we
have that $q_{i}=\widehat{S}^{*}v_{i}$ and
$\widehat{S}\widehat{W}_{\mu,\gamma}^{-1}\widehat{C}\widehat{S}^{*}v_{i}=\widehat{\sigma}_{i}^{2}v_{i}$.
Now, recalling that $\widehat{S}=[\mu^{-1/2}\,|\,0]\widehat{Z}_{\mu}$ we
obtain and that
$\widehat{Z}_{\mu}(\widehat{Z}_{\mu}^{*}\widehat{Z}_{\mu}+\mu\gamma
I)^{-1}{=}(\widehat{Z}_{\mu}\widehat{Z}_{\mu}^{*}+\mu\gamma
I)^{-1}\widehat{Z}_{\mu}$, we obtain

$[\mu^{-1/2}I\;|\;0]\textsc{F}_{\mu,\gamma}^{-1}\textsc{F}_{\mu}[\mu^{-1/2}I\;|\;0]^{\top}\textsc{K}v_{i}=\widehat{\sigma}_{i}^{2}v_{i},$

which after some algebra, using (39)

$\mu^{-1}(I-\gamma\textsc{J}_{\mu,\gamma}^{-1})\textsc{K}v_{i}=\widehat{\sigma}_{i}^{2}v_{i},$

i.e.
$\mu^{-1}(I-\gamma\textsc{J}_{\mu,\gamma}^{-1})\textsc{K}\textsc{V}_{r}=\textsc{V}_{r}\Sigma_{r}^{2}$.

By normalizing right singular value functions $q_{i}$ of
$\widehat{W}_{\mu,\gamma}^{-1/2}\widehat{C}$, that is by asking that ${\langle
q_{i},q_{i}\rangle}_{\mathcal{H}}=v_{i}^{\top}\textsc{K}v_{i}=1$, we obtain
that
$[\\![\widehat{W}_{\mu,\gamma}^{-1/2}\widehat{C}]\\!]_{r}=\widehat{W}_{\mu,\gamma}^{-1/2}\widehat{C}Q_{r}Q^{*}_{r}$,
for $Q_{r}=[q_{1}\;|\ldots|\;q_{r}]$. In other words, we have

$\widehat{G}_{\mu,\gamma}=\widehat{W}_{\mu,\gamma}^{-1}\widehat{Z}_{\mu}^{*}[\mu^{-1/2}I\;|\;0]^{\top}\widehat{S}\widehat{S}^{*}\textsc{V}_{r}\textsc{V}_{r}^{\top}\widehat{S}=\widehat{Z}_{\mu}^{*}\textsc{F}_{\mu,\gamma}^{-1}[\mu^{-1/2}I\;|\;0]^{\top}\textsc{K}\textsc{V}_{r}\textsc{V}_{r}^{\top}\widehat{S}=\widehat{Z}_{\mu}^{*}\textsc{U}_{r}\textsc{V}_{r}^{\top}\widehat{S},$

where we used that
$\textsc{J}_{\mu,\gamma}^{-1}\textsc{K}\textsc{V}_{r}=\gamma^{-1}[\textsc{K}\textsc{V}_{r}-\mu\textsc{V}_{r}\Sigma_{r}^{2}]$.
Once with this form, we apply [21, Theorem 2] to obtain the result. ∎

Finally, next result provides the reasoning for using empirical metric
distortion of Theorem 23.

######  Proposition 6.

Under the assumptions of Theorem 1, for every $i\in[r]$

$\widehat{\eta}_{i}=\frac{\lVert\widehat{h}_{i}\rVert}{\lVert\widehat{Z}_{\mu}\widehat{h}_{i}\rVert}=\sqrt{\frac{(w^{r}_{i})^{*}U_{r}^{\top}\textsc{F}_{\mu}U_{r}w^{r}_{i}}{\lVert\textsc{F}_{\mu}U_{r}w^{r}_{i}\rVert^{2}}},$
(40)

and

$\left|\widehat{\eta}_{i}-\eta(\widehat{h}_{i})\right|\leq\left(\eta(\widehat{h}_{i})\;\wedge\widehat{\eta}_{i}\right)\,\eta(\widehat{h}_{i})\,\widehat{\eta}_{i}\,\lVert\widehat{W}_{\mu}-W_{\mu}\rVert.$
(41)

###### Proof.

First, note that (40) follows directly from Theorem 1. Next, since for every
$i\in[r]$,

$(\widehat{\eta}_{i})^{-2}-(\eta(\widehat{h}_{i}))^{-2}=\frac{{\langle\widehat{h}_{i},(\widehat{W}_{\mu}-W_{\mu})\widehat{h}_{i})\rangle}}{\lVert\widehat{h}_{i}\rVert^{2}}\leq\lVert\widehat{W}_{\mu}-W_{\mu}\rVert,$

we obtain

$\left|\widehat{\eta}_{i}^{-1}-(\eta(\widehat{h}_{i}))^{-1}\right|\leq\frac{\left|\widehat{\eta}_{i}^{-2}-(\eta(\widehat{h}_{i}))^{-2}\right|}{(\eta(\widehat{h}_{i}))^{-1}\;\vee\;\widehat{\eta}_{i}^{-1}}\leq\left(\eta(\widehat{h}_{i})\;\wedge\widehat{\eta}_{i}\right)\lVert\widehat{W}_{\mu}-W_{\mu}\rVert.$

∎

##  Appendix E Learning bounds

###  E.1 Main assumptions

We start by observing that due to (BK), we have that operator
$Z_{\mu}\in{\rm{HS}}\left(\mathcal{H},\mathcal{W}^{\mu}_{\pi}(\mathcal{X})\right)$,
which according to the spectral theorem for positive self-adjoint operators,
has an SVD, i.e. there exists at most countable positive sequence
$(\sigma_{j})_{j\in J}$, where $J=\\{1,2,\ldots,\\}\subseteq\mathbb{N}$, and
ortho-normal systems $(z_{j})_{j\in J}$ and $(h_{j})_{j\in J}$ of
$\operatorname*{\text{\rm cl}}(\operatorname*{\text{\rm Im}}(Z_{\mu}))$ and
$\operatorname{\text{\rm Ker}}(Z_{\mu})^{\perp}$, respectively, such that
$Z_{\mu}h_{j}=\sigma_{j}z_{j}$ and $Z_{\mu}^{*}z_{j}=\sigma_{j}h_{j}$, $j\in
J$.

Now, given $\alpha\geq 0$, let us define scaled injection operator
$Z_{\mu,\alpha}\colon\mathcal{H}\to\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$ as

$Z_{\mu,\alpha}=\sum_{j\in J}\sigma_{j}^{\alpha}z_{j}\otimes h_{j}.$ (42)

Clearly, we have that $Z_{\mu}=Z_{\mu,1}$, while $\operatorname*{\text{\rm
Im}}{Z_{\mu,0}}=\operatorname*{\text{\rm cl}}(\operatorname*{\text{\rm
Im}}(Z_{\mu}))$. Next, we equip $\operatorname*{\text{\rm
Im}}(Z_{\mu,\alpha})$ with a norm $\lVert\cdot\rVert_{\mathcal{H},\alpha}$ to
build an interpolation space

$[\mathcal{H}]_{\alpha}=\left\\{f\in\operatorname*{\text{\rm
Im}}(Z_{\mu,\alpha})\;|\;\lVert f\rVert_{\mathcal{H},\alpha}^{2}=\sum_{j\in
J}\sigma_{j}^{-2\alpha}{\langle
f,z_{j}\rangle}_{\mathcal{W}}^{2}<\infty\right\\},$

noting that the inner product in $\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$ is
given by bilinear energy functional

${\langle f,g\rangle}_{\mathcal{W}}=\mu{\langle
f,g\rangle}_{\mathcal{L}^{2}_{\pi}}-{\langle
f,Lg\rangle}_{\mathcal{L}^{2}_{\pi}}.$

We remark that for $\alpha=1$ the space $[\mathcal{H}]_{\alpha}$ is just an
RKHS $\mathcal{H}$ seen as a subspace of
$\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$. Moreover, we have the following
injections

$[\mathcal{H}]_{\alpha_{1}}\hookrightarrow[\mathcal{H}]_{1}\hookrightarrow[\mathcal{H}]_{\alpha_{2}}\hookrightarrow[\mathcal{H}]_{0}=\mathcal{W}^{\mu}_{\pi}(\mathcal{X}),$

where $\alpha_{1}\geq 1\geq\alpha_{2}\geq 0$.

In addition, from (BK) we also have that RKHS $\mathcal{H}$ can be embedded
into

$W^{\mu,\infty}_{\pi}(\mathcal{X})=\\{f\in\mathcal{W}^{\mu}_{\pi}(\mathcal{X})\,|\,\lVert
f\rVert_{W^{\mu,\infty}_{\pi}}=\operatorname*{\text{\rm
ess\,sup}}_{x\sim\pi}[\lvert f(x)\rvert^{2}-f(x)[Lf](x)]<\infty\\}$

that is, for some $\tau\in(0,1]$

$[\mathcal{H}]_{1}\hookrightarrow[\mathcal{H}]_{\tau}\hookrightarrow
W^{\mu,\infty}_{\pi}(\mathcal{X})\hookrightarrow\mathcal{W}^{\mu}_{\pi}(\mathcal{X}).$

Now, according to [12], if
$Z_{\mu,\tau,\infty}\colon[\mathcal{H}]_{\tau}\hookrightarrow
L^{\infty}_{\pi}(\mathcal{X})$ denotes the injection operator, its boundedness
implies the polynomial decay of the singular values of $Z_{\mu}$, i.e.
$\sigma_{j}^{2}(Z_{\mu})\lesssim j^{-1/\tau}$, $j\in J$, and the condition
(KE) is assured.

Assumption (SD) allows one to quantify the effective dimension of
$\mathcal{H}$ in ambient space $\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$, while
the kernel embedding property (KE) allows one to estimate the norms of
whitened feature maps, in our generator setting vector-valued since they
define rank($1{+}p$) operators on $\mathcal{H}$,

$\xi(x):=W_{\mu,\gamma}^{-1/2}w\phi(x)\in\mathcal{H}^{1{+}p}.$ (43)

This object plays a key role in deriving the learning rates for regression
problems, [23] and the following result is bounding it.

######  Lemma 1.

Let (KE) hold for some $\tau\in[\beta,1]$ and $c_{\tau}\in(0,\infty)$. Then,

$\mathbb{E}_{x\sim\pi}\lVert\xi(x)\rVert_{\mathcal{H}^{1+p}}^{2}\leq\begin{cases}\frac{c_{\beta}^{\beta}}{1-\beta}(\mu\gamma)^{-\beta}&,\beta<1,\\\
c_{\tau}\,(\mu\gamma)^{-1}&,\beta=1.\end{cases}\text{ and
}\;\lVert\xi\rVert_{\infty}^{2}=\operatorname*{\text{\rm
ess\,sup}}_{x\sim\pi}\lVert\xi(x)\rVert_{\mathcal{H}^{1+p}}^{2}\leq{c_{\tau}}(\mu\gamma)^{-\tau}.$
(44)

###### Proof.

W.l.o.g. set $\mu=1$, observing that the only change in the proof is in
scaling $\gamma>0$. We first observe that for every $j\in J$ from definition
of $w\phi$ and fact that $h_{j}(x)=[Z_{\mu}h_{j}](x)$ $\pi$-a.e., it holds
that

$\sum_{i\in[1+p]}{\langle w_{i}\phi(x),h_{j}\rangle}^{2}{=}\mu\lvert
h_{j}(x)\rvert^{2}-h_{j}(x)[LZ_{\mu}h_{j}](x){=}\mu\lvert[Z_{\mu}h_{j}](x)\rvert^{2}{-}[Z_{\mu}h_{j}](x)[LZ_{\mu}h_{j}](x),\;\pi\text{-a.e.},$

implying that $\sum_{i\in[1+p]}{\langle
w_{i}\phi(x),h_{j}\rangle}^{2}\leq\sigma_{j}\mu\lvert
z_{j}(x)\rvert^{2}-z_{j}(x)[Lz_{j}](x)$. So, for every $\tau>0$

$\displaystyle\lVert\xi(x)\rVert_{\mathcal{H}^{1{+}p}}^{2}$
$\displaystyle{=}\sum_{j\in J}\sum_{i\in[1{+}p]}{\langle
W_{\mu,\gamma}^{-1/2}w_{i}\phi(x),h_{j}\rangle}^{2}=\sum_{j\in
J}\sum_{i\in[1{+}p]}\frac{1}{\sigma_{j}^{2}+\gamma}{\langle
w_{i}\phi(x),h_{j}\rangle}^{2}$ $\displaystyle{=}\sum_{j\in
J}\sum_{i\in[1{+}p]}\frac{\sigma_{j}^{2(1-\tau)}}{\sigma_{j}^{2}+\gamma}\frac{{\langle
w_{i}\phi(x),h_{j}\rangle}^{2}}{\sigma_{j}^{2}}\sigma_{j}^{2\tau}{=}\gamma^{-\tau}\sum_{j\in
J}\sum_{i\in[1{+}p]}\frac{(\sigma_{j}^{2}\gamma^{-1})^{1-\tau}}{\sigma_{j}^{2}\gamma^{-1}+1}\frac{{\langle
w_{i}\phi(x),h_{j}\rangle}^{2}}{\sigma_{j}^{2}}\sigma_{j}^{2\tau}$
$\displaystyle{\leq}\gamma^{-\tau}\sum_{j\in J}\frac{\mu\lvert
h_{j}(x)\rvert^{2}-h_{j}(x)[LZ_{\mu}h_{j}](x)}{\sigma_{j}^{2}}\sigma_{j}^{2\tau}{=}\gamma^{-\tau}\sum_{j\in
J}(\mu\lvert z_{j}(x)\rvert^{2}-z_{j}(x)[Lz_{j}](x))\sigma_{j}^{2\tau},$

and, due to (21), we obtain
$\lVert\xi\rVert_{\infty}^{2}\leq\gamma^{-\tau}c_{\tau}$. On the other hand,
we also have that

$\operatorname{\text{\rm
tr}}(\mathbb{E}_{x\sim\pi}[\xi(x)\otimes\xi(x)])=\operatorname{\text{\rm
tr}}(W_{\mu,\gamma}^{-1/2}W_{\mu}W_{\mu,\gamma}^{-1/2})=\operatorname{\text{\rm
tr}}(W_{\mu,\gamma}^{-1}W_{\mu}),$

which is an effective dimension of the RKHS $\mathcal{H}$ in
$\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$. Therefore, following the proof of
Fischer and Steinwart, [12, Lemma 11] for classical covariances in
$\mathcal{L}^{2}_{\pi}$, we show that the bound on the effective dimension is

$\operatorname{\text{\rm tr}}(W_{\mu,\gamma}^{-1}W_{\mu})=\sum_{j\in
N}\frac{\sigma_{j}^{2}}{\sigma_{j}^{2}+\gamma}\leq\begin{cases}\frac{c_{\beta}^{\beta}}{1-\beta}\gamma^{-\beta}&,\beta<1,\\\
c_{\tau}\,\gamma^{-1}&,\beta=1.\end{cases}$ (45)

For the case $\beta=1$, it suffices to see that

$\|z_{j}\|_{\mathcal{W}}=\mathfrak{E}_{\mu}[z_{j}]\leq\lVert
z_{j}\rVert_{W^{\mu,\infty}_{\pi}}=\operatorname*{\text{\rm
ess\,sup}}_{x\sim\pi}[\lvert z_{j}(x)\rvert^{2}-z_{j}(x)[LZ_{\mu}z_{j}](x)],$

and, hence

$\operatorname{\text{\rm
tr}}(W_{\mu,\gamma}^{-1}W_{\mu}){\leq}\gamma^{-1}\sum_{j\in
N}\sigma_{j}^{2}\|z_{j}\|_{\mathcal{W}}^{2}{=}\gamma^{-1}\sum_{j\in
N}\sigma_{j}^{2}\mathfrak{E}_{\mu}[z_{j}]{\leq}\gamma^{-1}c_{\tau}.$

For $\beta<1$ we can apply the same classical reasoning as in the proof of
Proposition 3 of [8]. ∎

######  Proposition 7.

If the eigenfunctions of $L$ belong to $[\mathcal{H}]_{\alpha}$, then
Condition (RC) is satisfied.

###### Proof.

Note first that the resolvent $(\mu I\,{-}\,L)^{-1}$ admits same
eigenfunctions as the generator $L$, meaning that $\operatorname*{\text{\rm
Im}}((\mu I\,{-}\,L)^{-1}Z_{\mu})\subseteq\operatorname*{\text{\rm
cl}}(\operatorname*{\text{\rm Im}}(Z_{\mu,\alpha}))$. But according to [47,
Theorem 2.2], this last condition is equivalent to Condition (RC). ∎

Finally we prove here Proposition 1

See 1

###### Proof.

Recall first that since
$Z_{\mu}\in{\rm{HS}}\left(\mathcal{H},\mathcal{W}^{\mu}_{\pi}(\mathcal{X})\right)$,
according to the spectral theorem for positive self-adjoint operators,
$Z_{\mu}$ admits an SVD. Its form is provided in (42) taking $\alpha=1$.

Now, recalling that $[\\![\cdot]\\!]_{r}$ denotes the $r$-truncated SVD, i.e.
$[\\![Z_{\mu}]\\!]_{r}=\sum_{j\in[r]}\sigma_{j}z_{j}\otimes h_{j}$, since
$\lVert
Z_{\mu}-[\\![Z_{\mu}]\\!]_{r}\rVert_{\rm{HS}}^{2}=\sum_{j>r}\sigma_{j}^{2}$,
for every $\delta>0$ there exists $r\in\mathbb{N}$ such that $\lVert
Z_{\mu}-[\\![Z_{\mu}]\\!]_{r}\rVert_{\rm{HS}}<\mu\delta/3$. Consequently since
all the eigenvalues of $L$ are non-positive, $\lVert(\mu
I\,{-}\,L)^{-1}(Z_{\mu}-[\\![Z_{\mu}]\\!]_{r})\rVert_{\rm{HS}}\leq\lVert
Z_{\mu}-[\\![Z_{\mu}]\\!]_{r}\rVert_{\rm{HS}}/\mu\leq\delta/3$. Next since
$\operatorname*{\text{\rm Im}}(P_{\mathcal{H}}(\mu
I\,{-}\,L)^{-1}Z_{\mu})\subseteq\operatorname*{\text{\rm
cl}}(\operatorname*{\text{\rm Im}}(Z_{\mu}))$, for any $j\in[r]$, there exists
$g_{j}\in\mathcal{H}$ s.t. $\lVert P_{\mathcal{H}}(\mu
I\,{-}\,L)^{-1}z_{j}-Z_{\mu}{g}_{j}\rVert\leq\frac{\delta}{3r}$, and, denoting
$B_{r}:=\sum_{j\in[r]}\sigma_{j}{g}_{j}\otimes{h}_{j}$ we conclude $\lVert
P_{\mathcal{H}}(\mu
I\,{-}\,L)^{-1}[\\![Z_{\mu}]\\!]_{r}-Z_{\mu}B_{r}\rVert_{\rm{HS}}\leq\delta/3$.
Finally we recall that the set of non-defective matrices is dense in the space
of matrices [40], implying that the set of non-defective rank-$r$ linear
operators is dense in the space of rank-$r$ linear operators on a Hilbert
space. Therefore, there exists a non-defective
$G\in{\rm{B}}_{r}({\mathcal{H}})$ such that $\lVert
G-B_{r}\rVert_{\rm{HS}}<\delta/(3\sigma_{1}(Z_{\mu}))$. So, we conclude

$\displaystyle\lVert(\mu I\,{-}\,L)^{-1}Z_{\mu}-Z_{\mu}G\rVert_{\rm{HS}}$
$\displaystyle\hskip 28.45274pt\leq\lVert(I\\!-\\!P_{\mathcal{H}})(\mu
I\,{-}\,L)^{-1}Z_{\mu}\rVert_{\rm{HS}}+\lVert(\mu
I\,{-}\,L)^{-1}Z_{\mu}-[\\![(\mu
I\,{-}\,L)^{-1}Z_{\mu}]\\!]_{r}\rVert_{\rm{HS}}$ $\displaystyle\hskip
56.9055pt+\lVert[\\![(\mu
I\,{-}\,L)^{-1}Z_{\mu}]\\!]_{r}-Z_{\mu}B_{r}\rVert_{\rm{HS}}+\lVert
Z_{\mu}(G-B_{r})\rVert_{\rm{HS}}$ $\displaystyle\hskip
28.45274pt\leq\lVert(I\\!-\\!P_{\mathcal{H}})(\mu
I\,{-}\,L)^{-1}Z_{\mu}\rVert_{\rm{HS}}+\delta.$

∎

######  Example 5.

These three conditions depend on the process $X$ (through its generator $L$)
as well as 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 process of (1)
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 $\mathcal{H}$ and
$\mathcal{L}^{2}_{\pi}$-norms, we have that

$\mathcal{E}(\widehat{G})\leq\underbrace{\lVert(\mu
I\,{-}\,L)^{-1}Z_{\mu}\\!-\\!Z_{\mu}G_{\mu,\gamma}\rVert_{\mathcal{H}\to\mathcal{W}}}_{\text{regularization
bias}}\\!+\\!\underbrace{\lVert
W_{\mu}^{1/2}(G_{\mu,\gamma}\\!\\!-\\!G^{r}_{\mu,\gamma})\rVert}_{\text{rank
reduction bias}}\\!+\\!\underbrace{\lVert
W_{\mu}^{1/2}(G^{r}_{\mu,\gamma}\\!\\!-\\!\widehat{G}^{r}_{\mu,\gamma})\rVert}_{\text{estimator's
variance}},$ (46)

and continue to prove the bound of the first term. Note that, while this proof
technique is standard for operator learning [23, 28], we present it here for
the sake of completeness.

######  Proposition 8.

Let $G_{\mu,\gamma}=W_{\mu,\gamma}^{-1}C$ for $\gamma>0$, and
$P_{\mathcal{H}}\colon\mathcal{W}^{\mu}_{\pi}(\mathcal{X})\to\mathcal{W}^{\mu}_{\pi}(\mathcal{X})$
be the orthogonal projector onto $\operatorname*{\text{\rm
cl}}(\operatorname*{\text{\rm Im}}(Z_{\mu}))$. If the assumptions (BK), (SD)
and (RC) hold, then $\lVert G_{\mu,\gamma}\rVert\leq
c_{\alpha}c_{\mathcal{W}}^{(\alpha-1)/2}$ for $\alpha\in[1,2]$, $\lVert
G_{\mu,\gamma}\rVert\leq c_{\alpha}\,(\mu\gamma)^{(\alpha-1)/2}$ for
$\alpha\in(0,1]$, and

$\|(\mu
I\,{-}\,L)^{-1}Z_{\mu}-Z_{\mu}G_{\mu,\gamma}\|_{\mathcal{H}{\to}\mathcal{W}}\leq
c_{\alpha}\,(\mu\gamma)^{\frac{\alpha}{2}}+\|(I-P_{\mathcal{H}})(\mu
I\,{-}\,L)^{-1}Z_{\mu}\|_{\mathcal{H}{\to}\mathcal{W}}.$ (47)

###### Proof.

Recalling that $P_{\mathcal{H}}:=\sum_{j\in J}z_{j}\otimes z_{j}$, start by
denoting the orthogonal projectors on the subspace of $k$ leading left
singular functions of $Z_{\mu}$ as $P_{k}:=\sum_{j\in[k]}z_{j}\otimes z_{j}$,
respectively. Next, observe that $C{=}Z_{\mu}^{*}(\mu I{-}L)^{{-1}}Z_{\mu}$
and $Z_{\mu}^{*}W_{\mu,\gamma}^{-1}{=}Z_{\mu}^{*}(Z_{\mu}^{*}Z_{\mu}+\mu\gamma
I)^{-1}{=}(Z_{\mu}Z_{\mu}^{*}+\mu\gamma I)^{-1}Z_{\mu}$. Therefore,

$\displaystyle(\mu I{-}L)^{{-1}}Z_{\mu}{-}Z_{\mu}G_{\mu,\gamma}$
$\displaystyle{=}(I{-}Z_{\mu}W_{\mu,\gamma}^{-1}Z_{\mu}^{*})(\mu
I{-}L)^{{-1}}Z_{\mu}$ $\displaystyle{=}(I{-}(Z_{\mu}Z_{\mu}^{*}{+}\mu\gamma
I)^{-1}Z_{\mu}Z_{\mu}^{*})(\mu I{-}L)^{{-1}}Z_{\mu}$
$\displaystyle{=}\mu\gamma(Z_{\mu}Z_{\mu}^{*}{+}\mu\gamma I)^{-1}(\mu
I{-}L)^{{-1}}Z_{\mu}$

and, hence

$\displaystyle(\mu I{-}L)^{{-1}}Z_{\mu}{-}Z_{\mu}G_{\mu,\gamma}$
$\displaystyle{=}\left(\sum_{j\in
J}\frac{\mu\gamma}{\sigma_{j}^{2}+\mu\gamma}z_{j}\otimes z_{j}\right)(\mu
I{-}L)^{{-1}}Z_{\mu}$ $\displaystyle=\left(\sum_{j\in
J}\frac{\mu\gamma}{(\sigma_{j}^{2}+\mu\gamma)\sigma_{j}}z_{j}\otimes(Z_{\mu}h_{j})\right)(\mu
I{-}L)^{{-1}}Z_{\mu}$ $\displaystyle{=}\left(\sum_{j\in
J}\frac{\mu\gamma}{(\sigma_{j}^{2}+\mu\gamma)\sigma_{j}}z_{j}\otimes
h_{j}\right)C.$

Therefore, for every $k\in J$, $\|P_{k}((\mu
I{-}L)^{{-1}}Z_{\mu}-Z_{\mu}G_{\mu,\gamma})\|_{\mathcal{H}{\to}\mathcal{W}}^{2}$
becomes

$\left\|\left(\sum_{j\in[k]}\frac{\mu\gamma}{(\sigma_{j}^{2}{+}\mu\gamma)\sigma_{j}}z_{j}\otimes
h_{j}\right)C^{2}\left(\sum_{j\in[k]}\frac{\gamma}{(\sigma_{j}^{2}{+}\mu\gamma)\sigma_{j}}h_{j}\otimes
z_{j}\right)\right\|_{\mathcal{H}\to\mathcal{W}},$

which, due to (RC), implies that

$\|P_{k}((\mu
I{-}L)^{{-1}}Z_{\mu}{-}Z_{\mu}G_{\mu,\gamma})\|_{\mathcal{H}{\to}\mathcal{W}}{\leq}c_{\alpha}\,\left\|\sum_{j\in[k]}\frac{\mu\gamma\,\sigma_{j}^{\alpha}}{\sigma_{j}^{2}+\mu\gamma}z_{j}\otimes
z_{j}\right\|_{\mathcal{W}\to\mathcal{W}}.$

On the other hand,

$\sum_{j\in[k]}\frac{\mu\gamma\,\sigma_{j}^{\alpha}}{\sigma_{j}^{2}+\mu\gamma}z_{j}\otimes
z_{j}{=}\gamma^{\frac{\alpha}{2}}\sum_{j\in[k]}\frac{(\sigma_{j}^{2}(\mu\gamma)^{-1})^{\frac{\alpha}{2}}}{\sigma_{j}^{2}(\mu\gamma)^{-1}+1}z_{j}\otimes
z_{j}\preceq(\mu\gamma)^{\frac{\alpha}{2}}\sum_{j\in[k]}z_{j}\otimes z_{j},$

where the inequality holds due to $x^{s}\leq x+1$ for all $x\geq 0$ and
$s\in[0,1]$. Since the norm of the projector equals one, we get $\lVert
P_{k}((\mu I{-}L)^{{-1}}Z_{\mu}-Z_{\mu}G_{\mu,\gamma})\rVert\leq
c_{\alpha}(\mu\gamma)^{\frac{\alpha}{2}}$.

Next, observe that

$\displaystyle\|(P_{\mathcal{H}}{-}P_{k})((\mu
I{-}L)^{{-1}}Z_{\mu}{-}Z_{\mu}G_{\mu,\gamma})\|_{\mathcal{H}{\to}\mathcal{W}}^{2}$
$\displaystyle{=}\left\|\sum_{j\in
J\setminus[k]}\frac{\mu^{2}\gamma^{2}}{(\sigma_{j}^{2}{+}\mu\gamma)^{2}}(Z_{\mu}^{*}z_{j})\otimes(Z_{\mu}^{*}z_{j})\right\|_{\mathcal{H}\to\mathcal{H}}$

which is bounded by

$\displaystyle\sum_{j\in
J\setminus[k]}\\!\frac{\mu^{2}\gamma^{2}}{(\sigma_{j}^{2}{+}\mu\gamma)^{2}}\lVert
Z_{\mu}^{*}z_{j}\rVert^{2}\leq\sum_{j\in
J\setminus[k]}\frac{\mu^{2}\gamma^{2}\,\sigma_{j}^{2\alpha}}{(\sigma_{j}^{2}{+}\mu\gamma)^{2}}\leq\sum_{j\in
J\setminus[k]}\\!\sigma_{j}^{2\alpha}.$

Using triangular inequality, for every $k\in J$ we have that
$\|P_{\mathcal{H}}((\mu
I{-}L)^{{-1}}Z_{\mu}{-}Z_{\mu}G_{\mu,\gamma})\|_{\mathcal{H}{\to}\mathcal{W}}$
is bounded by

$\|P_{k}((\mu
I{-}L)^{{-1}}Z_{\mu}{-}Z_{\mu}G_{\mu,\gamma})\|_{\mathcal{H}{\to}\mathcal{W}}{+}\|(P_{\mathcal{H}}{-}P_{k})((\mu
I{-}L)^{{-1}}Z_{\mu}{-}Z_{\mu}G_{\mu,\gamma})\|_{\mathcal{H}{\to}\mathcal{W}}$

and, therefore,

$\|P_{\mathcal{H}}((\mu
I{-}L)^{{-1}}Z_{\mu}{-}Z_{\mu}G_{\mu,\gamma})\|_{\mathcal{H}{\to}\mathcal{W}}\leq
c_{\alpha}(\mu\gamma)^{\frac{\alpha}{2}}+\sum_{j\in
J\setminus[k]}(\sigma_{j}^{2\beta})^{\frac{\alpha}{\beta}},$

and, hence, letting $k\to\infty$ we obtain $\lVert P_{\mathcal{H}}(\mu
I{-}L)^{{-1}}Z_{\mu}-Z_{\mu}G_{\mu,\gamma}\rVert\leq
c_{\alpha}(\mu\gamma)^{\frac{\alpha}{2}}$. Hence, (47) follows from triangular
inequality.

To estimate the $\lVert G_{\mu,\gamma}\rVert$, note that (RC) implies $\lVert
G_{\mu,\gamma}\rVert\leq\,c_{\alpha}\,\lVert
W_{\mu,\gamma}^{-1}W_{\mu}^{\frac{1+\alpha}{2}}\rVert$ and considering two
cases. First, if (RC) holds for some $\alpha\in[1,2]$, then, clearly $\lVert
G_{\mu,\gamma}\rVert\leq c_{\alpha}c_{\mathcal{W}}^{(\alpha-1)/2}$. On the
other hand, if $\alpha\in(0,1]$, then

$\frac{\sigma_{j}^{1+\alpha}}{\sigma_{j}^{2}+\mu\gamma}=(\mu\gamma)^{-1}\frac{\left(\sigma_{j}^{2}(\mu\gamma)^{-1}\right)^{\frac{1+\alpha}{2}}}{\sigma_{j}^{2}\,(\mu\gamma)^{-1}+1}\leq(\mu\gamma)^{\frac{\alpha-1}{2}},$

and, thus, $\lVert G_{\mu,\gamma}\rVert\leq
c_{\alpha}\,\gamma^{(\alpha-1)/2}$.

∎

######  Remark 2.

Inequality (47) says that the regularization bias is comprised of a term
depending on the choice of $\gamma$, and on a term depending on the
“alignment” between $\mathcal{H}$ and $\operatorname*{\text{\rm Im}}((\mu
I\,{-}\,L)^{-1}Z_{\mu})$. The term $\lVert(I-P_{\mathcal{H}})(\mu
I\,{-}\,L)^{-1}Z_{\mu}\rVert$ can be set to zero by two different approaches.
One is choose a kernel which in some way minimizes
$\lVert(I-P_{\mathcal{H}})(\mu I\,{-}\,L)^{-1}Z_{\mu}\rVert$. Another is to
choose a universal kernel [38, Chapter 4], for which $\operatorname*{\text{\rm
Im}}((\mu I\,{-}\,L)^{-1}S_{\pi})\subseteq\operatorname*{\text{\rm
cl}}(\operatorname*{\text{\rm Im}}(S_{\pi}))$. While we here develop theory
for universal kernels, deep learning approaches that leverage on 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 9.

Let $B:=W_{\mu,\gamma}^{-1/2}C$, let (RC) hold for some $\alpha\in(0,2]$, and
for $j\in\mathbb{N}$ denote

$\lambda_{j}^{\star}=\sigma_{j}^{2}((\mu
I{-}L)^{{-1}}Z_{\mu})=\lambda_{j}(S_{\pi}^{*}(\mu I\,{-}\,L)^{-1}S_{\pi})$
(48)

Then for every $j\in\mathbb{N}$

$\lambda_{j}^{\star}-c_{\alpha}^{2}\,c_{\mathcal{W}}^{\alpha/2}\,(\mu\gamma)^{\alpha/2}\leq\sigma_{j}^{2}(B)\leq\lambda_{j}^{\star}.$
(49)

###### Proof.

Start by observing that

$\displaystyle B^{*}B$ $\displaystyle=[(\mu
I{-}L)^{{-1}}Z_{\mu}]^{*}Z_{\mu}W_{\mu,\gamma}^{-1}Z_{\mu}^{*}(\mu
I{-}L)^{{-1}}Z_{\mu}$ $\displaystyle=[(\mu I{-}L)^{{-1}}Z_{\mu}]^{*}(\mu
I{-}L)^{{-1}}Z_{\mu}-\mu\gamma[(\mu
I{-}L)^{{-1}}Z_{\mu}]^{*}(Z_{\mu}Z_{\mu}^{*}+\mu\gamma I)^{-1}(\mu
I{-}L)^{{-1}}Z_{\mu},$

implies that, using $[(\mu I{-}L)^{{-1}}Z_{\mu}]^{*}{=}S_{\pi}$ and $[(\mu
I{-}L)^{{-1}}Z_{\mu}]^{*}(\mu I{-}L)^{{-1}}Z_{\mu}{=}S_{\pi}^{*}(\mu
I\,{-}\,L)^{-1}S_{\pi}$,

$S_{\pi}^{*}(\mu I\,{-}\,L)^{-1}S_{\pi}-\sum_{j\in
J}\frac{\mu\gamma}{\sigma_{j}^{2}+\mu\gamma}(S_{\pi}^{*}z_{j})\otimes(S_{\pi}^{*}z_{j})=B^{*}B\preceq
Z_{\mu}^{*}Z_{\mu}.$

Next, similarly to the above, for every $k\in J$, we have

$\displaystyle\left\|\sum_{j\in[k]}\frac{\mu\gamma}{\sigma_{j}^{2}+\mu\gamma}(S_{\pi}^{*}z_{j}){\otimes}(S_{\pi}^{*}z_{j})\right\|$
$\displaystyle{\leq}c_{\alpha}^{2}\left\|\sum_{j\in[k]}\frac{\sigma_{j}^{2\alpha}}{\sigma_{j}^{2}(\mu\gamma)^{-1}+1}z_{j}{\otimes}z_{j}\right\|$
$\displaystyle{=}c_{\alpha}^{2}\left\|\sum_{j\in[k]}\\!\\!\frac{(\sigma_{j}^{2}(\mu\gamma)^{-1})^{\alpha/2}\sigma_{j}^{\alpha}(\mu\gamma)^{\alpha/2}}{\sigma_{j}^{2}(\mu\gamma)^{-1}+1}z_{j}{\otimes}z_{j}\right\|{\leq}c_{\alpha}^{2}\gamma^{\alpha/2}\lVert
W_{\mu}\rVert^{\alpha/2},$

and

$\left\|\sum_{j\in
J\setminus[k]}\frac{\gamma}{\sigma_{j}^{2}+\mu\gamma}(S_{\pi}^{*}z_{j})\otimes(S_{\pi}^{*}z_{j})\right\|\leq
c_{\alpha}^{2}\sum_{j\in
J\setminus[k]}\frac{\gamma}{\sigma_{j}^{2}+\mu\gamma}\sigma_{j}^{2\alpha}\leq
c_{\alpha}^{2}\sum_{j\in J\setminus[k]}(\sigma_{j}^{2\beta})^{\alpha/\beta}.$

So, as before, letting $k\to\infty$ we get the result. ∎

As a consequence, we obtain the bound for the rank reduction bias of the RRR
method.

######  Proposition 10 (RRR).

Let (RC) hold for some $\alpha\in(0,2]$. Then the bias of $G^{r}_{\mu,\gamma}$
due to rank reduction, recalling (48) is bounded as

$\sqrt{\lambda^{\star}_{r+1}}-c_{\alpha}\,c_{\mathcal{W}}^{\alpha/4}(\mu\gamma)^{\alpha/4}-2\,c_{\alpha}\,(\mu\gamma)^{(1\wedge\alpha)/2}\leq\lVert
Z_{\mu}(G_{\mu,\gamma}-G^{r}_{\mu,\gamma})\rVert\leq\sqrt{\lambda^{\star}_{r+1}}.$
(50)

###### Proof.

Observe that

$\lVert Z_{\mu}(G_{\mu,\gamma}-G^{r}_{\mu,\gamma})\rVert\leq\lVert
W_{\mu,\gamma}^{1/2}(G_{\mu,\gamma}-G^{r}_{\mu,\gamma})\rVert=\lVert
B-[\\![B]\\!]_{r}\rVert=\sigma_{r+1}(B)\leq\sigma_{r+1}((\mu
I{-}L)^{{-1}}Z_{\mu})$

while

$\displaystyle\lVert Z_{\mu}(G_{\mu,\gamma}-G^{r}_{\mu,\gamma})\rVert$
$\displaystyle\geq\lVert
W_{\mu,\gamma}^{1/2}(G_{\mu,\gamma}-G^{r}_{\mu,\gamma})\rVert-(\mu\gamma)^{1/2}\lVert
G_{\mu,\gamma}-G^{r}_{\mu,\gamma}\rVert$ $\displaystyle\geq\sigma_{r+1}((\mu
I{-}L)^{{-1}}Z_{\mu})-c_{\alpha}\lVert
W_{\mu}\rVert^{\alpha/4}(\mu\gamma)^{\alpha/4}-2c_{\alpha}(\mu\gamma)^{(1\wedge\alpha)/2}.$

∎

###  E.3 Bounding the Variance

####  E.3.1 Concentration Inequalities

All the statistical bounds we present will relay on two versions of Bernstein
inequality. The first one is Pinelis and Sakhanenko inequality for random
variables in a separable Hilbert space, see [8, Proposition 2].

######  Proposition 11.

Let $A_{i}$, $i\in[n]$ be i.i.d copies of a random variable $A$ in a separable
Hilbert space with norm $\lVert\cdot\rVert$. If there exist constants
$\Lambda>0$ and $\sigma>0$ such that for every $m\geq 2$ $\mathbb{E}\lVert
A\rVert^{m}\leq\frac{1}{2}m!\Lambda^{m-2}\sigma^{2}$, then with probability at
least $1-\delta$

$\left\|\frac{1}{n}\sum_{i\in[n]}A_{i}-\mathbb{E}A\right\|\leq\frac{4\sqrt{2}}{\sqrt{n}}\log\frac{2}{\delta}\sqrt{\sigma^{2}+\frac{\Lambda^{2}}{n}}.$
(51)

On the other hand, we recall that in [32], a dimension-free version of the
non-commutative Bernstein inequality for finite-dimensional symmetric matrices
is proposed (see also Theorem 7.3.1 in [41] for an easier to read and slightly
improved version) as well as an extension to self-adjoint Hilbert-Schmidt
operators on a separable Hilbert spaces.

######  Proposition 12.

Let $A_{i}$, $i\in[n]$ be i.i.d copies of a Hilbert-Schmidt operator $A$ on
the separable Hilbert space. Let $\lVert A\rVert\leq c$ almost surely,
$\mathbb{E}A=0$ and let $\mathbb{E}[A^{2}]\preceq V$ for some trace class
operator $V$. Then with probability at least $1-\delta$

$\left\|\frac{1}{n}\sum_{i\in[n]}A_{i}\right\|\leq\frac{2c}{3n}\mathcal{L}_{A}(\delta)+\sqrt{\frac{2\lVert
V\rVert}{n}\mathcal{L}_{A}(\delta)},$ (52)

where

$\mathcal{L}_{A}(\delta)=\log\frac{4}{\delta}+\log\frac{\operatorname{\text{\rm
tr}}(V)}{\lVert V\rVert}.$

######  Proposition 13.

Given $\delta>0$, with probability in the i.i.d. draw of $(x_{i})_{i=1}^{n}$
from $\pi$, it holds that

$\mathbb{P}\\{\lVert\widehat{W}_{\mu}-W_{\mu}\rVert\leq\varepsilon_{n}(\delta)\\}\geq
1-\delta,$

where

$\varepsilon_{n}(\delta)=\frac{2c_{\mathcal{W}}}{3n}\mathcal{L}(\delta)+\sqrt{\frac{2\lVert
W_{\mu}\rVert}{n}\mathcal{L}(\delta)}\quad\text{ and
}\quad\mathcal{L}(\delta)=\log\frac{4\operatorname{\text{\rm
tr}}(W_{\mu})}{\delta\,\lVert W_{\mu}\rVert}.$ (53)

###### Proof.

Proof follows directly from Proposition 12 applied to rank-($1{+}p$) operators
$w\phi(x_{i})\otimes w\phi(x_{i})$ using the fact that
$W_{\mu}=\mathbb{E}[w\phi(x_{i})\otimes w\phi(x_{i})]$. ∎

######  Proposition 14.

Let (KE) hold for $\tau\in[\beta,1]$. Given $\delta>0$, with probability in
the i.i.d. draw of $(x_{i})_{i=1}^{n}$ from $\pi$, it holds that

$\mathbb{P}\left\\{\lVert
W_{\mu,\gamma}^{-1/2}(\widehat{W}_{\mu}-W_{\mu})W_{\mu,\gamma}^{-1/2}\rVert\leq\varepsilon^{1}_{n}(\gamma,\delta)\right\\}\geq
1-\delta,$ (54)

where

$\varepsilon_{n}^{1}(\gamma,\delta)=\frac{2c_{\tau}\mu^{-\tau}}{3n\gamma^{\tau}}\mathcal{L}^{1}(\gamma,\delta)+\sqrt{\frac{2\,c_{\tau}\mu^{-\tau}}{n\,\gamma^{\tau}}\mathcal{L}^{1}(\gamma,\delta)},$
(55)

and

$\mathcal{L}^{1}(\gamma,\delta)=\ln\frac{4}{\delta}+\ln\frac{\operatorname{\text{\rm
tr}}(W_{\mu,\gamma}^{-1}W_{\mu})}{\lVert W_{\mu,\gamma}^{-1}W_{\mu}\rVert}.$

Moreover,

$\mathbb{P}\left\\{\lVert
W_{\mu,\gamma}^{1/2}\widehat{W}_{\mu,\gamma}^{-1}W_{\mu,\gamma}^{1/2}\rVert\leq\frac{1}{1-\varepsilon_{n}^{1}(\gamma,\delta)}\right\\}\geq
1-\delta.$ (56)

###### Proof.

The idea is to apply Proposition 12 for operator $\xi(x)\otimes\xi(x)$, where
$\xi(x)$ is defined in (43). Due to (1), we have that $\lVert
A\rVert\leq\lVert\xi\rVert_{\infty}^{2}\leq(\mu\gamma)^{-\tau}c_{\tau}$. On
the other hand, we have that

$\mathbb{E}_{x\sim\pi}[\xi(x)\otimes\xi(x)]^{2}\preceq\lVert\xi\rVert_{\infty}^{2}\mathbb{E}_{x\sim\pi}[\xi(x)\otimes\xi(x)]=\lVert\xi\rVert_{\infty}^{2}W_{\mu,\gamma}^{-1/2}W_{\mu}W_{\mu,\gamma}^{-1/2},$

and, hence (54) follows. To complete the proof, observe that

$\lVert
I_{\mathcal{H}}-W_{\mu,\gamma}^{-1/2}\widehat{W}_{\mu,\gamma}W_{\mu,\gamma}^{-1/2}\rVert=\lVert
W_{\mu,\gamma}^{-1/2}(W_{\mu}-\widehat{W}_{\mu})W_{\mu,\gamma}^{-1/2}\rVert\leq\varepsilon_{n}^{1}(\gamma,\delta),$

and, hence for $\varepsilon_{n}^{1}(\gamma,\delta)$ smaller than one we obtain

$\lVert
W_{\mu,\gamma}^{1/2}\widehat{W}_{\mu,\gamma}^{-1}W_{\mu,\gamma}^{1/2}\rVert=\lVert(W_{\mu,\gamma}^{-1/2}\widehat{W}_{\mu,\gamma}W_{\mu,\gamma}^{-1/2})^{-1}\rVert\leq\frac{1}{1-\lVert
I_{\mathcal{H}}-W_{\mu,\gamma}^{-1/2}\widehat{W}_{\mu,\gamma}W_{\mu,\gamma}^{-1/2}\rVert}.$

∎

######  Proposition 15.

Let (RC), (SD) and (KE) hold for some $\alpha\in(0,2]$, $\beta\in(0,1]$ and
$\tau\in[\beta,1]$. Given $\delta>0$, with probability in the i.i.d. draw of
$(x_{i})_{i=1}^{n}$ from $\pi$, it holds

$\mathbb{P}\left\\{\lVert
W_{\mu,\gamma}^{-1/2}(\widehat{W}_{\mu}-W_{\mu})W_{\mu,\gamma}^{-1}C\rVert_{\rm{HS}}\leq\varepsilon^{2}_{n}(\gamma,\delta)\right\\}\geq
1-\delta,$

where

$\varepsilon_{n}^{2}(\gamma,\delta)=4\,\sqrt{2\,c_{\alpha}\,c_{\tau}}\ln\frac{2}{\delta}\,\sqrt{\frac{c_{\beta}\mu^{-\beta}}{n\gamma^{\beta}}+\frac{c_{\tau}\mu^{-\tau}}{n^{2}\gamma^{\tau}}}\;\begin{cases}(\mu\gamma)^{-(\tau-\alpha)/2}&,\alpha\leq\tau,\\\
c_{\mathcal{W}}^{(\alpha-\tau)/2}&,\alpha\geq\tau.\end{cases}$ (57)

###### Proof.

As before, w.l.o.g. set $\mu=1$. First, recall that
${\rm{HS}}\left(\mathcal{H}\right)$ equipped with
$\lVert\cdot\rVert_{\rm{HS}}$ is separable Hilbert space. Hence, we will apply
Proposition 11 for $A=\xi(x)\otimes\psi(x)$, where
$\psi(x)=CW_{\mu,\gamma}^{-1}w\phi(x)\in\mathcal{H}^{1{+}p}$. To that end,
observe that $\mathbb{E}\lVert
A\rVert_{\rm{HS}}^{m}{=}\mathbb{E}\,[\lVert\xi(x)\rVert^{m}_{\mathcal{H}^{1{+}p}}\,\lVert\psi(x)\rVert^{m}_{\mathcal{H}^{1{+}p}}]$,
and that

$\mathbb{E}[\lVert\xi(x)\rVert^{m}_{\mathcal{H}^{1{+}p}}]\leq\frac{1}{2}m!\left(\gamma^{-\tau/2}\,\sqrt{c_{\tau}}\right)^{m-2}\left(\sqrt{\operatorname{\text{\rm
tr}}(W_{\mu,\gamma}^{-1}W_{\mu})}\right)^{2}.$

Recalling Lemma 1, the task is to bound
$\lVert\psi(x)\rVert_{\mathcal{H}^{1{+}p}}^{2}{=}\sum_{i\in[1{+}p]}\|CW_{\mu,\gamma}^{-1}w_{i}\phi(x)\|_{\mathcal{H}}^{2}$.
Using (RC), we have that

$\lVert\psi(x)\rVert_{\mathcal{H}^{1{+}p}}^{2}{\leq}c_{\alpha}\sum_{i\in[1{+}p]}\|W_{\mu}^{(1{+}\alpha)/2}W_{\mu,\gamma}^{-1}w_{i}\phi(x)\|_{\mathcal{H}}^{2}{=}c_{\alpha}\sum_{j\in
J}\sum_{i\in[1{+}p]}{\langle
W_{\mu}^{(1{+}\alpha)/2}W_{\mu,\gamma}^{-1}w_{i}\phi(x),h_{j}\rangle}_{\mathcal{H}}^{2}.$

But, since

${\langle
W_{\mu}^{(1{+}\alpha)/2}W_{\mu,\gamma}^{-1}w_{i}\phi(x),h_{j}\rangle}_{\mathcal{H}}=\frac{\sigma_{j}^{1+\alpha}}{\sigma_{j}^{2}+\gamma}{\langle
w_{i}\phi(x),h_{j}\rangle}_{\mathcal{H}},$

expanding as in proof of Lemma 1 and using (KE), we have that

$\lVert\psi(x)\rVert_{\mathcal{H}^{1{+}p}}^{2}{\leq}c_{\alpha}\sum_{j\in
J}\sum_{i\in[1{+}p]}\left[\frac{\sigma_{j}^{(2+\alpha-\tau)}}{\sigma_{j}^{2}+\gamma}\right]^{2}\frac{{\langle
w_{i}\phi(x),h_{j}\rangle}^{2}}{\sigma_{j}^{2}}\sigma_{j}^{2\tau}\leq\begin{cases}\gamma^{-(\tau-\alpha)}&,\alpha\leq\tau,\\\
c_{\mathcal{W}}^{(\alpha-\tau)}&,\alpha\geq\tau.\end{cases}$

Therefore, we can set

$\Lambda^{2}=c_{\alpha}c_{\tau}^{2}\begin{cases}\gamma^{-(2\tau-\alpha)}&,\alpha\leq\tau,\\\
c_{\mathcal{W}}^{(\alpha-\tau)}\gamma^{-\tau}&,\alpha\geq\tau.\end{cases}\quad\text{
and
}\quad\sigma^{2}=c_{\alpha}c_{\tau}c_{\beta}\begin{cases}\gamma^{-(\beta+\tau-\alpha)}&,\alpha\leq\tau,\\\
c_{\mathcal{W}}^{(\alpha-\tau)}\gamma^{-\beta}&,\alpha\geq\tau.\end{cases}$

in Proposition 11 to obtain the bound ∎

######  Proposition 16.

Let (KE) hold for $\tau\in[\beta,1]$. Given $\delta>0$, with probability in
the i.i.d. draw of $(x_{i})_{i=1}^{n}$ from $\pi$, it holds

$\mathbb{P}\left\\{\lVert
W_{\mu,\gamma}^{-1/2}(\widehat{C}-C)\rVert_{\rm{HS}}\leq\varepsilon^{3}_{n}(\gamma,\delta)\right\\}\geq
1-\delta,$

where

$\varepsilon_{n}^{3}(\gamma,\delta)=\frac{4\,\sqrt{2\,c_{\mathcal{W}}}}{\mu}\,\ln\frac{2}{\delta}\,\sqrt{\frac{c_{\beta}\mu^{-\beta}}{n\gamma^{\beta}}+\frac{c_{\tau}\mu^{-\tau}}{n^{2}\gamma^{\tau}}}.$
(58)

###### Proof.

First, note that since
$\lVert\phi(x)\rVert^{2}_{\mathcal{H}}\leq\mu^{-1}\lVert
w\phi(x)\rVert^{2}_{\mathcal{H}^{1+p}}$ and $C\preceq\mu^{-1}W_{\mu}$, (KE)
and (SD) for $W_{\mu}$ imply analogous assumptions for $C$. Hence, we can
apply Proposition 13 from [23], which using the observation that $\lVert
W_{\mu,\gamma}^{-1/2}C_{\gamma}^{1/2}\rVert^{2}{=}\lVert C_{\gamma}^{1/2}(\mu
C_{\gamma}{-}T)^{-1}C_{\gamma}^{1/2}\rVert{\leq}\mu^{-1}\lVert
C_{\gamma}C_{\gamma}^{-1}\rVert{=}\mu^{-1}$, completes the proof. ∎

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 17.

Let (RC), (SD) and (KE) hold for some $\alpha\in(0,2]$, $\beta\in(0,1]$ and
$\tau\in[\beta,1]$. Given $\delta>0$ if
$\varepsilon^{1}_{n}(\gamma,\delta)<1$, then with probability at least
$1-\delta$ in the i.i.d. draw of $(x_{i},y_{i})_{i=1}^{n}$ from $\rho$

$\mathbb{P}\left\\{\lVert
W_{\mu,\gamma}^{1/2}(\widehat{G}_{\mu,\gamma}-G_{\mu,\gamma})\rVert\leq\frac{\varepsilon^{2}_{n}(\gamma,\delta/3)+\varepsilon^{3}_{n}(\gamma,\delta/3)}{1-\varepsilon^{1}_{n}(\gamma,\delta/3)}\right\\}\geq
1-\delta.$

###### Proof.

Note that
$W_{\mu,\gamma}^{1/2}(\widehat{G}_{\mu,\gamma}-G_{\mu,\gamma})=W_{\mu,\gamma}^{1/2}(\widehat{W}_{\mu,\gamma}^{-1}\widehat{C}-W_{\mu,\gamma}^{-1}C)$,
and, hence,

$\displaystyle W_{\mu,\gamma}^{1/2}(\widehat{G}_{\mu,\gamma}-G_{\mu,\gamma})$
$\displaystyle=W_{\mu,\gamma}^{1/2}\widehat{W}_{\mu,\gamma}^{-1}(\widehat{C}-\widehat{W}_{\mu,\gamma}W_{\mu,\gamma}^{-1}C\pm
C)$
$\displaystyle=W_{\mu,\gamma}^{1/2}\widehat{W}_{\mu,\gamma}^{-1}W_{\mu,\gamma}^{1/2}\left(W_{\mu,\gamma}^{-1/2}(\widehat{C}-C)-W_{\mu,\gamma}^{-1/2}(\widehat{W}_{\mu}-W_{\mu})W_{\mu,\gamma}^{-1}C\right).$
(59)

Thus, taking the norm and using $\lVert C_{\gamma}^{-1}T\rVert\leq
c_{\alpha}\,\sigma_{1}^{\alpha-1}(S_{\pi})$ with the Propositions 15 and 14 we
prove the first bound. ∎

####  E.3.3 Variance of Singular Values

In this section we prove concentration of singular values computed in Theorem
1, a necessary step to derive learining rates for RRR estimator.

######  Proposition 18.

Let (RC), (SD) and (KE) hold for some $\alpha\in(0,2]$, $\beta\in(0,1]$ and
$\tau\in[\beta,1]$. Let $B=W_{\mu,\gamma}^{-1/2}C$ and
$\widehat{B}=\widehat{W}_{\mu,\gamma}^{-1/2}\widehat{C}$. Given $\delta>0$ if
$\varepsilon^{1}_{n}(\gamma,\delta/5)<1$, then with probability at least
$1-\delta$ in the i.i.d. draw of $(x_{i})_{i=1}^{n}$ from $\pi$

$\lvert\sigma_{i}^{2}(\widehat{B}){-}\sigma_{i}^{2}(B)\rvert{\leq}\lVert\widehat{B}^{*}\widehat{B}-B^{*}B\rVert{\leq}\varepsilon^{4}_{n}(\gamma,\delta/3),$
(60)

where

$\varepsilon^{4}_{n}(\gamma,\delta/3)=(\varepsilon^{2}_{n}(\gamma,\delta/3){+}\varepsilon^{3}_{n}(\gamma,\delta/3))\Big{(}\frac{1}{\mu}{+}\frac{\varepsilon^{2}_{n}(\gamma,\delta/3){+}\varepsilon^{3}_{n}(\gamma,\delta/3)}{1-\varepsilon^{1}_{n}(\gamma,\delta/3)}\Big{)}.$
(61)

###### Proof.

We start from the Weyl’s inequalities for the square of singular values

$\lvert\sigma_{i}^{2}(\widehat{B})-\sigma_{i}^{2}(B)\rvert\leq\lVert\widehat{B}^{*}\widehat{B}-B^{*}B\rVert,\;i\in[n].$

But, since,

$\widehat{B}^{*}\widehat{B}-B^{*}B=\widehat{C}^{*}\widehat{W}_{\mu,\gamma}^{-1}\widehat{C}-C^{*}W_{\mu,\gamma}^{-1}C=(\widehat{C}-C)^{*}\widehat{W}_{\mu,\gamma}^{-1}\widehat{C}+C^{*}W_{\mu,\gamma}^{-1}(\widehat{C}-C)+C^{*}(\widehat{W}_{\mu,\gamma}^{-1}-W_{\mu,\gamma}^{-1})\widehat{C}$

denoting $M=W_{\mu,\gamma}^{-1/2}(\widehat{C}-C)$,
$N=W_{\mu,\gamma}^{-1/2}(\widehat{W}_{\mu}-W_{\mu})$ and
$R=W_{\mu,\gamma}^{1/2}(\widehat{G}_{\mu,\gamma}-G_{\mu,\gamma})$, we have

$\displaystyle\widehat{B}^{*}\widehat{B}-B^{*}B$
$\displaystyle=B^{*}M+M^{*}C_{\gamma}^{1/2}\widehat{G}_{\mu,\gamma}-B^{*}N\widehat{G}_{\mu,\gamma}=B^{*}M+(M^{*}C_{\gamma}^{1/2}-B^{*}N)(\widehat{G}_{\mu,\gamma}\pm
G_{\mu,\gamma})$
$\displaystyle=B^{*}M+M^{*}B-B^{*}NG_{\mu,\gamma}+(M^{*}-B^{*}NW_{\mu,\gamma}^{-1/2})R$
$\displaystyle=(G_{\mu,\gamma})^{*}(\widehat{C}-C){+}(\widehat{C}{-}C)G_{\mu,\gamma}{-}(G_{\mu,\gamma})^{*}(\widehat{W}_{\mu}{-}W_{\mu})G_{\mu,\gamma}{+}(M^{*}+(G_{\mu,\gamma})^{*}N^{*})R.$

Therefore, recalling that, due to (59),
$R=W_{\mu,\gamma}^{1/2}\widehat{W}_{\mu,\gamma}^{-1}W_{\mu,\gamma}^{1/2}(M-NG_{\mu,\gamma})$,
we conclude

$\displaystyle\widehat{B}^{*}\widehat{B}-B^{*}B=$
$\displaystyle(G_{\mu,\gamma})^{*}(\widehat{C}{-}C){+}(\widehat{C}{-}C)G_{\mu,\gamma}{-}(G_{\mu,\gamma})^{*}(\widehat{W}_{\mu}-W_{\mu})G_{\mu,\gamma}$
$\displaystyle+(M-NG_{\mu,\gamma})^{*}W_{\mu,\gamma}^{1/2}\widehat{W}_{\mu,\gamma}^{-1}W_{\mu,\gamma}^{1/2}(M-NG_{\mu,\gamma}).$
(62)

Next, observe that

  * •

$\lVert(\widehat{C}{-}C)G_{\mu,\gamma}\rVert{\leq}\lVert(\widehat{C}{-}C)W_{\mu,\gamma}^{1/2}\rVert\lVert
W_{\mu,\gamma}^{1/2}C\rVert{\leq}\mu^{-1}\lVert
W_{\mu,\gamma}^{1/2}(\widehat{C}{-}C)\rVert$ is bounded by Proposition 16

  * •

$\lVert M{-}NG_{\mu,\gamma}\rVert\leq\lVert
W_{\mu,\gamma}^{-1/2}(\widehat{C}{-}C)\rVert+\lVert