Regularization by spectral filtering

Spectral regularization is any of a class of regularization techniques used in machine learning to control the impact of noise and prevent overfitting. Spectral regularization can be used in a broad range of applications, from deblurring images to classifying emails into a spam folder and a non-spam folder. For instance, in the email classification example, spectral regularization can be used to reduce the impact of noise and prevent overfitting when a machine learning system is being trained on a labeled set of emails to learn how to tell a spam and a non-spam email apart.

Spectral regularization algorithms rely on methods that were originally defined and studied in the theory of ill-posed inverse problems (for instance, see) focusing on the inversion of a linear operator (or a matrix) that possibly has a bad condition number or an unbounded inverse. In this context, regularization amounts to substituting the original operator by a bounded operator called the "regularization operator" that has a condition number controlled by a regularization parameter, a classical example being Tikhonov regularization. To ensure stability, this regularization parameter is tuned based on the level of noise. The main idea behind spectral regularization is that each regularization operator can be described using spectral calculus as an appropriate filter on the eigenvalues of the operator that defines the problem, and the role of the filter is to "suppress the oscillatory behavior corresponding to small eigenvalues". Therefore, each algorithm in the class of spectral regularization algorithms is defined by a suitable filter function (which needs to be derived for that particular algorithm). Three of the most commonly used regularization algorithms for which spectral filtering is well-studied are Tikhonov regularization, Landweber iteration, and truncated singular value decomposition (TSVD). As for choosing the regularization parameter, examples of candidate methods to compute this parameter include the discrepancy principle, generalized cross validation, and the L-curve criterion.

It is of note that the notion of spectral filtering studied in the context of machine learning is closely connected to the literature on function approximation (in signal processing).

Notation

The training set is defined as $S = \{(x_1, y_1), \dots , (x_n, y_n)\}$, where $X$ is the $n \times d$ input matrix and $Y = (y_1,\dots,y_n)$ is the output vector. Where applicable, the kernel function is denoted by $k$, and the $n \times n$ kernel matrix is denoted by $K$ which has entries $K_{ij} = k(x_i,x_j)$ and $\mathcal{H}$ denotes the Reproducing Kernel Hilbert Space (RKHS) with kernel $k$. The regularization parameter is denoted by $\lambda$.

(Note: For $g \in G$ and $f \in F$, with $G$ and $F$ being Hilbert spaces, given a linear, continuous operator $L$, assume that $g = Lf$ holds. In this setting, the direct problem would be to solve for $g$ given $f$ and the inverse problem would be to solve for $f$ given $g$. If the solution exists, is unique and stable, the inverse problem (i.e. the problem of solving for $f$) is well-posed; otherwise, it is ill-posed.)

Relation to the theory of ill-posed inverse problems

The connection between the regularized least squares (RLS) estimation problem (Tikhonov regularization setting) and the theory of ill-posed inverse problems is an example of how spectral regularization algorithms are related to the theory of ill-posed inverse problems.

The RLS estimator solves

$\min_{f\in\mathcal{H}} \frac{1}{n}\sum_{i=1}^{n}(y_i-f(x_i))^2+\lambda \|f\|^2_\mathcal{H}$

and the RKHS allows for expressing this RLS estimator as $f_S^\lambda (X)=\sum_{i=1}^n c_i k(x,x_i)$ where $(K+n\lambda I)c=Y$ with $c=(c_1,\dots,c_n)$. The penalization term is used for controlling smoothness and preventing overfitting. Since the solution of empirical risk minimization $\min_{f\in\mathcal{H}} \frac{1}{n}\sum_{i=1}^{n}(y_i-f(x_i))^2$ can be written as $f_S^{\lambda}(X)=\sum_{i=1}^n c_i k(x,x_i)$ such that $Kc=Y$, adding the penalty function amounts to the following change in the system that needs to be solved:

$\bigg\{ \min_{f\in\mathcal{H}} \frac{1}{n}\sum_{i=1}^n (y_i-f(x_i))^2\rightarrow \min_{f\in\mathcal{H}} \frac{1}{n}\sum_{i=1}^n (y_i-f(x_i))^2+\lambda \|f\|^2_\mathcal{H}\bigg\}\equiv \bigg\{Kc=Y\rightarrow (K+n\lambda I)c=Y\bigg\}.$

In this learning setting, the kernel matrix can be decomposed as $K = Q\Sigma Q^T$, with

$\sigma = \operatorname{diag}(\sigma_1,\dots,\sigma_n),~\sigma_1 \geq \sigma_2 \geq \cdots \geq \sigma_n \geq 0$

and $q_1,\dots,q_n$ are the corresponding eigenvectors. Therefore, in the initial learning setting, the following holds:

$c=K^{-1}Y=Q\Sigma^{-1}Q^TY=\sum_{i=1}^n \frac{1}{\sigma_i} \langle q_i,Y \rangle q_i.$

