Early stopping and non-parametric regression: An optimal data-dependent stopping rule

Introduction

The phenomenon of overfitting is ubiquitous throughout statistics. It is especially problematic in nonparametric problems, where some form of regularization is essential to prevent overfitting. In the problem of nonparametric regression, the most classical form of regularization is that of Tikhonov regularization, where a quadratic smoothness penalty is added to the least-squares loss. An alternative and algorithmic approach to regularization is based on early stopping of an iterative algorithm, such as gradient descent applied to the unregularized loss function. The main advantage of early stopping for regularization, as compared to penalized forms, is lower computational complexity.

The idea of early stopping has a fairly lengthy history, dating back to the 1970’s in the context of the Landweber iteration. (For instance, see the paper by Strand , with follow-up work by Anderssen and Prenter as well as Wahba .) Early stopping has also been widely used in neural networks (e.g., ), in which stochastic gradient descent is used to estimate the network parameters. Past work provided intuitive arguments for the benefits of early stopping. It was argued that each step of an iterative algorithm will reduce bias but increase variance, so early stopping ensures the variance of the estimator is not too high. However, prior to the 1990s, there had been little theoretical justification for these claims. A more recent line of work has provided theoretical justification for various types of early stopping, including boosting algorithms (e.g., ), greedy methods , gradient descent over reproducing kernel Hilbert spaces (e.g. ), the conjugate gradient algorithm , and the power method for eigenvalue computation . Most relevant to our work is the work of Bühlmann and Yu , who derived optimal mean-squared error bounds for $L^{2}$-boosting with early stopping in the case of fixed design regression. However, these optimal rates are based on an “oracle” stopping rule, one that cannot be computed based on the data. Thus, their work left open the following natural question: is there a data-dependent and easily computable stopping rule that produces a minimax-optimal estimator?

The main contribution of this paper is to answer this question in the affirmative for a certain class of non-parametric regression problems, in which the underlying regression function belongs to a reproducing kernel Hilbert space (RKHS). In this setting, a standard estimator is the method of kernel ridge regression (e.g., ), which minimizes a weighted sum of the least-squares loss with a squared Hilbert norm penalty as a regularizer. Instead of a penalized form of regression, we analyze early stopping of an iterative update that is equivalent to gradient descent on the least-squares loss in an appropriately chosen coordinate system. By analyzing the mean-squared error of our iterative update, we derive a data-dependent stopping rule that provides the optimal trade-off between the estimated bias and variance at each iteration. In particular, our stopping rule is based on the first time that a running sum of step-sizes after $t$ steps increases above the critical trade-off between bias and variance. For Sobolev spaces and other types of kernel classes, we show that the function estimate obtained by this stopping rule achieves minimax-optimal estimation rates in both the empirical and population norms. Importantly, our stopping rule does not require the use of cross-validation or hold-out data.

Background and problem formulation

We begin by introducing some background on non-parametric regression and reproducing kernel Hilbert spaces, before turning to a precise formulation of the problem studied in this paper.

where $w_{i}:\,=y_{i}-f^{*}(x_{i})$ is a zero-mean noise random variable. Throughout this paper, we assume that the random variables $w_{i}$ are sub-Gaussian with parameter $\sigma$, meaning that

For instance, this sub-Gaussian condition is satisfied for normal variates $w_{i}\sim N(0,\sigma^{2})$, but it also holds for various non-Gaussian random variables. Parts of our analysis also apply to the fixed design setting, in which we condition on a particular realization $\{x_{i}\}_{i=1}^{n}$ of the covariates.

Since the eigenfunctions $\{\phi_{k}\}_{k=1}^{\infty}$ form an orthonormal basis, any function $f\in\mathcal{H}$ has an expansion of the form $f(x)=\sum_{k=1}^{\infty}\sqrt{\lambda_{k}}a_{k}\phi_{k}(x)$, where for all $k$ such that $\lambda_{k}>0$, the coefficients

The second inner product, denoted $\langle f,\,g\rangle_{\mathcal{H}}$, is the one that defines the Hilbert space; it can be written in terms of the expansion coefficients as

Using this definition, the Hilbert ball of radius $1$ for the Hilbert space $\mathcal{H}$ with eigenvalues $\{\lambda_{k}\}_{k=1}^{\infty}$ and eigenfunctions $\{\phi_{k}\}_{k=1}^{\infty}$ takes the form

The class of reproducing kernel Hilbert spaces contains many interesting classes that are widely used in practice, including polynomials of degree $d$, Sobolev spaces with smoothness $\nu$, and Gaussian kernels. For more background and examples on reproducing kernel Hilbert spaces, we refer the reader to various standard references .

Throughout this paper, we assume that any function $f$ in the unit ball of the Hilbert space is uniformly bounded, meaning that there is some constant $B<\infty$ such that

2 Gradient update equation

We now turn to the form of the gradient update that we study in this paper. Given the samples $\{(x_{i},y_{i})\}_{i=1}^{n}$, consider minimizing the least-squares loss function

