A self consistent theory of Gaussian Processes
captures feature learning effects in finite CNNs

Several previous papers derived leading order finite-DNN corrections to the GP results [30, 39, 12]. While these results are in principle extendable to any order in perturbation theory, such high order expansions have not been studied much, perhaps due to their complexity. In contrast, we develop an analytically tractable non-perturbative approach which we find crucial for obtaining non-negligible feature learning and associated performance enhancement effects.

Previous works [13, 14] studied how the behavior of infinite DNNs depends on the scaling of the top layer weights with its width. In [40] it is shown that the standard and NTK parameterizations of a neural network do not admit an infinite-width limit that can learn features, and instead suggest an alternative parameterization which can learn features in this limit. While unifying various viewpoints on infinite DNNs, this approach does not immediately lend itself to analytical analysis of the kind proposed here.

Several works [36, 25, 15, 2] show that finite width models can generalize either better or worse than their infinite width counterparts, and provide examples where the relative performance depends on the optimization details, the DNN architecture and the statistics of the data. Here we demonstrate analytically that finite DNNs outperform their GP counterparts when the latter have a prior that lacks some constraint found in the data (e.g. positive-definiteness [26], weight sharing [34]).

Deep linear networks (FCNs and CNNs) similar to our CNN toy example have been studied in the literature [4, 37, 20, 16]. These studies use different approaches and assumptions and do not discuss the target shift mechanism which applies also for non-linear CNNs. In addition, their analytical results hinge strongly on linearity whereas our approach could be useful whenever several leading cumulants of the DNN output are known.

A concurrent work [43] derived exact expressions for the output priors of finite FCNs induced by Gaussian priors over their weights. However, these results only apply to the limited case of a prior over a single training point and only for a FCN. In contrast, our approach applies to the setting of a large training set, it is not restricted to FCNs and yields results for the posterior predictions, not the prior. Focusing on deep linear fully connected DNNs, recent work [24] derived analytical finite-width renormalization results for the GP kernel, by sequentially integrating out the weights of the DNN, starting from the output layer and working backwards towards the input. Our analytical approach, its scope, and the models studied here differ substantially from that work.

Here we define a teacher-student toy model showing several qualitative real-world aspects of feature learning and analyze it via our self-consistent shifted target approach. Concretely, we consider the simplest student CNN $f({\mathbf{x}})$, having one hidden layer with linear activation, and a corresponding teacher CNN, $g({\mathbf{x}})$

This describes a CNN that performs 1-dimensional convolution where the convolutional weights for each channel are ${\mathbf{w}}_{c}\in\mathbb{R}^{S}$. These are dotted with a convolutional window of the input $\tilde{{\mathbf{x}}}_{i}=\left(x_{S\left(i-1\right)+1},\dots,x_{S\cdot i}\right)^{\mathsf{T}}\in\mathbb{R}^{S}$ and there are no overlaps between them so that ${\mathbf{x}}=\left(x_{1},\dots,x_{N\cdot S}\right)^{\mathsf{T}}=\left(\tilde{{\mathbf{x}}}_{1},\dots,\tilde{{\mathbf{x}}}_{N}\right)^{\mathsf{T}}\in\mathbb{R}^{N\cdot S}$. Namely, the input dimension is $d=NS$, where $N$ is the number of (non-overlapping) convolutional windows, $S$ is the stride of the conv-kernel and it is also the length of the conv-kernel, hence there is no overlap between the strides. The inputs ${\mathbf{x}}$ are sampled from $\mathcal{N}(0,I_{d})$.

Despite its simplicity, this model distils several key differences between feature learning models and lazy learning or GP models. Due to the lack of pooling layers, the GP associated with the student fails to take advantage of the weight sharing property of the underlying CNN [34]. In fact, here it coincides with a GP of a fully-connected DNN which is quite inappropriate for the task. We thus expect that the finite network will have good performance already for $n=C^{*}(N+S)$ whereas the GP will need $n$ of order of the dimension ($NS$) to learn well [10]. Thus, for $N+S\ll NS$ there should be a broad regime in the value of $n$ where the finite network substantially outperforms the corresponding GP. We later show (§IV.4) that this performance boost over GP is due to feature learning, as one may expect.

