End-to-End Kernel Learning with Supervised Convolutional Kernel Networks

Introduction

In the past years, deep neural networks such as convolutional or recurrent ones have become highly popular for solving various prediction problems, notably in computer vision and natural language processing. Conceptually close to approaches that were developed several decades ago (see, ), they greatly benefit from the large amounts of labeled data that have been made available recently, allowing to learn huge numbers of model parameters without worrying too much about overfitting. Among other reasons explaining their success, the engineering effort of the deep learning community and various methodological improvements have made it possible to learn in a day on a GPU complex models that would have required weeks of computations on a traditional CPU (see, e.g., ).

Before the resurgence of neural networks, non-parametric models based on positive definite kernels were one of the most dominant topics in machine learning . These approaches are still widely used today because of several attractive features. Kernel methods are indeed versatile; as long as a positive definite kernel is specified for the type of data considered—e.g., vectors, sequences, graphs, or sets—a large class of machine learning algorithms originally defined for linear models may be used. This family include supervised formulations such as support vector machines and unsupervised ones such as principal or canonical component analysis, or K-means and spectral clustering. The problem of data representation is thus decoupled from that of learning theory and algorithms. Kernel methods also admit natural mechanisms to control the learning capacity and reduce overfitting .

On the other hand, traditional kernel methods suffer from several drawbacks. The first one is their computational complexity, which grows quadratically with the sample size due to the computation of the Gram matrix. Fortunately, significant progress has been achieved to solve the scalability issue, either by exploiting low-rank approximations of the kernel matrix , or with random sampling techniques for shift-invariant kernels . The second disadvantage is more critical; by decoupling learning and data representation, kernel methods seem by nature incompatible with end-to-end learning—that is, the representation of data adapted to the task at hand, which is the cornerstone of deep neural networks and one of the main reason of their success. The main objective of this paper is precisely to tackle this issue in the context of image modeling.

Specifically, our approach is based on convolutional kernel networks, which have been recently introduced in . Similar to hierarchical kernel descriptors , local image neighborhoods are mapped to points in a reproducing kernel Hilbert space via the kernel trick. Then, hierarchical representations are built via kernel compositions, producing a sequence of “feature maps” akin to convolutional neural networks, but of infinite dimension. To make the image model computationally tractable, convolutional kernel networks provide an approximation scheme that can be interpreted as a particular type of convolutional neural network learned without supervision.

To perform end-to-end learning given labeled data, we use a simple but effective principle consisting of learning discriminative subspaces in RKHSs, where we project data. We implement this idea in the context of convolutional kernel networks, where linear subspaces, one per layer, are jointly optimized by minimizing a supervised loss function. The formulation turns out to be a new type of convolutional neural network with a non-standard parametrization. The network also admits simple principles to learn without supervision: learning the subspaces may be indeed achieved efficiently with classical kernel approximation techniques .

To demonstrate the effectiveness of our approach in various contexts, we consider image classification benchmarks such as CIFAR-10 and SVHN , which are often used to evaluate deep neural networks; then, we adapt our model to perform image super-resolution, which is a challenging inverse problem. On the SVHN and CIFAR-10 datasets, we obtain a competitive accuracy, with about $2\%$ and $10\%$ error rates, respectively, without model averaging or data augmentation. For image up-scaling, we outperform recent approaches based on classical convolutional neural networks .

We believe that these results are highly promising. Our image model achieves competitive performance in two different contexts, paving the way to many other applications. Moreover, our results are also subject to improvements. In particular, we did not use GPUs yet, which has limited our ability to exhaustively explore model hyper-parameters and evaluate the accuracy of large networks. We also did not investigate classical regularization/optimization techniques such as Dropout , batch normalization , or recent advances allowing to train very deep networks . To gain more scalability and start exploring these directions, we are currently working on a GPU implementation, which we plan to publicly release along with our current CPU implementation.

One of our goals is to make a bridge between kernel methods and deep networks, and ideally reach the best of both worlds. Given the potentially attractive features of such a combination, several attempts have been made in the past to unify these two schools of thought. A first proof of concept was introduced in with the arc-cosine kernel, which admits an integral representation that can be interpreted as a one-layer neural network with random weights and infinite number of rectified linear units. Besides, a multilayer kernel may be obtained by kernel compositions . Then, hierarchical kernel descriptors and convolutional kernel networks extend a similar idea in the context of images leading to unsupervised representations .