A direct approach would be to perform gradient descent on this form of the least-squares loss. For our purposes, it turns out to be more natural to perform gradient descent in the transformed co-ordinate system $\theta=\sqrt{K}\,\omega$. Some straightforward calculations (see Appendix A for details) yield that the gradient descent algorithm in this new co-ordinate system generates a sequence of vectors $\{\theta^{t}\}_{t=0}^{\infty}$ via the recursion

where $\{\alpha^{t}\}_{t=0}^{\infty}$ is a sequence of positive step sizes (to be chosen by the user). We assume throughout that the gradient descent procedure is initialized with $\theta^{0}=0$.

corresponds to the usual mean-squared error.

3 Overfitting and early stopping

Figure 1 provides plots of the squared prediction error $\|f_{t}-f^{*}\|_{n}^{2}$ as a function of the iteration number $t$. For both kernels, the prediction error decreases fairly rapidly, reaching a minimum before or around $T\approx 20$ iterations, before then beginning to increase.

As the analysis of this paper will clarify, too many iterations lead to fitting the noise in the data (i.e., the additive perturbations $w_{i}$), as opposed to the underlying function $f^{*}$. In a nutshell, the goal of this paper is to quantify precisely the meaning of “too many” iterations, and in a data-dependent and easily computable manner.

Main results and their consequences

In more detail, our main contribution is to formulate a data-dependent stopping rule, meaning a mapping from the data $\{(x_{i},y_{i})\}_{i=1}^{n}$ to a positive integer $\widehat{T}$, such that the two forms of prediction error $\|f_{\widehat{T}}-f^{*}\|_{n}$ and $\|f_{\widehat{T}}-f^{*}\|_{2}$ are minimal. In our formulation of such a stopping rule, two quantities play an important role: first, the running sum of the step sizes

and secondly, the eigenvalues $\widehat{\lambda}_{1}\geq\widehat{\lambda}_{2}\geq\cdots\geq\widehat{\lambda}_{n}\geq 0$ of the empirical kernel matrix $K$ previously defined (8). The kernel matrix and hence these eigenvalues are computable from the data. We also note that there is a large body of work on fast computation of kernel eigenvalues (e.g., see the paper and references therein).

Our stopping rule involves the use of a model complexity measure, familiar from past work on uniform laws over kernel classes , known as the local empirical Rademacher complexity. For the kernel classes studied in this paper, it takes the form

For a given noise variance $\sigma>0$, a closely related quantity—one of central importance to our analysis—is critical empirical radius $\widehat{\varepsilon}_{n}>0$, defined to be the smallest positive solution to the inequality

The existence and uniqueness of $\widehat{\varepsilon}_{n}$ is guaranteed for any reproducing kernel Hilbert space; see Appendix D for details. As clarified in our proof, this inequality plays a key role in trading off the bias and variance in a kernel regression estimate.

Our stopping rule is defined in terms of an analogous inequality that involves the running sum $\eta_{t}=\sum_{\tau=0}^{t-1}{\alpha^{\tau}}$ of the step sizes. Throughout this paper, we assume that the step sizes are chosen to satisfy the following properties:

Boundedness: $0\;\leq\;\alpha^{\tau}\;\leq\;\min\{1,1/\widehat{\lambda}_{1}\}$ for all $\tau=0,1,2,\ldots$.

Non-increasing: $\alpha^{\tau+1}\leq\alpha^{\tau}$ for all $\tau=0,1,2,\ldots$.

Infinite travel: the running sum $\eta_{t}=\sum_{\tau=0}^{t-1}\alpha^{\tau}$ diverges as $t\rightarrow+\infty$.

We refer to any sequence $\{\alpha^{\tau}\}_{\tau=0}^{\infty}$ that satisfies these conditions as a valid stepsize sequence. We then define the stopping time

As discussed in Appendix D, the integer $\widehat{T}$ belongs to the interval $[0,\infty)$ and is unique for any valid stepsize sequence. As will be clarified in our proof, the intuition underlying the stopping rule (15) is that the sum of the step-sizes $\eta_{t}$ acts as a tuning parameter that controls the bias-variance tradeoff. The stated choice of $\widehat{T}$ optimizes this trade-off.

The following result applies to any sequence $\{f_{t}\}_{t=0}^{\infty}$ of function estimates generated by the gradient iteration (9) with a valid stepsize sequence.

Given the stopping time $\widehat{T}$ defined by the rule (15), there are universal positive constants $(c_{1},c_{2})$ such that the following events both hold with probability at least $1-c_{1}\exp(-c_{2}n\widehat{\varepsilon}_{n}^{2})$:

For all iterations $t=1,2,...,\widehat{T}$:

At the iteration $\widehat{T}$ chosen according to the stopping rule (15), we have

Although the bounds (a) and (b) are stated as high probability claims, a simple integration argument can be used to show that the expected mean-squared error (over the noise variables, with the design fixed) satisfies a bound of the form