Conveniently, the cumulants of the student DNN of any order can be worked out exactly. Assuming $\gamma$ and $\sigma^{2}$ of the noisy GD training are chosen such that¹¹1Generically this requires $C$ dependent and layer dependent weight decay. $a_{i,c}\sim\mathcal{N}\left(0,\sigma_{a}^{2}/CN\right),\quad{\mathbf{w}}_{c}\sim\mathcal{N}\left(\bm{0},\frac{\sigma_{w}^{2}}{S}I_{S}\right)$ (and similarly for the teacher DNN) the covariance function for the associated GP reads

Denoting $\lambda:=\frac{\sigma_{a}^{2}}{N}\frac{\sigma_{w}^{2}}{S}$, the even cumulant of arbitrary order $2m$ is (see App. F):

while all odd cumulants vanish due to the sign flip symmetry of the last layer. In this notation, we mean that the $\bullet$’s stand for integers in $\left\{1,\dots,2m\right\}$ and e.g. $\left(1_{i_{1}},2_{i_{2}}\right)\equiv\left(\tilde{{\mathbf{x}}}_{i_{1}}^{1}\cdot\tilde{{\mathbf{x}}}_{i_{2}}^{2}\right)$ and the bracket notation $\left[\left(2m-1\right)!\right]$ stands for the number of ways to pair the integers $\left\{1,...,2m\right\}$ into the above form. This result can then be plugged in 7 to obtain the self consistent (saddle point) equations on the training set. See App. A.4 for a convergence criterion for the saddle point, supporting its application here.

In §III our description of the self consistent equations was for a finite and fixed training set. Further analytical insight can be gained if we consider the limit of a large training set, known in the GP literature as the Equivalent Kernel (EK) limit [35, 38]. For a short review of this topic, see App. D. In essence, in the EK limit we replace the discrete sums over a specific draw of training set, as in Eqs. 5, 7, 8, with integrals over the entire input distribution $\mu({\mathbf{x}})$. Given a kernel that admits a spectral decomposition in terms of its eigenvalues and eigenfunctions: $K\left({\mathbf{x}},{\mathbf{x}}^{\prime}\right)=\sum_{s}\lambda_{s}\psi_{s}\left({\mathbf{x}}\right)\psi_{s}\left({\mathbf{x}}^{\prime}\right)$, the standard result for the GP posterior mean at a test point is approximated by [35]

This has several advantages, already at the level of GP analysis. From a theoretical point of view, the integral expressions retain the symmetries of the kernel $K({\mathbf{x}},{\mathbf{x}}^{\prime})$ unlike the discrete sums that ruin these symmetries. Also, Eq. 15 does not involve computing the inverse matrix $\tilde{K}^{-1}$ which is costly for large matrices.

In the context of our theory, the EK limit allows for a derivation of a simple analytical form for the self consistent equations. As shown in App. E.1 in our toy CNN both $\Delta g$ and $\hat{\delta}g$ become linear in the target. Thus the self-consistent equations can be reduced to a single equation governing the proportionality factor ($\alpha$) between $\hat{\delta}g$ and $g$ ($\hat{\delta}g=\alpha g$). Thus starting from the general self consistent equations, 5, 7, 8, taking their EK limit, and plugging in the general cumulant for our toy model (14) we arrive at the following equation for $\alpha$

Setting for simplicity $\sigma_{a}^{2}=1=\sigma_{w}^{2}$ we have $\lambda=1/\left(NS\right)$ and we also introduced the constant $q\equiv\lambda^{-1}(1-\hat{\alpha}_{\rm{GP}})(\lambda+\sigma^{2}/n)$ where $\hat{\alpha}_{\mathrm{GP}}$ is computed using the empirical GP predictions on either the training set or test set: $\hat{\alpha}_{\rm{GP}}\equiv 1-\left(\sum_{\mu}f_{\mu}^{\rm{GP}}g_{\mu}\right)/\left(\sum_{\mu}g_{\mu}^{2}\right)$, or analytically in the perturbation theory approach developed in [10]. The quantity $q$ has an interpretation as $1/n$ corrections to the EK approximation [10] but here can be considered as a fitting parameter. It is non-negative and is typically $O(1)$; for more details and analytical estimates see App. E.2.