Multiple kernel learning is also related to our work since is it is a notable attempt to introduce supervision in the kernel design. It provides techniques to select a combination of kernels from a pre-defined collection, and typically requires to have already “good” kernels in the collection to perform well. More related to our work, the backpropagation algorithm for the Fisher kernel introduced in learns the parameters of a Gaussian mixture model with supervision. In comparison, our approach does not require a probabilistic model and learns parameters at several layers. Finally, we note that a concurrent effort to ours is conducted in the Bayesian community with deep Gaussian processes , complementing the Frequentist approach that we follow in our paper.

Learning Hierarchies of Subspaces with Convolutional Kernel Networks

In this section, we present the principles of convolutional kernel networks and a few generalizations and improvements of the original approach of . Essentially, the model builds upon four ideas that are detailed below and that are illustrated in Figure 1 for a model with a single layer.

Then, it is worth noting that the encoding function $\psi_{1}$ with kernel (2) is reminiscent of radial basis function networks (RBFNs) , whose hidden layer resembles (3) without the matrix $\kappa_{1}({\mathbf{Z}}^{\top}{\mathbf{Z}})^{-1/2}$ and with no normalization. The difference between RBFNs and our model is nevertheless significant. The RKHS mapping, which is absent from RBFNs, is indeed a key to the multilayer construction that will be presented shortly: a network layer takes points from the RKHS’s previous layer as input and use the corresponding RKHS inner-product. To the best of our knowledge, there is no similar multilayer and/or convolutional construction in the radial basis function network literature.

The multilayer scheme produces a sequence of maps $(I_{k})_{k\geq 0}$, where each vector $I_{k}(z)$ encodes a point—say $f_{k}(z)$—in the linear subspace $\mathcal{F}_{k}$ of $\mathcal{H}_{k}$. Thus, we implicitly represent an image at layer $k$ as a spatial map $f_{k}:\Omega_{k}\to\mathcal{H}_{k}$ such that $\langle I_{k}(z),I_{k}^{\prime}(z^{\prime})\rangle=\langle f_{k}(z),f_{k}^{\prime}(z^{\prime})\rangle_{\mathcal{H}_{k}}$ for all $z,z^{\prime}$. As mentioned previously, the mapping to the RKHS is a key to the multilayer construction. Given $I_{k}$, larger image neighborhoods are represented by patches of size $e_{k}\times e_{k}$ that can be mapped to a point in the Cartesian product space $\mathcal{H}_{k}^{e_{k}\times e_{k}}$ endowed with its natural inner-product; finally, the kernel $K_{k+1}$ defined on these patches can be seen as a kernel on larger image neighborhoods than $K_{k}$.

End-to-End Kernel Learning with Supervised CKNs

In the previous section, we have described a variant of convolutional kernel networks where linear subspaces are learned at every layer. This is achieved without supervision by a K-means algorithm leading to small projection residuals. It is thus natural to introduce also a discriminative approach.

Given a positive definite kernel $K$ on images, the classical empirical risk minimization formulation consists of finding a prediction function in the RKHS $\mathcal{H}$ associated to $K$ by minimizing the objective

where the parameter $\lambda$ controls the smoothness of the prediction function $f$ with respect to the geometry induced by the kernel, hence regularizing and reducing overfitting . After training a convolutional kernel network with $k$ layers, such a positive definite kernel may be defined as

where $I_{k},I_{k}^{\prime}$ are the $k$-th finite-dimensional feature maps of $I_{0},I_{0}^{\prime}$, respectively, and $f_{k},f_{k}^{\prime}$ the corresponding maps in $\Omega_{k}\to\mathcal{H}_{k}$, which have been defined in the previous section. The kernel is also indexed by $\mathcal{Z}$, which represents the network parameters—that is, the subspaces $\mathcal{F}_{1},\ldots,\mathcal{F}_{k}$, or equivalently the set of filters ${\mathbf{Z}}_{1},\ldots,{\mathbf{Z}}_{k}$ from Eq. (3). Then, formulation (5) becomes equivalent to

Since we consider a smooth loss function $L$, e.g., logistic, squared hinge, or square loss, optimizing (7) with respect to ${\mathbf{W}}$ can be achieved with any gradient-based method. Moreover, when $L$ is convex, we may also use fast dedicated solvers, (see, e.g., , and references therein). Optimizing with respect to the filters ${\mathbf{Z}}_{j}$, $j=1,\ldots,k$ is more involved because of the lack of convexity. Yet, the objective function is differentiable, and there is hope to find a “good” stationary point by using classical stochastic optimization techniques that have been successful for training deep networks.