Moreover, as will be clarified in corollaries to follow, Theorem 1 can be used to show that our stopping rule provides minimax-optimal rates for various function classes. The interpretation of Theorem 1 is as follows: if the sum of the step-sizes $\eta_{t}$ remains below the threshold defined by (15), applying the gradient update (9) reduces the prediction error. Moreover, note that for Hilbert spaces with a larger kernel complexity, the stopping time $\widehat{T}$ is smaller, since fitting functions in a larger class incurs a greater risk of overfitting.

Using this complexity measure, we define the critical population rate $\varepsilon_{n}$ to be the smallest positive solution to the inequality

(Our choice of the pre-factor $40$ is for later theoretical convenience.) In contrast to the critical empirical rate $\widehat{\varepsilon}_{n}$, this quantity is not data-dependent, since it depends on the population eigenvalues of the RKHS $\mathcal{H}$.

with probability at least $1-c_{1}\exp(-c_{2}n\widehat{\varepsilon}_{n}^{2})$.

Theorems 1 and 2 are general results that apply to any reproducing kernel Hilbert space. Their proofs involve combination of direct analysis of our iterative update (9) combined with techniques from empirical process theory and concentration of measure ; see Section 4 for the details.

To compare with the past work of Bühlmann and Yu , they also provide a theoretical analysis for gradient descent (referred to as $L^{2}$-boosting in their paper), focusing exclusively on the fixed design case. Our theory applies to random as well as fixed design, and a broader set of step-size choices. The most significant difference between Theorem 1 in our paper and Theorem 3 in the paper is that we provide a data-dependent stopping rule where as their analysis does not lead to a computable stopping rule.

2 Some consequences for specific kernel classes

Let us now illustrate some consequences of our general theory for special choices of kernels that are of interest in practice.

We begin with the class of RKHSs whose eigenvalues satisfy a polynomial decay condition, meaning that

For the uniform measure on $$, this class exhibits polynomial eigendecay (23) with$\nu=1$. For any class that satisfies the polynomial decay condition, we have the following corollary:

Suppose that in addition to the assumptions of Theorem 2, the kernel class $\mathcal{H}$ satisfies the polynomial eigenvalue decay (23) for some parameter $\nu>1/2$. Then there is a universal constant $c_{5}$ such that

Moreover, if $\lambda_{k}\geq c\,(1/k)^{2\nu}$ for all $k=1,2,\ldots$, then

The proof, provided in Section 4.3, involves showing that the population critical rate (20) is of the order $\mathcal{O}(n^{-\frac{2\nu}{2\nu+1}})$. By known results on non-parametric regression , the error bound (25) is minimax-optimal.

In the special case of the first-order spline family (24), Corollary 1 guarantees that

In order to test the accuracy of this prediction, we performed the following set of simulations. First, we generated samples from the observation model

where $x_{i}=i/n$, and $w_{i}\sim N(0,\sigma^{2})$ are i.i.d. noise terms. We present results for the function $f^{*}(x)=|x-1/2|-1/2$, a piecewise linear function belonging to the first-order Sobelev class. For all our experiments, the noise variance $\sigma^{2}$ was set to one, but so as to have a data-dependent method, this knowledge was not provided to the estimator. There is a large body of work on estimating the noise variance $\sigma^{2}$ in non-parametric regression (see e.g. Hall and Marron ). For our simulations, we use the simple estimator based on Hall and Marron . They proved that their estimator is ratio consistent, which is sufficient for our purposes.

For a range of sample sizes $n$ between $10$ and $300$, we performed the updates (9) with constant stepsize $\alpha=0.25$, stopping at the specified time $\widehat{T}$. For each sample size, we performed $10,000$ independent trials, and averaged the resulting prediction errors. In panel (a) of Figure 2, we plot the mean-squared error versus the sample size, which shows consistency of the method. We also plotted the mean-squared error raised to the power $-3/2$ versus the sample size. After this rescaling, the bound (27) predicts a linear relation, as is observed in panel (b) of Figure 2. We also performed the same experiments for the case of randomly drawn designs $x_{i}\sim\mbox{Unif}(0,1)$. In this case, we observed similar results but with more trials required to average out the additional randomness in the design.

If, in addition to the conditions of Theorem 2, the kernel has finite rank $m$, then

3 Comparison with other stopping rules

In this section, we provide a comparison of our stopping rule to two other stopping rules, as well as a oracle method (that involves knowledge of $f^{*}$, and so cannot be computed in practice).

First, we consider a simple hold-out method: it performs gradient descent using $50\%$ of the data, and uses the other $50\%$ of the data to estimate the risk (e.g. ). Assuming that the sample size is even for simplicity, we split the full data set $\{x_{i}\}_{i=1}^{n}$ into two equally sized subsets $S_{\mbox{\tiny{tr}}}$ and $S_{te}$. The data indexed by the training set $S_{\mbox{\tiny{tr}}}$ is used to estimate the function $f_{\mbox{\tiny{tr}}}^{t}$ using the gradient descent update (9). At each iteration $t=0,1,2,\ldots$, the data indexed by $S_{te}$ is used to estimate the risk via R_{\mbox{\tiny{HO}}}(f_{t})=\frac{1}{n}\sum_{i\in S_{te}}\big{(}y_{i}-f_{\mbox{\tiny{tr}}}^{t}(x_{i})\big{)}^{2}, which defines the stopping rule

