Non-parametric Stochastic Approximation with Large Step sizes

Introduction

Positive-definite-kernel-based methods such as the support vector machine or kernel ridge regression are now widely used in many areas of science of engineering. They were first developed within the statistics community for non-parametric regression using splines, Sobolev spaces, and more generally reproducing kernel Hilbert spaces (see, e.g., ). Within the machine learning community, they were extended in several interesting ways (see, e.g., ): (a) other problems were tackled using positive-definite kernels beyond regression problems, through the “kernelization” of classical unsupervised learning methods such as principal component analysis, canonical correlation analysis, or K-means, (b) efficient algorithms based on convex optimization have emerged, in particular for large sample sizes, and (c) kernels for non-vectorial data have been designed for objects like strings, graphs, measures, etc. A key feature is that they allow the separation of the representation problem (designing good kernels for non-vectorial data) and the algorithmic/theoretical problems (given a kernel, how to design, run efficiently and analyse estimation algorithms).

In practice, estimation algorithms based on regularized empirical risk minimization (e.g., penalized least-squares) face two challenges: (a) using the correct regularization parameter and (b) finding an approximate solution of the convex optimization problems. In this paper, we consider these two problems jointly by following a stochastic approximation framework formulated directly in the RKHS, in which each observation is used only once and overfitting is avoided by making only a single pass through the data–a form of early stopping, which has been considered in other statistical frameworks such as boosting . While this framework has been considered before , the algorithms that are considered either (a) require two sequences of hyperparameters (the step-size in stochastic gradient descent and a regularization parameter) or (b) do not always attain the optimal rates of convergence for estimating the regression function. In this paper, we aim to remove simultaneously these two limitations.

Traditional online stochastic approximation algorithms, as introduced by Robbins and Monro , lead in finite-dimensional learning problems (e.g., parametric least-squares regression) to stochastic gradient descent methods with step-sizes decreasing with the number of observations $n$, which are typically proportional to $n^{-\zeta}$, with $\zeta$ between $1/2$ and 1. Short step-sizes ($\zeta=1$) are adapted to well-conditioned problems (low dimension, low correlations between covariates), while longer step-sizes ($\zeta=1/2$) are adapted to ill-conditioned problems (high dimension, high correlations) but with a worse convergence rate—see, e.g., and references therein. More recently showed that constant step-sizes with averaging could lead to the best possible convergence rate in Euclidean spaces (i.e., in finite dimensions). In this paper, we show that using longer step-sizes with averaging also brings benefits to Hilbert space settings needed for non parametric regression.

With our analysis, based on positive definite kernels, under assumptions on both the objective function and the covariance operator of the RKHS, we derive improved rates of convergence , in both the finite horizon setting where the number of observations is known in advance and our bounds hold for the last iterate (with exact constants), and the online setting where our bounds hold for each iterate (asymptotic results only). It leads to an explicit choice of the step-sizes (which play the role of the regularization parameters) which may be used in stochastic gradient descent, depending on the number of training examples we want to use and on the assumptions we make.

In this paper, we make the following contributions:

We review in Section 2 a general though simple algebraic framework for least-squares regression in RKHS, which encompasses all commonly encountered situations. This framework however makes unnecessary topological assumptions, which we relax in Section 2.5 (with details in App. A).

We characterize in Section 3 the convergence rate of averaged least-mean-squares (LMS) and show how the proper set-up of the step-size leads to optimal convergence rates (as they were proved in ), extending results from finite-dimensional to infinite-dimensional settings. The problem we solve here was stated as an open problem in . Moreover, our results apply as well in the usual finite-dimensional setting of parametric least-squares regression, showing adaptivity of our estimator to the spectral decay of the covariance matrix of the covariates (see Section 4.1).

We compare our new results with existing work, both in terms of rates of convergence in Section 4, and with simulations on synthetic spline smoothing in Section 5.

Sketches of the proofs are given in Appendix B. Complete proofs are given in Appendices I, II.

Learning with positive-definite kernels

Throughout this section, we make the following assumption:

$\mathcal{X}$ is a compact topological space and $\mathcal{H}$ is an RKHS associated with a continuous kernel $K$ on the set $\mathcal{X}$.

RKHSs are well-studied Hilbert spaces which are particularly adapted to regression problems (see, e.g., ). They satisfy the following properties:

$\mathcal{H}$ contains all functions $K_{x}:t\mapsto K(x,t)$, for all $x$ in $\mathcal{X}$.

For any $x\in\mathcal{X}$ and $f\in\mathcal{H}$, the reproducing property holds:

2 Random variables

The space ${L}^{2}_{\rho_{{X}}}$ is then a Hilbert space with norm $\|\cdot\|_{{L}^{2}_{\rho_{{X}}}}$, which we will always assume separable (that is, with a countable orthonormal system).

Throughout this section, we make the following simple assumption regarding finiteness of moments:

Note that under these assumptions, any function in $\mathcal{H}$ in in ${L}^{2}_{\rho_{{X}}}$; however this inclusion is strict in most interesting situations.

3 Minimization problem

We are interested in minimizing the following quantity, which is the prediction error (or mean squared error) of a function $f$, defined for any function in ${L}^{2}_{\rho_{{X}}}$ as:

We are looking for a function with a low prediction error in the particular function space $\mathcal{H}$, that is we aim to minimize $\varepsilon(f)$ over $f\in\mathcal{H}$. We have for $f\in{L}^{2}_{\rho_{{X}}}$:

A key feature of our analysis is that we only considered $\|f-g_{\rho}\|^{2}_{{L}^{2}_{\rho_{{X}}}}$ as a measure of performance and do not consider convergences in stricter norms (which are not true in general). This allows us to neither assume that $g_{\rho}$ is in $\mathcal{H}$ nor that $\mathcal{H}$ is dense in ${L}^{2}_{\rho_{{X}}}$. We thus need to define a notion of the best estimator in $\mathcal{H}$. We first define the closure $\overline{{F}}$ (with respect to $\|\cdot\|_{{L}^{2}_{\rho_{{X}}}}$) of any set $F\subset{L}^{2}_{\rho_{{X}}}$ as the set of limits in ${L}^{2}_{\rho_{{X}}}$ of sequences in ${F}$. The space $\overline{\mathcal{H}}$ is a closed and convex subset in ${L}^{2}_{\rho_{{X}}}$. We can thus define $g_{\mathcal{H}}=\arg\min_{f\in\,\overline{\mathcal{H}}}\varepsilon(g)$, as the orthogonal projection of $g_{\rho}$ on $\overline{\mathcal{H}}$, using the existence of the projection on any closed convex set in a Hilbert space. See Proposition 8 in Appendix A for details. Of course we do not have $g_{\mathcal{H}}\in\mathcal{H}$, that is the minimum in $\mathcal{H}$ is in general not attained.

Seen as a function of $f\in\mathcal{H}$, our loss function $\varepsilon$ is not coercive (i.e., not strongly convex), as our covariance operator (see definition below) $\Sigma$ has no minimal strictly positive eigenvalue (the sequence of eigenvalues decreases to zero). As a consequence, even if $g_{\mathcal{H}}\in\mathcal{H}$, $g_{n}$ may not converge to $g_{\mathcal{H}}$ in $\mathcal{H}$, and when $g_{\mathcal{H}}\notin\mathcal{H}$, we shall even have $\|g_{n}\|_{\mathcal{H}}\rightarrow\infty$.

4 Covariance operator

We now define the covariance operator for the space $\mathcal{H}$ and probability distribution $\rho_{X}$. The spectral properties of such an operator have appeared to be a key point to characterize the convergence rates of estimators .

We implicitly define (via Riesz’ representation theorem) a linear operator $\Sigma:\mathcal{H}\rightarrow\mathcal{H}$ through

This operator is the covariance operator (defined on the Hilbert space $\mathcal{H}$). Using the reproducing property, we have:

where for any elements $g,h\in\mathcal{H}$, we denote by $g\otimes h$ the operator from $\mathcal{H}$ to $\mathcal{H}$ defined as:

Note that this expectation is formally defined as a Bochner expectation (an extension of Lebesgue integration theory to Banach spaces, see ) in $\mathcal{L}(\mathcal{H})$ the set of endomorphisms of $\mathcal{H}$.

We have defined the covariance operator on the Hilbert space $\mathcal{H}$. If $f\in\mathcal{H}$, we have for all $z\in\mathcal{X}$, using the reproducing property:

which shows that the operator $\Sigma$ may be extended to any square-integrable function $f\in{L}^{2}_{\rho_{{X}}}$. In the following, we extend such an operator as an endomorphism $T$ from ${L}^{2}_{\rho_{{X}}}$ to ${L}^{2}_{\rho_{{X}}}$.

Assume (A1-2). We define the operator $T$ as follows:

From the discussion above, if $f\in\mathcal{H}\subset{L}^{2}_{\rho_{{X}}}$, then $Tf=\Sigma f$. We give here some of the most important properties of $T$. The operator $T$ (which is an endomorphism of the separable Hilbert space ${L}^{2}_{\rho_{{X}}}$) may be reduced in some Hilbertian eigenbasis of ${L}^{2}_{\rho_{{X}}}$. It allows us to define the power of such an operator $T^{r}$, which will be used to quantify the regularity of the function $g_{\mathcal{H}}$. See proof in Appendix I.2, Proposition 19.

Assume (A1-2). $T$ is a bounded self-adjoint semi-definite positive operator on ${L}^{2}_{\rho_{{X}}}$, which is trace-class. There exists a Hilbertian eigenbasis $(\phi_{i})_{i\in I}$ of the orthogonal supplement $S$ of the null space ${\rm Ker}(T)$, with summable strictly positive eigenvalues $(\mu_{i})_{i\in I}$. That is:

$\forall i\in I,\ T\phi_{i}=\mu_{i}\phi_{i}$, $(\mu_{i})_{i\in I}$ strictly positive such that $\sum_{i\in I}\mu_{i}<\infty$.

${L}^{2}_{\rho_{{X}}}=\operatorname{Ker}(T)\overset{\perp}{\oplus}S$, that is, ${L}^{2}_{\rho_{{X}}}$ is the orthogonal direct sum of $\operatorname{Ker}(T)$ and $S$.

When the space $S$ has finite dimension, then $I$ has finite cardinality, while in general $I$ is countable. Moreover, the null space $\operatorname{Ker}(T)$ may be either reduced to $\{0\}$ (this is the more classical setting and such an assumption is often made), finite-dimensional (for example when the kernel has zero mean, thus constant functions are in $S$) or infinite-dimensional (e.g., when the kernel space only consists in even functions, the whole space of odd functions is in $S$).

Assume (A1-2). $\Sigma:\mathcal{H}\to\mathcal{H}$ is injective. The image of $T$ is included in $\mathcal{H}$: $\text{Im}(T)\subset\mathcal{H}$. Moreover, for any $i\in I$, $\phi_{i}=\frac{1}{\mu_{i}}T\phi_{i}\in\mathcal{H}$ , thus $\left(\mu_{i}^{1/2}\phi_{i}\right)_{i\in I}$ is an orthonormal eigen-system of $\Sigma$ and an Hilbertian basis of $\mathcal{H}$, i.e., for any $i$ in $I,\ \Sigma\phi_{i}=\mu_{i}\phi_{i}$.

This proposition will be generalized under relaxed assumptions (in particular as $\Sigma$ will no more be injective, see Section 2.5 and Appendix A).

We may now define all powers $T^{r}$ (they are always well defined because the sequence of eigenvalues is upper-bounded):

We define, for any $r\geqslant 0$, $T^{r}:{L}^{2}_{\rho_{{X}}}\rightarrow{L}^{2}_{\rho_{{X}}}$, for any $h\in\operatorname{Ker}(T)$ and $(a_{i})_{i\in I}$ such that $\sum_{i\in I}a_{i}^{2}<\infty$, through: $T^{r}\left(h+\sum_{i\in I}a_{i}\phi_{i}\right)=\sum_{i\in I}a_{i}\mu_{i}^{r}\phi_{i}.$ Moreover, for any $r>0$, $T^{r}$ may be defined as a bijection from $S$ into $\text{Im}(T^{r})$. We may thus define its unique inverse $T^{-r}:\ \ \text{Im}(T^{r})\rightarrow S.$

The following proposition is a consequence of Mercer’s theorem . It describes how the space $\mathcal{H}$ is related to the image of operator $T^{1/2}$.

Under assumptions (A1,2), $\mathcal{H}=T^{1/2}\left({L}^{2}_{\rho_{{X}}}\right)$ and $T^{1/2}:S\rightarrow\mathcal{H}$ is an isometrical isomorphism.

The proposition has the following consequences:

For any $r\geqslant 1/2,\ T^{r}(S)\subset\mathcal{H}$, because $T^{r}(S)\subset T^{1/2}(S)$, that is, with large enough powers $r$, the image of $T^{r}$ is in the Hilbert space.

$\forall r>0,\ \overline{T^{r}({L}^{2}_{\rho_{{X}}})}=S=\overline{T^{1/2}({L}^{2}_{\rho_{{X}}})}=\overline{\mathcal{H}}$, because (a) $T^{1/2}({L}^{2}_{\rho_{{X}}})={\mathcal{H}}$ and (b) for any $r>0$, $\overline{T^{r}({L}^{2}_{\rho_{{X}}})}=S$. In other words, elements of $\overline{\mathcal{H}}$ (on which our minimization problem attains its minimum), may seen as limits (in ${L}^{2}_{\rho_{{X}}}$) of elements of $T^{r}({L}^{2}_{\rho_{{X}}})$, for any $r>0$.

$\mathcal{H}$ is dense in ${L}^{2}_{\rho_{{X}}}$ if and only if $T$ is injective (which is equivalent to $\ker(T)=\{0\}$)

The sequence of spaces $\{T^{r}({L}^{2}_{\rho_{{X}}})\}_{r>0}$ is thus a decreasing (when $r$ is increasing) sequence of subspaces of ${L}^{2}_{\rho_{{X}}}$ such that any of them is dense in $\overline{\mathcal{H}}$, and $T^{r}({L}^{2}_{\rho_{{X}}})\subset\mathcal{H}$ if and only if $r\geqslant 1/2.$

In the following, the regularity of the function $g_{\mathcal{H}}$ will be characterized by the fact that $g_{\mathcal{H}}$ belongs to the space $T^{r}({L}^{2}_{\rho_{{X}}})$ (and not only to its closure), for a specific $r>0$ (see Section 2.7). This space may be described depending on the eigenvalues and eigenvectors as

We may thus see the spaces $T^{r}({L}^{2}_{\rho_{{X}}})$ as spaces of sequences with various decay conditions.

5 Minimal assumptions

$\mathcal{H}$ is a separable RKHS associated with kernel $K$ on the set $\mathcal{X}$.

In this section, we have to distinguish the set of square $\rho_{X}$-integrable functions $\mathcal{L}^{2}_{\rho_{{X}}}$ and its quotient ${L}^{2}_{\rho_{{X}}}$ that makes it a separable Hilbert space. We define $p$ the projection from $\mathcal{L}^{2}_{\rho_{{X}}}$ into ${L}^{2}_{\rho_{{X}}}$ (precise definitions are given in Appendix A). Indeed it is no more possible to identify the space $\mathcal{H}$, which is a subset of $\mathcal{L}^{2}_{\rho_{{X}}}$, and its canonical projection $p(\mathcal{H})$ in ${L}^{2}_{\rho_{{X}}}$.