For that, we need to compute the gradient by using the chain rule—also called “backpropagation” . We instantiate this rule in the next lemma, which we have found useful to simplify the calculation.

where the inner-product is the Frobenius’s one and $g_{j},h_{j}$ are linear functions. Then,

where $L^{\prime}$ denote the derivative of the smooth function $L$ with respect to its second argument.

The proof of this lemma is straightforward and follows from the definition of the Fréchet derivative. Nevertheless, it is useful to derive the closed form of the gradient in the next proposition.

Consider the quantities introduced in Lemma 1, but denote $I_{j}^{\mathcal{Z}}$ by $I_{j}$ for simplicity. By construction, we have for all $j\geq 1$,

where ${\mathbf{U}}$ is any matrix of the same size as $I_{j}$, $M_{j}={\mathbf{A}}_{j}\kappa_{j}({\mathbf{Z}}_{j}^{\top}{\mathbf{E}}_{j}(I_{{j-1}}){\mathbf{S}}_{j}^{-1}){\mathbf{S}}_{j}$ is the $j$-th feature map before the pooling step, $\odot$ is the Hadamart (elementwise) product, ${\mathbf{E}}_{j}^{\star}$ is the adjoint of ${\mathbf{E}}_{j}$, and

Computing the gradient requires a forward pass to obtain the maps $I_{j}$ through (11) and a backward pass that composes the functions $g_{j},h_{j}$ as in (10). The complexity of the forward step is dominated by the convolutions ${\mathbf{Z}}_{j}^{\top}{\mathbf{E}}_{j}(I_{{j-1}})$, as in convolutional neural networks. The cost of the backward pass is the same as the forward one up to a constant factor. Assuming $p_{j}\!\leq\!|\Omega_{{j-1}}|$, which is typical for lower layers that require more computation than upper ones, the most expensive cost is due to ${\mathbf{E}}_{j}(I_{{j-1}}){\mathbf{B}}_{j}^{\top}$ and ${\mathbf{Z}}_{j}{\mathbf{B}}_{j}$ which is the same as ${\mathbf{Z}}_{j}^{\top}{\mathbf{E}}_{j}(I_{{j-1}})$. We also pre-compute ${\mathbf{A}}_{j}^{1/2}$ and ${\mathbf{A}}_{j}^{3/2}$ by eigenvalue decompositions, whose cost is reasonable when performed only once per minibatch. Off-diagonal elements of $M_{j}^{\top}{\mathbf{U}}{\mathbf{P}}_{j}^{\top}-{\mathbf{E}}_{j}(I_{{j-1}})^{\top}{\mathbf{Z}}_{j}{\mathbf{B}}_{j}$ are also not computed since they are set to zero after elementwise multiplication with a diagonal matrix. In practice, we also replace ${\mathbf{A}}_{j}$ by $(\kappa_{j}({\mathbf{Z}}_{j}^{\top}{\mathbf{Z}}_{j})+\varepsilon{\mathbf{I}})^{-1/2}$ with $\varepsilon\!=\!0.001$, which corresponds to performing a regularized projection onto $\mathcal{F}_{j}$ (see Appendix A). Finally, a small offset of $0.00001$ is added to the diagonal entries of ${\mathbf{S}}_{j}$.

When using the kernel (2), the objective is differentiable with respect to the hyper-parameters $\alpha_{j}$. When large amounts of training data are available and overfitting is not a issue, optimizing the training loss by taking gradient steps with respect to these parameters seems appropriate instead of using a canonical parameter value. Otherwise, more involved techniques may be needed; we plan to investigate other strategies in future work.

2 Optimization and Practical Heuristics

The backpropagation rules of the previous section have set up the stage for using a stochastic gradient descent method (SGD). We now present a few strategies to accelerate it in our context.

