# A Theoretical Analysis of the Test Error of Finite-Rank Kernel Ridge
Regression

Tin Sum Cheng Aurelien Lucchi Ivan Dokmanić
Department of Mathematics and Computer Science
University of Basel
{tinsum.cheng, aurelien.lucchi<EMAIL_ADDRESS>
&Anastasis Kratsios
Department of Mathematics and Statistics,
McMaster University and Vector Institute,
<EMAIL_ADDRESS>&David Belius
Faculty of Mathematics and Computer Science
UniDistance Suisse
<EMAIL_ADDRESS>
###### Abstract

Existing statistical learning guarantees for general kernel regressors often
yield loose bounds when used with finite-rank kernels. Yet, finite-rank
kernels naturally appear in several machine learning problems, e.g. when fine-
tuning a pre-trained deep neural network’s last layer to adapt it to a novel
task when performing transfer learning. We address this gap for finite-rank
kernel ridge regression (KRR) by deriving sharp non-asymptotic upper and lower
bounds for the KRR test error of any finite-rank KRR. Our bounds are tighter
than previously derived bounds on finite-rank KRR, and unlike comparable
results, they also remain valid for any regularization parameters.

##  1 Introduction

Generalization is a central theme in statistical learning theory. The recent
renewed interest in kernel methods, especially in Kernel Ridge Regression
(KRR), is largely due to the fact that deep neural network (DNN) training can
be approximated using kernels under appropriate conditions Jacot et al.
(2018); Arora et al. (2019); Bordelon et al. (2020), in which the test error
is more tractable analytically and thus enjoys stronger theoretical
guarantees. However, many prior results have been derived under conditions
incompatible with practical settings. For instance Liang and Rakhlin (2020);
Liu et al. (2021a); Mei et al. (2021); Misiakiewicz (2022) give asymptotic
bounds on the KRR test error, which requires the input dimension $d$ to tend
to infinity. In reality, the input dimension of the data set and the target
function is typically finite. A technical difficulty to provide a sharp non-
asymptotic bound on the KRR test error comes from the infinite dimensionality
of the kernel 111For a fixed input dimension $d$ and sample size $N$, there
exists an (infinite-dimensional) subspace in the feature space in which we
cannot control the feature vector.. While this curse of dimensionality may be
unavoidable for infinite-rank kernels (at least without additional
restrictions), one may likely derive tighter bounds in a finite-rank setting.
Therefore, other works Bach (2023a); Amini et al. (2022) focused on a setting
where the kernel is of finite-rank, where non-asymptotic bounds and exact
formula of test error can be derived.

Since different generalization behaviours are observed depending on whether
the rank of the kernel rank $M$ is smaller than the sample size $N$, one
typically differentiates between the under-parameterized ($M<N$) and over-
parameterized regime ($M>N$). We focus on the former due to its relevance to
several applications in the field of machine learning, including random
feature models Rahimi and Recht (2007); Rudi and Rosasco (2017); Liu et al.
(2021b); Gonon et al. (2022); Herrera et al. (2021) such as reservoir
computers Gonon et al. (2022, 2020, 2023); Compagnoni et al. (2023) where all
the hidden layers of a deep learning model are randomly generated, and only
the final layer is trainable, or when fine-tuning the final layers of pre-
trained deep neural networks for transfer learning Azizpour et al. (2015) or
in few-shot learning Vinyals et al. (2016). The practice of only re-training
pre-trained deep neural networks final layer Yosinski et al. (2014); LI et al.
(2018) is justified by the fact that earlier layers encode general features
which are common to similar tasks thus fine-tuning can support’s a network’s
ability to generalize to new tasks Ju et al. (2022). In this case, the
network’s frozen hidden layers define a feature map into a finite-dimensional
RKHS, which induces a finite-rank kernel; however, similar finite-rank kernels
are often considered Donahue et al. (2014) in which features are extracted
directly from several hidden layers in a deep pre-trained neural network which
is then fed into a trainable linear regressor.

##### Contribution

Our main objective is to provide sharp non-asymptotic upper and lower bounds
for the finite-rank KRR test error in the under-parameterized regime. We make
the following contributions:

  1. (i)

Non-vacuous bound in ridgeless case: In contrast to prior work, our bounds
exhibit better accuracy when the ridge parameter $\lambda\to 0$, matching the
intuition that a smaller ridge yields a smaller test error;

  2. (ii)

Sharp non-asymptotic lower bound: We provide a sharp lower bound of test
error, which matches our derived upper bound as the sample size $N$ increases.
In this sense, our bounds are tight;

  3. (iii)

Empirical validation: We experimentally validate our results and show our
improvement in bounds over Bach (2023a).

As detailed in Section D, Table 1 compares our results to the available risk-
bounds for finite-rank KRR.

Table 1: Comparison of risk-bounds for finite-rank kernel ridge regressors. (*): The bound for the inconsistent case is implicit. See Section 5 for details. (**): By decay, we mean the difference between the upper bound and its limit as $N\to\infty$.  Assumptions / results | Mohri et al. (2018) | Amini et al. (2022) | Bach (2023a) | This paper
---|---|---|---|---
Include inconsistent case | ✓ | ✗ | (✓)* | ✓
Bias-variance decomposition | ✗ | ✓ | ✓ | ✓
Test error high probability upper bound | ✓ | ✗ | ✓ | ✓
Test error high probability lower bound | ✗ | ✗ | ✗ | ✓
Bounds improve with smaller ridge? | ✗ | - | ✗ | ✓
Upper bound decay rate** | ${\sqrt{\frac{\log N}{N}}}$ | - | ${\lambda+\frac{\log N}{\lambda N}}$ | ${(\lambda+\frac{1}{N})\sqrt{\frac{\log N}{N}}}$

##### Organization of Paper

Our paper is organized as follows: Section 2 motivates the importance of the
under-parameterized finite-rank KRR. Section 3 introduces the notation and
necessary background material required in our results’ formulation. Section 4
summarizes our main findings and illustrates their implications via a
numerical study. All notation is tabulated in Appendix A for the reader’s
convenience. All proofs can be found in Appendices B and C, and numerical
validation of our theoretical findings in Appendix D.

##  2 Applications

In this section, we motivate the importance of the under-parameterized finite-
rank KRR in practical settings. For readers more interested in our main
results, please start from Section 3.

##### Application: Fine-tuning Pre-Trained Deep Neural Networks For Transfer
Learning

Consider the transfer-learning problem of fine-tuning the final layer of a
pre-trained deep neural network model
$f:\mathbb{R}^{d}\rightarrow\mathbb{R}^{D}$ so that it can be adapted to a
task that is similar to what it was initially trained for. This procedure
defines a finite-rank kernel regressor because, for any $x\in\mathbb{R}^{d}$,
$f$ can be factored as

$\displaystyle f(x)$ $\displaystyle=A\phi(x),$ $\displaystyle\phi(x)$
$\displaystyle=\sigma\bullet W^{(L)}\circ\dots\sigma\bullet W^{(1)}(x),$

where $A$ is a $D\times d_{L}$-matrix, $\sigma\bullet$ denotes element-wise
application of a univariate non-linear function, for $l=1,\dots,L$, $L$ is a
positive integer, and
$W^{(l)}:\mathbb{R}^{d_{l-1}}\rightarrow\mathbb{R}^{d_{l}}$, $d=d_{0}$. In
pre-training, all parameters defining the affine maps $W^{(l)}$ in the hidden
layers $\sigma\bullet W^{(L)},\dots,\sigma\bullet W^{(1)}$ are simultaneously
optimized, while in fine-tuning, only the parameters in the final $A$ are
trained, and the others are frozen. This reduces $f$ to a finite-rank kernel
regressor with finite-rank kernel given for $x,\tilde{x}\in\mathbb{R}^{d}$ by

$K(x,\tilde{x})\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\phi(x)^{\top}\phi(\tilde{x}).$

Stably optimizing $f$ to the new task requires strict convexity of the KRR
problem

$\min_{A}\,\frac{1}{N}\sum_{n=1}^{N}(A\phi(x_{n})-y_{n})^{2}+\lambda\|A\|_{F}^{2},$
(1)

where the hyperparameter $\lambda>0$ ensures strong convexity of the KRR’s
loss function and where $\|A\|_{F}$ denotes the Frobenius norm of A. The
unique solution $\hat{A}$ to (1), determines the optimally trained KRR model
$\hat{f}(x)\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\hat{A}\phi(x)$
corresponding to the finite-rank kernel $K$.

##### Application: Random Feature Models

Popularized by Rahimi and Recht (2007) to enable more efficient kernel
computations, random feature model has recently seen a substantial spike in
popularity and has been the topic of many theoretical works Jacot et al.
(2020); Mei et al. (2021). Note that random feature models are finite-rank
kernel regressors once their features are randomly sampled Rahimi and Recht
(2007); Rudi and Rosasco (2017); Liu et al. (2021b); Gonon et al. (2022);
Compagnoni et al. (2023).

##### Application: General Use

We emphasize that, though fine-tuning and random feature models provide simple
typical examples of when finite-rank KRR arise in machine learning, there are
several other instances where our results apply, e.g. when deriving
generalization bounds for infinite-rank kernels by truncation, thereby
replacing them by finite-rank kernels; e.g. Mei et al. (2021).

##  3 Preliminary

We now formally define all concepts and notation used throughout this paper. A
complete glossary is found in Appendix A.

###  3.1 The Training and Testing Data

We fix a (non-empty) input space $\mathcal{X}\subset\mathbb{R}^{d}$ and a
target function $\tilde{f}:\mathcal{X}\to\mathbb{R}$ which is to be learned by
the KRR from a finite number of i.i.d. (possibly noisy) samples
$\mathbf{Z}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}(\mathbf{X},\mathbf{y})=\big{(}(x_{i})_{i=1}^{N},(y_{i})_{i=1}^{N}\big{)}\in\mathbb{R}^{d\times
N}\times\mathbb{R}^{N}$. The inputs $x_{i}$ are drawn from a sampling
distribution $\rho$ on $\mathcal{X}$ and outputs are modelled as
$y_{i}=\tilde{f}(x_{i})+\epsilon_{i}\in\mathbb{R}^{N\times 1}$ for some i.i.d.
independent random variable $\epsilon_{i}$ which is also independent of the
$\mathbf{X}$, satisfying $\mathbb{E}[\epsilon]=0$ and
$\mathbb{E}[\epsilon^{2}]\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\sigma^{2}\geq
0$. Our analysis is set in the space of all square-integrable “function” with
respect to the sampling distribution
$L_{2}({\rho})\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\\{f:\mathcal{X}\to\mathbb{R}:\,\mathbb{E}_{x\sim\rho}[f(x)^{2}]<\infty\\}$.

We abbreviate
$\mathbf{y}=\tilde{f}(\mathbf{X})+\bm{\epsilon}\in\mathbb{R}^{N\times 1}$,
where $\tilde{f}(\mathbf{X})=[\tilde{f}(x_{i})]_{i=1}^{N}$ and
$\bm{\epsilon}=[\epsilon_{i}]_{i=1}^{N}$.

###  3.2 The Assumption: The Finite-Rank Kernel Ridge Regressor

As in Amini et al. (2022); Bach (2023a), we fix a rank $M\in\mathbb{N}_{+}$
for the kernel $K$.

######  Definition 3.1 (Finite Rank Kernel).

Let $M$ be a positive integer. A (finite) rank-$M$ kernel $K$ is a map
$K:\mathcal{X}\times\mathcal{X}\rightarrow\mathbb{R}$ defined for any
$x,{x}^{\prime}\in\mathcal{X}$ by

$K(x,x^{\prime})=\sum_{k=1}^{M}\lambda_{k}\psi_{k}(x)\psi_{k}(x^{\prime}),$
(2)

where the positive numbers $\\{\lambda_{k}\\}_{k=1}^{M}$ are called the
eigenvalues, and orthonormal functions $\psi_{k}\in L^{2}_{\rho}$ are called
the eigenfunctions. 222This means that
$\int_{\mathcal{X}}\psi_{k}(x)\psi_{l}(x)d\rho(x)=\delta_{kl}$ for all $1\leq
k\leq l\leq M$.

######  Remark 3.2 (Eigenvalues are Ordered).

Without loss of generality, we will assume that the kernel $K$’s eigenvalues
are ordered $\lambda_{1}\geq\lambda_{2}\geq\dots\geq\lambda_{M}>0$.

We denote by
$\mathcal{H}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\text{span}\\{\lambda_{1}^{-1/2}\psi_{1},\dots,\lambda_{M}^{-1/2}\psi_{M}\\}$
the reproducing kernel Hilbert space (RKHS) associated to $K$. See Appendix B
for additional details on kernels and RKHSs.

Together, the kernel $K$ and a ridge $\lambda>0$ define an optimal regressor
for the training data $\mathbf{Z}$.

######  Definition 3.3 (Kernel Ridge Regressor (KRR)).

Fix a ridge $\lambda>0$, the regressor $f_{\mathbf{Z},\lambda}$ of the finite-
rank kernel $K$ is the (unique) minimizer of

$f_{\mathbf{Z},\lambda}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\min_{f\in\mathcal{H}}\frac{1}{N}\sum_{i=1}^{N}\left(f(x_{i})-y_{i}\right)^{2}+\lambda\|f\|_{\mathcal{H}}^{2}.$
(3)

###  3.3 The Target Function

The only regularity assumed of the target function $\tilde{f}$ is that it is
square-integrable with respect to the sampling distribution $\rho$; i.e.
$\tilde{f}\in L^{2}_{\rho}$. We typically assume that $\tilde{f}$ contains
strictly more features than can be expressed by the finite-rank kernel $K$.
Thus, $\tilde{f}$ can be arbitrarily complicated even if the kernel $K$ is
not. Since $\\{\psi_{k}\\}_{k=1}^{M}$ is an orthonormal set of $L_{\rho}^{2}$,
we decompose the target function $\tilde{f}$ as

$\tilde{f}=\underbrace{\sum_{k=1}^{M}\tilde{\gamma}_{k}\psi_{k}}_{\tilde{f}_{\leq
M}}+\tilde{\gamma}_{>M}\psi_{>M},$ (4)

for some real numbers $\tilde{\gamma}_{k}$’s and $\tilde{\gamma}_{>M}$ and for
some normal function $\psi_{>M}$ orthogonal to $\\{\psi_{k}\\}_{k=1}^{M}$. We
call $\psi_{>M}$ the orthonormal complement and $\tilde{\gamma}_{>M}$ the
complementary coefficient. The component $\tilde{f}_{\leq
M}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\sum_{k=1}^{M}\tilde{\gamma}_{k}\psi_{k}$
of $\tilde{f}$ is in $\mathcal{H}$. For the case $\tilde{\gamma}_{>M}=0$, we
call it a consistent case, as the target function $\tilde{f}$ lies in the
hypothesis set $\mathcal{H}$, else we call it an inconsistent case.

Alternatively, the orthonormal complement $\psi_{>M}$ can be understood as
some input-dependent noise. Assume we have chosen a suitable finite rank
kernel $K$ with corresponding RKHS $\mathcal{H}$ such that the target function
lies in $\mathcal{H}$. For this purpose, we can write the target function as
$\tilde{f}_{\leq M}=\sum_{k=1}^{M}\tilde{\gamma}_{k}\psi_{k}$ for some real
numbers $\tilde{\gamma}_{k}$’s. Suppose that we sample in a noisy environment;
then for each sample input $x_{i}$, the sample output $y_{i}$ can be written
as

$y_{i}=\underbrace{\tilde{f}_{\leq M}(x_{i})}_{\text{true
label}}+\underbrace{\tilde{\gamma}_{>M}\psi_{>M}(x_{i})}_{\text{input-
dependent noise}}+\underbrace{\epsilon_{i}}_{\text{input-independent noise}}.$
(5)

###  3.4 Test Error

Next, we introduce the subject of interest of this paper in more detail. Our
statistical analysis quantifies the deviation of the learned function from the
ground truth of the test error.

######  Definition 3.4 (KRR Test Error).

Fix a sample $\mathbf{Z}=(\mathbf{X},\tilde{f}(\mathbf{X})+\bm{\epsilon})$.
The finite-rank KRR (3)’s test error is

$\displaystyle\mathcal{R}_{\mathbf{Z},\lambda}$
$\displaystyle\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\mathbb{E}_{x,\epsilon}\left[(f_{\mathbf{Z},\lambda}(x)-\tilde{f}(x))^{2}\right]=\mathbb{E}_{\epsilon}\left[\int_{\mathcal{X}}\left(f_{\mathbf{Z},\lambda}(x)-\tilde{f}(x)\right)^{2}d\rho(x)\right].$
(6)

The analysis of this paper also follows the classical bias-variance
decomposition, thus we write

$\displaystyle\mathcal{R}_{\mathbf{Z},\lambda}$
$\displaystyle=\text{bias}+\text{variance},$

where bias measures the error when there is no noise in the label, that is,

$\text{bias}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\int_{\mathcal{X}}\left(f_{(\mathbf{X},\tilde{f}(\mathbf{X})),\lambda}(x)-\tilde{f}(x)\right)^{2}d\rho(x),$
(7)

and variance is defined to be the difference between test error and bias:
$\text{variance}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\mathcal{R}_{\mathbf{Z},\lambda}-\text{bias}$.

##  4 Main Result

We now present the informal version of our main result, which gives high-
probability upper and lower bounds on the test error. This result is obtained
by bounding both the bias and variance, and the probability that the bounds
hold is quantified with respect to the sampling distribution $\rho$. Here we
emphasize the condition that $N>M$ and hence our statement is valid only in
under-parametrized case.

We can assume the data-generating distribution $\rho$ and eigenfunctions
$\psi_{k}$ are well-behaved in the sense that:

######  Assumption 4.1 (Sub-Gaussian-ness).