Thus, for small eigenvalues, even small perturbations in the data can lead to considerable changes in the solution. Hence, the problem is ill-conditioned, and solving this RLS problem amounts to stabilizing a possibly ill-conditioned matrix inversion problem, which is studied in the theory of ill-posed inverse problems; in both problems, a main concern is to deal with the issue of numerical stability.

Implementation of algorithms

Each algorithm in the class of spectral regularization algorithms is defined by a suitable filter function, denoted here by $G_{\lambda}(\cdot)$. If the Kernel matrix is denoted by $K$, then $\lambda$ should control the magnitude of the smaller eigenvalues of $G_{\lambda}(K)$. In a filtering setup, the goal is to find estimators $f_S^{\lambda}(X):=\sum_{i=1}^n c_i k(x,x_i)$ where $c=G_{\lambda}(K)Y$. To do so, a scalar filter function $G_\lambda(\sigma)$ is defined using the eigen-decomposition of the kernel matrix:

$G_\lambda(K)=QG_{\lambda}(\Sigma)Q^T,$

which yields

$G_{\lambda}(K)Y~=~\sum_{i=1}^{n}G_{\lambda}(\sigma_i) \langle q_i,Y\rangle q_i.$

Typically, an appropriate filter function should have the following properties:

1. As $\lambda$ goes to zero, $G_\lambda (\sigma)~\rightarrow ~1/\sigma$.

2. The magnitude of the (smaller) eigenvalues of $G_{\lambda}$ is controlled by $\lambda$.

While the above items give a rough characterization of the general properties of filter functions for all spectral regularization algorithms, the derivation of the filter function (and hence its exact form) varies depending on the specific regularization method that spectral filtering is applied to.

Filter function for Tikhonov regularization

In the Tikhonov regularization setting, the filter function for RLS is described below. As shown in, in this setting, $c=(K+n\lambda I)^{-1}Y$. Thus,

$c=(K+n\lambda I)^{-1}Y=Q(\Sigma+n\lambda I)^{-1}Q^T Y=\sum_{i=1}^n \frac{1}{\sigma_i+n\lambda}<q_i,Y>q_i.$

The undesired components are filtered out using regularization:

• If $\sigma \gg \lambda n$, then $\frac{1}{\sigma_i+n\lambda} \sim \frac{1}{\sigma_i}$.
• If $\sigma \ll \lambda n$, then $\frac{1}{\sigma_i+n\lambda} \sim \frac{1}{\lambda n}$.

The filter function for Tikhonov regularization is therefore defined as:

$G_\lambda (\sigma)=\frac{1}{\sigma+n\lambda }.$

Filter function for Landweber iteration

The idea behind the Landweber iteration is gradient descent:

$c^0 =0 \,$

$\text{for } i=1,\dots,t-1$

$~~~~~ c^i=c^{i-1}+\eta(Y-Kc^{i-1})$

$\mathrm{end}$

In this setting, if $n$ is larger than $K$'s largest eigenvalue, the above iteration converges by choosing $\eta = 2/n$ as the step-size:. The above iteration is equivalent to minimizing $\frac{1}{n}||Y-Kc||_2^2$ (i.e. the empirical risk) via gradient descent; using induction, it can be proved that at the $t$-th iteration, the solution is given by

$c=\eta\sum_{i=0}^{t-1}(I-\eta K)^iY.$

Thus, the appropriate filter function is defined by:

$G_\lambda(\sigma)=\eta\sum_{i=0}^{t-1}(I-\eta\sigma)^i.$

It can be shown that this filter function corresponds to a truncated power expansion of $K^{-1}$; to see this, note that the relation $\sum_{i\geq0} x^i = 1/(1 - x)$, would still hold if $x$ is replaced by a matrix; thus, if $K$ (the kernel matrix), or rather $I - \eta K$, is considered, the following holds:

$K^{-1}=\eta\sum_{i=0}^\infty (I-\eta K)^i \sim \eta\sum_{i=0}^{t-1}(I-\eta K)^i.$

In this setting, the number of iterations gives the regularization parameter; roughly speaking, $t \sim 1/\lambda$. If $t$ is large, overfitting may be a concern. If $t$ is small, oversmoothing may be a concern. Thus, choosing an appropriate time for early stopping of the iterations provides a regularization effect.

Filter function for TSVD

In the TSVD setting, given the eigen-decomposition $K = Q\Sigma Q^T$ and using a prescribed threshold $\lambda n$, a regularized inverse can be formed for the kernel matrix by discarding all the eigenvalues that are smaller than this threshold. Thus, the filter function for TSVD can be defined as

$G_\lambda(\sigma) =\left\{\begin{array}[c]{lcll}1/\sigma & ,& \text{if }\sigma\geq\lambda n\\[0.05in]0&,& \text{otherwise}\\[0.05in]\end{array}\right..$

It can be shown that TSVD is equivalent to the (unsupervised) projection of the data using (kernel) Principal Component Analysis (PCA), and that it is also equivalent to minimizing the empirical risk on the projected data (without regularization). Note that the number of components kept for the projection is the only free parameter here.