A line of past work ) has analyzed stopping rules based on this type of hold-out rule. For instance, Caponetto analyzes a hold-out method, and shows that it yields rates that are optimal for Sobolev spaces with $\nu\leq 1$ but not in general. The major drawback of using hold-out as that a percentage of the data is lost which increases the risk.

which is easy to compute. This risk estimate defines the associated stopping rule

In contrast with hold-out, this approach makes use of all the data. However, we are not aware of any theoretical guarantees for early stopping using the stopping rule (32).

It can be shown for both stopping rules (30) and (32), a valid sequence of step-sizes guarantees existence and uniqueness of the stopping point. Note that our stopping rule $\widehat{T}$ based on (15) requires estimation of both the empirical eigenvalues, and the noise variance $\sigma^{2}$. In contrast, the SURE-based rule requires estimation of $\sigma^{2}$ but not the empirical eigenvalues, whereas the hold-out rule requires no parameters to be estimated, but a percentage of the data is used to estimate the risk.

As a third point of reference, we also plot the mean-squared error for an “oracle” method. It is allowed to base its stopping time on the exact prediction error $R_{\mbox{\tiny{OR}}}(f^{t})=\|f^{t}-f^{*}\|_{n}^{2}$, which defines the oracle stopping rule

Note that this stopping rule is not computable from the data, since it assumes exact knowledge of the function $f^{*}$ that we are trying to estimate.

In order to compare our stopping rule (15) with these alternatives, we generated i.i.d. samples from the previously described model (see equation (28) and the following discussion). We varied the sample size $n$ from $10$ to $300$, and for each sample size, we performed $10,000$ independent trials (randomizations of the noise variables $\{w_{i}\}_{i=1}^{n}$), and computed the average of squared prediction error.

Figure 3 plots the resulting mean-squared errors of our stopping rule, the hold-out stopping rule (30), the SURE-based stopping rule (32), and the oracle rule (33). Panel (a) shows the mean-squared error versus sample size, whereas panel (b) shows the same curves in terms of logarithm of mean-squared error. Our proposed rule exhibits better performance than the hold-out and SURE-based rules for sample sizes $n$ larger than $50$. On the flip side, since the construction of our stopping rule is based on the assumption that $f^{*}$ belongs to a known RKHS, it is unclear how robust it would be to model mis-specification. In contrast, the hold-out and SURE-based stopping rules are generic methods, not based directly on the RKHS structure, so might be more robust to model mis-specification. Thus, one interesting direction is to explore the robustness of our stopping rule. On the theoretical front, it would be interesting to determine whether the hold-out and/or SURE-based stopping rules can be proven to achieve minimax optimal rates for general kernels, as we have established for our stopping rule.

4 Connections to kernel ridge regression

We conclude by presenting an interesting link between our early stopping procedure and kernel ridge regression. The kernel ridge regression (KRR) estimate is defined as

where $\nu$ is the (inverse) regularization parameter. For any $\nu<\infty$, the objective is strongly convex, so that the KRR solution is unique.

Friedman and Popescu observed through simulations that the regularization paths for early stopping of gradient descent and ridge regression are similar, but did not provide any theoretical explanation of this fact. As an illustration of this empirical phenomenon, Figure 4 compares the prediction error $\|\widehat{f}_{\nu}-f^{*}\|_{n}^{2}$ of the kernel ridge regression estimate over the interval $\nu\in$ versus that of the gradient update (9) over the first $100$ iterations. Note that the curves, while not identical, are qualitatively very similar.

From past theoretical work [36, 26, e.g.,], kernel ridge regression, with the appropriate setting of the penalty parameter $\nu$, is known to achieve minimax-optimal error for various kernel classes, among them the Sobolev and finite-rank kernels for which stopping rule is provably optimal. In this section, we provide a theoretical basis for these connections, in particular by showing that if the inverse penalty parameter $\nu$ is chosen using the same criterion as our stopping rule, then the prediction error satisfies the same type of bounds (with $\nu$ now playing the role of the running sum $\eta_{t}$).

More precisely, suppose that we choose $\widehat{\nu}$ to be the smallest positive solution to the inequality

Note that this criterion is identical to the one underlying our stopping rule, except that the continuous parameter $\nu$ replaces the discrete parameter $\eta_{t}=\sum_{\tau=0}^{t-1}{\alpha^{\tau}}$.