Minimality: The separability assumption is necessary to be able to expand any element as an infinite sum, using a countable orthonormal family (this assumption is satisfied in almost all cases, for instance it is simple as soon as $\mathcal{X}$ admits a topology for which it is separable and functions in $\mathcal{H}$ are continuous, see for more details). Note that we do not make any topological assumptions regarding the set $\mathcal{X}$. We only assume that it is equipped with a probability measure.

Our assumptions are sufficient to analyze the minimization of $\varepsilon(f)$ with respect to $f\in\mathcal{H}$ and seem to allow the widest generality.

Comparison: These assumptions will include the previous setting, but also recover measures without full support (e.g., when the data lives in a small subspace of the whole space) and kernels on discrete objects (with non-finite cardinality).

Moreover, (A1’), (A2’) are stricly weeker than (A1), (A2). In previous work, (A2’) was sometimes replaced by the stronger assumptions $\sup_{x\in\mathcal{X}}K(x,x)<\infty$ and $|Y|$ bounded . Note that in functional analysis, the weaker hypothesis $\int_{\mathcal{X}\!\times\!\mathcal{X}}k(x,x^{\prime})^{2}d{\rho_{X}}(x)d\rho_{X}(x^{\prime})<\infty$ is often used , but it is not adapted to the statistical setting.

Main differences: The main difference here is that we cannot identify $\mathcal{H}$ and $p(\mathcal{H})$: there may exist functions $f\in\mathcal{H}\setminus\{0\}$ such that $\|f\|_{\mathcal{L}^{2}_{\rho_{{X}}}}=0$. This may for example occur if the support of $\rho_{X}$ is strictly included in $\mathcal{X}$, and $f$ is zero on this support, but not identically zero. See the Appendix I.5 for more details.

As a consequence, $\Sigma$ is no more injective and we do not have $\operatorname{Im}(T^{1/2})=\mathcal{H}$ any more. We thus denote $\mathscr{S}$ an orthogonal supplement of the null space $\text{Ker}(\Sigma)$. As we also need to be careful not to confuse $\mathcal{L}^{2}_{\rho_{{X}}}$ and ${L}^{2}_{\rho_{{X}}}$, we define an extension $\mathcal{T}$ of $\Sigma$ from $\mathcal{L}^{2}_{\rho_{{X}}}$ into $\mathcal{H}$, then $T=p\circ\mathcal{T}$. We can define for $r\geqslant 1/2$ the power operator $\mathcal{T}^{r}$ of $\mathcal{T}$ (from ${L}^{2}_{\rho_{{X}}}$ into $\mathcal{H}$), see App. A for details.

Conclusion: Our problem has the same behaviour under such assumptions. Proposition 1 remains unchanged. Decompositions in Prop. 2 and Corollary 1 must be slightly adapted (see Proposition 9 and Corollary 7 in Appendix A for details). Finally, Proposition 3 is generalized by the next proposition, which states that $p(\mathscr{S})=p(\mathcal{H})$ and thus $S$ and $p(\mathcal{H})$ are isomorphic (see proof in Appendix I.2, Proposition 19):

$\mathcal{T}^{1/2}:S\rightarrow\mathscr{S}$ is an isometry. Moreover, $\operatorname{Im}(T^{1/2})=p(\mathcal{H})$ and $T^{1/2}:S\rightarrow p(\mathcal{H})$ is an isomorphism.

We can also derive a version of Mercer’s theorem, which does not make any more assumptions that are required for defining RKHSs. As we will not use it in this article, this proposition is only given in Appendix A.

Convergence results: In all convergence results stated below, assumptions (A1, A2) may be replaced by assumptions (A1’, A2’).

6 Examples

The property $\overline{\mathcal{H}}=S$, stated after Proposition 3, is important to understand what the space $\overline{\mathcal{H}}$ is, as we are minimizing over this closed and convex set. As a consequence the space $\mathcal{H}$ is dense in ${L}^{2}_{\rho_{{X}}}$ if and only if $T$ is injective (or equivalently, ${\rm Ker}(T)=\{0\}\Leftrightarrow S={L}^{2}_{\rho_{{X}}}$). We detail below a few classical situations in which different configurations for the “inclusion” $\mathcal{H}\subset\overline{\mathcal{H}}\subset{L}^{2}_{\rho_{{X}}}$ appear:

Splines over the circle: When $X\sim\mathcal{U}[0;1]$ and $\mathcal{H}$ is the set of $m$-times periodic weakly differentiable functions (see Section 5), we have in general $\mathcal{H}\varsubsetneq\overline{\mathcal{H}}\varsubsetneq{L}^{2}_{\rho_{{X}}}$. In such a case, $\ker(T)=\operatorname{span}(x\mapsto 1)$, and $\overline{\mathcal{H}}\oplus\operatorname{span}(x\mapsto 1)={L}^{2}_{\rho_{{X}}}$, that is we can approximate any zero-mean function.

Many examples and more details may be found in . In particular, kernels on non-vectorial objects may be defined (e.g., sequences, graphs or measures).

7 Convergence rates

$g_{\mathcal{H}}\in T^{r}\left({L}^{2}_{\rho_{{X}}}\right)$ with $r\geqslant 0$, and as a consequence $\|T^{-r}(g_{\mathcal{H}})\|_{{L}^{2}_{\rho_{{X}}}}<\infty$.

We chose such assumptions in order to make the comparison with the existing literature as easy as possible, for example . However, some other assumptions may be found as in .

The two parameters $r$ and $\alpha$ intuitively parametrize the strengths of our assumptions:

In assumption (A3) a bigger $\alpha$ makes the assumption stronger: it means the reproducing kernel Hilbert space is smaller, that is if (A3) holds with some constant $\alpha$, then it also holds for any $\alpha^{\prime}<\alpha$. Moreover, if $T$ is reduced in the Hilbertian basis $(\phi_{i})_{i}$ of ${L}^{2}_{\rho_{{X}}}$, we have an effective search space S=\big{\{}\sum_{i=1}^{\infty}b_{i}\phi_{i}/\sum_{i=1}^{\infty}\frac{b_{i}^{2}}{\mu_{i}}<\infty\big{\}} : the smaller the eigenvalues, the smaller the space. Note that since $\operatorname{tr}(T)$ is finite, (A3) is always true for $\alpha=1$.

In assumption (A4), for a fixed $\alpha$, a bigger $r$ makes the assumption stronger, that is the function $g_{\mathcal{H}}$ is actually smoother. Indeed, considering that (A4) may be rewritten g_{\mathcal{H}}\in T^{r}\big{(}{L}^{2}_{\rho_{{X}}}\big{)} and for any r<r^{\prime},\ T^{r^{\prime}}\big{(}{L}^{2}_{\rho_{{X}}}\big{)}\subset T^{r}\big{(}{L}^{2}_{\rho_{{X}}}\big{)}. In other words, \big{\{}T^{r}\left({L}^{2}_{\rho_{{X}}}\right)\big{\}}_{r\geq 0} are decreasing ($r$ growing) subspaces of ${L}^{2}_{\rho_{{X}}}$.

For $r=1/2$, T^{1/2}\big{(}{L}^{2}_{\rho_{{X}}}\big{)}=\mathcal{H}; moreover, for $r\geqslant 1/2$, our best approximation function $g_{\mathcal{H}}\in\overline{\mathcal{H}}$ is in fact in $\mathcal{H}$, that is the optimization problem in the RKHS $\mathcal{H}$ is attained by a function of finite norm. However for $r<1/2$ it is not attained.

Furthermore, it is worth pointing the stronger assumption which is often used in the finite dimensional context, namely $\operatorname{tr}\left(\Sigma^{1/\alpha}\right)=\sum_{i\in I}\mu_{i}^{1/\alpha}$ finite. It turns out that this is a stronger assumption, indeed, since we have assumed that the eigenvalues $(\mu_{i})$ are arranged in non-increasing order, if $\operatorname{tr}\left(\Sigma^{1/\alpha}\right)$ is finite, then (A3) is satisfied for s^{2}=\big{[}2\operatorname{tr}\left(\Sigma^{1/\alpha}\right)\big{]}^{\alpha}. Such an assumption appears for example in Corollary 5.