We assume that the probability distribution of the random variable
$\psi_{k}(x)$, where $x\in\rho$, has sub-Gaussian norm bounded by a positive
constant $G>0$, for all $k\in\\{1,...,M\\}\cup\\{>M\\}$333it means the
orthonormal complement $\psi_{>M}$ is also mentioned in the assumption..

In particular, if the random variable $\psi_{k}(x)$ is bounded, the assumption
4.1 is fulfilled.

######  Theorem 4.2 (High Probability Bounds on Bias and Variance).

Suppose Assumption 4.1 holds, for $N>M$ sufficient large, there exists some
constants $C_{1},C_{2}$ independent to $N$ and $\lambda$ such that, with a
probability of at least $1-2/N$ w.r.t. random sampling, we have the following
results simultaneously:

  1. (i)

Upper-Bound on Bias: The bias is upper bounded by:

$\text{bias}\leq\tilde{\gamma}_{>M}^{2}+\lambda\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}+\left(\frac{1}{4}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}+2\lambda\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}\right)\sqrt{\frac{\log N}{N}}+C_{1}\frac{\log N}{N},$
(8)

where we denote $\tilde{f}_{\leq
M}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\sum_{k=1}^{M}\tilde{\gamma}_{k}\psi_{k}=\tilde{f}-\tilde{\gamma}_{>M}\psi_{>M}$;

  2. (ii)

Lower-Bound on Bias: The bias is lower bounded by an analogous result:

$\text{bias}\geq\tilde{\gamma}_{>M}^{2}+\frac{\lambda^{2}\lambda_{M}}{(\lambda_{M}+\lambda)^{2}}\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}-\left(\frac{1}{4}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}+\frac{2\lambda^{2}}{\lambda_{1}+\lambda}\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}\right)\sqrt{\frac{\log N}{N}}-C_{1}\frac{\log N}{N}.$
(9)

  3. (iii)

Upper-Bound on Variance: The variance is upper bounded by:

$\displaystyle\begin{split}\text{variance}&\leq\sigma^{2}\frac{M}{N}\left(1+\sqrt{\frac{\log
N}{N}}+C_{2}\frac{\log N}{N}\right);\\\ \end{split}$ (10)

  4. (iv)

Lower-Bound on Variance: The variance is lower bounded by an analogous result:

$\text{variance}\geq\frac{\lambda_{M}^{2}}{(\lambda_{M}+\lambda)^{2}}\sigma^{2}\frac{M}{N}\left(1-\sqrt{\frac{\log
N}{N}}\right)-C_{2}\sigma^{2}\frac{M}{N}\frac{\log N}{N}.$ (11)

For $\lambda\to 0$, we have a simpler bound on the bias: with a probability of
at least $1-2/N$, we have

$\displaystyle\begin{split}\lim_{\lambda\to
0}\text{bias}&\leq\tilde{\gamma}_{>M}^{2}\left(1+\frac{\log
N}{N}\right)+6\tilde{\gamma}_{>M}^{2}\left(\frac{\log
N}{N}\right)^{\frac{3}{2}};\\\ \lim_{\lambda\to
0}\text{bias}&\geq\tilde{\gamma}_{>M}^{2}\left(1-\frac{\log
N}{N}\right)-6\tilde{\gamma}_{>M}^{2}\left(\frac{\log
N}{N}\right)^{\frac{3}{2}}.\end{split}$ (12)

If $\tilde{\gamma}_{>M}^{2}=0$, then we are in the consistent case, meaning
that $\tilde{f}$ belongs to $\mathcal{H}$. In this case, we have a simpler
bound on the bias: with a probability of at least $1-2/N$, we have

$\text{bias}\leq\lambda\|\tilde{f}\|_{\mathcal{H}}^{2}\left(1+2\sqrt{\frac{\log
N}{N}}\right)+C_{1}\frac{\log N}{N}.$ (13)

Proof Sketch: The main technical tools are 1) more careful algebraic
manipulations when dealing with terms involving the regularizer $\lambda$ and
2) the use of a concentration result for a sub-Gaussian random covariance
matrix in Vershynin (2010) followed by the Neumann series expansion of a
matrix inverse. Hence, unlike previous work, our result holds for any
$\lambda$, which can be chosen independence to $N$. The complete proof of this
result, together with the explicit form of the lower bound, is presented in
Theorems C.19 and C.20 in the Appendix C. $\square$

######  Remark 4.3.

Note that one can also easily derive bounds on the expected value of the test
error to allow comparisons with prior works, such as Bach (2023a); Tsigler and
Bartlett (2023); Ronen et al. (2019).

######  Remark 4.4.

The main technical difference from prior works is that we consider the basis
$\\{\psi_{k}\\}_{k=1}^{M}$ in the function space $L_{\rho}^{2}$ instead of
$\\{\lambda_{k}^{-1/2}\psi_{k}\\}_{k=1}^{M}$ in the RKHS $\mathcal{H}$. This
way, we exploit decouple the effect of spectrum from the sampling randomness
to obtain a sharper bound. For further details, see Remark C.12 in Appendix C.

Combining Theorem 4.2 and the bias-variance decomposition in Definition 3.4,
we have both the upper and lower bounds of the test error on KRR.

######  Corollary 4.4.1.

Under mild conditions on the kernel $K$, for $N$ sufficiently large, there
exist some constants $C_{1},C_{2}$ independent to $N$ and $\lambda$ such that,
with a probability of at least $1-2/N$ w.r.t. random sampling, we have the
bounds on the test error $\mathcal{R}_{\mathbf{Z},\lambda}$:

$\displaystyle\mathcal{R}_{\mathbf{Z},\lambda}$
$\displaystyle\leq\tilde{\gamma}_{>M}^{2}+\lambda\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}+\left(\frac{1}{4}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}+2\lambda\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}\right)\sqrt{\frac{\log N}{N}}+C_{1}\frac{\log N}{N}$
$\displaystyle\quad+\sigma^{2}\frac{M}{N}\left(1+\sqrt{\frac{\log
N}{N}}+C_{2}\frac{\log N}{N}\right);$
$\displaystyle\mathcal{R}_{\mathbf{Z},\lambda}$
$\displaystyle\geq\tilde{\gamma}_{>M}^{2}+\frac{\lambda^{2}\lambda_{M}}{(\lambda_{k}+\lambda)^{2}}\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}-\left(\frac{1}{4}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}+\frac{2\lambda^{2}}{\lambda_{1}+\lambda}\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}\right)\sqrt{\frac{\log N}{N}}-C_{1}\frac{\log N}{N}$
$\displaystyle\quad+\sigma^{2}\frac{M}{N}\left(1-\sqrt{\frac{\log
N}{N}}-C_{2}\frac{\log N}{N}\right).$

In particular, we have:

$\displaystyle\lim_{\lambda\to 0}\mathcal{R}_{\mathbf{Z},\lambda}$
$\displaystyle\leq\left(1+\frac{\log
N}{N}\right)\tilde{\gamma}_{>M}^{2}+6\tilde{\gamma}_{>M}^{2}\left(\frac{\log
N}{N}\right)^{\frac{3}{2}}+\sigma^{2}\frac{M}{N}\left(1+\sqrt{\frac{\log
N}{N}}+C_{2}\frac{\log N}{N}\right).$

The corresponding lower bounds are given analogously:

$\lim_{\lambda\to 0}\mathcal{R}_{\mathbf{Z},\lambda}\geq\left(1-\frac{\log
N}{N}\right)\tilde{\gamma}_{>M}^{2}-6\tilde{\gamma}_{>M}^{2}\left(\frac{\log
N}{N}\right)^{\frac{3}{2}}+\sigma^{2}\frac{M}{N}\left(1-\sqrt{\frac{\log
N}{N}}-C_{2}\frac{\log N}{N}\right).$

See Section 5 for more details and validations of Theorem 4.2.

##  5 Discussion

In this section, we first elaborate on the result from Theorem 4.2, which we
then compare in detail to prior works, showcasing the improvements this paper
makes. Finally, we discuss several future research directions.

###  5.1 Elaboration on Main Result

##### Bias

From Eq. (8), we can draw the following observations: 1) The term
$\tilde{\gamma}_{>M}^{2}$ in the upper bound is the finite rank error due to
the inconsistency between RKHS and the orthogonal complement of the target
function, which cannot be improved no matter what sample size we have. We can
also view this term as the sample-dependent noise variance (see Eq. (5)).
Hence, unlike the sample-independent noise variance $\sigma$ in Eq. (10), the
sample-dependent noise variance does not vanish when $N\to\infty$. 2) Note
that the third term
$\frac{1}{4}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}\sqrt{\frac{\log N}{N}}$ is a
residue term proportional to $\|\tilde{f}\|_{L_{\rho}^{2}}^{2}$ and vanishes
when $N\to\infty$, which means we have better control of the sample-dependent
noise around its expected value for large $N$. Also, note that the factor
$\|\tilde{f}\|_{L_{\rho}^{2}}^{2}=\sum_{k=1}^{M}\tilde{\gamma}_{k}^{2}+\tilde{\gamma}_{>M}^{2}$
depends solely on the target function $\tilde{f}$ but not on the kernel $K$ or
ridge parameter $\lambda$. 3) The second plus the fourth terms
$\left(1+2\sqrt{\frac{\log N}{N}}\right)\lambda\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}$ depends strongly on the kernel training: the sum is
proportional to the ridge $\lambda$, and the RKHS norm square
$\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}=\sum_{k=1}^{M}\frac{\tilde{\gamma}_{k}^{2}}{\lambda_{k}}$
measures how well-aligned the target functions with the chosen kernel $K$ is.
Again, a larger $N$ favors the control as the residue $\sqrt{\frac{\log
N}{N}}\to 0$ as $N\to\infty$. 4) The fifth term $C_{1}\frac{\log N}{N}$ is a
small residue term with fast decay rate, which the other terms will overshadow
as $N\to\infty$.

##### Ridge parameter

The bounds in Eq. (12) demonstrate that in the ridgeless case, the bias can be
controlled solely by the finite rank error $\tilde{\gamma}_{k}^{2}$ with a
confidence level depending on $N$; also the upper and lower bounds coincide as
$N\to\infty$.

##### Variance

For the variance bounds in Eq. (10), we have similar results: the variance can
be controlled solely by the (sample-independent) noise variance $\sigma^{2}$
with a confidence level depending on $N$, also the upper and lower bounds
coincides as $N\to\infty$.

###  5.2 Comparison with Prior Works

We discuss how our results add to, and improve on, what is known about finite-
rank KRRs following the presentation in Table 1.

##### Classical tools to study generalization

The test error measures the average difference between the trained model and
the target function. It is one of the most common measures to study the
generalization performance of a machine learning model. Nevertheless,
generally applicable statistical learning-theoretic tools from VC-theory
(Bartlett and Mendelson, 2002), (local) Rademacher complexities (Kakade et
al., 2008; Bartlett et al., 2005; Mohri et al., 2018), PAC-Bayes methods
(Alquier et al., 2016; Mhammedi et al., 2019), or smooth optimal transport
theory (Hou et al., 2022) all yield pessimistic bounds on the KRR test error.
For example, Mohri et al. Mohri et al. (2018) bound the generalization error
using Rademacher Complexity: there exists some constant $C>0$ independent to
$N$ and $\lambda$ such that, with a probability of at least $1-\tau$:

$\text{test error}-\text{train
error}\leq\frac{C}{\sqrt{N}}\left(1+\frac{1}{2}\sqrt{\frac{\log\frac{1}{\tau}}{2}}\right).$
(14)

If we set $\tau$ to $2/N$, the decay in the generalization gap is
$\mathcal{O}\left(\sqrt{\frac{\log N}{N}}\right)$, which is too slow compared
to other kernel-specific analyses.

##### Truncated KRR

Amini et al. Amini et al. (2022) suggests an interesting type of finite-rank
kernels: for any (infinite-rank) kernel $K^{(\infty)}$, fix a sample
$\mathbf{Z}=(\mathbf{X},\mathbf{y})$ of size $N$ and define a rank-$N$ kernel
$K$ by the eigendecomposition of the kernel matrix $\mathbf{K}^{(\infty)}$.
Note that different random sample results in a completely different rank-$N$
kernel $K$. Assume that the target function $\tilde{f}$ lies in the
$N$-dimensional RKHS corresponding to $K$, then one can obtain an exact
formula for the expected value of the test error of the kernel $K$ (but not
only the original $K^{(\infty)}$). However, the formula obtained in Amini et
al. (2022) only takes into account the expectation over the noise $\epsilon$,
but not over the samples $x$. Hence, our paper yields a more general result
for the test error.

##### Upper bound Comparison

To the best of our knowledge, the following result is the closest and most
related to ours:

######  Theorem 5.1.

(Proposition 7.4 in Bach (2023a)) With notation as before, assume, in
addition, that: 1) $\tilde{\gamma}_{>M}=0$, that is,
$\tilde{f}\in\mathcal{H}$; 2) for all $k=1,...,M$, we have
$\sqrt{\sum_{k=1}^{M}\lambda_{k}}\leq R$ ; 3) and
$N\geq\left(\frac{4}{3}+\frac{R^{2}}{8\lambda}\log\frac{14R^{2}}{\lambda\tau}\right)$.
Then we have, with a probability of at least $1-\tau$,

bias $\displaystyle\leq 4\lambda\|\tilde{f}\|_{\mathcal{H}}^{2};$ variance
$\displaystyle\leq\frac{8\sigma^{2}R^{2}}{\lambda
N}\left(1+2\log\frac{2}{\tau}\right).$

In comparison to Bach (2023a), Theorem 4.2 makes the following significant
improvements:

  1. (i)

In Theorem 5.1, as the ridge parameter $\lambda\to 0$, the minimum requirement
of $N$ and the upper bound of the variance explode; while in our case, both
the minimum requirement of $N$ and the bounds on the variance are independent
of $\lambda$.

  2. (ii)

While Bach (2023a) also mentioned the case where $\tilde{\gamma}_{>M}\neq 0$,
however their bound for the inconsistent case is implicit; 444For inconsistent
cases (or misspecified model), Bach (2023a) only derives an upper bound on the
expected value of test error. Also, the bound involves a function defined as
an infimum of some set. while our work clearly states the impact of
$\tilde{\gamma}_{>M}$ on the test error.

  3. (iii)

Our upper bounds are sharper than Bach (2023a):

    1. (a)

For the bias, assume $\tilde{\gamma}_{>M}=0$ for comparison purposes. Then by
line (13) our upper bound would become:

$\text{bias}\leq\lambda\|\tilde{f}\|_{\mathcal{H}}^{2}\left(1+2\sqrt{\frac{\log
N}{N}}\right)+C_{1}\frac{\log N}{N},$

which improves Bach (2023a)’s upper bound $\text{bias}\leq
4\lambda\|\tilde{f}\|_{\mathcal{H}}^{2}$ significantly. First, we observe a
minor improvement in the constant (by a factor of $\frac{1}{4}$) when
$N\to\infty$. Second and more importantly, for finite $N$, we observe a non-
asymptotic decay on the upper bound, while Bach (2023a)’s upper bound is
insensitive to $N$.

    2. (b)

For the variance, if we replace $\tau$ by $2/N$ in Bach (2023a)’s upper bound:
with probability of at least $1-2/N$,

$\text{variance}\leq\frac{8\sigma^{2}R^{2}}{\lambda N}(1+2\log N),$ (15)

which has a decay rate of $\frac{\log N}{N}$, our work shows a much faster
decay of $\frac{1}{N}$. Moreover, we show that the "unscaled" variance, that
is $N\cdot(\text{variance})$, is actually bounded by
$\mathcal{O}\left(1+\sqrt{\frac{\log N}{N}}\right)$, while Bach (2023a)’s
upper bound on it would explode.

  4. (iv)

We also provide a lower bound for the test error that matches with the upper
bound when $N\to\infty$, while Bach (2023a) does not derive any lower bound.

While our work offers overall better bounds, there are some technical
differences between the bounds of Bach (2023a) and ours:

  1. (i)

We require some explicit mild condition on input distribution $\rho$ and on
kernel $K$, while Bach (2023a) rely on the assumption
$\sqrt{\sum_{k=1}^{M}\lambda_{k}}\leq R$ which gives an alternative implicit
requirement;

  2. (ii)

Our work eventually assumes under-parametrization $M<N$, while Bach (2023a)
does not requires this. Essentially, we exploit the assumption on the under-
parameterized regime to obtain a sharper bound than Bach (2023a). In short, we
select the $L_{\rho}^{2}$ basis $\psi_{k}$ for analysis while Bach (2023a)
selects the basis in RKHS, which eventually affects the assumptions and
conclusions. For more details, see Appendix C.

##### Nyström Subsampling and Random Feature Model

Nyström Subsampling and Random feature models are two setting closely related
to our work. Rudi et al. (2015); Rudi and Rosasco (2017) bound the test error
from above where the regularizer $\lambda=\lambda(N)$ decays as $N\to\infty$.
In comparison, one main contribution of our work is that we provide both
tighter upper bound and tighter lower bound than these prior works (they
derive the same convergence rate on the upper bound up to constants but they
do not derive a lower bound). Another major difference is that our bounds work
for any regularization independent to $N$.

##### Effective Dimension

In Rudi et al. (2015); Rudi and Rosasco (2017) and other related works on
kernel, the quantity
$\mathcal{N}(\lambda)=\sum_{k=1}^{M}\frac{\lambda_{k}}{\lambda_{k}+\lambda}=\operatorname{Tr}\left[\bar{\mathbf{P}}\right]$
is called the effective dimension and it appeared as a factor in the upper
bound. In our work, we can define a related quantity
$\mathcal{N}^{2}(\lambda)=\sum_{k=1}^{M}\frac{\lambda_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}=\operatorname{Tr}\left[\bar{\mathbf{P}}^{2}\right]$,
which appears in our bound (See Proposition C.14 for details.). Note that
$\mathcal{N}^{2}(\lambda)\leq\mathcal{N}(\lambda)\leq M$. Indeed, we can
sharpen the factor $\frac{M}{N}$ in line (10) to
$\frac{\mathcal{N}^{2}(\lambda)}{N}$.