Consider the kernel ridge regression estimator (34) applied to $n$ i.i.d. samples $\{(x_{i},y_{i})\}$ with $\sigma$-sub Gaussian noise. Then there are universal constants $(c_{1},c_{2},c_{3})$ such that for all $\delta>0$, the following claims hold with probability at least $1-c_{1}\exp(-c_{2}\,n\,\widehat{\varepsilon}_{n}^{2})$:

For all $0<\nu\leq\widehat{\nu}$, we have

With $\widehat{\nu}$ chosen according to the rule (35), we have

Moreover, for all $\nu>\widehat{\nu}$, we have

Note that (apart from a slightly different leading constant) the upper bound (36) is identical to the upper bound in equation (16) in Theorem 1. The only difference is that the inverse regularization parameter $\nu$ replaces the running sum $\eta_{t}=\sum_{\tau=0}^{t-1}{\alpha^{\tau}}$. Similarly, part (b) of Proposition 1 guarantees that the kernel ridge regression (34) has prediction error that is upper bounded by the empirical critical rate $\widehat{\varepsilon}_{n}^{2}$, as in part (b) of Theorem 1. Let us emphasize that bounds of this type on kernel ridge regression have been derived in past work [26, 36, e.g.,]. The novelty here is that the structure of our result reveals the intimate connection to early stopping, and in fact, the proofs follow a parallel thread.

In conjunction, Proposition 1 and Theorem 1 provide a theoretical explanation for why, as shown in Figure 4, the paths of the gradient descent update (9) and kernel ridge regression estimate (34) are so similar. However, it is important to emphasize that from a computational point of view, early stopping has certain advantages over kernel ridge regression. In general, solving a quadratic program of the form (34) requires on the order of $\mathcal{O}(n^{3})$ basic operations, and this must be done repeatedly at each new choice of $\nu$. On the other hand, by its very construction, the iterates of the gradient algorithm correspond to the desired path of solutions, and each gradient update involves multiplication by the kernel matrix, incurring $\mathcal{O}(n^{2})$ operations.

Proofs

We now turn to the proofs of our main results. The main steps in each proof are provided in the main text, with some of the more technical results deferred to the appendix.

corresponding to the $n$-vector obtained by evaluating the function $f^{t}$ at all design points, and the short-hand

corresponding to the vector of zero mean sub-Gaussian noise random variables. From equation (7), we have the relation

Consequently, by multiplying both sides of the gradient update (9) by $\sqrt{K}$, we find that the sequence $\{f^{t}(x_{1}^{n})\}_{t=0}^{\infty}$ evolves according to the recursion

Since $\theta^{0}=0$, the sequence is initialized with $f^{0}(x_{1}^{n})=0$. The recursion (41) lies at the heart of our analysis.

is the diagonal matrix of eigenvalues, augmented with $n-r$ zero eigenvalues as needed. We then define a sequence of diagonal shrinkage matrices $S^{t}$ as follows:

The matrix $S^{t}$ indicates the extent of shrinkage towards the origin; since $0\;\leq\;\alpha^{t}\;\leq\;\min\{1,1/\widehat{\lambda}_{1}\}$ for all iterations $t$, in the positive semodefinite ordering, we have the sandwich relation

See Appendix B.1 for the proof of this intermediate claim.

In order to complete the proof of the upper bound in Theorem 1, our next step is to obtain high probability upper bounds on these two terms. We summarize our conclusions in an additional lemma, and use it to complete the proof of Theorem 1(a) before returning to prove it.

For all iterations $t=1,2,\ldots$, the squared bias is upper bounded as

Moreover, there is a universal constant $c_{1}>0$ such that, for any iteration $t=1,2,\ldots,\widehat{T}$,

We can now complete the proof of Theorem 1(a). The bound (16) follows quickly: conditioned on the event V_{t}\leq 5\sigma^{2}\eta_{t}\mathcal{R}^{2}_{K}\big{(}1/\sqrt{\eta_{t}}\big{)}, we have

where inequality (i) follows from (42) in Lemma 1, and inequality (ii) follows from the bounds in Lemma 2 and (iii) follows since $t\leq\widehat{T}$. The lower bound (c) in equation (18) follows from (44).

Turning to the proof of part (b), using the upper bound from (a)

Based on the definition of $\widehat{T}$ and $\widehat{\varepsilon}_{n}$, we are guaranteed that $\frac{1}{\eta_{\widehat{T}+1}}\leq\widehat{\varepsilon}_{n}^{2}$, Moreover, by the non-decreasing nature of our step sizes, we have $\alpha^{\widehat{T}+1}\leq\alpha^{\widehat{T}}$, which implies that $\eta_{\widehat{T}+1}\leq 2\eta_{\widehat{T}}$, and hence

Putting together the pieces establishes the bound claimed in part (b).

It remains to establish the bias and variance bounds stated in Lemma 2, and we do so in the following subsections. The following auxiliary lemma plays a role in both proofs:

For all indices $j\in\{1,2,\ldots,r\}$, the shrinkage matrices $S^{t}$ satisfy the bounds

See Appendix B.2 for the proof of this result.