Recently, many incremental optimization techniques have been proposed for solving convex optimization problems of the form (7) when $n$ is large but finite (see and references therein). These methods usually provide a great speed-up over the stochastic gradient descent algorithm without suffering from the burden of choosing a learning rate. The price to pay is that they rely on convexity, and they require storing into memory the full training set. For solving (7) with fixed network parameters $\mathcal{Z}$, it means storing the $n$ maps $I_{k}^{i}$, which is often reasonable if we do not use data augmentation. To partially leverage these fast algorithms for our non-convex problem, we have adopted a minimization scheme that alternates between two steps: (i) fix $\mathcal{Z}$, then make a forward pass on the data to compute the $n$ maps $I_{k}^{i}$ and minimize the convex problem (7) with respect to ${\mathbf{W}}$ using the accelerated MISO algorithm ; (ii) fix ${\mathbf{W}}$, then make one pass of a projected stochastic gradient algorithm to update the $k$ set of filters ${\mathbf{Z}}_{j}$. The set of network parameters $\mathcal{Z}$ is initialized with the unsupervised learning method described in Section 2.

Choosing the right learning rate in stochastic optimization is still an important issue despite the large amount of work existing on the topic, see, e.g., and references therein. In our paper, we use the following basic heuristic: the initial learning rate $\eta_{t}$ is chosen “large enough”; then, the training loss is evaluated after each update of the weights ${\mathbf{W}}$. When the training loss increases between two epochs, we simply divide the learning rate by two, and perform “back-tracking” by replacing the current network parameters by the previous ones.

For classification tasks, “easy” samples have often negligible contribution to the gradient (see, e.g., ). For instance, for the squared hinge loss $L(y,\hat{y})=\max(0,1-y\hat{y})^{2}$, the gradient vanishes when the margin $y\hat{y}$ is greater than one. This motivates the following heuristic: we consider a set of active samples, initially all of them, and remove a sample from the active set as soon as we obtain zero when computing its gradient. In the subsequent optimization steps, only active samples are considered, and after each epoch, we randomly reactivate $10\%$ of the inactive ones.

Experiments

We now present experiments on image classification and super-resolution. All experiments were conducted on $8$-core and $10$-core $2.4$GHz Intel CPUs using C++ and Matlab.

We consider the datasets CIFAR-10 and SVHN , which contain $32\times 32$ images from $10$ classes. CIFAR-10 is medium-sized with $50\,000$ training samples and $10\,000$ test ones. SVHN is larger with $604\,388$ training examples and $26\,032$ test ones. We evaluate the performance of a $9$-layer network, designed with few hyper-parameters: for each layer, we learn $512$ filters and choose the RBF kernels $\kappa_{j}$ defined in (2) with initial parameters $\alpha_{j}\!=\!1/(0.5^{2})$. Layers $1,3,5,7,9$ use $3\!\times\!3$ patches and a subsampling pooling factor of $\sqrt{2}$ except for layer $9$ where the factor is $3$; Layers $2,4,6,8$ use simply $1\times 1$ patches and no subsampling. For CIFAR-10, the parameters $\alpha_{j}$ are kept fixed during training, and for SVHN, they are updated in the same way as the filters. We use the squared hinge loss in a one vs all setting to perform multi-class classification (with shared filters $\mathcal{Z}$ between classes). The input of the network is pre-processed with the local whitening procedure described in . We use the optimization heuristics from the previous section, notably the automatic learning rate scheme, and a gradient momentum with parameter $0.9$, following . The regularization parameter $\lambda$ and the number of epochs are set by first running the algorithm on a $80/20$ validation split of the training set. $\lambda$ is chosen near the canonical parameter $\lambda=1/n$, in the range $2^{i}/n$, with $i=-4,\ldots,4$, and the number of epochs is at most $100$. The initial learning rate is $10$ with a minibatch size of $128$.

We present our results in Table 1 along with the performance achieved by a few recent methods without data augmentation or model voting/averaging. In this context, the best published results are obtained by the generalized pooling scheme of . We achieve about $2\%$ test error on SVHN and about $10\%$ on CIFAR-10, which positions our method as a reasonably “competitive” one, in the same ballpark as the deeply supervised nets of or network in network of .

Due to lack of space, the results reported here only include a single supervised model. Preliminary experiments with no supervision show also that one may obtain competitive accuracy with wide shallow architectures. For instance, a two-layer network with (1024-16384) filters achieves $14.2\%$ error on CIFAR-10. Note also that our unsupervised model outperforms original CKNs . The best single model from gives indeed $21.7\%$. Training the same architecture with our approach is two orders of magnitude faster and gives $19.3\%$. Another aspect we did not study is model complexity. Here as well, preliminary experiments are encouraging. Reducing the number of filters to $128$ per layer yields indeed $11.95\%$ error on CIFAR-10 and $2.15\%$ on SVHN. A more precise comparison with no supervision and with various network complexities will be presented in another venue.

2 Image Super-Resolution from a Single Image

Image up-scaling is a challenging problem, where convolutional neural networks have obtained significant success . Here, we follow and replace traditional convolutional neural networks by our supervised kernel machine. Specifically, RGB images are converted to the YCbCr color space and the upscaling method is applied to the luminance channel only to make the comparison possible with previous work. Then, the problem is formulated as a multivariate regression one. We build a database of $200\,000$ patches of size $32\times 32$ randomly extracted from the BSD500 dataset after removing image 302003.jpg, which overlaps with one of the test images. $16\times 16$ versions of the patches are build using the Matlab function imresize, and upscaled back to $32\times 32$ by using bicubic interpolation; then, the goal is to predict high-resolution images from blurry bicubic interpolations.

The blurry estimates are processed by a $9$-layer network, with $3\times 3$ patches and $128$ filters at every layer without linear pooling and zero-padding. Pixel values are predicted with a linear model applied to the $128$-dimensional vectors present at every pixel location of the last layer, and we use the square loss to measure the fit. The optimization procedure and the kernels $\kappa_{j}$ are identical to the ones used for processing the SVHN dataset in the classification task. The pipeline also includes a pre-processing step, where we remove from input images a local mean component obtained by convolving the images with a $5\times 5$ averaging box filter; the mean component is added back after up-scaling.

For the evaluation, we consider three datasets: Set5 and Set14 are standard for super-resolution; Kodim is the Kodak Image database, available at http://r0k.us/graphics/kodak/, which contains high-quality images with no compression or demoisaicing artefacts. The evaluation procedure follows by using the code from the author’s web page. We present quantitative results in Table 2. For x3 upscaling, we simply used twice our model learned for x2 upscaling, followed by a 3/4 downsampling. This is clearly suboptimal since our model is not trained to up-scale by a factor 3, but this naive approach still outperforms other baselines that are trained end-to-end. Note that also proposes a data augmentation scheme at test time that slightly improves their results. In Appendix D, we also present a visual comparison between our approach and , whose pipeline is the closest to ours, up to the use of a supervised kernel machine instead of CNNs.

This work was supported by ANR (MACARON project ANR-14-CE23-0003-01).

References

First, we remark that the kernel $K_{1}$ is homogeneous such that for every patch ${\mathbf{x}}$ and scalar $\gamma>0$,

When $\kappa_{1}({\mathbf{Z}}^{\top}{\mathbf{Z}})$ is not invertible or simply badly conditioned, it is also common to use instead

where $\varepsilon>0$ is a small regularization that improves the condition number of $\kappa_{1}({\mathbf{Z}}^{\top}{\mathbf{Z}})$. Such a modification can be interpreted as performing a slightly regularized projection onto the finite-dimensional subspace $\mathcal{F}_{1}$.

Appendix B Computation of the Gradient with Respect to the Filters

To compute the gradient of the loss function, we use Lemma 1 and start by analyzing the effect of perturbing every quantity involved in (11) such that we may obtain the desired relations (8) and (9). Before proceeding, we recall the definition of the set $\mathcal{Z}+\mathcal{E}=\{{\mathbf{Z}}_{1}+\boldsymbol{\varepsilon}_{1},\ldots,{\mathbf{Z}}_{k}+\boldsymbol{\varepsilon}_{k}\}$ and the precise definition of the Landau notation $o(\|\mathcal{E}\|)$, which we use in (8). Here, it simply means a quantity that is negligible in front of the norm $\|\mathcal{E}\|=\sum_{j=1}^{k}\|\boldsymbol{\varepsilon}_{j}\|_{\text{F}}$—that is,

Then, we start by initializing a recursion: $I_{0}^{\mathcal{Z}}$ is unaffected by the perturbation and thus $\Delta I_{0}^{\mathcal{Z},\mathcal{E}}=0$. Consider now an index $j>0$ and assume that (8) holds for $j-1$ with $\Delta I_{{j-1}}^{\mathcal{Z},\mathcal{E}}=O(\|\mathcal{E}\|)$.

Then, the diagonal matrix ${\mathbf{S}}_{j}$ becomes after perturbation

The inverse diagonal matrix ${\mathbf{S}}_{j}^{-1}$ becomes

and the matrix ${\mathbf{A}}_{j}$ becomes

where we have used the relation $({\mathbf{I}}+{\mathbf{Q}})^{-1/2}={\mathbf{I}}-\frac{1}{2}{\mathbf{Q}}+o(\|{\mathbf{Q}}\|_{\text{F}})$. Note that the quantities $\Delta{\mathbf{A}}_{j},\Delta{\mathbf{S}}_{j},\Delta{\mathbf{S}}_{j}^{-1}$ that we have introduced are all $O(\|\mathcal{E}\|)$. Then, by replacing the quantities ${\mathbf{A}}_{j},{\mathbf{S}}_{j},{\mathbf{S}}_{j}^{-1},{\mathbf{I}}_{{j-1}}$ by their perturbed versions in the definition of $I_{j}$ given in (11), we obtain that $I_{j}^{\mathcal{Z}+\mathcal{E}}$ is equal to

Then, after short calculation, we obtain the desired relation $I_{j}^{\mathcal{Z}+\mathcal{E}}=I_{j}^{\mathcal{Z}+\mathcal{E}}+\Delta I_{j}^{\mathcal{Z},\mathcal{E}}+o(\|\mathcal{E}\|)$ with

First, we remark that $\Delta I_{j}^{\mathcal{Z},\mathcal{E}}=O(\|\mathcal{E}\|)$, which is one of the induction hypothesis we need. Then, after plugging in the values of $\Delta{\mathbf{A}}_{j},\Delta{\mathbf{S}}_{j},\Delta{\mathbf{S}}_{j}^{-1}$, and with further simplification, we obtain

where $M_{j}^{\mathcal{Z}}$ is the $j$-th feature map of $I_{0}$ before the $j$-th linear pooling step—that is, $I_{j}^{\mathcal{Z}}=M_{j}^{\mathcal{Z}}{\mathbf{P}}_{j}$. We now see that $\Delta I_{j}^{\mathcal{Z},\mathcal{E}}$ is linear in $\boldsymbol{\varepsilon}_{j}$ and $\Delta I_{{j-1}}^{\mathcal{Z},\mathcal{E}}$, which guarantees that there exist two linear functions $g_{j},h_{j}$ that satisfy (9). More precisely, we want for all matrix ${\mathbf{U}}$ of the same size as $\Delta I_{j}^{\mathcal{Z},\mathcal{E}}$

Then, it is easy to obtain the form of $g_{j},h_{j}$ given in (12), by using in the right order the following elementary calculus rules: (i) $\langle{\mathbf{U}}{\mathbf{V}},{\mathbf{W}}\rangle=\langle{\mathbf{U}},{\mathbf{W}}{\mathbf{V}}^{\top}\rangle=\langle{\mathbf{V}},{\mathbf{U}}^{\top}{\mathbf{W}}\rangle$, (ii) $\langle{\mathbf{U}},{\mathbf{V}}\rangle=\langle{\mathbf{U}}^{\top},{\mathbf{V}}^{\top}\rangle$, (iii) $\langle{\mathbf{U}}\odot{\mathbf{V}},{\mathbf{W}}\rangle=\langle{\mathbf{U}},{\mathbf{V}}\odot{\mathbf{W}}\rangle$ for any matrices ${\mathbf{U}},{\mathbf{V}},{\mathbf{W}}$ of appropriate sizes, and also (iv) $\langle{\mathbf{E}}_{j}({\mathbf{U}}),{\mathbf{V}}\rangle=\langle{\mathbf{U}},{\mathbf{E}}_{j}^{\star}({\mathbf{V}})\rangle$, by definition of the adjoint operator. We conclude by induction.

Appendix C Preconditioning Heuristic on the Sphere

With the change of variable ${\mathbf{z}}={\mathbf{Q}}^{1/2}{\mathbf{w}}$, this is equivalent to

This is exactly the update rule we have chosen in our paper, as a heuristic in a stochastic setting.

Appendix D Additional Results for Image Super-Resolution

We present a quantitative comparison in Table 3 using the structural similarity index measure (SSIM), which is known to better reflect the quality perceived by humans than the PSNR; it is commonly used to evaluate the quality of super-resolution methods, see . Then, we present a visual comparison between several approaches in Figures 2, 3, and 4. We focus notably on the classical convolutional neural network of since our pipeline essentially differs in the use of our supervised kernel machine instead of convolutional neural networks. After subjective evaluation, we observe that both methods perform equally well in textured areas. However, our approach recovers better thin high-frequency details, such as the eyelash of the baby in the first image. By zooming on various parts, it is easy to notice similar differences in other images. We also observed a few ghosting artefacts near object boundaries with the method of , which is not the case with our approach.