The assumptions (A3) and (A4) are adapted to our theoretical results, but some stricter assumptions are often used, that make comparison with existing work more direct. For comparison purposes, we will also use:

Assumption (a3) directly imposes that the eigenvalues of $T$ decay at rate $i^{-\alpha}$ (which imposes that there are infinitely many), and thus implies (A3). Together, assumptions (a3) and (a4), imply assumptions (A3) and (A4), with any $\delta>1+2\alpha r$. Indeed, we have

which is finite for $2\alpha r-\delta<-1$. Thus, the supremum element of the set of $r$ such that (A4) holds is such that $\delta=1+2\alpha r$. Thus, when comparing assumptions (A3-4) and (a3-4), we will often make the identification above, that is, $\delta=1+2\alpha r$.

The main advantage of the new assumptions is their interpretation when the basis $(\phi_{i})_{i\in I}$ is common for several RKHSs (such as the Fourier basis for splines, see Section 5): (a4) describes the decrease of the coordinates of the best function $g_{\mathcal{H}}\in{L}^{2}_{\rho_{{X}}}$ independently of the chosen RKHS. Thus, the parameter $\delta$ characterizes the prediction function, while the parameter $\alpha$ characterizes the RKHS.

Stochastic approximation in Hilbert spaces

In this section, we consider estimating a prediction function $g\in\mathcal{H}$ from observed data, and we make the following assumption:

Given $n$ observations, regularized empirical risk minimization corresponds to minimizing with respect to $g\in\mathcal{H}$ the following objective function:

In terms of convergence rates, assumptions (a3-4) allow to obtain convergence rates that decompose $\varepsilon(\hat{g})-\varepsilon(g_{\mathcal{H}})=\|\hat{g}-g_{\mathcal{H}}\|_{{L}^{2}_{\rho_{{X}}}}^{2}$ as the sum of two asymptotic terms :

Variance term: O\big{(}{\sigma^{2}}{n^{-1}\lambda^{-1/\alpha}}\big{)}, which is decreasing with $\lambda$, where $\sigma^{2}$ characterizes the noise variance, for example, in the homoscedastic case (i.i.d. additive noise), the marginal variance of the noise; see assumption (A6) for the detailed assumption that we need in our stochastic approximation context.

Bias term: O\big{(}\lambda^{\min\{(\delta-1)/\alpha,2\}}\big{)}, which is increasing with $\lambda$. Note that the corresponding $r$ from assumptions (A3-4) is $r=(\delta-1)/2\alpha$, and the bias term becomes proportional to $\lambda^{\min\{2r,2\}}$.

Optimal predictions: If $r<1$, then the optimal value of $\lambda$ (that minimizes the sum of two terms and makes them asymptotically equivalent) is proportional to $n^{-\alpha/(2r\alpha+1)}=n^{-\alpha/\delta}$ and the excess prediction error \|\hat{g}-g_{\mathcal{H}}\|_{{L}^{2}_{\rho_{{X}}}}^{2}=O\big{(}n^{-2\alpha r/(2\alpha r+1)}\big{)}=O\big{(}n^{-1+1/\delta}\big{)}, and the resulting procedure is then “optimal” in terms of estimation of $g_{\mathcal{H}}$ in ${L}^{2}_{\rho_{{X}}}$ (see Section 4 for details).

Saturation: If $r\geqslant 1$, where the optimal value of $\lambda$ (that minimizes the sum of two terms and makes them equivalent) is proportional to $n^{-\alpha/(2\alpha+1)}$, and the excess prediction error is less than O\big{(}n^{-2\alpha/(2\alpha+1)}\big{)}, which is suboptimal. Although assumption (A4) is valid for a larger $r$, the rate is the same than if $r=1$.

In this paper, we consider a stochastic approximation framework with improved running-time complexity and similar theoretical behavior than regularized empirical risk minimization, with the advantage of (a) needing a single pass through the data and (b) simple assumptions.

2 Stochastic approximation

We may also apply the recursion using representants. Indeed, if $g_{0}=0$, which we now assume, then for any $n\geqslant 1$,

with the following recursion on the sequence $(a_{n})_{n\geqslant 1}$:

We also output the averaged iterate defined as

The running time complexity is $O(i)$ for iteration $i$—if we assume that kernel evaluations are $O(1)$, and thus $O(n^{2})$ after $n$ steps. This is a serious limitation for practical applications. Several authors have considered expanding $g_{n}$ on a subset of all $(K_{x_{i}})$, which allows to bring down the complexity of each iteration and obtain an overall linear complexity is $n$ , but this comes at the expense of not obtaining the sharp generalization errors that we obtain in this paper. Note that when studying regularized least-squares problem (i.e., adding a penalisation term), one has to update every coefficient $(a_{i})_{1\leqslant i\leqslant n}$ at step $n$, while in our situation, only $a_{n}$ is computed at step $n$.

Similar algorithms have been studied before , under various forms. Especially, in a regularization term is added to the loss function (thus considering the following problem: $\arg\min_{f\in\mathcal{H}}\varepsilon(f)+\lambda||f||_{K}^{2}$). In , neither regularization nor averaging procedure are considered, but in the second case, multiple pass through the data are considered. In , a non-regularized averaged procedure equivalent to ours is considered. However, the step-sizes $\gamma_{n}$ which are proposed, as well as the corresponding analysis, are different. Our step-sizes are larger and our analysis uses more directly the underlying linear algebra to obtain better rates (while the proof of is applicable to all smooth losses).

We are mainly interested in two different types of step-sizes (a.k.a. learning rates): the sequence $(\gamma_{i})_{1\leqslant i\leqslant n}$ may be either:

a sequence of the type $\gamma_{i}=\Gamma(n)$ for $i\leqslant n$, which will be referred to as the “finite horizon setting”: in this situation the number of samples is assumed to be known and fixed and we chose a constant step-size which may depend on this number. Our bound then hold only for the last iterate.

In practice it is important to have an online procedure, to be able to deal with huge amounts of data (potentially infinite). However, the analysis is easier in the “finite horizon” setting. Some doubling tricks allow to pass to varying steps , but it is not fully satisfactory in practice as it creates jumps at every $n$ which is a power of two.

3 Extra regularity assumptions

We first present the results in the finite horizon setting in Section 3.4 before turning to the online setting in Section 3.5.

4 Main results (finite horizon)

We can first get some guarantee on the consistency of our estimator, for any small enough constant step-size:

Assume (A1-6), then for any constant choice $\gamma_{n}=\gamma_{0}<\frac{1}{2R^{2}}$, the prediction error of $\bar{g}_{n}$ converges to the one of $g_{\mathcal{H}}$, that is:

The expectation is considered with respect to the distribution of the sample $(x_{i},y_{i})_{1\leqslant i\leqslant n}$, as in all the following theorems (note that $\|\bar{g}_{n}-g_{\mathcal{H}}\|_{{L}^{2}_{\rho_{{X}}}}^{2}$ is itself a different expectation with respect to the law $\rho_{X}$).

Theorem 1 means that for the simplest choice of the learning rate as a constant, our estimator tends to the perform as well as the best estimator in the class $\mathcal{H}$. Note that in general, the convergence in $\mathcal{H}$ is meaningless if $r<1/2$. The following results will state some assertions on the speed of such a convergence; our main result, in terms of generality is the following:

Assume (A1-6) and $\gamma_{i}=\gamma=\Gamma(n)$, for $1\leqslant i\leqslant n$. If $\gamma R^{2}\leqslant 1/4$:

Where $q_{n,\gamma,s,r}:=(R^{2\alpha}\gamma^{1+\alpha}ns^{2})^{\frac{2r-1}{\alpha}}$ if $r\geqslant\frac{1}{2}$ and $q_{n,\gamma,s,r}:=0$ otherwise is a residual quantity.

Proof: Theorem 1 is directly derived from Theorem 2, which is proved in Appendix II.3: we derive for our algorithm a new error decomposition and bound the different sources of error via algebraic calculations. More precisely, following the proof in Euclidean space , we first analyze (in Appendix II.2) a closely related recursion (we replace $K_{x_{n}}\otimes K_{x_{n}}$ by its expectation $\Sigma$, and we thus refer to it as a semi-stochastic version of our algorithm):

It (a) leads to an easy computation of the main bias/variance terms of our result, (b) will be used to derive our main result by bounding the drifts between our algorithm and its semi-stochastic version. A more detailed sketch of the proof is given in Appendix B.

Bias/variance interpretation: The two main terms have a simple interpretation. The first one is a variance term, which shows the effect of the noise $\sigma^{2}$ on the error. It is bigger when $\sigma$ gets bigger, and moreover it also gets bigger when $\gamma$ is growing (bigger steps mean more variance). As for the second term, it is a bias term, which accounts for the distance of the initial choice (the null function in general) to the objective function. As a consequence, it is smaller when we make bigger steps.

Assumption (A4): Our assumption (A4) for $r>1$ is stronger than for $r=1$ but we do not improve the bound. Indeed the bias term (see comments below) cannot decrease faster than $O(n^{-2})$: this phenomenon in known as saturation . To improve our results with $r>1$ it may be interesting to consider another type of averaging. In the following, $r<1$ shall be considered as the main and most interesting case.

Relationship to regularized empirical risk minimization: Our bound ends up being very similar to bounds for regularized empirical risk minimization, with the identification $\lambda=\frac{1}{\gamma n}$. It is thus no surprise that once we optimize for the value of $\gamma$, we recover the same rates of convergence. Note that in order to obtain convergence, we require that the step-size $\gamma$ is bounded, which corresponds to an equivalent $\lambda$ which has to be lower-bounded by $1/n$.

Finite horizon: Once again, this theorem holds in the finite horizon setting. That is we first choose the number of samples we are going to use, then the learning rate as a constant. It allows us to chose $\gamma$ as a function of $n$, in order to balance the main terms in the error bound. The trade-off must be understood as follows: a bigger $\gamma$ increases the effect of the noise, but a smaller one makes it harder to forget the initial condition.

We may now deduce the following corollaries, with specific optimized values of $\gamma$:

Assume (A1-6) and a constant step-size $\gamma_{i}=\gamma=\Gamma(n)$, for $1\leqslant i\leqslant n$:

If $\frac{\alpha-1}{2\alpha}<r$ and $\Gamma(n)=\gamma_{0}\ n^{\frac{-2\alpha\min\{r,1\}-1+\alpha}{2\alpha\min\{r,1\}+1}}$, $\gamma_{0}R^{2}\leqslant 1/4$, we have:

with $A=4\left(1+{(\gamma_{0}s^{2})^{\frac{1}{\alpha}}}\right)\sigma^{2}+\frac{4(1+o(1))}{\gamma_{0}^{2r}}||L_{K}^{-r}g_{\mathcal{H}}||_{{L}^{2}_{\rho_{{X}}}}^{2}$.

If $0<r<\frac{\alpha-1}{2\alpha}$, with $\Gamma(n)=\gamma_{0}$ is constant, $\gamma_{0}R^{2}\leqslant 1/4$, we have:

Limit conditions: Assumption (A4), gives us some kind of “position” of the objective function with respect to our reproducing kernel Hilbert space. If $r\geqslant 1/2$ then $g_{\mathcal{H}}\in\mathcal{H}$. That means the regression function truly lies in the space in which we are looking for an approximation. However, it is not necessary neither to get the convergence result, which stands for any $r>0$, nor to get the optimal rate (see definition in Section 4.2), which is also true for $\frac{\alpha-1}{2\alpha}<r<1$ .

Evolution with $r$ and $\alpha$: As it has been noticed above, a bigger $\alpha$ or $r$ would be a stronger assumption. It is thus natural to get a rate which improves with a bigger $\alpha$ or $r$: the function $(\alpha,r)\mapsto\frac{2\alpha r}{2\alpha r+1}$ is increasing in both parameters.

The quantity $o(1)$ in Equation (3.3) stands for $(\gamma_{0}s^{2}n^{-2\alpha^{2}r+1})^{\frac{2r-1}{\alpha}}$ if $r\geqslant 1/2$ (0 otherwise) and is a quantity which decays to 0.

Different regions: in Figure 1(a), we plot in the plan of coordinates $\alpha,\delta$ (with $\delta=2\alpha r+1$) our limit conditions concerning our assumptions, that is, $r=1\Leftrightarrow\delta=2\alpha+1$ and $\frac{\alpha-1}{2\alpha}=r\Leftrightarrow\alpha=\delta$. The region between the two green lines is the region for which the optimal rate of estimation is reached. The magenta dashed lines stands for $r=1/2$, which has appeared to be meaningless in our context.

The region $\alpha\geqslant\delta\Leftrightarrow\frac{\alpha-1}{2\alpha}>r$ corresponds to a situation where regularized empirical risk minimization would still be optimal, but with a regularization parameter $\lambda$ that decays faster than $1/n$, and thus, our corresponding step-size $\gamma=1/(n\lambda)$ would not be bounded as a function of $n$. We thus saturate our step-size to a constant and the generalization error is dominated by the bias term.

The region $\alpha\leqslant(\delta-1)/2\Leftrightarrow r>1$ corresponds to a situation where regularized empirical risk minimization reaches a saturating behaviour. In our stochastic approximation context, the variance term dominates.

5 Online setting

We now consider the second case when the sequence of step-sizes does not depend on the number of samples we want to use (online setting).

The computation are more tedious in such a situation so that we will only state asymptotic theorems in order to understand the similarities and differences between the finite horizon setting and the online setting, especially in terms of limit conditions.

Assume (A1-6), assume for any $i$, $\gamma_{i}=\frac{\gamma_{0}}{i^{\zeta}}$, $\gamma_{0}R^{2}\leqslant 1/2$ :

If $0<r(1-\zeta)<1$, if $0<\zeta<\frac{1}{2}$ then

The constant in the $O(\cdot)$ notations only depend on $\gamma_{0}$ and $\alpha$.

Theorem 3 is proved in Appendix II.4. In the first case, the main bias and variance terms are the same as in the finite horizon setting, and so is the optimal choice of $\zeta$. However in the second case, the variance term behaviour changes: it does not decrease any more when $\zeta$ increases beyond $1/2$. Indeed, in such a case our constant averaging procedure puts to much weight on the first iterates, thus we do not improve the variance bound by making the learning rate decrease faster. Other type of averaging, as proposed for example in , could help to improve the bound.

Moreover, this extra condition thus changes a bit the regions where we get the optimal rate (see Figure 1(b)), and we have the following corollary:

Assume (A1-6) (in this corollary, $O(\cdot)$ stands for a constant depending on $\alpha,||L_{K}^{-r}g_{\mathcal{H}}||_{{L}^{2}_{\rho_{{X}}}},s,\sigma^{2},\gamma_{0}$ and universal constants):

If $\frac{\alpha-1}{2\alpha}<r<\frac{2\alpha-1}{2\alpha}$, with $\gamma_{n}=\gamma_{0}n^{\frac{-2\alpha r-1+\alpha}{2\alpha r+1}}$ for any $n\geqslant 1$ we get the rate:

If $\frac{2\alpha-1}{2\alpha}<r$, with $\gamma_{n}=\gamma_{0}n^{-1/2}$ for any $n\geqslant 1$, we get the rate:

If $0<r<\frac{\alpha-1}{2\alpha}$, with $\gamma_{n}=\gamma_{0}$ for any $n\geqslant 1$, we get the rate given in (3.4). Indeed the choice of a constant learning rate naturally results in an online procedure.

This corollary is directly derived from Theorem 3, balancing the two main terms. The only difference with the finite horizon setting is the shrinkage of the optimality region as the condition $r<1$ is replaced by $r<\frac{2\alpha-1}{2\alpha}<1$ (see Figure 1(b)). In the next section, we relate our results to existing work.

Links with existing results

In this section, we relate our results from the previous section to existing results.

Recently, Bach and Moulines showed in that for least squares regression, averaged stochastic gradient descent achieved a rate of $O(1/n)$, in a finite-dimensional Hilbert space (Euclidean space), under the same assumptions as above (except the first one of course), which is replaced by:

$\mathcal{H}$ is a $d$-dimensional Euclidean space.

Assume (A1-f), (A2-6). Then for $\gamma=\frac{1}{4R^{2}}$,

We show that we can deduce such a result from Theorem 2 (and even with comparable constants). Indeed under (A1-f) we have:

As we are in a finite-dimensional space (A4) is true for $r=1/2$ as $||T^{-1/2}g_{\mathcal{H}}||^{2}_{\mathcal{L}^{2}_{\rho_{{X}}}}=||g_{\mathcal{H}}||^{2}_{\mathcal{H}}$.

Under such remarks, the following corollary may be deduced from Theorem 2:

Assume (A1-f), (A2-6), then for any $\alpha>1$, with $\gamma R^{2}\leqslant 1/4$:

This bound is easily comparable to (4.1) and shows that our more general analysis has not lost too much. Moreover our learning rate is proportional to $n^{\frac{-1}{2\alpha+1}}$ with $r=1/2$, so tends to behave like a constant when $\alpha\rightarrow\infty$, which recovers the constant step set-up from .

Moreover, N. Flammarion proved (Personnal communication, 05/2014), using the same tpye of techniques, that their bound could be extended to:

a result that may be deduced of the following more general corollaries of our Theorem 2:

Assume (A1-f), (A2-6), and, for some $q\geqslant-1/2$, $||\Sigma^{-q}{g}_{\mathcal{H}}||^{2}_{\mathcal{H}}=||\Sigma^{-(q+1/2)}{g}_{\mathcal{H}}||^{2}_{{L}^{2}_{\rho_{{X}}}}<\infty$, then:

Such a result is derived from Theorem 2 and with the stronger assumption $\operatorname{tr}(\Sigma^{1/\alpha})<\infty$ clearly satisfied in finite dimension, and with $r=q+1/2$. Note that the result above is true for all values of $\alpha\geqslant 1$ and all $q\geqslant-1/2$ (for the ones with infinite $||\Sigma^{-(q+1/2)}{g}_{\mathcal{H}}||^{2}_{{L}^{2}_{\rho_{{X}}}}$, the statement is trivial). This shows that we may take the infimum over all possible $\alpha\leqslant 1$ and $q\geqslant 0$, showing adaptivity of the estimator to the spectral decay of $\Sigma$ and the smoothness of the optimal prediction function $g_{\mathcal{H}}$.

Thus with $\alpha\rightarrow\infty$, we obtain :

Assume (A1-f), (A2-6), and, for some $q\geqslant-1/2$, $||\Sigma^{-q}{g}_{\mathcal{H}}||^{2}_{\mathcal{H}}=||\Sigma^{-(q+1/2)}{g}_{\mathcal{H}}||^{2}_{{L}^{2}_{\rho_{{X}}}}<\infty$, then:

This final result bridges the gap between Proposition 5 ($q=0$), and its extension (4.2) ($q=1/2$). The constants 16 and 8 come from the upper bounds $(a+b)^{2}\leqslant a^{2}+b^{2}$ and $1+1/\sqrt{d}\leqslant 2$ and are thus non optimal.

Moreover, we can also derive from Corollary 5, with $\alpha=1$, $q=0$, and $\gamma\varpropto n^{-1/2}$, we recover the rate $O(n^{-1/2})$ (where the constant does not depend on the dimension $d$ of the Euclidean space). Such a rate was described, e.g., in .

Note that linking our work to the finite-dimensional setting is made using the fact that our assumption (A3) is true for any $\alpha>1$.

2 Optimal rates of estimation

In some situations, our stochastic approximation framework leads to “optimal” rates of prediction in the following sense. In [9, Theorem 2] a minimax lower bound was given: let $\mathcal{P}(\alpha,r)\ \ (\alpha>1,r\in[1/2,1])$ be the set of all probability measures $\rho$ on $\mathcal{X}\times\mathcal{Y}$, such that:

Then the following minimax lower rate stands:

for some constant $C>0$ where the infimum in the middle is taken over all algorithms as a map $((x_{i},y_{i})_{1\leqslant i\leqslant n})\mapsto g_{n}\in\mathcal{H}$.

When making assumptions (a3-4), the assumptions regarding the prediction problem (i.e., the optimal function $g_{\rho}$) are summarized in the decay of the components of $g_{\rho}$ in an orthonormal basis, characterized by the constant $\delta$. Here, the minimax rate of estimation (see, e.g., ) is $O(n^{-1+1/\delta})$ which is the same as O\big{(}n^{-2r\alpha/(2r\alpha+1)}\big{)} with the identification $\delta=2\alpha r+1$.

That means the rate we get is optimal for $\frac{\alpha-1}{2\alpha}<r<1$ in the finite horizon setting, and for $\frac{\alpha-1}{2\alpha}<r<\frac{2\alpha-1}{2\alpha}$ in the online setting. This is the region between the two green lines on Figure 1.

3 Regularized stochastic approximation

It is interesting to link our results to what has been done in and in the case of regularized least-mean-squares, so that the recursion is written:

Where $O_{\kappa}$ stands for a constant which depends on $\kappa$.

No assumption is made on the covariance operator beyond being trace class, but only on $\|T^{-r}g_{\rho}\|_{{L}^{2}_{\rho_{{X}}}}$ (thus no assumption (A3)). A few remarks may be made:

They get almost-sure convergence, when we only get convergence in expectation. We could perhaps derive a.s. convergence by considering moment bounds in order to be able to derive convergence in high probability and to use Borel-Cantelli lemma.

They only assume $\frac{1}{2}\leqslant r\leqslant 1$, which means that they assume the regression function to lie in the RKHS.

4 Unregularized stochastic approximation

In , Ying and Pontil studied the same unregularized problem as we consider, under assumption (A4). They obtain the same rates as above ($n^{-2r/(2r+1)}\log(n)$) in both online case (with $0\leqslant r\leqslant\frac{1}{2}$) and finite horizon setting ($0<r$).

They led as an open problem to improve bounds with some additional information on some decay of the eigenvalues of $T$, a question which is answered here.

Moreover, Zhang also studies stochastic gradient descent algorithms in an unregularized setting, also with averaging. As described in , his result is stated in the linear kernel setting but may be extended to kernels satisfying $\sup_{x\in\mathcal{X}}K(x,x)\leqslant R^{2}$. Ying and Pontil derive from Theorem 5.2 in the following proposition:

Assume we consider the algorithm defined in Section 3.2 and output $\overline{g}_{n}$ defined by equation (3.1). Assume the kernel $K$ satisfies $\sup_{x\in\mathcal{X}}K(x,x)\leqslant R^{2}$. Finally assume $g_{\rho}$ satisfies assumption (A4) with $0<r<1/2$. Then in the finite horizon setting, with $\Gamma(n)=\frac{1}{4R^{2}}n^{-\frac{2r}{2r+1}}$, we have:

Moreover, note that we may derive their result from Corollary 2. Indeed, using $\Gamma(n)=\gamma_{0}n^{\frac{-2r}{2r+1}}$, we get a bias term which is of order $n^{\frac{-2r}{2r+1}}$ and a variance term of order $n^{-1+\frac{1}{2r\alpha+\alpha}}$ which is smaller. Our analysis thus recovers their convergence rate with their step-size. Note that this step-size is significantly smaller than ours, and that the resulting bound is worse (but their result holds in more general settings than least-squares). See more details in Section 4.5.

5 Summary of results

All three algorithms are variants of the following:

But they are studied under different settings, concerning regularization, averaging, assumptions: we sum up in Table 1 the settings of each of these studies. For each of them, we consider the finite horizon settings, where results are generally better.

Worst-case result in $r$: in the setting of assumptions (a3,4), given $\delta$, the optimal rate of convergence is known to be $O(n^{-1+1/\delta})$, where $\delta=2\alpha r+1$. We thus get the optimal rate, as soon as $\alpha<\delta<2\alpha+1$, while the other algorithms get the suboptimal rate $n^{\frac{\delta-1}{\delta+\alpha-1}}$ under various conditions. Note that this sub-optimal rate becomes close to the optimal rate when $\alpha$ is close to one, that is, in the worst-case situation. Thus, in the worst-case ($\alpha$ arbitrarily close to one), all methods behave similarly, but for any particular instance where $\alpha>1$, our rates are better.

Choice of kernel: in the setting of assumptions (a3,4), given $\delta$, in order to get the optimal rate, we may choose the kernel (i.e., $\alpha$) such that $\alpha<\delta<2\alpha+1$ (that is neither too big, nor too small), while other methods need to choose a kernel for which $\alpha$ is as close to one as possible, which may not be possible in practice.

Improved bounds: Ying and Pontil only give asymptotic bounds, while we have exact constants for the finite horizon case. Moreover there are some logarithmic terms in which disappear in our analysis.

Saturation: our method does saturate for $r>1$, while the non-averaged framework of does not (but does not depend on the value of $\alpha$). We conjecture that a proper non-uniform averaging scheme (that puts more weight on the latest iterates), we should get the best of both worlds.

Experiments on artificial data

Following , we consider synthetic examples with smoothing splines on the circle, where our assumptions (A3-4) are easily satisfied.

The simplest example to match our assumptions may be found in . We consider $Y=g_{\rho}(X)+\varepsilon$, with $X\sim\mathcal{U}[\,0;1]$ is a uniform random variable in $$, and$g_{\rho}$ in a particular RKHS (which is actually a Sobolev space).

Let $\mathcal{H}$ be the collection of all zero-mean periodic functions on $[0;1]$ of the form

This means that the $m$-th derivative of $f$, $f^{(m)}$ is in $\mathcal{L}^{2}([0\,;1])$. We consider the inner product:

It is known that $\mathcal{H}$ is an RKHS and that the reproducing kernel $R_{m}(s,t)$ for $\mathcal{H}$ is

Moreover the study of Bernoulli polynomials gives a close formula for $R(s,t)$, that is:

with $B_{m}$ denoting the m-th Bernoulli polynomial and $\{s-t\}$ the fractional part of $s-t$ .

We can derive the following proposition for the covariance operator which means that our assumption (A3) is satisfied for our algorithm in $\mathcal{H}$ when $X\sim\mathcal{U}[0;\,1]$, with $\alpha=2m$, and $s=2(1/2\pi)^{m}$.

If $X\sim\mathcal{U}[0;\,1]$, then in $\mathcal{H}$:

the eigenvalues of $\Sigma$ are all of multiplicity 2 and are $\lambda_{i}=(2\pi i)^{-2m}$,

the eigenfunctions are $\phi_{i}^{c}:t\mapsto\sqrt{2}\cos(2\pi it)$ and $\phi_{i}^{s}:t\mapsto\sqrt{2}\sin(2\pi it)$.

For $\phi_{i}^{c}$ we have (a similar argument holds for $\phi_{i}^{s}$):

It is well known that $(\phi_{i}^{c},\phi_{i}^{s})_{i\geqslant 0}$ is an orthonormal system (the Fourier basis) of the functions in ${L}^{2}([\,0;1])$ with zero mean, and it is easy to check that $((2i\pi)^{-m}\phi_{i}^{c},(2i\pi)^{-m}\phi_{i}^{s})_{i\geqslant 1}$ is an orthonormal basis of our RKHS $\mathcal{H}$ (this may also be seen as a consequence of the fact that $T^{1/2}$ is an isometry). ∎

Here the eigenvectors do not depend on the kernel choice, only the re-normalisation constant depends on the choice of the kernel. Especially the eigenbasis of $T$ in ${L}^{2}_{\rho_{{X}}}$ does not depend on $m$. That can be linked with the previous remarks made in Section 4.

Assumption (A3) defines here the size of the RKHS: the smaller $\alpha=2m$ is, the bigger the space is, the harder it is to learn a function.

In the next section, we illustrate on such a toy model our main results and compare our learning algorithm to Ying and Pontil’s , Tarrès and Yao’s and Zhang’s algorithms.

2 Experimental set-up

We use $g_{\rho}(x)=B_{\delta/2}(x)$ with $\delta=2\alpha r+1$, as proposed above, with $B_{1}(x)=x-\frac{1}{2}$, $B_{2}(x)=x^{2}-x+\frac{1}{6}$ and $B_{3}(x)=x^{3}-\frac{3}{2}x^{2}+\frac{1}{2}x$.

We give in Figure 2 the functions used for simulations in a few cases that span our three regions. We also remind the choice of $\gamma$ proposed for the 4 algorithms. We always use the finite horizon setting.

3 Optimal learning rate for our algorithm

4 Comparison to competing algorithms

In this section, we compare the convergence rates of the four algorithms described in Section 4.5. We consider the different choices of $(r,\alpha)$ as described in Table 2 in order to go all over the different optimality situations. The main properties of each algorithm are described in Table 1. However we may note:

For our algorithm, $\Gamma(n)$ is chosen accordingly with Corollary 2, with $\gamma_{0}=\frac{1}{R^{2}}$.

For Ying and Pontil’s algorithm, accordingly to Theorem 6 in , we consider $\Gamma(n)=\gamma_{0}n^{-\frac{2r}{2r+1}}$. We choose $\gamma_{0}=\frac{1}{R^{2}}$ which behaves better than the proposed $\frac{r}{64(1+R^{4})(2r+1)}$.

For Tarrès and Yao’s algorithm, we refer to Theorem C in , and consider $\Gamma(n)=a\left({n_{0}+n}\right)^{-\frac{2r}{2r+1}}$ and $\Lambda(n)=\frac{1}{a}\left({n_{0}+n}\right)^{-\frac{1}{2r+1}}$. The theorem is stated for all $a\geqslant 4$: we choose $a=4$.

For Zhangl’s algorithm, we refer to Part 2.2 in , and choose $\Gamma(n)=\gamma_{0}n^{-\frac{2r}{2r+1}}$ with $\gamma_{0}=\frac{1}{R^{2}}$ which behaves better than the proposed choice $\frac{1}{4(1+R^{2})}$.

Finally, we sum up the rates that were both predicted and derived for the four algorithms in the four cases for $(\alpha,\delta)$ in Table 3. It appears that (a) we approximatively match the predicted rates in most cases (they would if $n$ was larger), (b) our rates improve on existing work.

Conclusion

In this paper, we have provided an analysis of averaged unregularized stochastic gradient methods for kernel-based least-squares regression. Our novel analysis allowed us to consider larger step-sizes, which in turn lead to optimal estimation rates for many settings of eigenvalue decay of the covariance operators and smoothness of the optimal prediction function. Moreover, we have worked on a more general setting than previous work, that includes most interesting cases of positive definite kernels.