Let us now prove the upper bound (43) on the squared bias. We bound each of the two terms in the definition (42) of $B_{t}^{2}$ in term. Applying the upper bound (45a) from Lemma 3, we see that

Here the final step follows from the fact that $\Psi$ is a unitary operator, so that $\|\Psi^{*}a\|_{2}^{2}\leq\|a\|_{2}^{2}=\|f^{*}\|_{\mathcal{H}}^{2}\leq 1$.

Turning to the second term in the definition (42), we have

where the final step uses the fact that $\Lambda^{1/2}_{jj}=0$ for all $j\in\{r+1,\ldots,n\}$ by construction. Combining the upper bounds (46) and (47) with the definition (42) of $B_{t}^{2}$ yields the claim (43).

1.2 Controlling the variance

Let us now prove the bounds (44) on the variance term $V_{t}$. (To simplify the proof, we assume throughout that $\sigma=1$; the general case can be recovered by a simple rescaling argument). By the definition of $V_{t}$, we have

where the final equality uses the definition of $\mathcal{R}_{K}$. Putting together the pieces, we see that

where $(|\!|\!|A|\!|\!|_{{\scriptsize{\mbox{op}}}},|\!|\!|A|\!|\!|_{{\scriptsize{\mbox{F}}}})$ are (respectively) the operator and Frobenius norms of the matrix $A=\{a_{ij}\}_{i,j=1}^{n}$.

If we apply this result with $A=\frac{2}{n}UQU^{T}$ and $Z_{i}=w_{i}$, then we have $Q=V_{t}$, and moreover

Consequently, the bound (49) implies that

Since $t\leq\widehat{T}$ setting \delta=3\sigma^{2}\eta_{t}\;\biggr{(}\mathcal{R}_{K}(1/\sqrt{\eta_{t}})\biggr{)}, the claim (44) follows.

2 Proof of Theorem 2

This proof is based on the following two steps:

second, showing the empirical critical radius $\widehat{\varepsilon}_{n}$ defined in equation (14) is upper bounded by the population critical radius $\varepsilon_{n}$ defined in equation (21).

Our proof is based on a number of more technical auxiliary lemmas, proved in the appendices. The first lemma provides a high probability bound on the Hilbert norm of the estimate $f_{\widehat{T}}$.

There exist universal constants $c_{1}$ and $c_{2}>0$ such that $\|f_{t}\|_{\mathcal{H}}\leq 2$ for all $t\leq\widehat{T}$ with probability greater than or equal to $1-c_{1}\exp(-c_{2}n\widehat{\varepsilon}_{n}^{2})$.

with probability at least $1-c_{2}\exp(-c_{3}nt^{2})$.

This claim follows from known results on reproducing kernel Hilbert spaces (e.g., Lemma 5.16 in the paper and Theorem 2.1 in the paper ). Our final lemma, proved in Appendix E.2, relates the critical empirical radius $\widehat{\varepsilon}_{n}$ to the population radius $\varepsilon_{n}$:

The inequality $\widehat{\varepsilon}_{n}\leq\varepsilon_{n}$ holds with probability at least $1-c_{1}\exp(-c_{2}n\widehat{\varepsilon}_{n}^{2})$.

With these lemmas in hand, the proof of the theorem is straightforward. First, from Lemma 4, we have $\|f_{\widehat{T}}\|_{\mathcal{H}}\leq 2$ and hence by triangle inequality, $\|f_{\widehat{T}}-f^{*}\|_{\mathcal{H}}\leq 3$ with high probability as well. Next, applying Lemma 5 with $t=\varepsilon_{n}$, we find that

with probability greater than $1-c_{2}\exp(-c_{3}n\varepsilon_{n}^{2})$. Finally, applying Lemma 6 yields that the bound $\|f_{\widehat{T}}-f^{*}\|_{2}^{2}\leq c\varepsilon_{n}^{2}$ holds with the claimed probability.

3 Proof of Corollaries

In each case, it suffices to upper bound the population critical rate $\varepsilon_{n}^{2}$ previously defined.

so that $\varepsilon_{n}^{2}=c^{\prime}\sigma^{2}\frac{m}{n}$.

4 Proof of Proposition 1

We now turn to the proof of our results on the kernel ridge regression estimate (34). The proof follows a very similar structure to that of Theorem 1. Recall the eigendecomposition $K=U\Lambda U^{T}$ of the empirical kernel matrix, and that we use $r$ to denote its rank. For each $\nu>0$, we define the ridge shrinkage matrix

We then have the following analog of Lemma 2 from the proof of Theorem 1:

For any $\nu>0$, the prediction error for the estimate $\widehat{f}_{\nu}$ is bounded as

Note that Lemma 7 is identical to Lemma 2 with the shrinkage matrices $S^{t}$ replaced by their analogues $R^{\nu}$. See Appendix C.1 for the proof of this claim.

Our next step is to show that the diagonal elements of the shrinkage matrices $R^{\nu}$ are bounded:

For all indices $j\in\{1,2,\ldots,r\}$, the diagonal entries $R^{\nu}$ satisfy the bounds

Note that this is the analog of Lemma 3 from Theorem 1, albeit with the constant $\frac{1}{4}$ in the bound (54a) instead of $\frac{1}{2e}$. See Appendix C.2 for the proof of this claim. With these lemmas in place, the remainder of the proof follows as in the proof of Theorem 1.

Discussion

Our analysis and stopping rule may be improved and extended in a number of ways. First, it would interesting to see how our stopping rule can be adapted to mis-specified models. As specified, our method relies on computation of the eigenvalues of the kernel matrix. A stopping rule based on approximate eigenvalue computations, for instance via some form of sub-sampling , would be interesting to study as well.

This work was partially supported by NSF grant DMS-1107000 to MJW and BY. In addition, BY was partially supported by the NSF grant SES-0835531 (CDI), ARO-W911NF-11-1-0114 and the Center for Science of Information (CSoI), an US NSF Science and Technology Center, under grant agreement CCF-0939370, and MJW was also partially supported ONR MURI grant N00014-11-1-086. During this work, GR received partial support from a Berkeley Graduate Fellowship.

Appendix A Derivation of gradient descent updates

In this appendix, we provide the details of how the gradient descent updates (9) are obtained. In terms of the transformed vector $\theta=\sqrt{K}\,\omega$, the least-squares objective takes the form

Given a sequence $\{\alpha^{t}\}_{t=0}^{\infty}$, the gradient descent algorithm operates via the recursion $\theta^{t+1}=\theta^{t}-\alpha^{t}\nabla\widetilde{\mathcal{L}}(\theta^{t})$. Taking the gradient of $\widetilde{\mathcal{L}}$ yields

Substituting into the gradient descent update yields the claim (9).

Appendix B Auxiliary lemmas for Theorem 1

In this appendix, we collect together the proofs of the lemmas for Theorem 1.

where step (i) uses the fact that $\widehat{\lambda}_{j}=0$ for all $j\in\{r+1,\ldots,n\}$. Finally, the orthogonality of $U$ implies that $\|\gamma^{t}-\gamma^{*}\|_{2}^{2}=\frac{1}{n}\|f^{t}(x_{1}^{n})-f^{*}(x_{1}^{n})\|_{2}^{2}$, from which the upper bound (42) follows.

B.2 Proof of Lemma 3

Using the definition of $S^{t}$ and the elementary inequality $1-u\leq\exp(-u)$, we have

Turning to the second set of inequalities, we have $1-[S^{t}]_{jj}=1-\prod_{\tau=0}^{t-1}{(1-\alpha^{\tau}\widehat{\lambda}_{j})}$. By induction, it can be shown that

where step (i) follows from the inequality $1-u\leq\exp(-u)$; and step (ii) follows from the inequality $\exp(-u)\leq(1+u)^{-1}$, valid for $u>0$.

Appendix C Auxiliary results for Proposition 1

In this appendix, we prove the auxiliary lemmas used in the proof of Proposition 1 on kernel ridge regression.

By definition of the KRR estimate, we have \biggr{(}K+\frac{1}{\nu}I\biggr{)}f_{\nu}(x_{1}^{n})=Ky_{1}^{n}. Consequently, some straightforward algebra yields the relation

where the shrinkage matrix $R^{\nu}$ was previously defined (52). The remainder of the proof follows using identical steps to the proof of Lemma 1 with $S^{t}$ replaced by $R^{\nu}$.

C.2 Proof of Lemma 8

By definition (52) of the shrinkage matrix, we have $[R^{\nu}]_{jj}^{2}=(1+\nu\widehat{\lambda}_{j})^{-2}\leq\frac{1}{4\nu\widehat{\lambda}_{j}}$. Moreover, we also have

Appendix D Properties of the empirical Rademacher complexity

In this section, we prove that the $\widehat{\varepsilon}_{n}$ lies in the interval $(0,\infty)$, and is unique. Recall that the stopping point $\widehat{T}$ is defined as \widehat{\varepsilon}_{n}:\,=\arg\min\biggr{\{}\epsilon>0\,\mid\widehat{\mathcal{R}}_{K}\big{(}\epsilon\big{)}\leq\epsilon^{2}/(2e\sigma)\biggr{\}}. Re-arranging and substituting for \widehat{\mathcal{R}}_{K}\big{(}\epsilon\big{)} yields the equivalent expression

Note that \sum_{i=1}^{n}\min\big{\{}\epsilon^{-2}\widehat{\lambda}_{i},1\big{\}} is non-increasing in $\epsilon$ while $n\epsilon^{2}$ is increasing in $\epsilon$. Furthermore when $\epsilon=0$, 0=n\epsilon^{2}<\sum_{i=1}^{n}\min\big{\{}\epsilon^{-2}\widehat{\lambda}_{i},1\big{\}}>0 while for $\epsilon=\infty$, \sum_{i=1}^{n}\min\big{\{}\eta_{t}\widehat{\lambda}_{i},1\big{\}}<n\epsilon^{2}. Hence $\widehat{\varepsilon}_{n}$ exists. Further, $\widehat{\mathcal{R}}_{K}(\epsilon)$ is a continuous function of $\epsilon$ since it is the sum of $n$ continuous functions, Therefore, the critical radius $\widehat{\varepsilon}_{n}$ exists, is unique and satisfies the fixed point equation