###  5.3 Experiments

We examine the test error and showcase the validity of the bounds we derived
for two different finite-rank kernels: the truncated neural tangent kernel
(tNTK) and the Legendre kernel (LK) (see details of their definitions in
Appendix D). Figure 1 shows the true test error as a function of the sample
size $N$. Figure 2 plots our upper bound compared to Bach (2023a), clearly
demonstrating significant improvements. We for instance note that the bound is
tighter for any value of $\lambda$. Finally, Figure 9 shows both the upper and
lower bound closely enclose the true test error. For more details on the
experiment setup, including evaluating the lower bound, please see Appendix D.

Figure 1: KRR on a finite-rank kernel. Left: KRR training; right: test error
for varying $N$, over 10 iterations. The upper and lower quartiles are used as
error bars.

###  5.4 Future Research Direction

An interesting extension would be to combine our results with Tsigler and
Bartlett (2023); Bordelon et al. (2020) in the over-parameterized regime to
give an accurate transition and observe double descent. However, since we are
interested in the sharpest possible bounds for the under-parameterized regime,
we treat these cases separately. We will do the same for the over-
parameterized regime in future work. A special case is where the
parameterization is at the threshold $M=N$ and $\lambda=0$. The blow-up of
variance as $M\to N$ is well-known, where Bach (2023b); Belkin et al. (2019)
gives empirical and theoretical reports on the blow-up for kernels. We expect
that we can formally prove this result for any finite-rank kernels using our
exact formula in Proposition C.7 and some anti-concentration inequalities on
random matrices.

Figure 2: Comparison of test error bounds. Left: upper bounds for varying $N$;
right: upper bounds for varying $\lambda$. In natural log scale. Figure 3:
Our bounds (u.b. = upper bound, l.w. = lower bound) comparing to the test
error with varying $N$, over 10 iterations, in natural log scale. The residues
with coefficients $C_{1}$ and $C_{2}$ are dropped by simplicity. As a result,
the bounds are not ‘bounding’ the averaged test error for small $N$. But for
large $N$, the residues become negligible.

## References

  * Jacot et al. [2018] Arthur Jacot, Franck Gabriel, and Clément Hongler.  Neural tangent kernel: Convergence and generalization in neural networks.  _Advances in neural information processing systems_ , 31, 2018.
  * Arora et al. [2019] Sanjeev Arora, Simon Du, Wei Hu, Zhiyuan Li, and Ruosong Wang.  Fine-grained analysis of optimization and generalization for overparameterized two-layer neural networks.  In _International Conference on Machine Learning_ , pages 322–332. PMLR, 2019.
  * Bordelon et al. [2020] Blake Bordelon, Abdulkadir Canatar, and Cengiz Pehlevan.  Spectrum dependent learning curves in kernel regression and wide neural networks.  In _International Conference on Machine Learning_ , pages 1024–1034. PMLR, 2020.
  * Liang and Rakhlin [2020] Tengyuan Liang and Alexander Rakhlin.  Just interpolate: Kernel “ridgeless” regression can generalize.  _The Annals of Statistics_ , 48(3), Jun 2020.  ISSN 0090-5364.  doi:10.1214/19-aos1849.  URL http://dx.doi.org/10.1214/19-AOS1849.
  * Liu et al. [2021a] Fanghui Liu, Zhenyu Liao, and Johan Suykens.  Kernel regression in high dimensions: Refined analysis beyond double descent.  In Arindam Banerjee and Kenji Fukumizu, editors, _Proceedings of The 24th International Conference on Artificial Intelligence and Statistics_ , volume 130 of _Proceedings of Machine Learning Research_ , pages 649–657. PMLR, 13–15 Apr 2021a.  URL https://proceedings.mlr.press/v130/liu21b.html.
  * Mei et al. [2021] Song Mei, Theodor Misiakiewicz, and Andrea Montanari.  Generalization error of random feature and kernel methods: hypercontractivity and kernel matrix concentration.  _Applied and Computational Harmonic Analysis_ , 2021.
  * Misiakiewicz [2022] Theodor Misiakiewicz.  Spectrum of inner-product kernel matrices in the polynomial regime and multiple descent phenomenon in kernel ridge regression, 2022.  URL https://arxiv.org/abs/2204.10425.
  * Bach [2023a] Francis Bach.  Learning theory from first principles.  _Online version_ , 2023a.  URL https://www.di.ens.fr/~fbach/ltfp_book.pdf.
  * Amini et al. [2022] Arash Amini, Richard Baumgartner, and Dai Feng.  Target alignment in truncated kernel ridge regression.  _Advances in Neural Information Processing Systems_ , 35:21948–21960, 2022.
  * Rahimi and Recht [2007] Ali Rahimi and Benjamin Recht.  Random features for large-scale kernel machines.  _Advances in neural information processing systems_ , 20, 2007.
  * Rudi and Rosasco [2017] Alessandro Rudi and Lorenzo Rosasco.  Generalization properties of learning with random features.  _Advances in neural information processing systems_ , 30, 2017.
  * Liu et al. [2021b] Fanghui Liu, Xiaolin Huang, Yudong Chen, and Johan AK Suykens.  Random features for kernel approximation: A survey on algorithms, theory, and beyond.  _IEEE Transactions on Pattern Analysis and Machine Intelligence_ , 44(10):7128–7148, 2021b.
  * Gonon et al. [2022] L. Gonon, L. Grigoryeva, and J.P Ortega.  Approximation bounds for random neural networks and reservoir systems.  _Annals of Applied Probability_ , 2022.
  * Herrera et al. [2021] Calypso Herrera, Florian Krach, Pierre Ruyssen, and Josef Teichmann.  Optimal stopping via randomized neural networks.  _arXiv preprint arXiv:2104.13669_ , 2021.
  * Gonon et al. [2020] Lukas Gonon, Lyudmila Grigoryeva, and Juan-Pablo Ortega.  Risk bounds for reservoir computing.  _Journal of Machine Learning Research (JMLR)_ , 21, 2020.
  * Gonon et al. [2023] Lukas Gonon, Lyudmila Grigoryeva, and Juan-Pablo Ortega.  Infinite-dimensional reservoir computing.  _arXiv preprint arXiv:2304.00490_ , 2023.
  * Compagnoni et al. [2023] Enea Monzio Compagnoni, Anna Scampicchio, Luca Biggio, Antonio Orvieto, Thomas Hofmann, and Josef Teichmann.  On the effectiveness of randomized signatures as reservoir for learning rough dynamics.  _IEEE IJCNN 2023_ , 2023.
  * Azizpour et al. [2015] Hossein Azizpour, Ali Sharif Razavian, Josephine Sullivan, Atsuto Maki, and Stefan Carlsson.  Factors of transferability for a generic convnet representation.  _IEEE transactions on pattern analysis and machine intelligence_ , 38(9):1790–1802, 2015.
  * Vinyals et al. [2016] Oriol Vinyals, Charles Blundell, Timothy Lillicrap, Daan Wierstra, et al.  Matching networks for one shot learning.  _Advances in neural information processing systems_ , 29, 2016.
  * Yosinski et al. [2014] Jason Yosinski, Jeff Clune, Yoshua Bengio, and Hod Lipson.  How transferable are features in deep neural networks?  _Advances in neural information processing systems_ , 27, 2014.
  * LI et al. [2018] Xuhong LI, Yves Grandvalet, and Franck Davoine.  Explicit inductive bias for transfer learning with convolutional networks.  In Jennifer Dy and Andreas Krause, editors, _Proceedings of the 35th International Conference on Machine Learning_ , volume 80 of _Proceedings of Machine Learning Research_ , pages 2825–2834. PMLR, 10–15 Jul 2018.  URL https://proceedings.mlr.press/v80/li18a.html.
  * Ju et al. [2022] Haotian Ju, Dongyue Li, and Hongyang R Zhang.  Robust fine-tuning of deep neural networks with hessian-based generalization guarantees.  In Kamalika Chaudhuri, Stefanie Jegelka, Le Song, Csaba Szepesvari, Gang Niu, and Sivan Sabato, editors, _Proceedings of the 39th International Conference on Machine Learning_ , volume 162 of _Proceedings of Machine Learning Research_ , pages 10431–10461. PMLR, 17–23 Jul 2022.  URL https://proceedings.mlr.press/v162/ju22a.html.
  * Donahue et al. [2014] Jeff Donahue, Yangqing Jia, Oriol Vinyals, Judy Hoffman, Ning Zhang, Eric Tzeng, and Trevor Darrell.  Decaf: A deep convolutional activation feature for generic visual recognition.  In Eric P. Xing and Tony Jebara, editors, _Proceedings of the 31st International Conference on Machine Learning_ , volume 32 of _Proceedings of Machine Learning Research_ , pages 647–655, Bejing, China, 22–24 Jun 2014. PMLR.
  * Mohri et al. [2018] Mehryar Mohri, Afshin Rostamizadeh, and Ameet Talwalkar.  _Foundations of machine learning_.  MIT press, 2018.
  * Jacot et al. [2020] Arthur Jacot, Berfin Simsek, Francesco Spadaro, Clément Hongler, and Franck Gabriel.  Implicit regularization of random feature models.  In _International Conference on Machine Learning_ , pages 4631–4640. PMLR, 2020.
  * Vershynin [2010] Roman Vershynin.  Introduction to the non-asymptotic analysis of random matrices.  _arXiv preprint arXiv:1011.3027_ , 2010.
  * Tsigler and Bartlett [2023] Alexander Tsigler and Peter L. Bartlett.  Benign overfitting in ridge regression.  _Journal of Machine Learning Research_ , 24(123):1–76, 2023.  URL http://jmlr.org/papers/v24/22-1398.html.
  * Ronen et al. [2019] Basri Ronen, David Jacobs, Yoni Kasten, and Shira Kritchman.  The convergence rate of neural networks for learned functions of different frequencies.  _Advances in Neural Information Processing Systems_ , 32, 2019.
  * Bartlett and Mendelson [2002] Peter L. Bartlett and Shahar Mendelson.  Rademacher and Gaussian complexities: risk bounds and structural results.  _J. Mach. Learn. Res._ , 3(Spec. Issue Comput. Learn. Theory):463–482, 2002.  ISSN 1532-4435.  doi:10.1162/153244303321897690.  URL https://doi.org/10.1162/153244303321897690.
  * Kakade et al. [2008] Sham M Kakade, Karthik Sridharan, and Ambuj Tewari.  On the complexity of linear prediction: Risk bounds, margin bounds, and regularization.  _Advances in neural information processing systems_ , 21, 2008.
  * Bartlett et al. [2005] Peter L. Bartlett, Olivier Bousquet, and Shahar Mendelson.  Local Rademacher complexities.  _Ann. Statist._ , 33(4):1497–1537, 2005.  ISSN 0090-5364.  doi:10.1214/009053605000000282.  URL https://doi.org/10.1214/009053605000000282.
  * Alquier et al. [2016] Pierre Alquier, James Ridgway, and Nicolas Chopin.  On the properties of variational approximations of Gibbs posteriors.  _J. Mach. Learn. Res._ , 17:Paper No. 239, 41, 2016.  ISSN 1532-4435.
  * Mhammedi et al. [2019] Zakaria Mhammedi, Peter Grünwald, and Benjamin Guedj.  Pac-bayes un-expected bernstein inequality.  _Advances in Neural Information Processing Systems_ , 32, 2019.
  * Hou et al. [2022] Songyan Hou, Parnian Kassraie, Anastasis Kratsios, Jonas Rothfuss, and Andreas Krause.  Instance-dependent generalization bounds via optimal transport.  _arXiv preprint arXiv:2211.01258_ , 2022.
  * Rudi et al. [2015] Alessandro Rudi, Raffaello Camoriano, and Lorenzo Rosasco.  Less is more: Nyström computational regularization.  _Advances in Neural Information Processing Systems_ , 28, 2015.
  * Bach [2023b] Francis Bach.  High-dimensional analysis of double descent for linear regression with random projections.  _arXiv preprint arXiv:2303.01372_ , 2023b.
  * Belkin et al. [2019] Mikhail Belkin, Daniel Hsu, Siyuan Ma, and Soumik Mandal.  Reconciling modern machine-learning practice and the classical bias–variance trade-off.  _Proceedings of the National Academy of Sciences_ , 116(32):15849–15854, 2019.
  * De Vito et al. [2006] Eernesto De Vito, Lorenzo Rosasco, and Andrea Caponnetto.  Discretization error analysis for tikonov regularization.  _Analysis and Applications_ , 04(01):81–99, 2006\.  doi:10.1142/S0219530506000711.  URL https://doi.org/10.1142/S0219530506000711.
  * Wainwright [2019] Martin J. Wainwright.  _High-Dimensional Statistics: A Non-Asymptotic Viewpoint_.  Cambridge Series in Statistical and Probabilistic Mathematics. Cambridge University Press, 2019.  doi:10.1017/9781108627771.
  * Petersen and Pedersen [2012] K. B. Petersen and M. S. Pedersen.  The matrix cookbook, nov 2012.  URL http://www2.compute.dtu.dk/pubdb/pubs/3274-full.html.  Version 20121115.

Appendix

##  Appendix A Glossary

Table 2: Glossary Name | Notation | Expression | Dimension
---|---|---|---
sampling distribution | $\rho$ | - | $\mathcal{X}\to\mathbb{R}^{+}$
sampling size | $N$ | - | integer
input matrix | $\mathbf{X}$ | $(x_{i})_{i=1}^{N}\underset{iid}{\sim}\rho_{\mathcal{X}}$ | $N\times d$
output vector | $\mathbf{y}$ | $(y_{i})_{i=1}^{N}$ | $N\times 1$
sample | $\mathbf{Z}$ | $(\mathbf{X},\mathbf{y})$ | $N\times(d+1)$
noise | $\varepsilon$ | - | random scalar
noise variance | $\sigma^{2}$ | $\mathbb{E}[\varepsilon^{2}]$ | scalar
ridge | $\lambda$ | - | scalar
finite-rank kernel | $K$ | $\sum_{k=1}^{M}\lambda_{k}\psi_{k}(\cdot)\psi_{k}(\cdot)$ | $\mathcal{X}\times\mathcal{X}\to\mathbb{R}$
kernel rank | $M$ | - | integer
$k$th eigenfunction | $\psi_{k}$ | - | $\mathcal{X}\to\mathbb{R}$
$k$th value | $\lambda_{k}$ | - | scalar
- | $\bm{\psi}(x)$ | $[\psi_{k}(x)]_{k=1}^{M}$ | $M\times 1$
- | $\bm{\Psi}$ | $[\psi_{k}(x_{i})]_{k,i}$ | $M\times N$
- | $\bm{\Lambda}$ | $\operatorname{diag}\big{[}\lambda_{k}\big{]}$ | $M\times M$
kernel matrix | $\mathbf{K}$ | $[K(x_{i},x_{j})]_{i,j}=\bm{\Psi}^{\top}\bm{\Lambda}\bm{\Psi}$ | $N\times N$
resolvent | $\mathbf{R}$ | $(\mathbf{K}+\lambda N\mathbf{I}_{N})^{-1}$ | $N\times N$
target function | $\tilde{f}$ | $\sum_{k=1}^{M}\tilde{\gamma}_{k}\psi_{k}+\tilde{\gamma}_{>M}\psi_{>M}$ | $\mathcal{X}\to\mathbb{R}$
- | $\tilde{f}_{\leq M}$ | $\sum_{k=1}^{M}\tilde{\gamma}_{k}\psi_{k}$ | $\mathcal{X}\to\mathbb{R}$
$k$th target coefficient | $\tilde{\gamma}_{k}$ | $\int_{\mathcal{X}}\tilde{f}(x)\psi_{k}(x)d\rho_{\mathcal{X}}(x)$ | scalar
- | $\bm{\gamma}$ | $[\gamma_{k}]$ | $M\times 1$
orthonormal complement | $\psi_{>M}$ | - | $\mathcal{X}\to\mathbb{R}$
complementary coefficient | $\tilde{\gamma}_{>M}$ | - | scalar
- | $\bm{\Psi}_{>M}$ | $[\psi_{>M}(x_{i})]$ | $1\times N$
test error | $\mathcal{R}_{\mathbf{Z},\lambda}$ | $\mathbb{E}_{x,\epsilon}\left[(f_{\mathbf{Z},\lambda}(x)-\tilde{f}(x))^{2}\right]$ | scalar
bias | - | $\int_{\mathcal{X}}\left(f_{(\mathbf{X},\tilde{f}(\mathbf{X})),\lambda}(x)-\tilde{f}(x)\right)^{2}d\rho(x)$ | scalar
variance | - | $\mathcal{R}_{\mathbf{Z},\lambda}-\text{bias}=\mathbb{E}_{x,\varepsilon}\left(\mathbf{K}_{x}^{\top}\mathbf{R}\bm{\varepsilon}\right)^{2}$ | scalar
fluctuation matrix | $\bm{\Delta}$ | $\frac{1}{N}\bm{\Psi}\bm{\Psi}^{\top}-\mathbf{I}_{M}$ | $M\times M$
fluctuation | $\delta$ | $\|\bm{\Delta}\|_{\text{op}}$ | scalar
error vector | $\bm{E}$ | $[\eta_{k}]$ | $M\times 1$
- | $\eta_{k}$ | $\frac{1}{N}\sum_{i=1}^{N}\psi_{k}(x_{i})\psi_{>M}(x_{i})$ | scalar
- | $\mathbf{B}$ | $(\mathbf{I}_{M}+\bm{\Delta}+\lambda\bm{\Lambda}^{-1})^{-1}$ | $M\times M$
- | $\bar{\mathbf{P}}$ | $\operatorname{diag}\left[\frac{\lambda_{k}}{\lambda_{k}+\lambda}\right]$ | $M\times M$

##  Appendix B Classical KRR Theory

In an effort to keep our manuscript as self-contained as possible, we recall
the Mercer decomposition, representer theorem for kernel ridge regression as
well as the form of the bias-variance tradeoff in the KRR context.

###  B.1 Mercer Decomposition

We begin with a general kernel
$K^{(\infty)}:\mathcal{X}\times\mathcal{X}\to\mathbb{R}$.

######  Proposition B.1.

De Vito et al. [2006] Fix a sample distribution $\rho$. Let
$K^{(\infty)}:\mathcal{X}\times\mathcal{X}\to\mathbb{R}$ be a reproducing
kernel with corresponding RKHS $\mathcal{H}^{(\infty)}$. There exists a
decreasing sequence of real numbers $\lambda_{1}\geq\lambda_{2}\geq...$,
called the eigenvalues of the kernel $K^{(\infty)}$; and a sequence of
pairwise-orthonormal functions $\\{\psi_{k}\\}_{k=1}^{\infty}\subset
L_{\rho}^{2}$, called the eigenfunctions of $K^{(\infty)}$, such that for all
$x,x^{\prime}\in\mathcal{X}$, we have

$\displaystyle K^{(\infty)}(x,x^{\prime})$
$\displaystyle=\sum_{k=1}^{\infty}\lambda_{k}\psi_{k}(x)\psi_{k}(x^{\prime})$
(16)

In particular, we assume $\lambda_{k}=0,\ \forall k>M$. In this case, we say
the kernel
$K(x,x^{\prime})=\sum_{k=1}^{M}\lambda_{k}\psi_{k}(x)\psi_{k}(x^{\prime})$ is
of finite rank $M$ with corresponding (finite-dimensional) RKHS $\mathcal{H}$,
recovering equation (2).

The first of these results, allows us to explicitly express the finite-rank
kernel ridge regressor $f_{\mathbf{Z},\lambda}$.

######  Proposition B.2 (Representer Theorem - [Wainwright, 2019, Chapter
12]).

Let $\mathbf{R}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}(\mathbf{K}+\lambda
N\mathbf{I}_{N})^{-1}\in\mathbb{R}^{N\times N}$ be the resolvent matrix and
recall the kernel ridge regressor $f_{\mathbf{Z},\lambda}$ given by equation
(3):

$f_{\mathbf{Z},\lambda}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\operatorname*{arg\,min}_{f\in\mathcal{H}}\frac{1}{N}\sum_{i=1}^{N}\left(f(x_{i})-y_{i}\right)^{2}+\lambda\|f\|_{\mathcal{H}}^{2}$

Then, for every $x\in\mathcal{X}$, we have the expression

$f_{\mathbf{Z},\lambda}(x)=\mathbf{y}^{\top}\mathbf{R}\mathbf{K}_{x},\ \forall
x\in\mathcal{X},$ (17)

where
$\mathbf{K}_{x}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}[K(x_{i},x)]_{i=1}^{N}\in\mathbb{R}^{N\times
1}$.

###  B.2 Compact Matrix Expression

First, let
$\bm{\Psi}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}(\psi_{k}(x_{i}))_{k=1,i=1}^{M,N}$
be the random $M\times N$ matrix defined by evaluating the $M$ eigenfunctions
on all input training instances
$\bm{X}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}(x_{i})_{i=1}^{N}$,
$\bm{\Lambda}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\operatorname{diag}[\lambda_{k}]\in\mathbb{R}^{M\times
M}$, and
$\bm{\psi}(x)\stackrel{{\scriptstyle\mbox{def.}}}{{=}}[\psi_{k}(x)]_{k=1}^{M}\in\mathbb{R}^{M\times
1}$. The advantage of this notation is that we can rewrite the equations in a
more compact form. For equation (17):

$f_{\mathbf{Z},\lambda}(x)=\mathbf{y}^{\top}\underbrace{(\bm{\Psi}^{\top}\bm{\Lambda}\bm{\Psi}+\lambda
N\mathbf{I}_{M})^{-1}}_{\mathbf{R}}\bm{\Psi}^{\top}\bm{\Lambda}\bm{\psi}(x);$
(18)

for equation (4):

$\tilde{f}(x)=\tilde{\bm{\gamma}}^{\top}\bm{\psi}(x)+\tilde{\gamma}_{>M}\psi_{>M}(x).$
(19)

Last but not least, we define some important quantities for later analysis.

######  Definition B.3 (Fluctuation matrix).

The fluctuation matrix is the random $M\times M$-matrix given by
$\bm{\Delta}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\frac{1}{N}\bm{\Psi}\bm{\Psi}^{\top}-\mathbf{I}_{M}.$
Our analysis will often involve the operator norm of $\bm{\Delta}$, which we
denote by
$\delta\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\|\bm{\Delta}\|_{\text{op}}$.

The fluctuation matrix $\bm{\Delta}$ measures the first source of randomness
in the KRR’s test error. Namely it encodes the degree of non-orthonormality
between the vectors obtained by evaluating of the $M$ eigenfunctions
$\psi_{1},\dots,\psi_{M}$ on the input $\mathbf{X}$.

The second source of randomness in the KRR’s test error comes from the
empirical evaluation of the dot product of the eigenfunction $\psi_{k}$’s and
the orthogonal complement $\psi_{>M}$:

######  Definition B.4 (Error Vector).

$\bm{E}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\frac{1}{N}\bm{\Psi}\psi_{>M}(\mathbf{X})$
is called the error vector.

The random matrix $\bm{\Delta}$ and the random vector $\bm{E}$ are centered;
i.e. $\mathbb{E}_{\mathbf{X}}[\bm{\Delta}]=0$ and
$\mathbb{E}_{\mathbf{X}}[\bm{E}]=0$.

###  B.3 Bias-Variance Decomposition

The derivation of several contemporary KRR generalization bounds [Bach, 2023a,
Liu et al., 2021a, Mei et al., 2021] involves the classical Bias-Variance
Trade-off:

######  Proposition B.5 (Bias-Variance Trade-off).

Fix a sample $\mathbf{Z}$. Recall the definition 3.4 of test error
$\mathcal{R}_{\mathbf{Z},\lambda}$, bias, and variance:

$\displaystyle\mathcal{R}_{\mathbf{Z},\lambda}$
$\displaystyle\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\mathbb{E}_{x,\epsilon}\left[(f_{\mathbf{Z},\lambda}(x)-\tilde{f}(x))^{2}\right]=\mathbb{E}_{\epsilon}\left[\int_{\mathcal{X}}\left(f_{\mathbf{Z},\lambda}(x)-\tilde{f}(x)\right)^{2}d\rho(x)\right];$
bias
$\displaystyle\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\int_{\mathcal{X}}\left(f_{(\mathbf{X},\tilde{f}(\mathbf{X})),\lambda}(x)-\tilde{f}(x)\right)^{2}d\rho(x);$
variance
$\displaystyle\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\mathcal{R}_{\mathbf{Z},\lambda}-\text{bias}.$

Then, we can write
$\text{variance}_{\text{test}}=\mathbb{E}_{x,\varepsilon}\left(\mathbf{K}_{x}^{\top}\mathbf{R}\bm{\varepsilon}\right)^{2}$
and hence the test error $\mathcal{R}_{\mathbf{Z},\lambda}$ admits a
decomposition:

$R_{\mathbf{Z},\lambda}=\text{bias}+\mathbb{E}_{x,\varepsilon}\left(\mathbf{K}_{x}^{\top}\mathbf{R}\bm{\varepsilon}\right)^{2}.$

###### Proof.

See the proof of Theorem C.8. ∎

##  Appendix C Proofs

In this section, we will derive the essential lemmata and propositions for
proving the main theorems.

###  C.1 Formula Derivation

We begin with writing the test error in convenient forms.

####  C.1.1 Bias

We first derive, from the definition of the bias, a convenient expression to
proceed:

######  Proposition C.1 (Bias Expression).

Let
$\bm{\Psi}_{>M}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}[\psi_{>M}(x_{i})]_{i=1}^{N}$
as an $1\times N$\- row vector, $\binom{\bm{\Psi}}{\bm{\Psi}_{>M}}$ as an
$(M+1)\times N$ matrix. Denote
$\mathbf{P}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\begin{pmatrix}\mathbf{P}_{\leq
M}&\mathbf{P}_{>M}\end{pmatrix}=\bm{\Lambda}\bm{\Psi}\mathbf{R}\begin{pmatrix}\bm{\Psi}^{\top}&\bm{\Psi}_{>M}^{\top}\end{pmatrix}\in\mathbb{R}^{M\times(M+1)}$
, $\mathbf{P}_{\leq M}\in\mathbb{R}^{M\times M}$ and
$\mathbf{P}_{>M}\in\mathbb{R}^{M\times 1}$. Then the bias admits the following
expression:

bias $\displaystyle=\underbrace{\tilde{\gamma}_{>M}^{2}}_{\mbox{Finite Rank
Error}}+\underbrace{\|\tilde{\bm{\gamma}}-\mathbf{P}_{\leq
M}\tilde{\bm{\gamma}}-\tilde{{\gamma}}_{>M}\mathbf{P}_{>M}\|_{2}^{2}}_{\mbox{Fitting
Error}}.$

###### Proof.

Recall that, by equations (18) and (19), we can write

$\displaystyle\tilde{f}(x)$
$\displaystyle=\tilde{\bm{\gamma}}^{\top}\bm{\psi}(x)+\tilde{\gamma}_{>M}\psi_{>M}(x),$
$\displaystyle f_{(\mathbf{X},\tilde{f}(\mathbf{X})),\lambda}(x)$
$\displaystyle=(\tilde{\bm{\gamma}}^{\top}\bm{\Psi}+\tilde{\gamma}_{>M}\bm{\Psi}_{>M}^{\top})\mathbf{R}\bm{\Psi}^{\top}\bm{\Lambda}\bm{\Psi}(x).$

Hence

bias
$\displaystyle=\mathbb{E}_{x}\left[\big{(}\tilde{\bm{\gamma}}^{\top}\bm{\psi}(x)+\tilde{\gamma}_{>M}\psi_{>M}(x)-(\tilde{\bm{\gamma}}^{\top}\bm{\Psi}+\tilde{\gamma}_{>M}\bm{\Psi}_{>M}^{\top})\mathbf{R}\bm{\Psi}^{\top}\bm{\Lambda}\bm{\Psi}(x)\big{)}^{2}\right]$
$\displaystyle=\Big{\|}\binom{\tilde{\bm{\gamma}}}{\tilde{\gamma}_{>M}}-\begin{pmatrix}\mathbf{P}\binom{\tilde{\bm{\gamma}}}{\tilde{\gamma}_{>M}}\\\
0\end{pmatrix}\Big{\|}_{2}^{2}$ (20)
$\displaystyle=\underbrace{\tilde{\gamma}_{>M}^{2}}_{\mbox{Finite Rank
Error}}+\underbrace{\|\tilde{\bm{\gamma}}-\mathbf{P}_{\leq
M}\tilde{\bm{\gamma}}-\tilde{{\gamma}}_{>M}\mathbf{P}_{>M}\|_{2}^{2}}_{\mbox{Fitting
Error}},$

in line (20), we use Parseval’s identity. ∎

We proceed by reformulating the projection matrix $\mathbf{P}$, first with the
left matrix $\mathbf{P}_{\leq M}$:

######  Lemma C.2.

Recall the following notations

$\displaystyle\mathbf{K}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\bm{\Psi}^{\top}\bm{\Lambda}\bm{\Psi},\quad\mathbf{R}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}(\mathbf{K}+\lambda
N\mathbf{I}_{M})^{-1},\quad\mathbf{P}_{\leq
M}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\bm{\Lambda}\bm{\Psi}\mathbf{R}\bm{\Psi}^{\top}.$

Define the symmetric random matrix
$\mathbf{B}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}(\mathbf{I}_{M}+\bm{\Delta}+\lambda\bm{\Lambda}^{-1})^{-1}$.
It holds that

$\displaystyle\mathbf{P}_{\leq M}$
$\displaystyle=\mathbf{I}_{M}-\lambda\mathbf{B}\bm{\Lambda}^{-1}.$

###### Proof.

We first observe that

$\displaystyle\bm{\Psi}\bm{\Psi}^{\top}\mathbf{P}_{\leq M}$
$\displaystyle=\bm{\Psi}\bm{\Psi}^{\top}\bm{\Lambda}\bm{\Psi}\mathbf{R}\bm{\Psi}^{\top}$
(21) $\displaystyle=\bm{\Psi}\mathbf{K}(\mathbf{K}+\lambda
N\mathbf{I}_{M})^{-1}\bm{\Psi}^{\top}$
$\displaystyle=\bm{\Psi}\left(\mathbf{I}_{M}-\lambda N(\mathbf{K}+\lambda
N\mathbf{I}_{M})^{-1}\right)\bm{\Psi}^{\top}$
$\displaystyle=\bm{\Psi}\bm{\Psi}^{\top}-\lambda N\bm{\Psi}(\mathbf{K}+\lambda
N\mathbf{I}_{M})^{-1}\bm{\Psi}^{\top}.$ (22)

From lines (21)- (22) and the definition of the fluctuation matrix
$\bm{\Delta}$ we deduce

$\displaystyle\frac{1}{N}\bm{\Psi}\bm{\Psi}^{\top}(\mathbf{I}_{M}-\mathbf{P}_{\leq
M})$ $\displaystyle=\lambda\bm{\Psi}(\mathbf{K}+\lambda
N\mathbf{I}_{M})^{-1}\bm{\Psi}^{\top}$
$\displaystyle(\mathbf{I}_{M}+\bm{\Delta})(\mathbf{I}_{M}-\mathbf{P}_{\leq
M})$ $\displaystyle=\lambda\bm{\Psi}\mathbf{R}\bm{\Psi}^{\top}$
$\displaystyle(\mathbf{I}_{M}+\bm{\Delta})(\mathbf{I}_{M}-\mathbf{P}_{\leq
M})$ $\displaystyle=\lambda\bm{\Lambda}^{-1}\mathbf{P}_{\leq M}$
$\displaystyle(\bm{\Lambda}+\bm{\Lambda}\bm{\Delta})(\mathbf{I}_{M}-\mathbf{P}_{\leq
M})$ $\displaystyle=\lambda\mathbf{P}_{\leq M}$
$\displaystyle\bm{\Lambda}+\bm{\Lambda}\bm{\Delta}$
$\displaystyle=(\bm{\Lambda}+\bm{\Lambda}\bm{\Delta}+\lambda\mathbf{I}_{M})\mathbf{P}_{\leq
M}.$ (23)

Rearranging (23) and applying the definition of $\bm{B}$ we find that

$\displaystyle\mathbf{P}_{\leq M}$
$\displaystyle=(\bm{\Lambda}+\bm{\Lambda}\bm{\Delta}+\lambda\mathbf{I}_{M})^{-1}(\bm{\Lambda}+\bm{\Lambda}\bm{\Delta})$
(24)
$\displaystyle=\mathbf{I}_{M}-\lambda(\bm{\Lambda}+\bm{\Lambda}\bm{\Delta}+\lambda\mathbf{I}_{M})^{-1}$
$\displaystyle=\mathbf{I}_{M}-\lambda(\bm{\Lambda}+\bm{\Lambda}\bm{\Delta}+\lambda\mathbf{I}_{M})^{-1}\bm{\Lambda}\bm{\Lambda}^{-1}$
$\displaystyle=\mathbf{I}_{M}-\lambda\mathbf{B}\bm{\Lambda}^{-1}.$ (25)

∎

Arguing analogously for the right matrix $\mathbf{P}_{>M}$, we draw the
subsequent similar conclusion.

######  Lemma C.3.

Recall the following notations

$\displaystyle\mathbf{K}$
$\displaystyle\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\bm{\Psi}^{\top}\bm{\Lambda}\bm{\Psi},\quad\mathbf{R}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}(\mathbf{K}+\lambda
N\mathbf{I}_{M})^{-1},\quad\mathbf{P}_{>M}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\bm{\Lambda}\bm{\Psi}\mathbf{R}\bm{\Psi}_{>M}^{\top},$
$\displaystyle\bm{E}$
$\displaystyle\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\frac{1}{N}\bm{\Psi}\bm{\Psi}_{>M}^{\top},\quad\mathbf{B}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}(\mathbf{I}_{M}+\bm{\Delta}+\lambda\bm{\Lambda}^{-1})^{-1}.$

We have that $\mathbf{P}_{>M}=\mathbf{B}\bm{E}.$

###### Proof.

Similarly to (21)- (22) we note that

$\displaystyle\bm{\Psi}\bm{\Psi}^{\top}\mathbf{P}_{>M}$
$\displaystyle=\bm{\Psi}\bm{\Psi}^{\top}\bm{\Lambda}\bm{\Psi}\mathbf{R}\bm{\Psi}_{>M}^{\top}$
$\displaystyle=\bm{\Psi}\mathbf{K}(\mathbf{K}+\lambda
N\mathbf{I}_{M})^{-1}\bm{\Psi}_{>M}^{\top}$
$\displaystyle=\bm{\Psi}\left(\mathbf{I}_{M}-\lambda N(\mathbf{K}+\lambda
N\mathbf{I}_{M})^{-1}\right)\bm{\Psi}_{>M}^{\top}$
$\displaystyle=\bm{\Psi}\bm{\Psi}_{>M}^{\top}-\lambda
N\bm{\Psi}(\mathbf{K}+\lambda N\mathbf{I}_{M})^{-1}\bm{\Psi}_{>M}^{\top}.$

Analogously to the computations in (24)-(25)

$\displaystyle(\mathbf{I}_{M}+\bm{\Delta})\mathbf{P}_{>M}$
$\displaystyle=\bm{E}-\lambda\bm{\Psi}(\mathbf{K}+\lambda
N\mathbf{I}_{M})^{-1}\bm{\Psi}_{>M}^{\top}$
$\displaystyle(\mathbf{I}_{M}+\bm{\Delta})\mathbf{P}_{>M}$
$\displaystyle=\bm{E}-\lambda\bm{\Lambda}^{-1}\bm{\Lambda}\bm{\Psi}(\mathbf{K}+\lambda
N\mathbf{I}_{M})^{-1}\bm{\Psi}_{>M}^{\top}$
$\displaystyle(\mathbf{I}_{M}+\bm{\Delta})\mathbf{P}_{>M}$
$\displaystyle=\bm{E}-\lambda\bm{\Lambda}^{-1}\mathbf{P}_{>M}$
$\displaystyle(\bm{\Lambda}+\bm{\Lambda}\bm{\Delta})\mathbf{P}_{>M}$
$\displaystyle=\bm{\Lambda}\bm{E}-\lambda\mathbf{P}_{>M}$
$\displaystyle(\bm{\Lambda}+\bm{\Lambda}\bm{\Delta}+\lambda\mathbf{I}_{M})\mathbf{P}_{>M}$
$\displaystyle=\bm{\Lambda}\bm{E}$ $\displaystyle\mathbf{P}_{>M}$
$\displaystyle=(\bm{\Lambda}+\bm{\Lambda}\bm{\Delta}+\lambda\mathbf{I}_{M})^{-1}\bm{\Lambda}\bm{E}$
$\displaystyle\mathbf{P}_{>M}$ $\displaystyle=\mathbf{B}\bm{E}.$

∎

######  Lemma C.4 (Fitting Error).

Recall the notation

fitting error $\displaystyle=\|\tilde{\bm{\gamma}}-\mathbf{P}_{\leq
M}\tilde{\bm{\gamma}}-\tilde{\gamma}_{>M}\mathbf{P}_{>M}\|_{2}^{2},$
$\displaystyle\mathbf{B}$
$\displaystyle\stackrel{{\scriptstyle\mbox{def.}}}{{=}}(\mathbf{I}_{M}+\bm{\Delta}+\lambda\bm{\Lambda}^{-1})^{-1}.$

We have $\text{fitting
error}=\left\|\mathbf{B}\left(\lambda\bm{\Lambda}^{-1}\tilde{\bm{\gamma}}-\bm{E}\tilde{\bm{\gamma}}_{>M}\right)\right\|_{2}^{2}.$

###### Proof.

By lemmata C.2 and C.3,

$\displaystyle\|\tilde{\bm{\gamma}}-\mathbf{P}_{\leq
M}\tilde{\bm{\gamma}}-\tilde{\gamma}_{>M}\mathbf{P}_{>M}\|_{2}^{2}$
$\displaystyle=\left\|\tilde{\bm{\gamma}}-(\mathbf{I}_{M}-\lambda\mathbf{B}\bm{\Lambda}^{-1})\tilde{\bm{\gamma}}\tilde{\bm{\gamma}}-\mathbf{B}\bm{E}\tilde{\bm{\gamma}}_{>M}\right\|_{2}^{2}$
$\displaystyle=\left\|\mathbf{B}\left(\lambda\bm{\Lambda}^{-1}\tilde{\bm{\gamma}}-\bm{E}\tilde{\bm{\gamma}}_{>M}\right)\right\|_{2}^{2}.$

∎

Hence we come up with a new expression of the bias:

######  Proposition C.5 (Bias).

Recall that
$\mathbf{B}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}(\mathbf{I}_{M}+\bm{\Delta}+\lambda\bm{\Lambda}^{-1})^{-1}$.
The bias
$\mathbb{E}_{x}\big{(}f^{\lambda}_{\mathbf{X}}(x)-\tilde{f}(x)\big{)}^{2}$ has
the following expression:

bias
$\displaystyle=\tilde{\gamma}_{>M}^{2}+\left\|\mathbf{B}\left(\lambda\bm{\Lambda}^{-1}\tilde{\bm{\gamma}}-\tilde{\gamma}_{>M}\bm{E}\right)\right\|_{2}^{2}.$

###### Proof.

We apply Proposition C.1 and Lemma C.4 to obtain the result. ∎

####  C.1.2 Variance

If we consider noise in the label, we have to compute the variance part of the
test error.

######  Proposition C.6 (Variance Expression).

Define

$\displaystyle\mathbf{M}$
$\displaystyle\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\mathbb{E}_{x}[\mathbf{K}_{x}\mathbf{K}_{x}^{\top}]$
$\displaystyle=\mathbb{E}_{x}[\bm{\Psi}^{\top}\bm{\Lambda}\bm{\psi}(x)\bm{\psi}(x)^{\top}\bm{\Lambda}\bm{\Psi}]$
$\displaystyle=\bm{\Psi}^{\top}\bm{\Lambda}\mathbb{E}_{x}[\bm{\psi}(x)\bm{\psi}(x)^{\top}]\bm{\Lambda}\bm{\Psi}$
$\displaystyle=\bm{\Psi}^{\top}\bm{\Lambda}\mathbf{I}_{M}\bm{\Lambda}\bm{\Psi}$
$\displaystyle=\bm{\Psi}^{\top}\bm{\Lambda}^{2}\bm{\Psi}.$

We can further simplify the variance part:

variance
$\displaystyle\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\mathbb{E}_{x,\varepsilon}\left[\left(\mathbf{K}_{x}^{\top}\mathbf{R}\bm{\varepsilon}\right)^{2}\right]$
$\displaystyle=\mathbb{E}_{x,\varepsilon}\left[\bm{\varepsilon}^{\top}\mathbf{R}\mathbf{K}_{x}\mathbf{K}_{x}^{\top}\mathbf{R}\bm{\varepsilon}\right]$
$\displaystyle=\mathbb{E}_{\varepsilon}\left[\bm{\varepsilon}^{\top}\mathbf{R}\mathbf{M}\mathbf{R}\bm{\varepsilon}\right]$
$\displaystyle=\sigma^{2}\cdot\operatorname{Tr}[\mathbf{R}\mathbf{M}\mathbf{R}].$

######  Theorem C.7 (Variance).

Recall that
$\mathbf{B}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}(\mathbf{I}_{M}+\bm{\Delta}+\lambda\bm{\Lambda}^{-1})^{-1}$.
The variance part, variance, can be expressed as:

variance
$\displaystyle=\frac{\sigma^{2}}{N}\operatorname{Tr}\left[\mathbf{B}^{2}(\mathbf{I}_{M}+\bm{\Delta})\right].$

###### Proof.

We argue similarly as in lemma C.2. Since

$\displaystyle\bm{\Psi}\bm{\Psi}^{\top}\bm{\Lambda}\bm{\Psi}\mathbf{R}$
$\displaystyle=\bm{\Psi}\mathbf{K}(\mathbf{K}+\lambda N\mathbf{I}_{M})^{-1}$
$\displaystyle=\bm{\Psi}(\mathbf{I}_{M}-\lambda N\mathbf{R})$
$\displaystyle=\bm{\Psi}-\lambda N\bm{\Psi}\mathbf{R},$

therefore, we deduce that

$\displaystyle(\mathbf{I}_{M}+\bm{\Delta})\bm{\Lambda}\bm{\Psi}\mathbf{R}$
$\displaystyle=\frac{1}{N}\bm{\Psi}-\lambda\bm{\Psi}\mathbf{R}$ (26)
$\displaystyle(\mathbf{I}_{M}+\bm{\Delta})\bm{\Lambda}\bm{\Psi}\mathbf{R}$
$\displaystyle=\frac{1}{N}\bm{\Psi}-\lambda\bm{\Lambda}^{-1}\bm{\Lambda}\bm{\Psi}\mathbf{R}$
$\displaystyle(\mathbf{I}_{M}+\bm{\Delta}+\lambda\bm{\Lambda}^{-1})\bm{\Lambda}\bm{\Psi}\mathbf{R}$
$\displaystyle=\frac{1}{N}\bm{\Psi}$
$\displaystyle\bm{\Lambda}\bm{\Psi}\mathbf{R}$
$\displaystyle=\frac{1}{N}(\mathbf{I}_{M}+\bm{\Delta}+\lambda\bm{\Lambda}^{-1})^{-1}\bm{\Psi}$
$\displaystyle\bm{\Lambda}\bm{\Psi}\mathbf{R}$
$\displaystyle=\frac{1}{N}\mathbf{B}\bm{\Psi}.$ (27)

By leveraging the identity
$\mathbf{M}=\bm{\Psi}^{\top}\bm{\Lambda}^{2}\bm{\Psi}$ and elementary
properties of the trace map, the computations in (26)-(27) imply that

$\displaystyle\operatorname{Tr}[\mathbf{R}\mathbf{M}\mathbf{R}]$
$\displaystyle=\operatorname{Tr}[\mathbf{R}\bm{\Psi}^{\top}\bm{\Lambda}^{2}\bm{\Psi}\mathbf{R}]$
(28)
$\displaystyle=\operatorname{Tr}\left[\left(\bm{\Lambda}\bm{\Psi}\mathbf{R}\right)^{\top}\left(\bm{\Lambda}\bm{\Psi}\mathbf{R}\right)\right]$
(29)
$\displaystyle=\operatorname{Tr}\left[\left(\frac{1}{N}\mathbf{B}\bm{\Psi}\right)^{\top}\left(\frac{1}{N}\mathbf{B}\bm{\Psi}\right)\right]$
(30)
$\displaystyle=\frac{1}{N}\operatorname{Tr}\left[\frac{1}{N}\bm{\Psi}^{\top}\mathbf{B}^{\top}\mathbf{B}\bm{\Psi}\right]$
$\displaystyle=\frac{1}{N}\operatorname{Tr}\left[\mathbf{B}^{\top}\mathbf{B}\cdot\frac{1}{N}\bm{\Psi}\bm{\Psi}^{\top}\right]$
(31)
$\displaystyle=\frac{1}{N}\operatorname{Tr}\left[\mathbf{B}^{\top}\mathbf{B}(\mathbf{I}_{M}+\bm{\Delta})\right]$
(32)
$\displaystyle=\frac{1}{N}\operatorname{Tr}\left[\mathbf{B}^{2}(\mathbf{I}_{M}+\bm{\Delta})\right];$
(33)

in more detail: in line (28), we use the definition of $\mathbf{M}$; in line
(29), we use the fact that both $\bm{\Lambda}$ and $\mathbf{R}$ are symmetric;
in line (30), we use line (27); in line (31), we use the cyclicity of the
trace; in line (32), we use the definition of $\bm{\Delta}$; in line (33), we
use the symmetry of $\mathbf{B}$. We obtain the result upon applying Lemma
C.6. ∎

####  C.1.3 Test Error

The Bias-Variance trade-off (see Proposition B.5) decomposed the KRR’s test
error into two terms, the bias and variance. Since Propositions C.5 and C.7
give us exact expressions for the bias and variance, respectively, we deduce
the following exact expression for the KRR’s test error.

######  Theorem C.8 (Exact Formula for KRR’s Test Error).

The test error $\mathcal{R}_{\mathbf{Z},\lambda}$ of KRR equals

$\mathcal{R}_{\mathbf{Z},\lambda}=\underbrace{\overbrace{\left\|\mathbf{B}\left(\lambda\bm{\Lambda}^{-1}\tilde{\bm{\gamma}}-\tilde{{\gamma}}_{>M}\bm{E}_{M}\right)\right\|_{2}^{2}}^{\text{fitting
error}}+\overbrace{\tilde{{\gamma}}_{>M}^{2}}^{\text{finite rank
error}}}_{\mbox{bias}}+\underbrace{\frac{\sigma_{\text{noise}}^{2}}{N}\operatorname{Tr}\left[\mathbf{B}^{2}(\mathbf{I}_{M}+\bm{\Delta})\right]}_{\mbox{variance}},$

where
$\mathbf{B}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}(\mathbf{I}_{M}+\bm{\Delta}+\lambda\bm{\Lambda}^{-1})^{-1}.$

###### Proof.

We begin with the bias/variance decomposition:

$\displaystyle R_{\mathbf{Z}}^{\lambda}$
$\displaystyle\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\mathbb{E}_{y}\|f^{\lambda}_{\mathbf{Z}}-\tilde{f}\|^{2}_{L^{2}_{\rho_{\mathcal{X}}}}$
$\displaystyle=\mathbb{E}_{x,y}\left(\mathbf{K}_{x}^{\top}\mathbf{R}\mathbf{y}-\tilde{f}(x)\right)^{2}$
$\displaystyle=\mathbb{E}_{\varepsilon,x}\left(\mathbf{K}_{x}^{\top}\mathbf{R}(\tilde{f}(\mathbf{X})+\bm{\varepsilon})-\tilde{f}(x)\right)^{2}$
$\displaystyle=\mathbb{E}_{x}\left(f_{\mathbf{X}}^{\lambda}(x)-\tilde{f}(x)\right)^{2}+\mathbb{E}_{x,\varepsilon}\left[\left(\mathbf{K}_{x}^{\top}\mathbf{R}\bm{\varepsilon}\right)^{2}\right]$
$\displaystyle=\text{bias}+\text{variance},$

then we apply Propositions C.5 and C.7. ∎

For the validation of the Theorem C.8, please see Appendix D for details.

The matrix $\mathbf{B}$ plays an important role in the expression since it
encodes most information of the KRR. Therefore, the following subsection will
discuss the approximation of the matrix $\mathbf{B}$.

###  C.2 Matrix Approximation

Recall that the matrix
$\mathbf{B}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}(\mathbf{I}_{M}+\bm{\Delta}+\lambda\bm{\Lambda}^{-1})^{-1}$
is the inverse of a random matrix. The following lemma helps to approximate
$\mathbf{B}$. Informally, it says that: given that
$\delta\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\|\bm{\Delta}\|_{\text{op}}<1$.
We have

$\mathbf{B}=\sum_{s=0}^{\infty}(-\bar{\mathbf{P}}\bm{\Delta})^{s}\bar{\mathbf{P}}$

in operator norm $\|\cdot\|_{op}$ for an $M\times M$ matrix $\bar{P}$
depending only on the $M$ eigenvalues $\\{\lambda_{k}\\}_{k=1}^{M}$ and on the
ridge $\lambda>0$. More precisely we have the following.

######  Lemma C.9 ($\mathbf{B}$-Expansion).

Given that
$\delta\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\|\bm{\Delta}\|_{\text{op}}<1$.
It holds that

$\lim\limits_{n\uparrow\infty}\,\Big{\|}\mathbf{B}-\sum_{s=0}^{n}(-\bar{\mathbf{P}}\bm{\Delta})^{s}\bar{\mathbf{P}}\Big{\|}_{op}=0$

where
$\bar{\mathbf{P}}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\operatorname{diag}\left[\frac{\lambda_{k}}{\lambda_{k}+\lambda}\right]_{k}=\bm{\Lambda}(\bm{\Lambda}+\lambda\mathbf{I}_{M})^{-1}\in\mathbb{R}^{M\times
M}$.

###### Proof.

Set $\mathbf{A}=\mathbf{I}_{M}+\lambda\bm{\Lambda}^{-1}$ and repeatedly use
the formula
$(\mathbf{A}+\bm{\Delta})^{-1}=\mathbf{A}^{-1}-\mathbf{A}^{-1}\bm{\Delta}(\mathbf{A}+\bm{\Delta})^{-1}$
from Petersen and Pedersen [2012], we have

$\displaystyle\mathbf{B}$
$\displaystyle\stackrel{{\scriptstyle\mbox{def.}}}{{=}}(\mathbf{I}_{M}+\bm{\Delta}+\lambda\bm{\Lambda}^{-1})^{-1}$
$\displaystyle=(\mathbf{A}+\bm{\Delta})^{-1}$
$\displaystyle=\mathbf{A}^{-1}-\mathbf{A}^{-1}\bm{\Delta}(\mathbf{A}+\bm{\Delta})^{-1}$
$\displaystyle=\mathbf{A}^{-1}-\mathbf{A}^{-1}\bm{\Delta}\left(\mathbf{A}^{-1}-\mathbf{A}^{-1}\bm{\Delta}(\mathbf{A}+\bm{\Delta})^{-1}\right)$
$\displaystyle=\mathbf{A}^{-1}-\mathbf{A}^{-1}\bm{\Delta}\mathbf{A}^{-1}+(\mathbf{A}^{-1}\bm{\Delta})^{2}(\mathbf{A}+\bm{\Delta})^{-1}$
$\displaystyle=\sum_{s=0}^{n}(-\mathbf{A}^{-1}\bm{\Delta})^{s}\mathbf{A}^{-1}+(-\mathbf{A}^{-1}\bm{\Delta})^{n+1}(\mathbf{A}+\bm{\Delta})^{-1}$

Note that
$\mathbf{A}^{-1}=(\mathbf{I}_{M}+\lambda\bm{\Lambda}^{-1})^{-1}=\bm{\Lambda}(\bm{\Lambda}+\lambda\mathbf{I}_{M})^{-1}=\bar{\mathbf{P}}$
with operator norm $\frac{\lambda_{1}}{\lambda_{1}+\lambda}<1$, hence we have
$(\mathbf{A}^{-1}\bm{\Delta})^{n+1}=(-\bar{\mathbf{P}}\bm{\Delta})^{n+1}\to 0$
in operator norm as $n\to\infty$. Hence

$\displaystyle\mathbf{B}$
$\displaystyle=\sum_{s=0}^{\infty}(-\mathbf{A}^{-1}\bm{\Delta})^{s}\mathbf{A}^{-1}$
$\displaystyle=\sum_{s=0}^{\infty}(-\bar{\mathbf{P}}\bm{\Delta})^{s}\bar{\mathbf{P}}$

in operator norm. ∎

Due to the convergence result in lemma C.9, it is natural to define:

######  Definition C.10.

For any $n\in\mathbb{N}\cup\\{\infty\\}$, write
$\mathbf{B}^{(n)}=\sum_{s=0}^{n}(-\bar{\mathbf{P}}\bm{\Delta})^{s}\bar{\mathbf{P}}$.
For example, We have

$\displaystyle\mathbf{B}^{(0)}$ $\displaystyle=\bar{\mathbf{P}}$
$\displaystyle\mathbf{B}^{(1)}$
$\displaystyle=\bar{\mathbf{P}}-\bar{\mathbf{P}}\bm{\Delta}\bar{\mathbf{P}}$
$\displaystyle\mathbf{B}^{(\infty)}$ $\displaystyle=\mathbf{B}$

Although lemma C.9 is valid when $\delta<1$, we need a slightly stronger
condition that $\delta$ is upper bounded by an arbitrary constant strictly
small than 1. For simplicity, we assume this constant to be $\frac{1}{2}$ in
the following lemma:

######  Lemma C.11 ($\mathbf{B}$-Approximation).

Assume that
$\delta\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\left\|\bm{\Delta}\right\|_{\text{op}}<\frac{1}{2}$.
Let
$\mathbf{B}^{(n)}=\sum_{s=0}^{n}(-\bar{\mathbf{P}}\bm{\Delta})^{s}\bar{\mathbf{P}}$
be the $n$th-order approximation of the matrix $\mathbf{B}$ as in definition
C.10. Then we have

$\left\|\mathbf{B}-\mathbf{B}^{(n)}\right\|_{\text{op}}<2\delta^{n+1}.$

###### Proof.

We first bound the operator norm of the matrix $\mathbf{B}$: since the minimum
singular value of the matrix $\bar{\mathbf{P}}^{-1}+\bm{\Delta}$ is at least

$\frac{\lambda_{k}+\lambda}{\lambda_{k}}-\left\|\Delta\right\|_{\text{op}}\geq
1+\frac{\lambda}{\lambda_{k}}-\frac{1}{2}>\frac{1}{2},$

and hence

$\left\|\mathbf{B}\right\|_{\text{op}}=\left\|\left(\bar{\mathbf{P}}^{-1}+\bm{\Delta}\right)^{-1}\right\|_{\text{op}}<2.$

Also, we have

$\displaystyle\mathbf{B}-\mathbf{B}^{(n)}$
$\displaystyle=\sum_{s=n+1}^{\infty}(-\bar{\mathbf{P}}\bm{\Delta})^{s}\bar{\mathbf{P}}$
$\displaystyle=(-\bar{\mathbf{P}}\bm{\Delta})^{n+1}\sum_{s=0}^{\infty}(-\bar{\mathbf{P}}\bm{\Delta})^{s}\bar{\mathbf{P}}$
$\displaystyle=(-\bar{\mathbf{P}}\bm{\Delta})^{n+1}\mathbf{B}.$

Hence
$\left\|\mathbf{B}-\mathbf{B}^{(n)}\right\|_{\text{op}}\leq\left\|\bar{\mathbf{P}}\bm{\Delta}\right\|_{\text{op}}^{n+1}\left\|\mathbf{B}\right\|_{\text{op}}<\left\|\bm{\Delta}\right\|_{\text{op}}^{n+1}\cdot
2=2\delta^{n+1}$, since we have
$\left\|\bar{\mathbf{P}}\right\|_{\text{op}}=\frac{\lambda_{1}}{\lambda_{1}+\lambda}<1$.
∎

Note that the upper bound $\frac{1}{2}$ of $\delta$ can be replaced by any
constant strictly small than 1 to get a similar conclusion.

######  Remark C.12.

Using the concentration result from random matrix theory, for $M<N$, one can
show with high probability that the operator norm $\delta$ of the fluctuation
matrix $\bm{\Delta}$ is less than 1. 555From there, we differentiate the
approach from Bach Bach [2023a]: From Propositions C.5 and C.7, it is
inevitable to approximate the matrix $\mathbf{B}$, and we have
$\mathbf{I}_{M}$ as support of the inverse. Bach instead uses RHKS basis to
express the fluctuation matrix and is hence forced to use
$\lambda\mathbf{I}_{M}$ as the support. As a result, he would need to require
that the fluctuation is less than $\lambda$ and hence his requirement on $N$
is antiproportional to $\lambda$ in Theorem 5.1.

See subsection C.3 for details. Then we can use the the above lemmata C.9 and
C.11 to approximate the test error of KRR:

######  Proposition C.13 (Bias Approximation).

Fix a sample $\mathbf{Z}$ of $\rho$ such that
$\delta\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\|\bm{\Delta}\|_{\text{op}}<\frac{1}{2}$.
Then the $\text{bias}_{\text{test}}$ term is bounded above and below by

$\left|\text{bias}-\left(\left\|\bar{\mathbf{P}}w\right\|_{2}^{2}+\tilde{\gamma}_{>M}^{2}\right)\right|\leq
2\delta\left\|\bar{\mathbf{P}}w\right\|_{2}^{2}+\left\|w\right\|_{2}^{2}\delta^{2}p(\delta),$

where
$\bar{\mathbf{P}}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\bm{\Lambda}(\bm{\Lambda}+\lambda\mathbf{I}_{M})^{-1}$,
$w=\lambda\bm{\Lambda}^{-1}\tilde{\bm{\gamma}}-\tilde{{\gamma}}_{>M}\bm{E}$,
and $p(\delta)\stackrel{{\scriptstyle\mbox{def.}}}{{=}}5+4\delta+4\delta^{2}$.
By writing $\bm{E}=(\eta_{k})_{k=1}^{M}$, the bounds simplify to

bias
$\displaystyle\leq\tilde{\gamma}_{>M}^{2}+(1+2\delta)\,\sum_{k=1}^{M}\frac{(\lambda\tilde{\gamma}_{k}-\tilde{\gamma}_{>M}\eta_{k}\lambda_{k})^{2}}{(\lambda_{k}+\lambda)^{2}}+\left\|w\right\|_{2}^{2}\delta^{2}p(\delta);$
bias
$\displaystyle\geq\tilde{\gamma}_{>M}^{2}+(1-2\delta)\,\sum_{k=1}^{M}\frac{(\lambda\tilde{\gamma}_{k}-\tilde{\gamma}_{>M}\eta_{k}\lambda_{k})^{2}}{(\lambda_{k}+\lambda)^{2}}-\left\|w\right\|_{2}^{2}\delta^{2}p(\delta).$

###### Proof.

Let
$w=\lambda\bm{\Lambda}^{-1}\tilde{\bm{\gamma}}-\tilde{{\gamma}}_{>M}\bm{E}$.
We apply lemma C.4 followed by the 1st-order approximation $\mathbf{B}^{(1)}$
of the matrix $\mathbf{B}$ in lemma C.11:

fitting error
$\displaystyle=\|\mathbf{B}w\|_{2}^{2}=\left\|\mathbf{B}^{(1)}w+\left(\mathbf{B}-\mathbf{B}^{(1)}\right)w\right\|_{2}^{2}$
$\displaystyle=\left\|\mathbf{B}^{(1)}w\right\|_{2}^{2}+w^{\top}\left(\mathbf{B}-\mathbf{B}^{(1)}\right)\mathbf{B}^{(1)}w+w^{\top}\mathbf{B}^{(1)}\left(\mathbf{B}-\mathbf{B}^{(1)}\right)w+\left\|\left(\mathbf{B}-\mathbf{B}^{(1)}\right)w\right\|_{2}^{2}$
$\displaystyle\leq\left\|\mathbf{B}^{(1)}w\right\|_{2}^{2}+2\left\|\mathbf{B}^{(1)}\right\|_{\text{op}}\left\|\mathbf{B}-\mathbf{B}^{(1)}\right\|_{\text{op}}\left\|w\right\|_{2}^{2}+\left\|\mathbf{B}-\mathbf{B}^{(1)}\right\|_{\text{op}}^{2}\left\|w\right\|_{2}^{2}$
$\displaystyle\leq\left\|\mathbf{B}^{(1)}w\right\|_{2}^{2}+2\cdot(1+\delta)\cdot
2\delta^{2}\left\|w\right\|_{2}^{2}+4\delta^{4}\left\|w\right\|_{2}^{2}$
$\displaystyle\leq\left\|\mathbf{B}^{(1)}w\right\|_{2}^{2}+4\left\|w\right\|_{2}^{2}\delta^{2}(1+\delta+\delta^{2})$
$\displaystyle\leq\left\|\left(\bar{\mathbf{P}}-\bar{\mathbf{P}}\bm{\Delta}\bar{\mathbf{P}}\right)w\right\|_{2}^{2}+4\left\|w\right\|_{2}^{2}\delta^{2}(1+\delta+\delta^{2})$
$\displaystyle\leq\left\|\mathbf{I}_{M}-\bar{\mathbf{P}}\bm{\Delta}\right\|_{\text{op}}^{2}\left\|\bar{\mathbf{P}}w\right\|_{2}^{2}+4\left\|w\right\|_{2}^{2}\delta^{2}(1+\delta+\delta^{2})$
$\displaystyle\leq\left(1+2\left\|\bar{\mathbf{P}}\bm{\Delta}\right\|_{\text{op}}+\left\|\bar{\mathbf{P}}\bm{\Delta}\right\|_{\text{op}}^{2}\right)\left\|\bar{\mathbf{P}}w\right\|_{2}^{2}+4\left\|w\right\|_{2}^{2}\delta^{2}(1+\delta+\delta^{2})$
$\displaystyle\leq\left(1+2\left\|\bar{\mathbf{P}}\bm{\Delta}\right\|_{\text{op}}\right)\left\|\bar{\mathbf{P}}w\right\|_{2}^{2}+\left\|w\right\|_{2}^{2}\delta^{2}(5+4\delta+4\delta^{2})$
$\displaystyle\leq(1+2\delta)\left\|\bar{\mathbf{P}}w\right\|_{2}^{2}+\left\|w\right\|_{2}^{2}\delta^{2}(5+4\delta+4\delta^{2}).$

Hence we have the upper bound:

$\text{bias}\leq\tilde{\gamma}_{>M}^{2}+(1+2\delta)\left\|\bar{\mathbf{P}}w\right\|_{2}^{2}+\left\|w\right\|_{2}^{2}\delta^{2}p(\delta).$

We argue similarly for the lower bound using:
$\left\|\mathbf{A}\right\|_{\text{op}}\left\|v\right\|_{2}^{2}\geq
v^{\top}\mathbf{A}v\geq-\left\|\mathbf{A}\right\|_{\text{op}}\left\|v\right\|_{2}^{2}$
for any $\mathbf{A}\in\mathbb{R}^{M\times M},\ v\in\mathbb{R}^{M\times 1}$. ∎

We argue similarly for variance.

######  Proposition C.14 (Variance Approximation).

Fix a sampling $\mathbf{Z}$ such that
$\delta\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\|\bm{\Delta}\|_{\text{op}}<\frac{1}{2}$.
Then we have

$\displaystyle\left|\text{variance}-\frac{\sigma^{2}}{N}\sum_{k=1}^{M}\frac{\lambda_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}\right|\leq\delta\frac{\sigma^{2}}{N}\sum_{k=1}^{M}\frac{\lambda_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}+M\frac{\sigma^{2}}{N}(1+\delta)\delta^{2}p(\delta),$

where
$p(\delta)\stackrel{{\scriptstyle\mbox{def.}}}{{=}}5+4\delta+4\delta^{2}$, and
$\sigma^{2}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\mathbb{E}[\epsilon^{2}]$
is the noise variance.

###### Proof.

Note that $\operatorname{Tr}{\mathbf{A}}\leq
M\left\|\mathbf{A}\right\|_{\text{op}}$ for any matrix
$\mathbf{A}\in\mathbb{R}^{M\times M}$. Since
$\mathbf{B}^{2}(\mathbf{I}_{M}+\bm{\Delta})=(\mathbf{B}^{(1)})^{2}(\mathbf{I}_{M}+\bm{\Delta})+2\mathbf{B}^{(1)}\left(\mathbf{B}-\mathbf{B}^{(1)}\right)(\mathbf{I}_{M}+\bm{\Delta})+\left(\mathbf{B}-\mathbf{B}^{(1)}\right)^{2}(\mathbf{I}_{M}+\bm{\Delta})$,
we can bound the residue term by $\delta$:

$\displaystyle\operatorname{Tr}\left[2\mathbf{B}^{(1)}\left(\mathbf{B}-\mathbf{B}^{(1)}\right)(\mathbf{I}_{M}+\bm{\Delta})+\left(\mathbf{B}-\mathbf{B}^{(1)}\right)^{2}(\mathbf{I}_{M}+\bm{\Delta})\right]$
$\displaystyle\leq
M(1+\delta)\left\|\mathbf{B}-\mathbf{B}^{(1)}\right\|_{\text{op}}(2\left\|\mathbf{B}^{(1)}\right\|_{\text{op}}+\left\|\mathbf{B}-\mathbf{B}^{(1)}\right\|_{\text{op}})$
$\displaystyle\leq M(1+\delta)\cdot 2\delta^{2}(2(1+\delta)+2\delta^{2})$
$\displaystyle\leq 4M\delta^{2}(1+\delta)(1+\delta+\delta^{2}),$

For the main terms, we have

$\displaystyle\operatorname{Tr}[(\mathbf{B}^{(1)})^{2}(\mathbf{I}_{M}+\bm{\Delta})]$
$\displaystyle\leq\operatorname{Tr}[\bar{\mathbf{P}}^{2}]\cdot\left\|(\mathbf{I}_{M}-\bm{\Delta}\bar{\mathbf{P}})^{2}(\mathbf{I}_{M}+\bm{\Delta})\right\|_{\text{op}}$
$\displaystyle=\operatorname{Tr}[\bar{\mathbf{P}}^{2}]\left\|\mathbf{I}_{M}+\bm{\Delta}(\mathbf{I}_{M}-2\bar{\mathbf{P}})+(\bm{\Delta}\bar{\mathbf{P}})^{2}-2\bm{\Delta}\bar{\mathbf{P}}\bm{\Delta}+(\bm{\Delta}\bar{\mathbf{P}})^{2}\bm{\Delta}\right\|_{\text{op}}$
$\displaystyle\leq\operatorname{Tr}[\bar{\mathbf{P}}^{2}]\left\|\mathbf{I}_{M}+\bm{\Delta}(\mathbf{I}_{M}-2\bar{\mathbf{P}})\right\|_{\text{op}}+M\left\|(\bm{\Delta}\bar{\mathbf{P}})^{2}-2\bm{\Delta}\bar{\mathbf{P}}\bm{\Delta}+(\bm{\Delta}\bar{\mathbf{P}})^{2}\bm{\Delta}\right\|_{\text{op}}$
$\displaystyle\leq\operatorname{Tr}[\bar{\mathbf{P}}^{2}](1+\delta)+M\delta^{2}(1+\delta).$

We apply Theorem C.7 to yield a bound on variance:

$\displaystyle\left|\text{variance}-\frac{\sigma^{2}}{N}\sum_{k=1}^{M}\frac{\lambda_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}\right|$
$\displaystyle\leq\delta\frac{\sigma^{2}}{N}\sum_{k=1}^{M}\frac{\lambda_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}+M\frac{\sigma^{2}}{N}(1+\delta)\delta^{2}p(\delta).$

∎

Note that the above propositions C.13 and C.14 give absolute (non-
probabilistic) bounds on the test error, once $\delta$ is controlled.

###  C.3 Concentration Results

In this subsection, we focus on bounding the operator norm $\delta$ of the
fluctuation matrix $\bm{\Delta}$.

First, we establish some concentration results.

######  Lemma C.15 (Theorem 3.59 in Vershynin [2010]).

Let $\mathbf{A}$ be an $n\times N$ matrix with independent isotropic sub-
Gaussian columns in $\mathbb{R}^{n}$ which sub-gaussian norm is bounded by a
positive constant $G>0$. Then for all $t\geq 0$, with probability at least
$1-2\exp(-\frac{1}{3}t^{2})$, we have

$\displaystyle\left\|\frac{1}{N}\mathbf{A}\mathbf{A}^{\top}-\mathbf{I}_{n}\right\|_{\text{op}}\leq\max(a,a^{2}),$
(34)

where
$a\stackrel{{\scriptstyle\mbox{def.}}}{{=}}C\sqrt{\frac{n}{N}}+\frac{t}{\sqrt{N}}$,
for all constant $C\geq 12G^{2}$ .

###### Proof.

Let
$a\stackrel{{\scriptstyle\mbox{def.}}}{{=}}C\sqrt{\frac{n}{N}}+\frac{t}{\sqrt{N}}$
with $C>0$ to be determined, and
$\epsilon\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\max\\{a,a^{2}\\}$. The
first step to show that :

$\displaystyle\max_{x\in\mathcal{N}}\left|\frac{1}{N}\left\|\mathbf{A}^{\top}x\right\|_{2}^{2}-1\right|\leq\epsilon$

for some $\frac{1}{4}$-net $\mathcal{N}$ on the sphere
$\mathbb{S}^{n-1}\subset\mathbb{R}^{n}$. Choose such a net $\mathcal{N}$ with
$|\mathcal{N}|<\left(1+\frac{2}{1/4}\right)^{n}=9^{n}$. Let $\mathbf{A}_{i}$
be the $i$th column of the matrix $\mathbf{A}$ and let
$Z_{i}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\mathbf{A}_{i}^{\top}x$ be a
random variable. By definition of $\mathbf{A}$, $Z_{i}$ is centered with unit
variance with sub-Gaussian norm upper bounded by $G$. Note that
$G\geq\frac{1}{\sqrt{2}}\mathbb{E}[Z_{i}^{2}]^{1/2}=\frac{1}{\sqrt{2}}$, and
the random variable $Z_{i}^{2}-1$ is centered and has sub-exponential norm
upper bounded by $4G^{2}$. Hence by an exponential deviation inquality 666This
inequality is Corollary 5.17 from Vershynin [2010]., we have, for any
$x\in\mathbb{S}^{n-1}$:

$\displaystyle\mathbb{P}\left\\{\left|\frac{1}{N}\left\|\mathbf{A}^{\top}x\right\|_{2}^{2}-1\right|\geq\frac{\epsilon}{2}\right\\}$
$\displaystyle=\mathbb{P}\left\\{\left|\frac{1}{N}\sum_{i=1}^{N}Z_{i}^{2}-1\right|\geq\frac{\epsilon}{2}\right\\}$
$\displaystyle\leq
2\exp\left(-\frac{1}{2}e^{-1}G^{-4}\min\\{\epsilon,\epsilon^{2}\\}\right)$
$\displaystyle=\leq 2\exp\left(-\frac{1}{2}e^{-1}G^{-4}a^{2}\right)$
$\displaystyle\leq 2\exp\left(-\frac{1}{2}e^{-1}G^{-4}(C^{2}n+t^{2})\right).$

Then by union bound, we have

$\displaystyle\mathbb{P}\left\\{\max_{x\in\mathcal{N}}\left|\frac{1}{N}\left\|\mathbf{A}^{\top}x\right\|_{2}^{2}-1\right|\geq\frac{\epsilon}{2}\right\\}$
$\displaystyle\leq 9^{n}\cdot
2\exp\left(-\frac{1}{2}e^{-1}G^{-4}(C^{2}n+t^{2})\right)$ $\displaystyle\leq
2\exp\left(-\frac{1}{2}e^{-1}G^{-4}t^{2}\right),$

for $C\geq\sqrt{2e\log 9}G^{2}$. Since $12>\sqrt{2e\log 9}$, for simplicity,
we assume $C>12G^{2}$. Moreover, since $G\geq\frac{1}{\sqrt{2}}$, we have
$\frac{1}{2}e^{-1}G^{-4}\leq\frac{1}{3}$, we have

$\displaystyle\mathbb{P}\left\\{\max_{x\in\mathcal{N}}\left|\frac{1}{N}\left\|\mathbf{A}^{\top}x\right\|_{2}^{2}-1\right|\geq\frac{\epsilon}{2}\right\\}$
$\displaystyle\leq 2\exp\left(-\frac{1}{3}t^{2}\right).$

Then by the $\frac{1}{4}$-net argument, with probability at least
$1-2\exp\left(-\frac{1}{3}t^{2}\right)$, we have

$\displaystyle\left\|\frac{1}{N}\mathbf{A}\mathbf{A}^{\top}-\mathbf{I}_{n}\right\|_{\text{op}}$
$\displaystyle\leq\frac{4}{2}\max_{x\in\mathcal{N}}\left|\frac{1}{N}\left\|\mathbf{A}^{\top}x\right\|_{2}^{2}-1\right|$
$\displaystyle\leq\epsilon=\max\\{a,a^{2}\\}.$

∎

######  Lemma C.16.

Assume Assumption 4.1 holds and that $N>\exp(4(12G^{2})^{2}(M+1))$. Then with
a probability of at least $1-2/N$, we have

$\displaystyle\max\left\\{\delta,\left\|\bm{E}_{M}\right\|_{2}\right\\}\leq\sqrt{\frac{\log
N}{N}}.$

###### Proof.

Set $n=M+1$, $\mathbf{A}=\binom{\bm{\Psi}_{\leq
M}}{\psi_{>M}(\mathbf{X})^{\top}}\in\mathbb{R}^{(M+1)\times N}$. Then

$\displaystyle\frac{1}{N}\mathbf{A}\mathbf{A}^{\top}-\mathbf{I}_{n}=\begin{pmatrix}\frac{1}{N}\bm{\Psi}_{\leq
M}\bm{\Psi}_{\leq M}^{\top}&\bm{E}_{M}\\\
\bm{E}_{M}^{\top}&\eta_{>M}+1\end{pmatrix}-\mathbf{I}_{n}=\begin{pmatrix}\bm{\Delta}_{M}&\bm{E}_{M}\\\
\bm{E}_{M}^{\top}&\eta_{>M}\end{pmatrix}.$

where
$\eta_{>M}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\frac{1}{N}\sum_{i=1}^{N}\psi_{>M}(x_{i})^{2}-1$.
On one hand, the operator norm of the above matrix bounds $\delta$ and
$\left\|\bm{E}_{M}\right\|_{2}$ from above:

$\displaystyle\left\|\begin{pmatrix}\bm{\Delta}_{M}&\bm{E}_{M}\\\
\bm{E}_{M}^{\top}&\eta_{>M}\end{pmatrix}\right\|_{\text{op}}$
$\displaystyle=\max_{\left\|\mathbf{u}\right\|_{2}^{2}+v^{2}=1}\left\|\begin{pmatrix}\bm{\Delta}_{M}&\bm{E}_{M}\\\
\bm{E}_{M}^{\top}&\eta_{>M}\end{pmatrix}\begin{pmatrix}\mathbf{u}\\\
v\end{pmatrix}\right\|_{2}$
$\displaystyle=\max_{\left\|\mathbf{u}\right\|_{2}^{2}+v^{2}=1}\left\|\begin{pmatrix}\bm{\Delta}_{M}\mathbf{u}+v\bm{E}_{M}\\\
\bm{E}_{M}^{\top}\mathbf{u}+\eta_{>M}v\end{pmatrix}\right\|_{2}$
$\displaystyle\geq\max_{\left\|\mathbf{u}\right\|_{2}^{2}+v^{2}=1}\left\|\bm{\Delta}_{M}\mathbf{u}+v\bm{E}_{M}\right\|_{2}$
$\displaystyle\geq\max_{\left\|\mathbf{u}\right\|_{2}^{2}=1,v=0}\left\|\bm{\Delta}_{M}\mathbf{u}+v\bm{E}_{M}\right\|_{2}$
$\displaystyle\geq\max_{\left\|\mathbf{u}\right\|_{2}^{2}=1}\left\|\bm{\Delta}_{M}\mathbf{u}\right\|_{2}=\delta,$

and

$\displaystyle\left\|\begin{pmatrix}\bm{\Delta}_{M}&\bm{E}_{M}\\\
\bm{E}_{M}^{\top}&\eta_{>M}\end{pmatrix}\right\|_{\text{op}}\geq\max_{\left\|\mathbf{u}\right\|_{2}^{2}+v^{2}=1}\left\|\bm{\Delta}_{M}\mathbf{u}+v\bm{E}_{M}\right\|_{2}\geq\max_{\left\|\mathbf{u}\right\|_{2}^{2}=0,|v|=1}\left\|\bm{\Delta}_{M}\mathbf{u}+v\bm{E}_{M}\right\|_{2}=\left\|\bm{E}_{M}\right\|_{2}.$

On the other hand, set $t=\frac{1}{2}\sqrt{\log N},\ C=12G^{2}$, since
$N>\exp(4C^{2}(M+1))$, we have

$\displaystyle
a=C\sqrt{\frac{n}{N}}+\frac{t}{\sqrt{N}}=12G^{2}\sqrt{\frac{M+1}{N}}+\frac{1}{2}\sqrt{\frac{\log
N}{N}}\leq\sqrt{\frac{\log N}{N}}<1.$

By Lemma C.16, then with probability of at least
$1-2\exp(-\frac{1}{3}t^{2})=1-2\exp(-\frac{1}{12})/N>1-2/N$, we have

$\displaystyle\left\|\begin{pmatrix}\bm{\Delta}_{M}&\bm{E}_{M}\\\
\bm{E}_{M}^{\top}&\eta_{>M}\end{pmatrix}\right\|_{\text{op}}\leq\max\\{a,a^{2}\\}=a\leq\sqrt{\frac{\log
N}{N}}.$

Combine the both results and we conclude the upper bounds. ∎

In particular, as $N\rightarrow\infty$, $\delta$ vanishes almost surely. In
empirical calculation, if the requirement $N>\exp(4(12G^{2})^{2}(M+1))$
exponential in $M$ is too demanding for a large integer $M$, we can take
$t=N^{s}$ for any positive number $s\in\left(0,\frac{1}{2}\right)$ instead of
$t=\frac{1}{2}{\log N}$. In this way, we decrease the requirement to $N$
polynomial in $M$ in sacrificing the decay from
$\mathcal{O}\left(\sqrt{\frac{logN}{N}}\right)$ to
$\mathcal{O}\left(N^{s-1/2}\right)$. For simplicity purpose, we do not list
out the result with this decay in this paper.

###  C.4 Refined Test Error Analysis

We can apply the above concentration results to refine the following bounds on
the finite-rank KRR test error. First of all, we realize the decay of target
function coefficient comparable to the spectral decay:

######  Definition C.17 (Comparable Decay).

Denote
$\underline{r}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\min_{k}\\{|\tilde{\gamma}_{k}/\lambda_{k}|\\}$
and
$\overline{r}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\max_{k}\\{|\tilde{\gamma}_{k}/\lambda_{k}|\\}$.

####  C.4.1 Refined Bounds on Bias

Recall that Proposition C.13 bounding the bias in terms of $\delta$ and
$\eta_{k}$. For the former one, we can choose: for
$N>\max\\{\exp(4(12G^{2})^{2}(M+1)),9\\}$, by Lemma C.16, with probability of
at least $1-2/N$, we have $\delta\leq\sqrt{\frac{\log N}{N}}<\sqrt{\frac{\log
9}{9}}<\frac{1}{2}.$ For the latter one, we have to control the vector $w$:

######  Lemma C.18.

Let
$w=\lambda\bm{\Lambda}^{-1}\tilde{\bm{\gamma}}-\tilde{{\gamma}}_{>M}\bm{E}$.
We have

$\displaystyle\left\|w\right\|_{2}^{2}$
$\displaystyle\leq\left(\lambda\overline{r}\sqrt{M}+|\tilde{\gamma}_{>M}|\left\|\bm{E}\right\|_{2}\right)^{2};$
$\displaystyle\frac{\lambda^{2}\lambda_{M}}{(\lambda_{M}+\lambda)^{2}}\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}-\frac{1}{2}|\tilde{\gamma}_{>M}|\|\tilde{f}_{\leq
M}\|_{L_{\rho}^{2}}\left\|\bm{E}\right\|_{2}$
$\displaystyle\leq\left\|\bar{\mathbf{P}}w\right\|_{2}^{2}\leq\lambda\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}+\frac{1}{2}|\tilde{\gamma}_{>M}|\|\tilde{f}_{\leq
M}\|_{L_{\rho}^{2}}\left\|\bm{E}\right\|_{2}+\tilde{\gamma}_{>M}^{2}\left\|\bm{E}\right\|_{2}^{2}.$

###### Proof.

Since
$\lambda^{2}\underline{r}^{2}M\leq\left\|\lambda\bm{\Lambda}^{-1}\tilde{\bm{\gamma}}\right\|_{2}^{2}\leq\lambda^{2}\overline{r}^{2}M$
and
$\left\|\tilde{{\gamma}}_{>M}\bm{E}\right\|_{2}^{2}=\tilde{{\gamma}}_{>M}^{2}\left\|\bm{E}\right\|_{2}^{2}$,
we have

$\left\|w\right\|_{2}^{2}\leq\left(\lambda\overline{r}\sqrt{M}+|\tilde{\gamma}_{>M}|\left\|E\right\|_{2}\right)^{2}.$

Similarly, we can bound $\left\|\bar{\mathbf{P}}w\right\|_{2}$. Observe that:

$\displaystyle\left\|\bar{\mathbf{P}}w\right\|_{2}^{2}$
$\displaystyle=\underbrace{\lambda^{2}\sum_{k=1}^{M}\frac{\tilde{\gamma}_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}}_{I}\underbrace{-2\lambda\tilde{\gamma}_{>M}\sum_{k=1}^{M}\frac{\tilde{\gamma}_{k}\lambda_{k}\eta_{k}}{(\lambda_{k}+\lambda)^{2}}}_{II}+\underbrace{\tilde{\gamma}_{>M}^{2}\sum_{k=1}^{M}\frac{\lambda_{k}^{2}\eta_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}}_{III}.$

Since
$1\geq\frac{\lambda}{\lambda_{k}+\lambda}\geq\frac{\lambda}{\lambda_{M}+\lambda}$,
we have the upper bound:

$I=\lambda^{2}\sum_{k=1}^{M}\frac{\tilde{\gamma}_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}\leq\lambda\sum_{k=1}^{M}\frac{\lambda}{\lambda_{k}+\lambda}\frac{\tilde{\gamma}_{k}^{2}}{\lambda_{k}+\lambda}\leq\lambda\sum_{k=1}^{M}\frac{\tilde{\gamma}^{2}}{\lambda_{k}}=\lambda\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}.$ (35)

where $\tilde{f}_{\leq
M}\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\sum_{k=1}^{M}\tilde{\gamma}_{k}\psi_{k}=\tilde{f}-\tilde{\gamma}_{>M}\psi_{>M}$.
For the lower bound, we have:

$I=\lambda^{2}\sum_{k=1}^{M}\frac{\tilde{\gamma}_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}\geq\lambda^{2}\sum_{k=1}^{M}\frac{\lambda_{k}}{(\lambda_{k}+\lambda)^{2}}\frac{\tilde{\gamma}_{k}^{2}}{\lambda_{k}}\geq\lambda^{2}\frac{\lambda_{M}}{(\lambda_{M}+\lambda)^{2}}\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}$ (36)

Similarly, since $4\lambda\lambda_{k}\leq(\lambda_{k}+\lambda)^{2}$,

$\displaystyle|II|$
$\displaystyle=2\lambda|\tilde{\gamma}_{>M}|\sum_{k=1}^{M}\frac{|\tilde{\gamma}_{k}|\lambda_{k}|\eta_{k}|}{(\lambda_{k}+\lambda)^{2}}\leq\frac{1}{2}|\tilde{\gamma}_{>M}|\sum_{k+1}^{M}|\tilde{\gamma}_{k}\eta_{k}|\leq\frac{1}{2}|\tilde{\gamma}_{>M}|\sqrt{\sum_{k+1}^{M}\tilde{\gamma}_{k}^{2}\sum_{k=1}^{M}\eta_{k}^{2}}\leq\frac{1}{2}|\tilde{\gamma}_{>M}|\|\tilde{f}_{\leq
M}\|_{L_{\rho}^{2}}\left\|\bm{E}\right\|_{2}.$

And

$\displaystyle III$
$\displaystyle=\tilde{\gamma}_{>M}^{2}\sum_{k=1}^{M}\frac{\lambda_{k}^{2}\eta_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}\leq\tilde{\gamma}_{>M}^{2}\sum_{k=1}^{M}\eta_{k}^{2}=\tilde{\gamma}_{>M}^{2}\left\|\bm{E}\right\|_{2}^{2}.$

∎

Combining the above result, we state the following theorem:

######  Theorem C.19.

For $N>\max\left\\{\exp(4(12G^{2})^{2}(M+1)),9\right\\}$ and for any constant
$C_{1}>8\left(\lambda\overline{r}\sqrt{M}+\frac{1}{2}|\tilde{\gamma}_{>M}|\right)^{2}+\frac{5}{2}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}$
(independent to $N$), with a probability of at least $1-2/N$, we have the
upper and lower bounds of bias:

bias $\displaystyle\leq\tilde{\gamma}_{>M}^{2}+\lambda\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}+\left(\frac{1}{4}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}+2\lambda\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}\right)\sqrt{\frac{\log N}{N}}+C_{1}\frac{\log N}{N};$
bias
$\displaystyle\geq\tilde{\gamma}_{>M}^{2}+\frac{\lambda^{2}\lambda_{M}}{(\lambda_{M}+\lambda)^{2}}\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}-\left(\frac{1}{4}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}+\frac{2\lambda^{2}}{\lambda_{1}+\lambda}\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}\right)\sqrt{\frac{\log N}{N}}-C_{1}\frac{\log N}{N}.$

For $\lambda\to 0$, we have a simpler bound: with a probability of at least
$1-2/N$, we have

$\displaystyle\begin{split}\lim_{\lambda\to
0}\text{bias}&\leq\left(1+\frac{\log
N}{N}\right)\tilde{\gamma}_{>M}^{2}+6\tilde{\gamma}_{>M}^{2}\left(\frac{\log
N}{N}\right)^{\frac{3}{2}};\\\ \lim_{\lambda\to
0}\text{bias}&\geq\left(1-\frac{\log
N}{N}\right)\tilde{\gamma}_{>M}^{2}-6\tilde{\gamma}_{>M}^{2}\left(\frac{\log
N}{N}\right)^{\frac{3}{2}}.\end{split}$ (37)

For $\tilde{\gamma}_{>M}^{2}=0$, that is $\tilde{f}\in\mathcal{H}$, we have a
simpler upper bound on bias: with a probability of at least $1-2/N$, we have

$\displaystyle\begin{split}\text{bias}&\leq\lambda\|\tilde{f}\|_{\mathcal{H}}^{2}\left(1+2\sqrt{\frac{\log
N}{N}}\right)+C_{1}\frac{\log N}{N};\\\
\text{bias}&\geq\frac{\lambda^{2}\lambda_{M}}{(\lambda_{M}+\lambda)^{2}}\|\tilde{f}\|_{\mathcal{H}}^{2}\left(1-2\sqrt{\frac{\log
N}{N}}\right)-C_{1}\frac{\log N}{N}.\end{split}$ (38)

###### Proof.

By Proposition C.13 and Lemma C.18,

fitting error
$\displaystyle\leq(1+2\delta)\left\|\bar{\mathbf{P}}w\right\|_{2}^{2}+\left\|w\right\|_{2}^{2}\delta^{2}p(\delta)$
(39) $\displaystyle\leq(1+2\delta)(\lambda\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}+\frac{1}{2}|\tilde{\gamma}_{>M}|\|\tilde{f}_{\leq
M}\|_{L_{\rho}^{2}}\left\|\bm{E}\right\|_{2}+\tilde{\gamma}_{>M}^{2}\left\|\bm{E}\right\|_{2}^{2})+\left\|w\right\|_{2}^{2}\delta^{2}p(\delta)$
(40) $\displaystyle\leq(1+2\delta)(\lambda\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}+\frac{1}{4}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}\left\|\bm{E}\right\|_{2}+\tilde{\gamma}_{>M}^{2}\left\|\bm{E}\right\|_{2}^{2})+\left\|w\right\|_{2}^{2}\delta^{2}p(\delta).$
(41)

where in line (39), we use Proposition C.13; in line (40), we use Lemma C.18;
in line (40), we use the fact that $2ab\leq a^{2}+b^{2}$ where
$a=|\tilde{\gamma}_{>M}|,b=\|\tilde{f}_{\leq M}\|_{L_{\rho}^{2}}$.

Now we apply the concentration result in Lemma C.16: with a probability of at
least $1-2/N$:

fitting error $\displaystyle\leq\left(1+2\sqrt{\frac{\log
N}{N}}\right)\lambda\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}+\frac{1}{4}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}\sqrt{\frac{\log
N}{N}}$ $\displaystyle\quad+\frac{\log
N}{N}\left(\left\|w\right\|_{2}^{2}p(\delta)+(1+2\delta)\tilde{\gamma}_{>M}^{2}+\frac{1}{2}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}\right)$
$\displaystyle\leq\lambda\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}+\left(\frac{1}{4}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}+2\lambda\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}\right)\sqrt{\frac{\log N}{N}}+C_{1}\frac{\log N}{N},$

where we choose $C_{1}>0$ to be such that:

$\displaystyle\left\|w\right\|_{2}^{2}p(\delta)+(1+2\delta)\tilde{\gamma}_{>M}^{2}+\frac{1}{2}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}$
$\displaystyle\leq\left\|w\right\|_{2}^{2}p(\delta)+\left(1+2\delta+\frac{1}{2}\right)\|\tilde{f}\|_{L_{\rho}^{2}}^{2}$
$\displaystyle\leq\left\|w\right\|_{2}^{2}p\left(\frac{1}{2}\right)+\left(1+2\cdot\frac{1}{2}+\frac{1}{2}\right)\|\tilde{f}\|_{L_{\rho}^{2}}^{2}$
$\displaystyle\leq
8\left(\lambda\overline{r}\sqrt{M}+|\tilde{\gamma}_{>M}|\left\|E\right\|_{2}\right)^{2}+\frac{5}{2}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}$
$\displaystyle\leq
8\left(\lambda\overline{r}\sqrt{M}+\frac{1}{2}|\tilde{\gamma}_{>M}|\right)^{2}+\frac{5}{2}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}<C_{1}.$

Hence we have an upper bound for the bias. We argue similarly for the lower
bound:

fitting error $\displaystyle\geq\left(1-2\sqrt{\frac{\log
N}{N}}\right)\frac{\lambda^{2}\lambda_{M}}{(\lambda_{M}+\lambda)^{2}}\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}-\frac{1}{4}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}\sqrt{\frac{\log
N}{N}}$ $\displaystyle\quad-\frac{\log
N}{N}\left(\left\|w\right\|_{2}^{2}p(\delta)+(1+2\delta)\tilde{\gamma}_{>M}^{2}+|\tilde{\gamma}_{>M}|\|\tilde{f}_{\leq
M}\|_{L_{\rho}^{2}}\right)$
$\displaystyle\geq\frac{\lambda^{2}\lambda_{M}}{(\lambda_{M}+\lambda)^{2}}\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}-\left(\frac{1}{4}\|\tilde{f}\|_{L_{\rho}^{2}}^{2}+\frac{2\lambda^{2}}{\lambda_{1}+\lambda}\|\tilde{f}_{\leq
M}\|_{\mathcal{H}}^{2}\right)\sqrt{\frac{\log N}{N}}-C_{1}\frac{\log N}{N}.$

For $\lambda\to 0$, note that $w\to-\tilde{\gamma}_{>M}\bm{E}$. This yields

$\displaystyle\lim_{\lambda\to 0}\text{fitting error}$
$\displaystyle\leq\lim_{\lambda\to
0}\left\\{(1+2\delta)\left\|\bar{\mathbf{P}}w\right\|_{2}^{2}+\left\|w\right\|_{2}^{2}\delta^{2}p(\delta)\right\\}$
$\displaystyle=(1+2\delta)\left\|-\tilde{\gamma}_{>M}\bm{E}\right\|_{2}^{2}+\left\|-\tilde{\gamma}_{>M}\bm{E}\right\|_{2}^{2}\delta^{2}p(\delta)$
$\displaystyle=\tilde{\gamma}_{>M}^{2}\left\|\bm{E}\right\|_{2}^{2}(1+\delta(2+\delta
p(\delta))).$

Hence, by plugging in $\delta<\frac{1}{2}$, with probability of at least
$1-2/N$,

$\displaystyle\lim_{\lambda\to 0}\text{fitting error}$
$\displaystyle\leq\tilde{\gamma}_{>M}^{2}\frac{\log
N}{N}\left(1+6\sqrt{\frac{\log N}{N}}\right)$ $\displaystyle\lim_{\lambda\to
0}\text{bias}$ $\displaystyle\leq\left(1+\frac{\log
N}{N}\right)\tilde{\gamma}_{>M}^{2}+6\tilde{\gamma}_{>M}^{2}\left(\frac{\log
N}{N}\right)^{\frac{3}{2}}.$

For lower bound, it follows similarly:

$\displaystyle\lim_{\lambda\to 0}\text{bias}$
$\displaystyle\geq\left(1-\frac{\log
N}{N}\right)\tilde{\gamma}_{>M}^{2}-6\tilde{\gamma}_{>M}^{2}\left(\frac{\log
N}{N}\right)^{\frac{3}{2}},$

and we obtain line (37). For the case where $\tilde{\gamma}_{>M}=0$,
recalculate and simplify line (40) to obtain line (38). ∎

####  C.4.2 Refined Bounds on Variance

Similarly, we can refine Theorem C.14 to get a bound on the variance:

######  Theorem C.20.

For $N>\max\left\\{(12G)^{4}(M+1)^{2},9\right\\}$, and set $C_{2}=12$
(independent to $N$), with a probability of at least $1-2/N$, we have the
upper and lower bounds of variance:

variance $\displaystyle\leq\sigma^{2}\frac{M}{N}\left(1+\sqrt{\frac{\log
N}{N}}+C_{2}\frac{\log N}{N}\right);$ variance
$\displaystyle\geq\frac{\lambda_{M}^{2}}{(\lambda_{M}+\lambda)^{2}}\sigma^{2}\frac{M}{N}\left(1-\sqrt{\frac{\log
N}{N}}\right)-C_{2}\sigma^{2}\frac{M}{N}\frac{\log N}{N}.$

###### Proof.

We argue analogously as in Theorem C.19: by Proposition C.14 and Lemma C.16,
we have

variance
$\displaystyle\leq(1+\delta)\frac{\sigma^{2}}{N}\sum_{k=1}^{M}\frac{\lambda_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}+M\frac{\sigma^{2}}{N}(1+\delta)\delta^{2}p(\delta)$
$\displaystyle\leq(1+\delta)\sigma^{2}\frac{M}{N}+\sigma^{2}\frac{M}{N}(1+\delta)\delta^{2}p(\delta)$
$\displaystyle\leq\left(1+\sqrt{\frac{\log
N}{N}}\right)\sigma^{2}\frac{M}{N}+\sigma^{2}\frac{M}{N}\frac{\log
N}{N}\left(1+\frac{1}{2}\right)p\left(\frac{1}{2}\right)$
$\displaystyle\leq\left(1+\sqrt{\frac{\log
N}{N}}\right)\sigma^{2}\frac{M}{N}+12\sigma^{2}\frac{M}{N}\frac{\log N}{N}$

Hence we can choose $C_{2}=12$. For the lower bound, since
$\frac{\lambda_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}>\frac{\lambda_{M}^{2}}{(\lambda_{M}+\lambda)^{2}}$,
we have

variance
$\displaystyle\geq(1-\delta)\frac{\sigma^{2}}{N}\sum_{k=1}^{M}\frac{\lambda_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}-M\frac{\sigma^{2}}{N}(1+\delta)\delta^{2}p(\delta)$
$\displaystyle\geq\frac{\lambda_{M}^{2}}{(\lambda_{M}+\lambda)^{2}}\sigma^{2}\frac{M}{N}\left(1-\sqrt{\frac{\log
N}{N}}\right)-12\sigma^{2}\frac{M}{N}\frac{\log N}{N}.$

∎

Note that in both Theorems C.19 and C.20, the constants $C_{1},C_{2}>0$ is not
optimized.

##  Appendix D Numerical Validation

In this section, we illustrate our result for KRR with two different finite
rank kernels.

###  D.1 Truncated NTK

First, we need to define a finite-rank kernel
$K:\mathcal{X}\times\mathcal{X}\to\mathbb{R}$. We set
$\mathcal{X}=\mathbb{S}^{1}\subset\mathbb{R}^{2}$. By reparametrization, we
write $\mathbb{S}^{1}\cong[0,2\pi]/_{0\sim 2\pi}$. We assume the data are
drawn uniformly on the circle, that is
$\rho_{\mathcal{X}}=\operatorname{unif}[\mathbb{S}^{1}]$. We can use the
Fourier functions $\cos(k\cdot),\sin(k\cdot)$ as the orthogonal eigenfunctions
of the kernel. We define the NTK

$K^{(\infty)}(\theta,\theta^{\prime})\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\frac{\cos(\theta-\theta^{\prime})\left(\pi-|\theta-\theta^{\prime}|\right)}{2\pi}$

for all $\theta,\theta^{\prime}\in[0,2\pi]$. 3) We choose a rank-$M$
truncation
$K(\theta,\theta^{\prime})=\sum_{k=1}^{M}\lambda_{k}\psi_{k}(\theta)\psi_{k}(\theta^{\prime})$
for all $\theta,\theta^{\prime}\in[0,2\pi]$. For the first few eigenvalues of
the kernel, please see Table 3 for example.

$k$ | 1 | 2 | 3 | 4 | 5 | 6 | 7 | $\infty$
---|---|---|---|---|---|---|---|---
$\lambda_{k}$ | $\frac{1}{\pi^{2}}$ | $\frac{1}{8}$ | $\frac{1}{8}$ | $\frac{5}{9\pi^{2}}$ | $\frac{5}{9\pi^{2}}$ | $\frac{17}{225\pi^{2}}$ | $\frac{17}{225\pi^{2}}$ | -
$\psi_{k}(\theta)$ | 1 | $\sqrt{2}\cos(\theta)$ | $\sqrt{2}\sin(\theta)$ | $\sqrt{2}\cos(2\theta)$ | $\sqrt{2}\sin(2\theta)$ | $\sqrt{2}\cos(4\theta)$ | $\sqrt{2}\sin(4\theta)$ | -
$\sum_{k^{\prime}=0}^{k}\lambda_{k^{\prime}}$ | 0.1013 | 0.2263 | 0.3513 | 0.4076 | 0.4639 | 0.4716 | 0.4792 | 0.5

Table 3: The first few eigenvalues of the NTK

Before proceeding to test error computation, we present a training example,
Figure 4, to give readers an intuition on the truncated NTK (tNTK).

Figure 4:  (left): NTK training; (right): tNTK training where $N=50,M=7$.
$\sigma^{2}=0.05,\lambda=\sigma^{2}/N$.

###  D.2 Test Error Computations

In the following tNTK training, we set the hyperparameters as follows:

##### Target function

We choose a simple target function $\tilde{f}(x)=\cos
x=\frac{1}{\sqrt{2}}\psi_{2}(x)$. Throughout the experiment, we set the noise
variance $\sigma^{2}=0.05$.

##### Ridge

We choose $\lambda=\frac{\sigma^{2}}{N}$. In Figure 5 (left), we set
$N=50,\lambda=0.05/50$ for tNTK training; (right) we set set $\lambda=0.05/50$
for varying $N$ from 10 to 200.

##### Error bars

In Figure 7 (right), for each value of $N$, we run over 10 iterations of
random samples and compute the test error. The error bars are shown as the
difference between the upper and the lower quartiles.

Figure 5:  (left): tNTK training; (right): the decay of test error as $N$
varies.

##### Lower bound

See the subsection below.

###  D.3 Bound Comparison

We continue with the experiment on the tNTK this time with varying $N$ and
compare our upper bound with Bach [2023a].

##### Upper bounds

In Figure 6, the expression of Bach’s and our upper bounds are directly
computed:

Bach’s upper bound
$\displaystyle=4\lambda\|\tilde{f}\|_{\mathcal{H}}^{2}+\frac{8\sigma^{2}R^{2}}{\lambda
N}(1+2\log N)$ Our upper bound without residue
$\displaystyle=\lambda\|\tilde{f}\|_{\mathcal{H}}^{2}\left(1+2\sqrt{\frac{\log
N}{N}}\right),$

where the constants $\|\tilde{f}\|_{\mathcal{H}}^{2}$ and $R^{2}$ can be
computed directed from the choice of kernel and target function. For
simplicity reason, we drop the residue term $C_{1}\frac{\log N}{N}$ since it
is overshadowed by the other terms and the constant $C_{1}$ is not optimized.

Figure 6:  Test error bound improvement on tNTK. Same as Figure 2.

###  D.4 Legendre Kernel

To illustrate the bounds with another finite-rank, we choose a simple legendre
kernel (LK):

$K(x,z)=\sum_{k=0}^{M}\lambda_{k}P_{k}(x)P_{k}(z)$

where $P_{k}$ is the Legendre polynomial of degree $k$, and $\lambda_{k}>0$
are the eigenvalues.

##### Eigenvalues

To better compare the Legendre kernel $K$ with the NTK, we choose
$\lambda_{k}=C\cdot(k+1)^{-2}$ of quadratic decay such that the spectral sums
are the same: $\sum_{k=0}^{\infty}\lambda_{k}=0.5$. Hence we choose
$C=0.5/\sum_{k=1}^{\infty}k^{-2}=\frac{3}{\pi^{2}}$.

##### Target function

We choose a simple target function
$\tilde{f}(x)=x^{2}=\frac{1}{3}P_{0}(x)+\frac{2}{3}P_{2}(x)$. Throughout the
experiment, we set the noise variance $\sigma^{2}=0.05$.

###  D.5 Test Error Computation

##### Ridge

As before, our bound suggests that, to balance the bias and the variance with
a fixed $N$, we can choose $\lambda=\frac{\sigma^{2}}{N}$. In Figure 7 (left),
we set $N=50,\lambda=0.05/50$ for KRR training; (right) we set set
$\lambda=0.05/50$ for varying $N$ from 10 to 200.

##### Error bars

In Figure 7 (right), for each value of $N$, we run over 10 iterations of
random samples and compute the test error. The error bars are shown as the
different between the upper and the lower quartiles. The median is taken as
average.

Figure 7:  (left): LK training; (right): the decay of test error as $N$
varies. Same as Figure 1.

##### Upper bounds

In Figure 8, the expression of Bach’s and our upper bounds are directly
computed:

Bach’s upper bound
$\displaystyle=4\lambda\|\tilde{f}\|_{\mathcal{H}}^{2}+\frac{8\sigma^{2}R^{2}}{\lambda
N}(1+2\log N)$ Our upper bound without residue
$\displaystyle=\lambda\|\tilde{f}\|_{\mathcal{H}}^{2}\left(1+2\sqrt{\frac{\log
N}{N}}\right),$

where the constants $\|\tilde{f}\|_{\mathcal{H}}^{2}$ and $R^{2}$ can be
computed directed from the choice of kernel and target function.

Figure 8:  Test error bound improvement on LK. Same as Figure 2.

##### Lower Bound

Last but not least, we need to show our lower bound is valid. To see this
clearly, we need to write the bound in exact sums instead of in HKRS norm
square $\|\tilde{f}\|_{\mathcal{H}}^{2}$: namely, we compute $I$

$\frac{\lambda^{2}\lambda_{M}}{(\lambda_{M}+\lambda)^{2}}\|\tilde{f}\|_{\mathcal{H}}^{2}\leq
I=\lambda^{2}\sum_{k=1}^{M}\frac{\tilde{\gamma}_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}\leq\lambda\|\tilde{f}\|_{\mathcal{H}}^{2},$
(42)

instead of using the inequality (42) in Lemma C.18; and

$M\frac{\lambda_{M}^{2}}{(\lambda_{M}+\lambda)^{2}}\leq\sum_{k=1}^{M}\frac{\lambda_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}\leq
M,$ (43)

instead of using the inequality (42) in Theorem C.20. Then we can compute our
bounds as:

Our upper bound without residue
$\displaystyle=\lambda^{2}I\left(1+2\sqrt{\frac{\log
N}{N}}\right)+\frac{\sigma^{2}}{N}\sum_{k=1}^{M}\frac{\lambda_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}\left(1+\sqrt{\frac{\log
N}{N}}\right),$ Our lower bound without residue
$\displaystyle=\lambda^{2}I\left(1-2\sqrt{\frac{\log
N}{N}}\right)+\frac{\sigma^{2}}{N}\sum_{k=1}^{M}\frac{\lambda_{k}^{2}}{(\lambda_{k}+\lambda)^{2}}\left(1-\sqrt{\frac{\log
N}{N}}\right),$

and we drop the residue terms $C_{1}\frac{\log N}{N}$ and
$C_{2}\frac{\sigma^{2}}{N}M\frac{\log N}{N}$ by the same reason as before.
From Figure 3, we can see that our bounds precisely describe the decay of the
test error. Our bounds are not ‘bounding’ the test errors in smaller instances
due to the absence of the residue terms, which increases the interval of
confidence of our approximation. But for larger instances, say $N>100$, all
upper and lower bounds, and the averaged test error converge to the same
limit.

Figure 9:  Our bounds comparing to the averaged test error with varying $N$,
over 10 iterations. Same as Figure 3.