Equation 16 is the second main analytical result of this work. It simplifies the highly non-linear inference problem to a single equation that embodies strong non-linear finite-DNN effect and feature learning (see also §IV.4). In practice, to compute $\alpha_{\mathrm{test}}$ we numerically solve 16 using $q_{\mathrm{train}}$ for the training set to get $\alpha_{\mathrm{train}}$, and then set $\alpha=\alpha_{\mathrm{train}}$ in the r.h.s. of 16 but use $q=q_{\mathrm{test}}$. Equation 16 can also be used to bound $\alpha$ analytically on both the training set and test point, given the reasonable assumption that $\alpha$ changes continuously with $C$. Indeed, at large $C$ the pole in this equation lays at $\alpha_{\rm{pole}}=(\sigma^{2}/n)(C/\lambda)^{1/2}\gg 1$ whereas $\alpha\approx\alpha_{\rm{GP}}<\alpha_{\rm{pole}}$. As $C$ diminishes, continuity implies that $\alpha$ must remain smaller than $\alpha_{\rm{pole}}$. The latter decays as $\sigma^{2}\sqrt{CNS}/n$ implying that the amount of data required for good performance scales as $\sqrt{CNS}$ rather than as $NS$ in the GP case.

In this section we numerically verify the predictions of the self consistent theory of Sec. §IV.2, by training linear shallow student CNNs on a teacher with $C^{*}=1$ as in Eq. 12, using noisy gradients as in Eq. 1, and averaging their outputs across noise realizations and across dynamics after reaching equilibrium.

For simplicity we used $N=S$ and $n\in\left\{62,200,650\right\},S\in\left\{15,30,60\right\}$ so that $n\propto S^{1.7}$. The latter scaling places us in the poorly performing regime of the associated GP while allowing good performance of the CNN. Indeed, as aforementioned, the GP here requires $n$ of the scale of $\lambda^{-1}=NS=O(S^{2})$ for good performance [10], while the CNN requires $n$ of scale of the number of parameters ($C(N+S)=O(S)$).

The results are shown in Fig 1 where we compare the theoretical predictions given by the solutions of the self consistent equation (16) to the empirical values of $\alpha$ obtained by training actual CNNs and averaging their outputs across the ensemble.

At this point there is evidence that our self-consistent shifted target approach works well within the feature learning regime of the toy model. Indeed GP is sub-optimal here, since it does not represent the CNN’s weight sharing present in the teacher network. Weight sharing is intimately tied with feature learning in the first layer, since it aggregates the information coming from all convolutional windows to refine a single set of repeating convolution-filters. Empirically, we observed a large performance gap of finite $C$ CNNs over the infinite-$C$ (GP) limit, which was also observed previously in more realistic settings [23, 13, 34]. Taken together with the existence of a clear feature in the teacher, a natural explanation for this performance gap is that feature learning, which is completely absent in GPs, plays a major role in the behavior of finite $C$ CNNs.

To analyze this we wish to track how the feature of the teacher ${\mathbf{w}}^{*}$ are reflected in the student network’s first layer weights ${\mathbf{w}}_{c}$ across training time (after reaching equilibrium) and across training realizations. However, as our formalism deals with ensembles of DNNs, computing averages of ${\mathbf{w}}_{c}$ with respect to these ensembles would simply give zero. Indeed, the chance of a DNN with specific parameters $\theta=\left\{a_{i,c},\,{\mathbf{w}}_{c}\right\}$ appearing is the same as that of $-\theta$. Consequently, to detect feature learning the first reasonable object to examine is the empirical covariance matrix $\Sigma_{W}=\frac{S}{C}WW^{\mathsf{T}}$, where the matrix $W\in\mathbb{R}^{S\times C}$ has ${\mathbf{w}}_{c}$ as its $c$’th column. This $\Sigma_{W}$ is invariant under such a change of signs and provides important information on the statistics of ${\mathbf{w}}_{c}$.

As shown in App. G, using our field-theory or function-space formulation, we find that to leading order in $1/C$ the ensemble average of the empirical covariance matrix, for a teacher with a single feature ${\mathbf{w}}^{*}$, is

A first conclusion that could be drawn here, is that given access to an ensemble of such trained CNNs, feature learning happens for any finite $C$ as a statistical property. We turn to discuss the more common setting where one wishes to use the features learned by a specific randomly chosen CNN from this ensemble.

To this end, we follow Ref. [27] and model $\Sigma_{W}$ as a Wishart matrix with a rank-one perturbation. The variance of the matrix and details of the rank one perturbation are then determined by the above equation. Consequently the eigenvalue distribution is expected to follow a spiked Marchenko-Pastur (MP), which was studied extensively in [6]. To test this modeling assumption, for each snapshot of training time (after reaching equilibrium) and noise realization we compute $\Sigma_{W}$’s eigenvalues and aggregate these across the ensemble. In Fig. 2 we plot the resulting empirical spectral distribution for varying values of $C$ while keeping $S$ fixed. Note that, differently from the usual spiked-MP model, varying $C$ here changes both the distribution of the MP bulk (which is determined by the ratio $S/C$) as well as the strength of the low-rank perturbation.

Our main finding is a phase transition between two regimes which becomes sharp as one takes $n,S\rightarrow\infty$. In the regime of large $C$ the eigenvalue distribution of $\Sigma_{W}$ is indistinguishable from the MP distribution, whereas in the regime of small $C$ an outlier eigenvalue $\lambda_{m}$ departs from the support of the bulk MP distribution and the associated top eigenvector has a non-zero overlap with ${\mathbf{w}}^{*}$, see Fig. 2. We refer to the latter as the feature-learning regime, since the feature ${\mathbf{w}}^{*}$ is manifested in the spectrum of the students weights, whereas the former is the non-feature learning regime. We use the quantity $\mathcal{Q}\equiv{\mathbf{w}}^{*\mathsf{T}}\Sigma_{W}{\mathbf{w}}^{*}$ as a surrogate for $\lambda_{m}$, as it is valid on both sides of the transition. Having established the correspondence to the MP plus low rank model, we can use the results of [6] to find the exact location of the phase transition, which occurs at the critical value $C_{\rm{crit}}$ given by

where we assumed for simplicity $N=S$ so that $\lambda=S^{-2}$.

Another setting where GPs are expected to under-perform finite-DNNs is the case of quadratic fully connected teacher and student DNNs where the teacher is rank-1, also known as the phase retrieval problem [26]. Here we consider some positive target of the form $g({\mathbf{x}})=({\mathbf{w}}_{*}\cdot{\mathbf{x}})^{2}-\sigma_{w}^{2}||{\mathbf{x}}||^{2}$ where ${\mathbf{w}}_{*},{\mathbf{x}}\in{\mathbb{R}}^{d}$ and a student DNN given by $f({\mathbf{x}})=\sum_{m=1}^{M}({\mathbf{w}}_{m}\cdot{\mathbf{x}})^{2}-\sigma_{w}^{2}||{\mathbf{x}}||^{2}$. We consider training this DNN on $n$ train points $\left\{{\mathbf{x}}_{\mu}\right\}_{\mu=1}^{n}$ using noisy GD training with weight decay $\gamma=2M\sigma^{2}/\sigma_{w}^{2}$.

Similarly to the previous toy model, here too the GP associated with the student at large $M$ (and finite $\sigma^{2}$) overlooks a qualitative feature of the finite DNN — the fact that the first term in $f({\mathbf{x}})$ is non-negative. Interestingly, this feature provides a strong performance boost [26] in the $\sigma^{2}\rightarrow 0$ limit compared to the associated GP. Namely the DNN, even at large $M$, performs well for $n>2d$ [26] whereas the associated GP is expected to work well only for $n=O(d^{2})$ [10].

We wish to solve for the predictions of this model with our self consistent GP based approach. As shown in App. I, the cumulants of this model can be obtained from the following cumulant generating function

The associated GP kernel is given by $K({\mathbf{x}}_{\mu},{\mathbf{x}}_{\nu})=2M^{-1}\sigma_{w}^{4}({\mathbf{x}}_{\mu}\cdot{\mathbf{x}}_{\nu})^{2}$. Following this, the target shift equation, at the saddle point level, appears as

In App. I, we solve these equations numerically for $\sigma^{2}=10^{-5}$ and show that our approach captures the correct $n=2d$ threshold value. An analytic solution of these equations at low $\sigma^{2}$ using EK or other continuum approximations is left for future work (see Refs. [10, 7, 8] for potential approaches). As a first step towards this goal, in App. I we consider the simpler case of $\sigma^{2}=1$ and derive the asymptotics of the learning curves which deviate strongly from those of GP for $M\ll d$.