Finally, we show that the integer $\widehat{T}$ belongs to the interval $[0,\infty)$ and is unique for any valid sequence of step-sizes. Be the definition of $\widehat{T}$ given by the stopping rule (15) and $\widehat{\varepsilon}_{n}$, we have $\frac{1}{\eta_{\widehat{T}+1}}\leq\widehat{\varepsilon}_{n}^{2}\leq\frac{1}{\eta_{\widehat{T}}}$. Since $\eta_{0}=0$ and $\eta_{t}\rightarrow\infty$ as $t\rightarrow\infty$ and $\widehat{\varepsilon}_{n}\in(0,\infty)$, there exists a unique stopping point $\widehat{T}$ in the interval $[0,\infty)$.

Appendix E Auxiliary results for Theorem 2

This appendix is devoted to the proofs of auxiliary lemmas used in the proof for Theorem 2.

Recall the eigendecomposition $K=U\Lambda U^{T}$ with $\Lambda=\mbox{diag}(\widehat{\lambda}_{1},\widehat{\lambda}_{2},\ldots\widehat{\lambda}_{r})$, and the relation $U^{T}f^{t}(x_{1}^{n})=(I-S^{t})U^{T}y_{1}^{n}$. Substituting into equation (56) yields

where equality (i) follows from the observation equation $y_{1}^{n}=f^{*}(x_{1}^{n})+w$. From Lemma 3, we have $1-S^{t}_{jj}\leq 1$, and hence $C_{t}\leq\frac{1}{n}f^{*}(x_{1}^{n})^{T}U\Lambda^{-1}U^{T}f^{*}(x_{1}^{n})\,\stackrel{{\scriptstyle(i)}}{{\leq}}\,1$, where the last step follows from the analysis in Section 4.1.1.

It remains to derive upper bounds on the random variables $A_{t}$ and $B_{t}$.

where the final inequality follows from the analysis in Section 4.1.1.

where $Q=\mbox{diag}\{\frac{(1-S^{t}_{jj})^{2}}{\widehat{\lambda_{j}}},\;j=1,2,\ldots r\}$. Consequently, $B_{t}$ is a quadratic form in zero-mean sub-Gaussian variables, and using the tail bound (49), we have

But by the definition (15) of the stopping rule and the fact that $t\leq\widehat{T}$, we have

Using the definition of the empirical kernel complexity, we have

where the final inequality holds for $t\leq{\widehat{T}}$, using the definition of the stopping rule.

Putting together the pieces, we have shown that

for all $t\leq\widehat{T}$. Since $\frac{1}{\eta_{t}}\geq\widehat{\varepsilon}_{n}^{2}$ for any $t\leq\widehat{T}$, the claim follows.

E.2 Proof of Lemma 6

In this section, we need to show that $\widehat{\varepsilon}_{n}\leq\varepsilon_{n}$. Recall that $\widehat{\varepsilon}_{n}$ and $\varepsilon_{n}$ satisfy

It suffices to prove that $\widehat{\mathcal{R}}_{K}(\varepsilon_{n})\leq\frac{\varepsilon_{n}^{2}}{2e\sigma}$ using the definition of $\widehat{\varepsilon}_{n}$.

In order to prove the claim, we define the random variables

where $w_{i}\sim N(0,1)$ are i.i.d. standard normal, as well as the associated (deterministic) functions

so that it is also Lipschitz $\frac{t}{\sqrt{n}}$. Consequently, standard concentration results imply that

Let us condition on the two events $\mathcal{A}(t,t_{0}):\,=\{|\widehat{Z}_{n}(w,t)-\widehat{\mathcal{Q}}_{n}(t)|\leq t_{0}\}$ and $\mathcal{A}^{\prime}(t,t_{0}):\,=\{|Z_{n}(w,t)-\mathcal{Q}_{n}(t)|\leq t_{0}\}$. We then have

where inequality (a) follows the first bound in equation (59) with $t_{0}=\frac{\varepsilon_{n}^{2}}{4e\sigma}$ and $t=\varepsilon_{n}^{2}$, inequality (b) follows from Lemma 5 with $t=\varepsilon_{n}$, inequality (c) follows from the second bound (59) with $t_{0}=\frac{\varepsilon_{n}^{2}}{8e\sigma}$ and $t=\varepsilon_{n}^{2}$, and inequality (d) follows from the definition of $\varepsilon_{n}$. Since the events $\mathcal{A}(t,t_{0})$ and $\mathcal{A}^{\prime}(t,t_{0})$ hold with the stated probability, the claim follows.

References