Our work can be extended in a number of interesting ways: First, (a) we have considered results in expectation; following the higher-order moment bounds from in the Euclidean case, we could consider higher-order moments, which in turn could lead to high-probability results or almost-sure convergence. Moreover, (b) while we obtain optimal convergence rates for a particular regime of kernels/objective functions, using different types of averaging (i.e., non uniform) may lead to optimal rates in other regimes. Besides, (c) following , we could extend our results for infinite-dimensional least-squares regression to other smooth loss functions, such as for logistic regression, where an online Newton algorithm with the same running-time complexity would also lead to optimal convergence rates. Also, (d) the running-time complexity of our stochastic approximation procedures is still quadratic in the number of samples $n$, which is unsatisfactory when $n$ is large; by considering reduced set-methods , we hope to able to obtain a complexity of $O(d_{n}n)$, where $d_{n}$ is such that the convergence rate is $O(d_{n}/n)$, which would extend the Euclidean space result, where $d_{n}$ is constant equal to the dimension. Finally, (e) in order to obtain the optimal rates when the bias term dominates our generalization bounds, it would be interesting to combine our spectral analysis with recent accelerated versions of stochastic gradient descent which have been analyzed in the finite-dimensional setting .

Acknowledgements

This work was partially supported by the European Research Council (SIERRA Project). We thank Nicolas Flammarion for helpful discussions.

Appendix A Minimal assumptions

We first define the set of square $\rho_{X}$-integrable functions $\mathcal{L}^{2}_{\rho_{{X}}}$:

we will always make the assumptions that this space is separable (this is the case in most interesting situations. See for more details.) ${L}^{2}_{\rho_{{X}}}$ is its quotient under the equivalence relation given by

which makes it a separable Hilbert space (see, e.g., ).

A.2 Isomorphism

As it has been explained in the main text, the minimization problem will appear to be an approximation problem in $\mathcal{L}^{2}_{\rho_{{X}}}$, for which we will build estimates in $\mathcal{H}$. However, to derive theoretical results, it is easier to consider it as an approximation problem in the Hilbert space ${L}^{2}_{\rho_{{X}}}$, building estimates in $p(\mathcal{H})$.

We thus need to define a notion of the best estimation in $p(\mathcal{H})$. We first define the closure $\overline{{F}}$ (with respect to $\|\cdot\|_{{L}^{2}_{\rho_{{X}}}}$) of any set $F\subset{L}^{2}_{\rho_{{X}}}$ as the set of limits of sequences in ${F}$. The space $\overline{p(\mathcal{H})}$ is a closed and convex subset in ${L}^{2}_{\rho_{{X}}}$. We can thus define $g_{\mathcal{H}}=\arg\min_{f\in\,\overline{p(\mathcal{H})}}\varepsilon(g)$, as the orthogonal projection of $g_{\rho}$ on $\overline{p(\mathcal{H})}$, using the existence of the projection on any closed convex set in a Hilbert space. See Proposition 8 in Appendix A for details.

Assume (A1-2). The minimum of $\varepsilon(f)$ in $\overline{p(\mathcal{H})}$ is attained at a certain $g_{\mathcal{H}}$ (which is unique and well defined in ${L}^{2}_{\rho_{{X}}}$).

Where \overline{{p(\mathcal{H})}}=\Big{\{}f\in{L}^{2}_{\rho_{{X}}}\ /\ \exists(f_{n})\subset{{p(\mathcal{H})}},\|f_{n}-f\|_{{L}^{2}_{\rho_{{X}}}}\rightarrow 0\Big{\}} is the set of functions $f$ for which we can hope for consistency, i.e., having a sequence $(f_{n})_{n}$ of estimators in $\mathcal{H}$ such that $\varepsilon(f_{n})\rightarrow\varepsilon(f)$.

The properties of our estimator, especially its rate of convergence will strongly depend on some properties of both the kernel, the objective function and the distributions, which may be seen through the properties of the covariance operator which is defined in the main text. We have defined the covariance operator, $\Sigma:\mathcal{H}\rightarrow\mathcal{H}$. In the following, we extend such an operator as an endomorphism $\mathcal{T}$ from ${L}^{2}_{\rho_{{X}}}$ to $\mathcal{L}^{2}_{\rho_{{X}}}$ and by projection as an endomorphism $T=p\circ\mathcal{T}$ from ${L}^{2}_{\rho_{{X}}}$ to ${L}^{2}_{\rho_{{X}}}$. Note that $\mathcal{T}$ is well defined as $\int_{\mathcal{X}}g(t)\ K_{t}\ d\rho_{\mathcal{X}}(t)$ does not depend on the function $g$ chosen in the class of equivalence of $g$.

Assume (A1-2). We define the operator $\mathcal{T}$ as follows (this expectation is formally defined as a Bochner expectation in $\mathcal{H}$.):

A first important remark is that $\Sigma f=0$ implies $\langle f,\Sigma f\rangle=\|f\|_{{L}^{2}_{\rho_{{X}}}}^{2}=0$, that is $p(\operatorname{Ker}(\Sigma))=\{0\}$. However, $\Sigma$ may not be injective (unless $\|f\|_{{L}^{2}_{\rho_{{X}}}}^{2}\Rightarrow f=0$, which is true when $f$ is continuous and $\rho_{X}$ has full support). $\Sigma$ and $\mathcal{T}$ may independently be injective or not.

The operator $T$ (which is an endomorphism of the separable Hilbert space ${L}^{2}_{\rho_{{X}}}$) can be reduced in some Hilbertian eigenbasis of ${L}^{2}_{\rho_{{X}}}$. The linear operator $\mathcal{T}$ happens to have an image included in $\mathcal{H}$, and the eigenbasis of $T$ in ${L}^{2}_{\rho_{{X}}}$ may also be seen as eigenbasis of $\Sigma$ in $\mathcal{H}$ (See proof in Appendix I.2, Proposition 18):

Assume (A1-2). The image of $\mathcal{T}$ is included in $\mathcal{H}$: $\text{Im}(\mathcal{T})\subset\mathcal{H}$, that is, for any $f\in{L}^{2}_{\rho_{{X}}}$, $\mathcal{T}f\in\mathcal{H}$. Moreover, for any $i\in I$, $\phi_{i}^{H}=\frac{1}{\mu_{i}}\mathcal{T}\phi_{i}\in\mathcal{H}\subset\mathcal{L}^{2}_{\rho_{{X}}}$ is a representant for the equivalence class $\phi_{i}$, that is $p(\phi_{i}^{H})=\phi_{i}$. Moreover $\mu_{i}^{1/2}\phi_{i}^{H}$ is an orthonormal eigen-system of the orthogonal supplement $\mathscr{S}$ of the null space $\text{Ker}(\Sigma)$. That is:

$\forall i\in I,\ \Sigma\phi_{i}^{H}=\mu_{i}\phi_{i}^{H}$.

$\mathcal{H}=\operatorname{Ker}(\Sigma)\overset{\perp}{\oplus}\mathscr{S}$.

Such decompositions allow to define $\mathcal{T}^{r}:{L}^{2}_{\rho_{{X}}}\rightarrow\mathcal{H}$ for $r\geqslant 1/2$. Indeed , completeness allows to define infinite sums which satisfy a Cauchy criterion. See proof in Appendix I.2, Proposition 19. Note the different condition concerning $r$ in the definitions. For $r\geqslant 1/2$, $T^{r}=p\circ\mathcal{T}^{r}$. We need $r\geqslant 1/2$, because $(\mu_{i}^{1/2}\phi^{H})$ is an orthonormal system of $\mathscr{S}$.

We define, for any $r\geqslant 1/2$, $\mathcal{T}^{r}:{L}^{2}_{\rho_{{X}}}\rightarrow\mathcal{H}$, for any $h\in\operatorname{Ker}(T)$ and $(a_{i})_{i\in I}$ such that $\sum_{i\in I}a_{i}^{2}<\infty$, through: