Approximation and Learning with Deep Convolutional Models: a Kernel Perspective

Alberto Bietti

Introduction

Deep convolutional models have been at the heart of the recent successes of deep learning in problems where the data consists of high-dimensional signals, such as image classification or speech recognition. Convolution and pooling operations have notably contributed to the practical success of these models, yet our theoretical understanding of how they enable efficient learning is still limited.

One key difficulty for understanding such models is the curse of dimensionality: due to the high-dimensionality of the input data, it is hopeless to learn arbitrary functions from samples. For instance, classical non-parametric regression techniques for learning generic target functions typically require either low dimension or very high degrees of smoothness in order to obtain good generalization (e.g., Wainwright, 2019), which makes them impractical for dealing with high-dimensional signals. Thus, further assumptions on the target function are needed to make the problem tractable, and we seek assumptions that make convolutions a useful modeling tool. Various works have studied approximation benefits of depth with models that resemble deep convolutional architectures (Cohen & Shashua, 2017; Mhaskar & Poggio, 2016; Schmidt-Hieber et al., 2020). Nevertheless, while such function classes may provide improved statistical efficiency in theory, it is unclear if there exist efficient algorithms to learn such models, and hence, whether they might correspond to what convolutional networks learn in practice. To overcome this issue, we consider instead function classes based on kernel methods (Schölkopf & Smola, 2001; Wahba, 1990), which are known to be learnable with efficient (polynomial-time) algorithms, such as kernel ridge regression or gradient descent.

We consider “deep” structured kernels known as convolutional kernels, which yield good empirical performance on standard computer vision benchmarks (Li et al., 2019; Mairal, 2016; Mairal et al., 2014; Shankar et al., 2020), and are related to over-parameterized convolutional networks (CNNs) in so-called “kernel regimes” (Arora et al., 2019; Bietti & Mairal, 2019b; Daniely et al., 2016; Garriga-Alonso et al., 2019; Jacot et al., 2018; Novak et al., 2019; Yang, 2019). Such regimes may be seen as providing a first-order description of what common deep models trained with gradient methods may learn. Studying the corresponding function spaces (reproducing kernel Hilbert spaces, or RKHS) may then provide insight into the benefits of various architectural choices. For fully-connected architectures, such kernels are rotation-invariant, and the corresponding RKHSs are well understood in terms of regularity properties on the sphere (Bach, 2017a; Smola et al., 2001), but do not show any major differences between deep and shallow kernels (Bietti & Bach, 2021; Chen & Xu, 2021; Geifman et al., 2020). In contrast, in this work we show that even in the kernel setting, multiple layers of convolution and pooling operations can be crucial for efficient learning of functions with specific structures that are well-suited for natural signals. Our work paves the way for further studies of the inductive bias of optimization algorithms on deep convolutional networks beyond kernel regimes, for instance by incorporating adaptivity to low-dimensional structure (Bach, 2017a; Chizat & Bach, 2020; Wei et al., 2019) or hierarchical learning (Allen-Zhu & Li, 2020).

We revisit convolutional kernel networks (Mairal, 2016), finding that simple two or three layers models with Gaussian pooling and polynomial kernels of degree 2-4 at higher layers provide competitive performance with state-of-the-art convolutional kernels such as Myrtle kernels (Shankar et al., 2020) on Cifar10.

For such kernels, we provide an exact description of the RKHS functions and their norm, illustrating representation benefits of multiple convolutional and pooling layers for capturing additive and interaction models on patches with certain spatial regularities among interaction terms.

We provide generalization bounds that illustrate the benefits of architectural choices such as pooling and patches for learning additive interaction models with spatial invariance in the interaction terms, namely, improvements in sample complexity by polynomial factors in the size of the input signal.

Convolutional kernel networks were introduced by Mairal et al. (2014); Mairal (2016). Empirically, they used kernel approximations to improve computational efficiency, while we evaluate the exact kernels in order to assess their best performance, as in (Arora et al., 2019; Li et al., 2019; Shankar et al., 2020). Bietti & Mairal (2019a; b) show invariance and stability properties of its RKHS functions, and provide upper bounds on the RKHS norm for some specific functions (see also Zhang et al., 2017); in contrast, we provide exact characterizations of the RKHS norm, and study generalization benefits of certain architectures. Scetbon & Harchaoui (2020) study statistical properties of simple convolutional kernels without pooling, while we focus on the role of architecture choices with an emphasis on pooling. Cohen & Shashua (2016; 2017); Mhaskar & Poggio (2016); Poggio et al. (2017) study expressivity and approximation with models that resemble CNNs, showing benefits thanks to hierarchy or local interactions, but such models are not known to be learnable with tractable algorithms, while we focus on (tractable) kernels. Regularization properties of convolutional models were also considered in (Gunasekar et al., 2018; Heckel & Soltanolkotabi, 2020), but in different regimes or architectures than ours. Li et al. (2021); Malach & Shalev-Shwartz (2021) study benefits of convolutional networks with efficient algorithms, but do not study the gains of pooling. Du et al. (2018) study sample complexity of learning CNNs, focusing on parametric rather than non-parametric models. Mei et al. (2021) study statistical benefits of global pooling for learning invariant functions, but only consider one layer with full-size patches. Concurrently to our work, Favero et al. (2021); Misiakiewicz & Mei (2021) study benefits of local patches, but focus on one-layer models.

Deep Convolutional Kernels

In this section, we recall the construction of multi-layer convolutional kernels on discrete signals, following most closely the convolutional kernel network (CKN) architectures studied by Mairal (2016); Bietti & Mairal (2019a). These architectures rely crucially on pooling layers, typically with Gaussian filters, which make them empirically effective even with just two convolutional layers. These kernels define function spaces that will be the main focus of our theoretical study of approximation and generalization in the next sections. In particular, when learning a target function of the form f∗(x)=∑ifi(x)f^{*}(x)=\sum_{i}f_{i}(x), we will show that they are able to efficiently exploit two useful properties of f∗f^{*}: (locality) each fif_{i} may depend on only one or a few small localized patches of the signal; (invariance) many different terms fif_{i} may involve the same function applied to different input patches. We provide further background and motivation in Appendix A.

We note that our construction closely resembles kernels derived from infinitely wide convolutional networks, known as conjugate or NNGP kernels (Garriga-Alonso et al., 2019; Novak et al., 2019), and is also related to convolutional neural tangent kernels (Arora et al., 2019; Bietti & Mairal, 2019b; Yang, 2019). The Myrtle family of kernels (Shankar et al., 2020) also resembles our models, but they use small average pooling filters instead of Gaussian filters, which leads to deeper architectures due to smaller receptive fields.

Approximation with (Deep) Convolutional Kernels

In this section, we present our main results on the approximation properties of convolutional kernels, by characterizing functions in the RKHS as well as their norms. We begin with the one-layer case, which does not capture interactions between patches but highlights the role of pooling, before moving multiple layers, where interaction terms play an important role. Proofs are given in Appendix E.

We begin by considering the case of a single convolutional layer, which can already help us illustrate the role of patches and pooling. Here, the kernel is given by

with Φ(x)[u]=φ(xu)\Phi(x)[u]=\varphi(x_{u}), where we use the shorthand xu=Px[u]x_{u}=Px[u] for the patch at position uu. We now characterize the RKHS of K1K_{1}, showing that it consists of additive models of functions in H{\mathcal{H}} defined on patches, with spatial regularities among the terms, induced by the pooling operator AA. (Notation: A∗A^{*} and A†A^{\dagger} denote the adjoint and pseudo-inverse of an operator AA, respectively.)

The RKHS of K1K_{1} consists of functions f(x)=⟨G,Φ(x)⟩L2(Ω,H)=∑u∈ΩG[u](xu)f(x)=\langle G,\Phi(x)\rangle_{L^{2}(\Omega,{\mathcal{H}})}=\sum_{u\in\Omega}G[u](x_{u}), with G∈Range(A∗)G\in\text{Range}(A^{*}), and with RKHS norm

Note that if A∗A^{*} is not invertible (for instance in the presence of downsampling), the constraint G∈Range⁡(A∗)G\in\operatorname{Range}(A^{*}) is active and A†∗A^{\dagger*} is its pseudo-inverse. In the extreme case of global average pooling, we have A=(1,…,1)⊗Id:L2(Ω,H)→HA=(1,\ldots,1)\otimes Id:L^{2}(\Omega,{\mathcal{H}})\to{\mathcal{H}}, so that G∈Range⁡(A∗)G\in\operatorname{Range}(A^{*}) is equivalent to G[u]=gG[u]=g for all uu, for some fixed g∈Hg\in{\mathcal{H}}. In this case, the penalty in \eqrefeq:onelayernorm\eqref{eq:one_layer_norm} is simply the squared RKHS norm ∥g∥H2\|g\|_{{\mathcal{H}}}^{2}.

2 The Multi-Layer Case

We now study the case of convolutional kernels with more than one convolutional layer. While the patch kernels used at higher layers are typically similar to the ones from the first layer, we show empirically on Cifar10 that they may be replaced by simple polynomial kernels with little loss in accuracy. We then proceed by studying the RKHS of such simplified models, highlighting the role of depth for capturing interactions between different patches via kernel tensor products.

Table 1 shows the performance of a given 2-layer convolutional kernel architecture, with different choices of patch kernels κ1\kappa_{1} and κ2\kappa_{2}. The reference model uses exponential kernels in both layers, following the construction in Mairal (2016). We find that replacing the second layer kernel by a simple polynomial kernel of degree 3, κ2(u)=u3\kappa_{2}(u)=u^{3}, leads to roughly the same test accuracy. By changing κ2\kappa_{2} to κ2(u)=u2\kappa_{2}(u)=u^{2}, the test accuracy is only about 1% lower, while doing the same for the first layer decreases it by about 3%. The shallow kernel with a single non-linear convolutional layer (shown in the last line of Table 1) performs significantly worse. This suggests that the approximation properties described in Section 3.1 may not be sufficient for this task, while even a simple polynomial kernel of order 2 at the second layer may substantially improve things by capturing interactions, in a way that we describe below.

Let Φ(x)=(φ1(xu)⊗φ1(xv))u,v∈Ω∈L2(Ω2,H⊗H)\Phi(x)=(\varphi_{1}(x_{u})\otimes\varphi_{1}(x_{v}))_{u,v\in\Omega}\in L^{2}(\Omega^{2},{\mathcal{H}}\otimes{\mathcal{H}}). The RKHS of K2K_{2} when k2(z,z′)=(⟨z,z′⟩)2k_{2}(z,z^{\prime})=(\langle z,z^{\prime}\rangle)^{2} consists of functions of the form

where Gpq∈L2(Ω2,H⊗H)G_{pq}\in L^{2}(\Omega^{2},{\mathcal{H}}\otimes{\mathcal{H}}) obeys the constraints Gpq∈Range⁡(Epq)G_{pq}\in\operatorname{Range}(E_{pq}) and diag⁡((LpA1⊗LqA1)†∗Gpq)∈Range⁡(A2∗)\operatorname{diag}((L_{p}A_{1}\otimes L_{q}A_{1})^{\dagger*}G_{pq})\in\operatorname{Range}(A_{2}^{*}). Here, Epq:L2(Ω1)→L2(Ω2)E_{pq}:L^{2}(\Omega_{1})\to L^{2}(\Omega^{2}) is a linear operator given by

The squared RKHS norm ∥f∥HK22\|f\|_{{\mathcal{H}}_{K_{2}}}^{2} is then equal to the minimum over such decompositions of the quantity

As discussed in the one-layer case, the inverses should be replaced by pseudo-inverses if needed, e.g., when using downsampling. In particular, if A2∗A_{2}^{*} is singular, the second constraint plays a similar role to the one-layer case. In order to understand the first constraint, we show in Figure 2 the outputs of EpqxE_{pq}x for Dirac delta signals x[v]=δu[v]x[v]=\delta_{u}[v]. We can see that if the pooling filter h1h_{1} has a small support of size mm, then Gpq[u−p,v−q]G_{pq}[u-p,v-q] must be zero when ∣u−v∣>m|u-v|>m, which highlights that the functions in GpqG_{pq} may only capture interactions between pairs of patches where the (signed) distance between the first and the second is close to p−qp-q.

where F2=F⊗F\mathcal{F}_{2}=\mathcal{F}\otimes\mathcal{F} is the 2D discrete Fourier transform. Thus, this penalizes the variations of gz,z′g_{z,z^{\prime}} in both dimensions, encouraging the interaction functions to not rely too strongly on the specific positions of the two patches. This regularization is stronger when the spatial bandwidth of h1h_{1} is large, since this leads to a more localized filter in the frequency domain, with stronger penalties on high frequencies. In addition to this 2D smoothness, the penalty in Proposition 2 also encourages smoothness along the p−qp-q diagonal of this resulting 2D image using the pooling operator A2A_{2}. This has a similar behavior to the one-layer case, where the penalty prevents the functions from relying too much on the absolute position of the patches. Since A2A_{2} typically has a larger bandwidth than A1A_{1}, interaction functions Gpq[u,u+r]G_{pq}[u,u+r] are allowed to vary with rr more rapidly than with uu. The regularity of the resulting “smoothed” interaction terms as a function of the input patches is controlled by the RKHS norm of the tensor product kernel k1⊗k1k_{1}\otimes k_{1} as described in Appendix A.2.

When using a polynomial kernel k2(z,z′)=(⟨z,z′⟩)αk_{2}(z,z^{\prime})=(\langle z,z^{\prime}\rangle)^{\alpha} with α>2\alpha>2, we obtain a similar picture as above, with higher-order interaction terms. For example, if α=3\alpha=3, the RKHS contains functions with interaction terms of the form Gpqr[u,v,w](xu,xv,xw)G_{pqr}[u,v,w](x_{u},x_{v},x_{w}), with a penalty

where A1c=LcA1A_{1c}=L_{c}A_{1}. Similarly to the quadratic case, the first-layer pooling operator encourages smoothness with respect to relative positions between patches, while the second-layer pooling penalizes dependence on the global location. One may extend this further to higher orders to capture more complex interactions, and our experiments suggest that a two-layer kernel of this form with a degree-4 polynomial at the second layer may achieve state-of-the-art accuracy for kernel methods on Cifar10 (see Table 2). We note that such fixed-order choices for κ2\kappa_{2} lead to convolutional kernels that lower-bound richer kernels with, e.g., an exponential kernel at the second layer, in the Loewner order on positive-definite kernels. This imples in particular that the RKHS of these “richer” kernels also contains the functions described above. For more than two layers with polynomial kernels, one similarly obtains higher-order interactions, but with different regularization properties (see Appendix D).

Generalization Properties

In this section, we study generalization properties of the convolutional kernels studied in Section 3, and show improved sample complexity guarantees for architectures with pooling and small patches when the problem exhibits certain invariance properties.

As discussed in Section 3, the RKHS of 1-layer CKNs consists of sums of functions that are localized on patches, each belonging to the RKHS H{\mathcal{H}} of the patch kernel k1k_{1}. The next result illustrates the benefits of pooling when f∗f^{*} is translation invariant.

When using two layers with polynomial kernels at the second layer, we saw in Section 3 that the RKHS of CKNs consists of additive models of interaction terms of the order of the polynomial kernel used. The next proposition illustrates how pooling filters and patch sizes at the second layer may affect generalization on a simple target function consisting of order-2 interactions.

As an example, consider f∗(x)=∑u,vg(xu,xv)f^{*}(x)=\sum_{u,v}g(x_{u},x_{v}) for g∈H⊗Hg\in{\mathcal{H}}\otimes{\mathcal{H}} of minimal norm. The following table illustrates the obtained generalization bounds R(f^n)−R(f∗)R(\hat{f}_{n})-R(f^{*}) for KRR with various two-layer architectures (δ\delta: Dirac filter; 1\mathbf{1}: global average pooling):

Numerical Experiments

In this section, we provide additional experiments illustrating numerical properties of the convolutional kernels considered in this paper. We focus here on the Cifar10 dataset, and on CKN architectures based on the exponential kernel. Additional results are given in Appendix B.

We consider classification on Cifar10 dataset, which consists of 50k training images and 10k test images with 10 different output categories. We pre-process the images using a whitening/ZCA step at the patch level, which is commonly used for such kernels on images (Mairal, 2016; Shankar et al., 2020; Thiry et al., 2021). This may help reduce the effective dimensionality of patches, and better align the dominant eigen directions to the target function, a property which may help kernel methods (Ghorbani et al., 2020). Our convolutional kernel evaluation code is written in C++ and leverages the Eigen library for hardware-accelerated numerical computations. The computation of kernel matrices is distributed on up to 1000 cores on a cluster consisting of Intel Xeon processors. Computing the full Cifar10 kernel matrix typically takes around 10 hours when running on all 1000 cores. Our results use kernel ridge regression in a one-versus-all approach, where each class uses labels 0.90.9 for the correct label and −0.1-0.1 for the other labels. We report the test accuracy for a fixed regularization parameter λ=10−8\lambda=10^{-8} (we note that the performance typically remains the same for smaller values of λ\lambda). The exponential kernel always refers to κ(u)=e1σ2(u−1)\kappa(u)=e^{\frac{1}{\sigma^{2}}(u-1)} with σ=0.6\sigma=0.6. Code is available at https://github.com/albietz/ckn_kernel.

Table 2 shows test accuracies for different architectures compared to Table 1, including 3-layer models and 2-layer models with larger patches. In both cases, the full models with exponential kernels outperform the 2-layer architecture of Table 1, and provide comparable accuracy to the Myrtle10 kernel of Shankar et al. (2020), with an arguably simpler architecture. We also see that using degree-3 or 4 polynomial kernels at the second second layer of the two-layer model essentially provides the same performance to the exponential kernel, and that degree-2 at the second and third layer of the 3-layer model only results in a 0.3% accuracy drop. The two-layer model with degree-2 at the second layer loses about 1% accuracy, suggesting that certain Cifar10 images may require capturing interactions between at least 3 different patches in the image for good classification, though even with only second-order interactions, these models significantly outperform single-layer models. While these results are encouraging, computing such kernels is prohibitively costly, and we found that applying the Nyström approach of Mairal (2016) to these kernels with more layers or larger patches requires larger models than for the architecture of Table 1 for a similar accuracy. Figure 3(left) shows learning curves for different architectures, with slightly better convergence rates for more expressive models involving higher-order kernels or more layers; this suggests that their approximation properties may be better suited for these datasets.

Figure 3 shows the spectral decays of the empirical kernel matrix on 1000 Cifar images, which may help assess the “effective dimensionality” of the data, and are related to generalization properties (Caponnetto & De Vito, 2007). While multi-layer architectures with pooling seem to provide comparable decays for various depths, removing pooling leads to significantly slower decays, and hence much larger RKHSs. In particular, the “strided pooling” architecture (i.e., with Dirac pooling filters and downsampling) shown in Figure 3(right), which resembles the kernel considered in (Scetbon & Harchaoui, 2020), obtains less than 40% accuracy on 10k examples. This suggests that the regularization properties induced by pooling, studied in Section 3, are crucial for efficient learning on these problems, as shown in Section 4. Appendix B provides more empirics on different pooling configurations.

Discussion and Concluding Remarks

In this paper, we studied approximation and generalization properties of convolutional kernels, showing how multi-layer models with convolutional architectures may effectively break the curse of dimensionality on problems where the input consists of high-dimensional natural signals, by modeling localized functions on patches and interactions thereof. We also show how pooling induces additional smoothness constraints on how interaction terms may or may not vary with global and relative spatial locations. An important question for future work is how optimization of deep convolutional networks may further improve approximation properties compared to what is captured by the kernel regime presented here, for instance by selecting well-chosen convolution filters at the first layer, or interaction patterns in subsequent layers, perhaps in a hierarchical manner.

The author would like to thank Francis Bach, Alessandro Rudi, Joan Bruna, and Julien Mairal for helpful discussions.

References

Appendix A Further Background

This section provides further background on the problem of approximation of functions defined on signals, as well as on the kernels considered in the paper. We begin by introducing and motivating the problem of learning functions defined on signals such as images, which captures tasks such as image classification where deep convolutional networks are predominant. We then recall properties of dot-product kernels and kernel tensor products, which are key to our study of approximation.

using samples from the data distribution ρ\rho. If f∗f^{*} is only assumed to be Lipschitz, learning requires a number of samples that scales exponentially in the dimension (see, e.g., von Luxburg & Bousquet (2004); Wainwright (2019)), a phenomenon known as the curse of dimensionality. In the case of natural signals, the dimension d=p∣Ω∣d=p|\Omega| scales with the size of the domain ∣Ω∣|\Omega| (e.g., the number of pixels), which is typically very large and thus makes this intractable. One common way to alleviate this is to assume that f∗f^{*} is smooth, however the order of smoothness typically needs to be of the order of the dimension in order for the problem to become tractable, which is a very strong assumption here when dd is very large. This highlights the need for more structured assumptions on f∗f^{*} which may help overcome the curse of dimensionality.

where kk is a “simple” kernel such as a dot-product kernel, as discussed in Appendix C. In contrast, if KK is a dot-product kernel on the entire image, corresponding to an infinite-width limit of a fully-connected network, then approximation is more difficult and is generally cursed by the full dimension (see Appendix C). While some models of wide fully-connected networks provide some adaptivity to low-dimensional structures such as the variables in a patch (Bach, 2017a), no tractable algorithms are currently known to achieve such behavior provably, and it is reasonable to instead encode such prior information in a convolutional architecture.

Modeling interactions between elements of a system at different scales, possibly hierarchically, is important in physics and complex systems, in order to efficiently handle systems with large numbers of variables (Beylkin & Mohlenkamp, 2002; Hackbusch & Kühn, 2009). As an example, one may consider target functions f∗(x)f^{*}(x) that consist of interaction functions of the form g(xp,xq)g(x_{p},x_{q}), where p,qp,q denote locations of the corresponding patches, and higher-order interactions may also be considered. In the context of image recognition, while functions of a single patch may capture local texture information such as edges or color, such an interaction function may also respond to specific spatial configurations of relevant patches, which could perhaps help identify properties related to the “shape” of an object, for instance. If such functions gg are too general, then the curse of dimensionality may kick in again when one considers more than a handful of patches. Certain idealized models of approximation may model such interactions more efficiently through hierarchical compositions (e.g., Poggio et al. (2017)) or tensor decompositions (Cohen & Shashua, 2016; 2017), though no tractable algorithms are known to find such models. In this work, we tackle this in a tractable way using multi-layer convolutional kernels. We show that they can model interactions through kernel tensor products, which define functional spaces that are typically much smaller and more structured than for a generic kernel on the full vector (xp,xq)(x_{p},x_{q}).

A.2 Dot-Product Kernels and their Tensor Products

In this section, we review some properties of dot-product kernels, their induced RKHS and regularization properties. We then recall the notion of tensor product of kernels, which allows us to describe the RKHS of products of kernels in terms of that of individual kernels.

For more than one layer, the convolutional kernels we study in Section 3 can be expressed in terms of products of kernels on patches, of the form

is a feature map for KK, and the corresponding RKHS, denoted H⊗m=H⊗⋯⊗H{\mathcal{H}}^{\otimes m}={\mathcal{H}}\otimes\cdots\otimes{\mathcal{H}}, contains all functions

Appendix B Additional Experiments

In this section, we provide additional experiments to those presented in Section 5, using different patch kernels, patch sizes, pooling filters, preprocessings, and datasets.

Table 3 provides more results on 3-layer architectures compared to Table 2, including different changes in the degrees of polynomial kernels at the second and third layer. In particular we see that the architecture with degree-2 kernels at both layers, which captures interactions of order 4, also outperforms the simpler ones using degree-4 kernels at either layer, suggesting that a deeper architecture may better model relevant interactions terms on this problem.

Table 4 shows the variations in test performance when changing the size of the second patches at the second layer. We see that intermediate sizes between 3x3 and 9x9 work best, but that performance degrades when using patches that are too large or too small. For very large patches, this may be due to the large variance in (10), or perhaps instability (Bietti & Mairal, 2019a). For ∣S2∣=|S_{2}|=1x1, note that while pooling after the first layer allows even 1x1 patches to capture interactions across different input image patches, these may be limited to short range interactions when the pooling filter is localized (see Proposition 2), which may limit the expressivity of the model.

In Table 5, we consider 2-layer convolutional kernels with a similar architecture to those considered in Table 1, but where we use arc-cosine kernels arising from ReLU activations instead of the exponential kernel used in Section 5, given by

Table 6 shows the accuracy for one-layer convolutional kernels with 6x6 patchesNote that in this case the ZCA/whitening step is applied on these larger 6x6 patches. and various pooling sizes, with a highest accuracy of 75.8% for a pooling size of 88. While this improves on the accuracy obtained with 3x3 patches (slightly above 74% for the architectures in Tables 1 and 3 with a single non-linear kernel at the first layer), these accuracies remain much lower than those achieved by two-layer architectures with even quadratic kernels at the second layer. While using larger patches may allow capturing patterns that are less localized compared to small 3x3 patches, the neighborhoods that they model need to remain small in order to avoid the curse of dimensionality when using dot-product kernels, as discussed in Section 2. Instead, the multi-layer architecture may model information at larger scales with a much milder dependence on the size of the neighborhood, thanks to the structure imposed by tensor product kernels (see Section A.2) and the additional regularities induced by pooling.

We also found that larger patches at the first layer may hurt performance in multi-layer models: when considering the architecture of Table 1 with exponential kernels, using 5x5 patches instead of 3x3 at the first layer yields an accuracy of 79.6% instead of 80.5% on Cifar10 when training on the same 10k images. This again reflects the benefits of using small patches at the first layer for allowing better approximation on small neighborhoods, while modeling larger scales using interaction models according to the structure of the architecture. We note nevertheless that for standard deep networks, larger patches are often used at the first layer (e.g., He et al., 2016), as the feature selection capabilities of SGD may alleviate the dependence on dimension, e.g., by finding Gabor-like filters.

Table 7 shows the differences in performance between two or three layer architectures considered in Table 2, when Gaussian pooling filters are replaced by average pooling filters. For both architectures considered, average pooling leads to a significant performance drop. This suggests that one may need deeper architectures in order for such average pooling filters to work well, as in (Shankar et al., 2020), either with multiple 3x3 convolutional layers before applying pooling, or by applying multiple average pooling layers in a row as in certain Myrtle kernels. Note that iterating multiple average pooling layers in a row is equivalent to using a larger and more smooth pooling filter (with one more order of smoothness at each layer), which may then be more comparable to our Gaussian pooling filters.

We now consider the SVHN dataset, which consists of 32x32 images of digits from Google Street View images, 73 257 for training and 26 032 for testing. Due to the larger dataset size, we only consider the kernel approximation approach of Mairal (2016) based on the Nyström method, which projects the patch kernel feature maps at each layer to finite-dimensional subspaces generated by a set of anchor points (playing the role of convolutional filters), themselves computed via a K-means clustering of patches.We use the PyTorch implementation available at https://github.com/claying/CKN-Pytorch-image. We train one-versus-all classifiers on the resulting finite-dimensional representations using regularized ERM with the squared hinge loss, and simply report the best test accuracy over a logarithmic grid of choices for the regularization parameter, ignoring model selection issues in order to assess approximation properties. We use the same ZCA preprocessing as on Cifar10 and the same architecture as in Table 1, with a relatively small number of filters (256 at the first layer, 4096 at the second layer, leading to representations of dimension 65 536), noting that the accuracy can further improve when increasing this number. Our observations are similar to those for the Cifar10 dataset: using a degree-3 polynomial kernel at the second layer reaches very similar accuracy to the exponential kernel; using a degree-2 polynomial leads to a slight drop, but a smaller drop than when making this same change at the first layer; using a linear kernel at the second layer leads to a much larger drop. This again highlights the importance of using non-linear kernels on top of the first layer in order to capture interactions at larger scales than the scale of a single patch.

Recall that our pre-processing is based on a patch-level whitening or ZCA on each image, following Mairal (2016). In practice, this is achieved by whitening extracted patches from each image, and reconstructing the image from whitened patches via averaging. In contrast, other approaches use global whitening of the entire image Lee et al. (2020); Shankar et al. (2020). For the 2-layer model shown in Table 2 with 5x5 patches at the second layer, we found global ZCA to provide significantly worse performance, with a drop from 88.3% to about 80%.

The work Shankar et al. (2020) introduces Myrtle kernels but also consider similar architectures for usual CNNs with finite-width, trained with stochastic gradient descent. Obtaining competitive architectures for the finite-width case is not the goal of our work, which focuses on good architectures for the kernel setup, yet it remains interesting to consider this question. In the case of Shankar et al. (2020), training the finite-width networks yields better accuracy compared to their “infinite-width” kernel counterparts, a commonly observed phenomenon which may be due to better “adaptivity” of optimization algorithms compared to kernel methods, which have a fixed representation and thus may not learn representations adapted to the data (see, e.g., Allen-Zhu & Li, 2020; Bach, 2017a; Chizat et al., 2019). Nevertheless, we found that for the two-layer architecture considered in Table 1, which has many fewer layers compared to the Myrtle architectures of Shankar et al. (2020), using a finite-width ReLU network yields poorer performance compared to the kernel (around 83% at best, compared to 87.9%). This may suggest that for convolutional networks, deeper networks may have additional advantages when using optimization algorithms, in terms of adapting to possibly relevant structure of the problem, such as hierarchical representations (see, e.g., Allen-Zhu & Li (2020); Chen et al. (2020); Poggio et al. (2017) for theoretical justifications of the benefits of depth in non-kernel regimes).

Appendix C Complexity of Spatially Localized Functions

If we define the kernel Ku(x,x′)=k(xu,xu′)K_{u}(x,x^{\prime})=k(x_{u},x^{\prime}_{u}), where kk is a dot-product kernel arising from a one-hidden layer network with positively-homogeneous activation such as the ReLU, and further assume patches to be bounded and g∗g^{*} to be bounded, then the uniform approximation error bound of Bach (2017a, Proposition 6) together with a simple O(1/n)O(1/\sqrt{n}) Rademacher complexity bound on estimation error shows that we may achieve a generalization bound with a rate that only depends on the patch dimension p∣S∣p|S| rather than p∣Ω∣p|\Omega| in this setup (i.e., a sample complexity that is exponential in p∣S∣p|S|, which is much smaller than p∣Ω∣p|\Omega|).

If we consider the kernel K(x,x′)=∑u∈Ωk(xu,xu′)K(x,x^{\prime})=\sum_{u\in\Omega}k(x_{u},x^{\prime}_{u}), the RKHS contains all functions in the RKHS of KuK_{u} for all u∈Ωu\in\Omega, with the same norm (this may be seen as an application of Theorem 6 with a feature map given by concatenating the kernel maps of each KuK_{u}), so that we may achieve the same approximation error as above, and thus a similar generalization bound that is not cursed by dimension. This kernel also allows us to obtain similar generalization guarantees when f∗f^{*} consists of linear combinations of such spatially localized functions on different patches within the image.

In contrast, when using a similar dot-product kernel on the full signal, corresponding to using a fully-connected network in a kernel regime, one may construct functions f∗(x)=g∗(xu)f^{*}(x)=g^{*}(x_{u}) with g∗g^{*} Lipschitz where an RKHS norm that is exponentially large in the (full) dimension p∣Ω∣p|\Omega| is needed for a small approximation error (see Bach, 2017a, Appendix D.5).

Related to this, Malach & Shalev-Shwartz (2021) show a separation in the different setting of learning certain parity functions on the hypercube using gradient methods; their upper bound for convolutional networks is based on a similar kernel regime as above. We note that kernels that exploit such a localized structure have also been considered in the context of structured prediction for improved statistical guarantees (Ciliberto et al., 2019).

Appendix D Extensions to More Layers

In this section, we study the RKHS for convolutional kernels with more than 2 convolutional layers, by considering the simple example of a 3-layer convolutional kernel K3K_{3} defined by the feature map

We may then describe the RKHS as follows.

The RKHS of K3K_{3} when k2k_{2} and k3k_{3} are quadratic kernels (⟨⋅,⋅⟩)2(\langle\cdot,\cdot\rangle)^{2} consists of functions of the form

where Gα∈L2(Ω4,H⊗4)G_{\alpha}\in L^{2}(\Omega^{4},{\mathcal{H}}^{\otimes 4}) obeys the constraint

where the linear operator Eα:L2(Ω3)→L2(Ω4)E_{\alpha}:L^{2}(\Omega_{3})\to L^{2}(\Omega^{4}) for α=(p,q,r,p′,q′,r′)\alpha=(p,q,r,p^{\prime},q^{\prime},r^{\prime}) (with p,p′∈S3p,p^{\prime}\in S_{3} and q,r,q′,r′∈S2q,r,q^{\prime},r^{\prime}\in S_{2}) is defined by

The operators A1,αA_{1,\alpha} and A2,αA_{2,\alpha} denote:

The squared RKHS norm ∥f∥HK32\|f\|_{{\mathcal{H}}_{K_{3}}}^{2} is then equal to the minimum over decompositions \eqrefeq:fdecompthreelayers\eqref{eq:f_decomp_three_layers} of the quantity

Appendix E Proofs

We recall the following result about reproducing kernel Hilbert spaces, which characterizes the RKHS of kernels defined by explicit Hilbert space features maps (see, e.g., Saitoh, 1997, §2.1).

Let HH be some Hilbert space, ψ:X→H\psi:{\mathcal{X}}\to H a feature map, and K(x,x′)=⟨ψ(x),ψ(x′)⟩HK(x,x^{\prime})=\langle\psi(x),\psi(x^{\prime})\rangle_{H} a kernel on X{\mathcal{X}}. The RKHS H{\mathcal{H}} of KK consists of functions f=⟨g,ψ(⋅)⟩Hf=\langle g,\psi(\cdot)\rangle_{H}, with norm

We also state here the generalization bound for kernel ridge regression used in Section 4, adapted from Bach (2021, Proposition 7.1).

for i.i.d. data (xi,yi)∼ρ(x_{i},y_{i})\sim\rho, i=1,…,ni=1,\ldots,n. Let

Under the conditions of the theorem, we may apply (Bach, 2021, Proposition 7.1), which states that for λ≤1\lambda\leq 1 and n≥5λ(1+log⁡(1/λ))n\geq\frac{5}{\lambda}(1+\log(1/\lambda)), we have

From Theorem 6, the RKHS contains functions of the form

with RKHS norm equal to the minimum of ∥F∥L2(Ω1,H)\|F\|_{L^{2}(\Omega_{1},{\mathcal{H}})} over such decompositions.

We may alternatively write f(x)=⟨G,Φ(x)⟩L2(Ω,H)f(x)=\langle G,\Phi(x)\rangle_{L^{2}(\Omega,{\mathcal{H}})} with G=A∗FG=A^{*}F. The mapping from FF to GG is one-to-one if G∈Range⁡(A∗)G\in\operatorname{Range}(A^{*}). Then, we obtain that equivalently, the RKHS contains functions of this form, with G∈Range⁡(A∗)G\in\operatorname{Range}(A^{*}), and with RKHS norm equal to the minimum of ∥A†∗G∥L2(Ω1,H)\|A^{\dagger*}G\|_{L^{2}(\Omega_{1},{\mathcal{H}})} over such decompositions. ∎

From Theorem 6, the RKHS contains functions of the form

with RKHS norm equal to the minimum of ∥F∥L2(Ω2,H2)\|F\|_{L^{2}(\Omega_{2},{\mathcal{H}}_{2})} over such decompositions. Here, Φ1(x)∈L2(Ω,H)\Phi_{1}(x)\in L^{2}(\Omega,{\mathcal{H}}) is given by Φ1(x)[u]=φ1(xu)\Phi_{1}(x)[u]=\varphi_{1}(x_{u}), so that Φ(x)\Phi(x) in the statement is given by Φ(x)=Φ1(x)⊗Φ1(x)\Phi(x)=\Phi_{1}(x)\otimes\Phi_{1}(x). We also have that H2=(H⊗H)∣S2∣×∣S2∣{\mathcal{H}}_{2}=({\mathcal{H}}\otimes{\mathcal{H}})^{|S_{2}|\times|S_{2}|}, so that we may write F=(Fpq)p,q∈S2F=(F_{pq})_{p,q\in S_{2}} with Fpq∈L2(Ω2,H⊗H)F_{pq}\in L^{2}(\Omega_{2},{\mathcal{H}}\otimes{\mathcal{H}}).

For p,q∈S2p,q\in S_{2}, denoting by LcL_{c} the translation operator Lcx[u]=x[u−c]L_{c}x[u]=x[u-c], we have

We may then write this as ⟨Gpq,Φ(x)⟩L2(Ω2,H⊗H)\langle G_{pq},\Phi(x)\rangle_{L^{2}(\Omega^{2},{\mathcal{H}}\otimes{\mathcal{H}})} with

and the mapping between Fpq∈L2(Ω2,H⊗H)F_{pq}\in L^{2}(\Omega_{2},{\mathcal{H}}\otimes{\mathcal{H}}) and Gpq∈L2(Ω2,H⊗H)G_{pq}\in L^{2}(\Omega^{2},{\mathcal{H}}\otimes{\mathcal{H}}) is one-to-one if Gpq∈Range⁡((LpA1⊗LqA1)∗)G_{pq}\in\operatorname{Range}((L_{p}A_{1}\otimes L_{q}A_{1})^{*}), and diag⁡((LpA1⊗LqA1)†∗Gpq)∈Range⁡(A2∗)\operatorname{diag}((L_{p}A_{1}\otimes L_{q}A_{1})^{\dagger*}G_{pq})\in\operatorname{Range}(A_{2}^{*}). We may then equivalently write the RKHS norm as the minimum over GpqG_{pq} satisfying such constraints for all p,q∈S2p,q\in S_{2}, of the quantity

Let Φ(x)=(φ1(xu))u∈L2(Ω,H)\Phi(x)=(\varphi_{1}(x_{u}))_{u}\in L^{2}(\Omega,{\mathcal{H}}), so that we may write

for some G∈L2(Ω4,H⊗4)G\in L^{2}(\Omega^{4},{\mathcal{H}}^{\otimes 4}).

From Theorem 6, the RKHS contains functions of the form

with RKHS norm equal to the minimum of ∥F∥L2(Ω3,H3)\|F\|_{L^{2}(\Omega_{3},{\mathcal{H}}_{3})} over such decompositions. Here, Φ2(x)∈L2(Ω1,H2)=L2(Ω1,(H⊗H)∣S2∣×∣S2∣)\Phi_{2}(x)\in L^{2}(\Omega_{1},{\mathcal{H}}_{2})=L^{2}(\Omega_{1},({\mathcal{H}}\otimes{\mathcal{H}})^{|S_{2}|\times|S_{2}|}) is given as in the proof of Proposition 2, by

for q,r∈S2q,r\in S_{2}. A patch P3A2Φ2(x)[u]P_{3}A_{2}\Phi_{2}(x)[u] is then given by

Applying the quadratic feature map given by φ3(z)=z⊗z∈(H⊗4)(∣S3∣×∣S2∣×∣S2∣)2\varphi_{3}(z)=z\otimes z\in({\mathcal{H}}^{\otimes 4})^{(|S_{3}|\times|S_{2}|\times|S_{2}|)^{2}} for z∈(H⊗H)∣S3∣×∣S2∣×∣S2∣z\in({\mathcal{H}}\otimes{\mathcal{H}})^{|S_{3}|\times|S_{2}|\times|S_{2}|}, we obtain for α=(p,q,r,p′,q′,r′)∈(S3×S2×S2)2\alpha=(p,q,r,p^{\prime},q^{\prime},r^{\prime})\in(S_{3}\times S_{2}\times S_{2})^{2},

Now, one can check that we have the following relation:

Since H3=(H⊗4)(∣S3∣×∣S2∣×∣S2∣)2{\mathcal{H}}_{3}=({\mathcal{H}}^{\otimes 4})^{(|S_{3}|\times|S_{2}|\times|S_{2}|)^{2}}, we may write F=(Fα)α∈(S3×S2×S2)2F=(F_{\alpha})_{\alpha\in(S_{3}\times S_{2}\times S_{2})^{2}}, with each Fα∈L2(Ω3,H⊗4)F_{\alpha}\in L^{2}(\Omega_{3},{\mathcal{H}}^{\otimes 4}). We then have

We may write this as ⟨Gα,Φ(x)⊗4⟩L2(Ω4,H⊗4)\langle G_{\alpha},\Phi(x)^{\otimes 4}\rangle_{L^{2}(\Omega^{4},{\mathcal{H}}^{\otimes 4})}, with

The mapping from FαF_{\alpha} to GαG_{\alpha} is bijective if GαG_{\alpha} is constrained to the lie in the range of the operator EαE_{\alpha}. If GαG_{\alpha} satisfies this constraint, we may write

Then, the resulting penalty on GαG_{\alpha} is as desired.

E.4 Proof of Proposition 3 (generalization for one-layer CKN)

It remains to verify that ∥f∗∥HK1=∣Ω∣∥g∥H\|f^{*}\|_{{\mathcal{H}}_{K_{1}}}=\sqrt{|\Omega|}\|g\|_{\mathcal{H}}. Note that if we denote G=(g)u∈Ω∈L2(Ω,H)G=(g)_{u\in\Omega}\in L^{2}(\Omega,{\mathcal{H}}), then we have A∗G=GA^{*}G=G, since ∑vh[v−u]G[v]=(∑vh[v−u])g=g\sum_{v}h[v-u]G[v]=(\sum_{v}h[v-u])g=g. This implies that A∗†G=GA^{*\dagger}G=G, regardless of which pooling filter is used. Then we have, by (5) that ∥f∗∥2≤∣Ω∣∥g∥2\|f^{*}\|^{2}\leq|\Omega|\|g\|^{2}. Further, since gg is of minimal norm, no other G∈L2(Ω,H)G\in L^{2}(\Omega,{\mathcal{H}}) may lead to a smaller norm, so that we can conclude ∥f∗∥2=∣Ω∣∥g∥2\|f^{*}\|^{2}=|\Omega|\|g\|^{2}. ∎

while the sum of coefficients bounded by ϵ\epsilon is upper bounded by ∣Ω∣∣S2∣2|\Omega||S_{2}|^{2}. Overall, this yields

When using a Dirac filter at the first layer, we need ∣S2∣=∣Ω∣|S_{2}|=|\Omega| in order to capture log-range interaction terms, and we represent f∗f^{*} as a sum of GpG_{p} that are non-zero and equal to gg only on the pp-th diagonal, i.e., Gp[u,v]=gG_{p}[u,v]=g when v=u+pv=u+p, and zero otherwise. We then verify ∑p∑u,vGp[u,v](xu,xv)=∑p∑ug(xu,xu+p)=f∗(x)\sum_{p}\sum_{u,v}G_{p}[u,v](x_{u},x_{v})=\sum_{p}\sum_{u}g(x_{u},x_{u+p})=f^{*}(x). Then, the expression (7) on this decomposition yields ∣S2∣∣Ω∣∥g∥H⊗H2=∣Ω∣2∥g∥H⊗H2|S_{2}||\Omega|\|g\|^{2}_{{\mathcal{H}}\otimes{\mathcal{H}}}=|\Omega|^{2}\|g\|^{2}_{{\mathcal{H}}\otimes{\mathcal{H}}}, for any choice of pooling h2h_{2} such that ∥h2∥1=1\|h_{2}\|_{1}=1 (using similar arguments to the proof of Proposition 3). This is then equal to the squared norm of f∗f^{*} due to the minimality of ∥g∥H⊗H\|g\|_{{\mathcal{H}}\otimes{\mathcal{H}}}.

When using average pooling at the first layer and ∣S2∣=∣Ω∣|S_{2}|=|\Omega|, we may use a single term G[u,v]G[u,v] with all entries equal to gg, i.e., a decomposition f∗(x)=∑u,vG[u,v](xu,xv)f^{*}(x)=\sum_{u,v}G[u,v](x_{u},x_{v}). Using (7), we obtain an upper bound ∣Ω∣∥g∥2|\Omega|\|g\|^{2} on the squared norm. The same decomposition can be used when ∣S2∣=1|S_{2}|=1, leading to the same bound.

Appendix F Generalization Gains under Specific Data Models

In this section, we consider simple models of architectures and data distribution where we may quantify more precisely the improvements in sample complexity guarantees thanks to pooling.

We note that when patches are in high dimension, overlapping patches may become near-orthogonal, which could allow extensions of our arguments below to the case with overlap, yet this may require different tools similar to Mei et al. (2021). We leave these questions to future work.

where II denotes the modified Bessel function of the first kind. For arc-cosine kernels, it may be obtained by leveraging the random feature expansion of the kernel (Bach, 2017a). More generally, we also note that when the patch dimension dd is large, we have

since ωd−2ωd−1(1−t2)d−32\frac{\omega_{d-2}}{\omega_{d-1}}(1-t^{2})^{\frac{d-3}{2}} is a probability density that converges weakly to a Dirac mass at 0. In particular, for the exponential kernel with σ=0.6\sigma=0.6, we have κ(0)=e−1/σ2≈0.06\kappa(0)=e^{-1/\sigma^{2}}\approx 0.06. When learning a translation-invariant function, the bound in Prop. 3 then shows that global average pooling yields an improvement w.r.t. no pooling of order ∣Ω∣/(1+0.06∣Ω∣)|\Omega|/(1+0.06|\Omega|). Note that removing the constant component of κ\kappa, i.e., using the kernel κ(u)−κ(0)κ(1)−κ(0)\frac{\kappa(u)-\kappa(0)}{\kappa(1)-\kappa(0)}, may further improve this bound, leading to a denominator very close to 11 when dd is large, and hence an improvement in sample complexity of order ∣Ω∣|\Omega|. We also remark that the dependence on κ(0)\kappa(0) may be removed by using a finer generalization analysis beyond uniform convergence that leverages spectral properties of the kernel (see Section F.2).

If ∣{u,v,u′,v′}∣=4|\{u,v,u^{\prime},v^{\prime}\}|=4, we have

If u=vu=v and ∣{u,u′,v′}∣=3|\{u,u^{\prime},v^{\prime}\}|=3, then we have

F.2 Fast rates

In this section, we derive spectral decompositions of 1-layer CKN architectures with non-overlapping patches under the product of spheres distribution described in the previous section. This allows us to derive fast rates that depend on the complexity of the target functions on patches, and shows similar improvement factors to those derived in Section F.1, without the κ(0)\kappa(0) term, which in fact turns out to only be due to a single eigenspace, namely constant functions. We note that our derivation extends (Favero et al., 2021) to the case of generic pooling filters, and considers a different data distribution.

Before studying the 1-layer case, we remark that while it may seem natural to extend such decompositions to the 2-layer case using tensor products of spherical harmonics, as done by Scetbon & Harchaoui (2020) in the case without pooling, it appears that pooling may make it more challenging to find an eigenbasis since subspaces consisting of tensor products of spherical harmonics with fixed total degree are no longer left stable by the kernel.For instance, the term k1(xw,yu)k1(xw,yv)k_{1}(x_{w},y_{u})k_{1}(x_{w},y_{v}), which may only appear in the presence of pooling, maps the polynomial Yk(xu)Yk(xv)Y_{k}(x_{u})Y_{k}(x_{v}) of degree 2k2k to a polynomial μk2Yk(xw)2\mu_{k}^{2}Y_{k}(x_{w})^{2} which is not necessarily orthogonal to all spherical harmonics tensor products of degree smaller than 2k. We thus leave such a study to future work.

We begin by considering the following Mercer decomposition of the patch kernel

where Yk,jY_{k,j} for k≥0k\geq 0 and j=1,…,N(d,k)j=1,\ldots,N(d,k) are spherical harmonic polynomials of degree kk forming an orthonormal basis of L2(dτ)L^{2}(d\tau). The μk\mu_{k} here are Legendre/Gegenbauer coefficients of the function κ\kappa (see, e.g., Bach, 2017a; Smola et al., 2001).

Note that the 1-layer kernel may be written as

and let λw:=h⊛hˉ^[w]=∣h^[w]∣2\lambda_{w}:=\widehat{h\circledast\bar{h}}[w]=|\hat{h}[w]|^{2}. Note that when the filter is normalized s.t. ∥h∥1=1\|h\|_{1}=1, we have λ0=h^=1\lambda_{0}=\hat{h}=1. We will also use the Parseval identity ∥h∥22=(∑wλw)/∣Ω∣\|h\|_{2}^{2}=(\sum_{w}\lambda_{w})/|\Omega|. Using the inverse DFT, it holds

Note that when k=0k=0, we have Y0,1(xu)=1Y_{0,1}(x_{u})=1 for all uu, hence

since e0=∣Ω∣−1/2(1,…,1)e_{0}=|\Omega|^{-1/2}(1,\ldots,1) and ∑uew[u]=0\sum_{u}e_{w}[u]=0 for w>0w>0. Then, we may write

with ϕ0(x)=1\phi_{0}(x)=1, ϕw,k,j(x)=∑uew[u]Yk,j(xu)\phi_{w,k,j}(x)=\sum_{u}e_{w}[u]Y_{k,j}(x_{u}), and ϕ∗\phi^{*} denotes the complex conjugate of ϕ\phi.

It is then easy to check that the ϕ0\phi_{0} and ϕw,k,j\phi_{w,k,j} form an orthonormal basis of L2(dτ⊗∣Ω∣)L^{2}(d\tau^{\otimes|\Omega|}). We thus have obtained a Mercer decomposition of the kernel KhK_{h} w.r.t. the data distribution, so that its eigenvalues are also the eigenvalues of the covariance operator (Caponnetto & De Vito, 2007), and control generalization performance, typically through the degrees of freedom

where Σ\Sigma is the covariance operator, and (ξm)m(\xi_{m})_{m} is the collection of its eigenvalues. In particular, if we have a decay ξm≍m−α\xi_{m}\asymp m^{-\alpha} (with α<1\alpha<1), then we have N(λ)≤O(λ−1/α)\mathcal{N}(\lambda)\leq O(\lambda^{-1/\alpha}), which then leads to a fast rate of n−α/(α+1)n^{-\alpha/(\alpha+1)} on the excess risk when optimizing for λ\lambda in kernel ridge regression (Bach, 2021; Caponnetto & De Vito, 2007). In our case, the degrees of freedom takes the form

Given that the number of spatial frequencies ww is fixed, the asymptotic decay rate of eigenvalues associated to the ϕw,k,j\phi_{w,k,j} (that is, λwμk\lambda_{w}\mu_{k}, each with multiplicity N(d,k)N(d,k)) is the same as that of the eigenvalues associated to ϕw0,k,j\phi_{w_{0},k,j} for some fixed w0w_{0}, which in turn corresponds to the decay for the corresponding dot-product kernel on the sphere. For instance, if κ\kappa is the arc-cosine kernel, we have ξm≍m−α\xi_{m}\asymp m^{-\alpha} with α=d+2d−1\alpha=\frac{d+2}{d-1}, and more generally α=2sd−1\alpha=\frac{2s}{d-1} for a kernel resembling a Sobolev space with ss bounded derivatives. This then leads to a rate n2s/(2s+(d−1))n^{2s/(2s+(d-1))}, which only depends on the dimension dd of patches, rather than the full dimension d∣Ω∣d|\Omega|. One may also add more general assumption on the smoothness of the localized components of f∗f^{*} (such as gg in Proposition 3) in order to get rates that depend explicitly on the order of smoothness ss of such components (as in Caponnetto & De Vito, 2007).

The eigenvalue associated to ϕ0\phi_{0} plays a minor role as by itself as it only contributes at most τρ2/n\tau_{\rho}^{2}/n to the excess risk, which is negligible compared to the rest of the eigenvalues which lead to a slower n−α/(α+1)n^{-\alpha/(\alpha+1)} rate.

With no pooling (hh is a Dirac delta), we have λw=∣h^[w]∣2=1\lambda_{w}=|\hat{h}[w]|^{2}=1 for all ww. We then have

where we defined Nκ(λ):=∑k≥1N(d,k)μk/(λ+μk)\mathcal{N}_{\kappa}(\lambda):=\sum_{k\geq 1}N(d,k)\mu_{k}/(\lambda+\mu_{k}).

With global pooling (h=1/∣Ω∣h=1/|\Omega| is constant), we have λ0=1\lambda_{0}=1, and λw=0\lambda_{w}=0 for w>0w>0. This yields

This then yields an improvement by a factor ∣Ω∣|\Omega| in sample complexity guarantees compared to the scenario above with no pooling, namely the dominant term in the excess risk bound will be C(1/n)α/(α+1)C(1/n)^{\alpha/(\alpha+1)} compared to C(∣Ω∣/n)α/(α+1)C(|\Omega|/n)^{\alpha/(\alpha+1)}.

For more general pooling, one may exploit specific decays of λw\lambda_{w} to obtain finer bounds. We may also obtain the following bound by Jensen’s inequality

where we used that λˉ=(∑wλw)/∣Ω∣=∥h∥22\bar{\lambda}=(\sum_{w}\lambda_{w})/|\Omega|=\|h\|_{2}^{2}, by Parseval’s identity. When Nκ(λ)≤Cκλ−1/α\mathcal{N}_{\kappa}(\lambda)\leq C_{\kappa}\lambda^{-1/\alpha}, we get a bound

instead of a bound N(λ)≤1+C∣Ω∣λ−1/α\mathcal{N}(\lambda)\leq 1+C|\Omega|\lambda^{-1/\alpha} for the case of no pooling, that is, the improvement is again controlled by ∥h∥22\|h\|_{2}^{2}, which goes from 1/∣Ω∣1/|\Omega| for global pooling, to 11 for no pooling.

The following result provides an example of a generalization bound for invariant target functions, which illustrates that there is no curse of dimensionality in the rate if the patch dimension is much smaller than the full dimension (i.e., d≪d∣Ω∣d\ll d|\Omega|), as well as the benefits of pooling.

(capacity condition) Nκ(λ)≤Cκλ−1/α\mathcal{N}_{\kappa}(\lambda)\leq C_{\kappa}\lambda^{-1/\alpha},

(source condition) g∗=Tκrg0g^{*}=T_{\kappa}^{r}g_{0} and ∥g0∥L2(dτ)≤C∗\|g_{0}\|_{L^{2}(d\tau)}\leq C_{*}, where TκT_{\kappa} is the integral operator of the kernel κ\kappa on L2(dτ)L^{2}(d\tau), with r>α−12αr>\frac{\alpha-1}{2\alpha}.

Then, kernel ridge regression with the one-layer CKN kernel KhK_{h} with pooling filter hh (with ∥h∥1=1\|h\|_{1}=1) satisfies, for nn large enough,

where CC is independent of ∣Ω∣|\Omega| and hh. For global pooling, the factor ∥h∥22/α=∣Ω∣−1/α\|h\|_{2}^{2/\alpha}=|\Omega|^{-1/\alpha} can be improved to ∣Ω∣−1|\Omega|^{-1}. In contrast, with no pooling we have ∥h∥22/α=1\|h\|_{2}^{2/\alpha}=1, i.e., nn needs to be ∣Ω∣|\Omega| times larger for the same guarantee. Note that for α→1\alpha\to 1 and r=1/2r=1/2, the resulting bound resembles that of Proposition 3.

We note that if g∗g^{*} is assumed to be ss-smooth on the sphere, the source condition with 2αr=2sd−12\alpha r=\frac{2s}{d-1} corresponds to a Sobolev condition of order ss, and leads to the bound

which highlights that the rate only depends on the patch dimension dd instead of the full dimension d∣Ω∣d|\Omega|.

Under the conditions of the theorem, we may apply (Bach, 2021, Proposition 7.2), which states that for λ≤1\lambda\leq 1 and n≥5λ(1+log⁡(1/λ))n\geq\frac{5}{\lambda}(1+\log(1/\lambda)), we have

where the degrees of freedom Nh(λ){\mathcal{N}}_{h}(\lambda) is given in (22) and satisfies the upper bound (24), and the approximation error is given by

where HKh{\mathcal{H}}_{K_{h}} is the RKHS of KhK_{h}. Denote by

the decompositions of ff and f∗f^{*} in the orthonormal basis defined above. If g∗=∑k,jgk,jYk,jg^{*}=\sum_{k,j}g_{k,j}Y_{k,j} is the spherical harmonic decomposition of g∗g^{*}, then we have

where the second line uses a0∗=0a_{0}^{*}=0 and considers bk,j=a0,k,j/∣Ω∣b_{k,j}=a_{0,k,j}/\sqrt{|\Omega|}, while the third line uses g0,1=0g_{0,1}=0 and the fact that λ0=1\lambda_{0}=1 regardless of the choice of pooling filter hh. Thus, Ah(λ,f∗)A_{h}(\lambda,f^{*}) does not depend on hh, and corresponds to the approximation error Aκ(λ,g∗)A_{\kappa}(\lambda,g^{*}) of the patch kernel κ\kappa on the sphere (up to a factor ∣Ω∣|\Omega|). Under the source condition, we then have, by (Cucker & Smale, 2002, Theorem 3, p.33),

with a0,k,j∗=a0∗a^{*}_{0,k,j}=a^{*}_{0}. Combining with (24) and plugging this into (27), we obtain

In the case of global pooling, the term in parentheses scales as (∣Ω∣−1/α/n)2αr2αr+1(|\Omega|^{-1/\alpha}/n)^{\frac{2\alpha r}{2\alpha r+1}}, but can be improved to (1/∣Ω∣n)2αr2αr+1(1/|\Omega|n)^{\frac{2\alpha r}{2\alpha r+1}} by using the bound (23) on Nh(λ)\mathcal{N}_{h}(\lambda) instead of (24).

When r>α−12αr>\frac{\alpha-1}{2\alpha}, we can choose nn large enough so that n≥5λn(1+log⁡(1/λn))n\geq\frac{5}{\lambda_{n}}(1+\log(1/\lambda_{n})) is satisfied, and the higher order terms in nn are negligible. ∎