Learning Curves for Deep Neural Networks: A Gaussian Field Theory Perspective

TL;DR

This paper introduces a physics-inspired field-theory approach to analyze deep neural networks, accurately predicting learning curves and explaining their generalization despite over-parameterization.

Contribution

It develops a novel field-theory formalism for supervised deep learning, providing precise learning curve predictions and insights into generalization mechanisms.

Findings

Accurate predictions of deep DNN learning curves on polynomial regression tasks

Efficient hyper-parameter optimization based on the theoretical model

Explanation of generalization in over-parameterized DNNs via entropic bias

Abstract

In the past decade, deep neural networks (DNNs) came to the fore as the leading machine learning algorithms for a variety of tasks. Their raise was founded on market needs and engineering craftsmanship, the latter based more on trial and error than on theory. While still far behind the application forefront, the theoretical study of DNNs has recently made important advancements in analyzing the highly over-parameterized regime where some exact results have been obtained. Leveraging these ideas and adopting a more physics-like approach, here we construct a versatile field-theory formalism for supervised deep learning, involving renormalization group, Feynman diagrams and replicas. In particular we show that our approach leads to highly accurate predictions of learning curves of truly deep DNNs trained on polynomial regression tasks and that these predictions can be used for efficient hyper-parameter optimization. In addition, they explain how DNNs generalize well despite being highly over-parameterized, this due to an entropic bias to simple functions which, for the case of fully-connected DNNs with data sampled on the hypersphere, are low order polynomials in the input vector. Being a complex interacting system of artificial neurons, we believe that such tools and methodologies borrowed from condensed matter physics would prove essential for obtaining an accurate quantitative understanding of deep learning.

Tables

Table 1. Table 1: Comparison of the performance of networks trained in the NTK regime using different hyper-parameters ( 𝝈 𝝈 \bm{\sigma} ). The Test column shows an estimate of the DNN test loss, the Prediction column the loss as predicted by our learning-curve, and the GPR column an estimate of the dataset averaged expected loss of the corresponding GP. Worst, Median and Best refer to one of 21 networks with random hyper-parameters ranked by test loss. Typical and Optimized refer to networks with 𝝈 Typical subscript 𝝈 Typical \bm{\sigma}_{\textrm{Typical}} (defined in the text) and the optimal hyper-parameters following our optimization scheme. For more experiment results see App. L

	Test	Prediction	GPR
Worst	0.413	0.406	0.382
Median	0.313	0.316	0.319
Best	0.175	0.198	0.214
Typical	0.307	0.307	0.317
Optimized	0.078	0.110	0.141

Table 2. Table 2: Hyper-parameter performance comparison

	Test	Prediction	GPR	${\sigma }_r^2$	Train	${\sigma }_{w_1}$	${\sigma }_{b_1}$	${\sigma }_{w_2}$	${\sigma }_{b_2}$
Random 1	0.389	0.400	0.364	0.068	1.29e-04	1.555	0.032	2.028	0.026
Random 2	0.316	0.287	0.331	0.004	4.17e-04	0.914	0.029	0.880	0.062
Random 3	0.191	0.219	0.250	0.003	3.36e-04	0.922	0.049	0.760	0.046
Random 4	0.300	0.324	0.306	0.045	1.39e-04	1.552	0.058	1.644	0.054
Random 5	0.268	0.338	0.273	0.070	1.24e-04	1.585	0.072	2.020	0.029
Random 6	0.413	0.406	0.382	0.065	1.28e-04	2.037	0.032	1.512	0.053
Random 7	0.332	0.297	0.335	0.010	8.76e-04	0.994	0.030	1.190	0.071
Random 8	0.228	0.245	0.262	0.015	6.41e-04	1.165	0.059	1.271	0.068
Random 9	0.337	0.308	0.332	0.018	1.37e-03	0.909	0.027	1.758	0.030
Random 10	0.371	0.392	0.355	0.069	1.19e-04	1.658	0.041	1.908	0.057
Random 11	0.313	0.316	0.319	0.032	1.31e-04	1.440	0.050	1.502	0.045
Random 12	0.335	0.340	0.329	0.042	1.33e-04	2.106	0.065	1.175	0.040
Random 13	0.397	0.336	0.373	0.018	1.71e-03	1.546	0.029	1.037	0.044
Random 14	0.236	0.253	0.271	0.017	4.18e-04	1.388	0.067	1.133	0.050
Random 15	0.293	0.288	0.306	0.018	6.37e-04	1.534	0.058	1.068	0.065
Random 16	0.175	0.198	0.214	0.011	1.04e-03	1.132	0.074	1.117	0.053
Random 17	0.206	0.226	0.246	0.013	7.66e-04	1.095	0.061	1.275	0.072
Random 18	0.264	0.264	0.292	0.011	1.00e-03	1.035	0.043	1.230	0.066
Random 19	0.239	0.248	0.282	0.006	5.34e-04	0.835	0.037	1.107	0.050
Random 20	0.385	0.378	0.364	0.054	1.43e-04	1.810	0.039	1.540	0.046
Random 21	0.318	0.300	0.323	0.018	9.58e-04	1.282	0.042	1.266	0.055
Worst	0.413	0.406	0.382	0.065	1.28e-04	2.037	0.032	1.512	0.053
Median	0.313	0.316	0.319	0.032	1.31e-04	1.440	0.050	1.502	0.045
Best	0.175	0.198	0.214	0.011	1.04e-03	1.132	0.074	1.117	0.053
Typical	0.307	0.307	0.317	0.028	1.41e-04	1.414	0.050	1.414	0.050
Optimized	0.078	0.110	0.141	0.002	1.66e-04	0.707	0.075	0.707	0.027

Equations

$$\displaystyle\left\{\begin{array}[]{l}h^{l+1}=x^{l}W^{l+1}+b^{l+1}\\ x^{l+1}=\phi\left(h^{l+1}\right)\\ f\left(x;W^{1},\ldots,W^{l},b^{1},\ldots,b^{l}\right)=x^{L+1}\end{array}\right.$$

$$g^{\ast }(x_{\ast })$$

$${\Vert g-g^{\ast }\Vert }^2=\int d{\mu }_x{\left(g(x)-g^{\ast }(x)\right)}^2$$

$$⟨g^{\ast }(x){⟩}_{D_N}\approx g_{EK,N}^{\ast }(x)=\sum _n\frac{{\lambda }_n}{{\lambda }_n+\frac{{\sigma }^2}{N}}{g}_n{\varphi }_n(x)$$

$$\frac{d{w}_i}{dt}=-{\partial }_{w_i}\left(\mathcal{L}[z_W]+\sum _j\frac{T{w}_j^2}{2{\sigma }_w^2}\right)+\sqrt{2T}{\xi }_i(t)$$

$$P[f]$$

$$\propto e^{-\frac{1}{2{\sigma }^2}\mathcal{L}\left[f\right]}\int dW{e}^{-\frac{1}{2{\sigma }_w^2}{\sum }_i{w}_i^2}\delta [f-z_W]$$

$$\propto P_{nd}[f]e^{-\frac{1}{2{\sigma }^2}\mathcal{L}\left[f\right]}$$

$$P_0\left[f\right]\propto e^{-\frac{1}{2}{\Vert f\Vert }_K^2}$$

$${\Vert f\Vert }_K^2=\int d{\mu }_{x}d{\mu }_{x^{\prime }}f(x)K^{-1}(x,x^{\prime })f(x^{\prime })$$

$$P_0[f]=\int dW{P}_0(W)\delta [f-z_W]$$

$$P\left[f\right]\propto e^{-\frac{1}{2}{\Vert f\Vert }_K^2-\frac{1}{2{\sigma }^2}{\sum }_{i=1}^N{\left(f(x_i)-g(x_i)\right)}^2}$$

$$P\left[f\right]\propto e^{-\frac{1}{2}{\Vert f\Vert }_K^2-\frac{1}{2{\sigma }^2}\mathcal{L}\left[f\right]}$$

$$Z[\alpha ]=\int \mathcal{D}f{e}^{-S[f]+\int dx\alpha (x)f(x)}$$

$$g^{\ast }(x_{\ast })={\frac{\delta \log (Z[\alpha ])}{\delta \alpha (x_{\ast })}|}_{\alpha =0}=\frac{1}{Z[0]}\int \mathcal{D}f\cdot f(x_{\ast })e^{-S[f]}$$

$$⟨...{⟩}_{\eta }$$

$$⟨g^{\ast }(x_{\ast }){⟩}_{\eta }=\underset{M\to 0}{\lim }{\frac{1}{M}\frac{\delta ⟨Z^M{⟩}_{\eta }}{\delta \alpha (x_{\ast })}|}_{\alpha =0}$$

$$⟨Z^M{⟩}_{\eta }=$$

$$e^{-{\sum }_{m=1}^M\left(\frac{1}{2}{\Vert f_m\Vert }_K^2-\int dx\alpha f_m\right)+\eta \int d{\mu }_x{e}^{-\frac{{\sum }_{m=1}^M{\mathcal{L}}_{f_m}(x)}{2{\sigma }^2}}}$$

$$e^{-\frac{{\sum }_{m=1}^M{\mathcal{L}}_{f_m}(x)}{2{\sigma }^2}}$$

$$+\frac{1}{2}{\left[\frac{{\sum }_{m=1}^M{\mathcal{L}}_{f_m}(x)}{2{\sigma }^2}\right]}^2+...$$

$$⟨Z^M{⟩}_{\eta }=Z_{EK}^M+O\left(1/{\eta }^2\right)$$

$$Z_{EK}\left[\alpha \right]=\int \mathcal{D}f{e}^{-S_{EK}[f]+\int dx\alpha (x)f(x)}$$

$${⟨g^{\ast }(x_{\ast })⟩}_{\eta }$$

$$=\arg \min \left[S_{EK}[f]\right]+O\left(1/{\eta }^2\right)$$

$$=g_{EK,\eta }^{\ast }(x_{\ast })+O\left(1/{\eta }^2\right)$$

$${\mathrm{Cov}}_0\left[f(x),f(y)\right]=\sum _i{\left(\frac{1}{{\lambda }_i}+\frac{\eta }{{\sigma }^2}\right)}^{-1}{\varphi }_i\left(x\right){\varphi }_i\left(y\right)$$

$$⟨g^{\ast }(x_{\ast }){⟩}_{\eta }$$

$$g_{SL,\eta }^{\ast }(x_{\ast })=$$

$${\mathrm{Cov}}_0\left[f\left(x\right),f\left(x\right)\right]{\mathrm{Cov}}_0\left[f\left(x\right),f\left(x_{\ast }\right)\right]$$

$$g_{SL,\eta }^{\ast }(x_{\ast })$$

$$-\frac{\eta }{{\sigma }^4}$$

$${\Lambda }_{i,j,k}$$

$${\lambda }_l$$

$$g_{SL,\eta }^{\ast }(x_{\ast })$$

$$K(x,y)$$

Full text

I Introduction

Deep artificial neural networks (DNNs) have been rapidly advancing the state-of-the-art in machine learning, showing human and sometimes super-human performance in image recognition [1], speech recognition [2], reinforcement learning [3] and natural language processing tasks [4]. Their raise to prominence was largely results-driven, with little theoretical support or guarantee [5]. Such mode of invention is very different from how, say, the transistor was discovered, and more akin to how new materials, such as lithium-ion batteries, are discovered. Indeed, being huge interacting systems of artificial neurons, DNNs are more analogous to a complex meta-material than to an electronic component [6]. Due to this complexity, a general theory of deep learning with predictive power is still lacking.

Notwithstanding, recently several results were obtained in the highly over-parameterized regime [7, 8] where the role played by any specific DNN weight is small. This facilitated the proofs of various bounds [9, 10, 11] on generalization for shallow networks and, more relevant for this work, two correspondences between fully-trained DNNs and a different type of inference models called Gaussian Processes (GPs) [12]. As shown below, these can be thought of as non-interacting scalar field-theories with disorder and a non-local action.

The first such correspondence [8] between GPs and trained DNNs is known as the Neural Tangent Kernel result, which we would refer to here as the NTK correspondence. It holds when highly over-parameterized DNNs are initialized according to standard practice and trained with Mean-Square-Error (MSE) loss at vanishing learning rate and without weight decay.

The second correspondence [13] (the NNSP, Neural Networks Stochastic Process correspondence) applies when DNNs are trained using a similar protocol which involves random noise, roughly mimicking the Stochastic Gradient Descent (SGD) optimization, as well as weight-decay. It relates the outputs of the trained DNN to a Stochastic Process (SP) which, in the highly over-parameterized limit, tends to a GP. It thus yields an additional training protocol, complementary in some ways to the previous one, which is analytically tractable.

How much of deep learning can be explained through such correspondences remains to be seen. On the one hand, some aspects such as learning sharp filters (features) in the first DNN layers, seem out of reach as specific DNN weights change only infinitesimally in the NTK case and remain largely random, apart from a small bias, in the NNSP case. In addition, learning in the NTK regime, sometimes dubbed ”lazy-learning”, often lags behind state-of-the-art training protocols (see 14 and Refs. therein) where finite learning-rates, widths, and mini-batches are used. On the other hand, lazy-learning or more generally GP methods are being extended and improved [15, 16, 17, 18] by importing technologies such as pooling and data-augmentation. Currently GP models corresponding to DNNs are competitive with deep learning on the UCI datasets [15] as well as Fashion-MNIST [16], whereas on the CIFAR-10 dataset the performance of the best GPs currently lags $5\%$ behind celebrated DNNs such as AlexNet [16, 18], while surpassing pre-AlexNet non-deep methods by $8\%$ [19]. In addition, there is the prospect of extending these correspondences by include non-linearities coming from finite-width [13] and finite-learning rates [20, 21] settings. These results and prospects invite further study of how such DNNs trained in the NTK and NNSP regimes make predictions.

In this work we introduce a versatile field-theory formalism for analyzing deep neural networks, which involves replicas, Feynman diagrams, and renormalization group techniques. In its most basic version, studied in depth below, it applies to DNNs trained using the protocols for which the NTK and NNSP correspondence hold exactly and lead to GP models. For these cases we provide expressions for the generalization power of fully-connected DNNs in the form of learning curves. These learning curves depend on the dataset distribution and the target which we learn. For uniform datasets on the hyper-sphere and any target function, our learning curves become fully explicit and provide a clear picture of how such DNNs generalize. This includes the more challenging case of the NTK correspondence where certain infinities in the action are removed by our renormalization group transformation. To the best of our knowledge, the accuracy at which our learning-curves capture the empirical ones far exceeds the current theoretical state of the art.

In addition our formalism can also accommodate various extension of these correspondences. For the case of the NNSP correspondence, we can work with loss functions different than MSE as well as corrections to the infinite over-parameterization limit. Furthermore, recent results on extensions of the NTK correspondence [20], suggest that high-learning-rate leads to a renormalized NTK correspondence whose performance can again be analyzed using our approach. Such extensions may prove useful in addressing the gap [22], between GPs and their DNN counterparts.

We hope that the results and formalism introduced here would aid in developing a more physics-like paradigm for studying DNNs, complementary to the proof-based approach common in theoretical computer science (see also [23, 24, 25, 26, 27]). Such a paradigm should fill in the gap, typically large in complex systems, between what can be predicted following some reasonable assumptions and what can be proven rigorously.

This paper is structured as followed. In section III we provide the necessary background on Deep neural networks, Gaussian processes and the correspondences between the two. Section IV describes our novel field theory approach and analytical results. Section V considers the case of uniformly distributed data on the hypersphere, where further analytical simplifications can be carried. Section VI introduces the RG approach used to tackle the noiseless NTK case. Section VII applies our results to concrete examples and compares them with empirical results. Section VIII shows how our results can be used to perform efficient hyper-parameters optimization on actual DNNs, and Section IX summarizes the results and discusses possible directions for future work.

II Prior works

Learning curves for GPs have been analyzed using a variety of techniques (see [12] for a review) most of which focus on a GP-teacher averaged case where the target/teacher is drawn from the same GP used for inference (matched priors) and is furthermore averaged over. Fixed-teacher or fixed-target learning curves have been analyzed using a grand-canonical/Poisson-averaged approach [28] similar to the one we used. However, their treatment of the resulting partition function was variational whereas we take a different, perturbation-theory based, approach. In addition, previous cited results for MSE-loss break in the noiseless limit [28]. To the best of our knowledge, noiseless GPs learning-curves have been analyzed analytically only in the teacher-averaged case and limited to the following settings: For matched priors, exact results are known for one dimensional data [29, 12] and two dimensional data with some limitations of how one samples the inputs (in the context of optimal design) [30, 31]. In addition [32] derived a lower bound on generalization. For noiseless inference with partially mismatched-priors (matching features, mismatching eigenvalues) and at large input dimension the teacher and dataset averaging involved in obtaining learning curves, has been performed analytically and the resulting matrix traces analyzed numerically [33]. Notably none of these cited results apply in any straightforward manner in the NTK-regime.

Considering kernel eigenvalues, explicit expression for the features and eigenvalues of dot-product kernels ($K=K(x\cdot x^{\prime})$) were given in [34]. The fact that the $l$-th eigenvalue of such kernels scales as $d^{-l}$ ($d$ being the input dimension), which we used in our derivation of the bound, has been noticed in [33]. Kernels with a trimmed spectrum where the spectrum is trimmed after the first $r$’s leading eigenvalues, has previously been suggested as a way of reducing the computational cost of GP inference [35]. In contrast we trim the Taylor expansion of the kernel function rather than the spectrum (which has a very different effect) and show that an effective observation noise compensates for our trimming/renormalization procedure.

Several interesting recent works give bounds on generalization [9, 10, 11] which show $O(1/\sqrt{N})$ asymptotic decay of the learning-curve (at best). In contrast our predictions are typically well below this bound.

III Theoretical background

III.1 DNNs, expected error, and learning curves

We begin with the standard definitions of DNNs as they apply to this work. While the majority of this work is applicable to many network architectures, we will focus on a simple feed forward network for the sake of simplicity. A fully connected feed forward DNN with $L$ hidden layers of width $n_{l}$ for $l=1,\ldots,L$ and readout layer $n_{L+1}=k$ is a function $f$ defined recursively by:

[TABLE]

where $\phi$ is a point-wise activation function, $x^{0}\in\mathbb{R}^{d}$ is the input of the network and $W^{l+1}\in\mathbb{R}^{n_{l}\times n_{l+1}}$, $b^{l+1}\in\mathbb{R}^{n_{l+1}}$ are trainable weights and biases, which will be collectively referred to as weights from here on. Each component of the weights is usually initialized randomly from a normal distribution ${\cal N}\left(0,\sigma_{w}^{2}\right)$ for the weights and ${\cal N}\left(0,\sigma_{b}^{2}\right)$ for the biases.

In the usual setting one starts with a training set – a set of input points $D=\left\{x_{n}\right\}^{N}_{n=1}$ where $x_{n}\in\mathbb{R}^{d}$ along with their labels $\left\{l_{n}\right\}^{N}_{n=1}$ where $l_{n}\in\mathbb{R}^{k}$. One then picks weights for the network by minimizing a loss function ${\cal L}\left(f\left(D\right),\left\{l_{n}\right\}\right)$ which compares the values of network function over $D$ to the labels $\left\{l_{n}\right\}^{N}_{n=1}$, assigning a smaller value to points where the network function and labels are similar. One then finds weights which minimize the loss by some variation of gradient descent , usually stochastic gradient descent (SGD) wherein one approximates the gradient at each iteration using a random batch of the training set (see [36] for details). The performance of the network is then evaluated by computing the loss function over a set of labeled points, different from the training set, known as the test set. This is known as the test error of the network, and is used as a proxy for the expected error – the loss averaged over draws from the dataset distribution.

One of the most detailed objects quantifying the performance of a machine learning algorithm, and the main focus of this work, is its learning-curve – a graph of how the expected error diminishes with the number of data points ($N$). There are currently no analytical predictions or bounds we are aware of for DNN learning-curves which are tight even just in terms of their scaling with $N$, let alone tight in an absolute sense (see App. II)

III.2 Gaussian processes regression

In this work we will investigate the properties of DNNs by their correspondence with GPs. We supply here some standard definitions of GPs and their usage in regression tasks. Regression here simply means approximating a function ($g(x)$) based on discrete samples ($\{g(x_{n})\}_{n=1}^{N}$). A GP is commonly defined as a stochastic process of which any finite subset of random variables follow a multivariate normal distribution [12]. In a similar fashion to multivariate normal variables, GPs are also determined by their first and second moments. The first is typically taken to be zero, and second is known as the covariance function or the kernel $K_{xx^{\prime}}=\mathbb{E}\left[f(x)f(x^{\prime})\right]$, where $\mathbb{E}[\cdot]$ here denotes expectation with respect to the GP distribution. The main appeal of GPs is that Bayesian Inference with GP priors is tractable [12]. In GP inference we use the mean of the GP distribution conditioned on the data (posterior) as the predictor $g^{*}$, and it is given by:

[TABLE]

where $x_{*}$ is a new data point, $l_{m}$ are the training targets, $x_{n}$ are the training data-points, $[K(D)]_{nm}=K_{x_{n},x_{m}}$ is the covariance-matrix (the covariance-function projected on the training dataset $D$), and $\sigma^{2}$ is the variance of the assumed Gaussian noise of the labels, which also acts as a regulator of the prediction. Some intuition for this formula can be gained by verifying that in the noiseless case ($\sigma^{2}=0$) the prediction at some training point $x_{*}=x_{q}$ coincides with that point’s label $g^{*}=l_{q}$.

The quantity of interest in this paper, which we define now, is the expected error averaged over all the possible datasets. Throughout this paper we will assume that both train and test points are drawn from a probability measure $d\mu_{x}=P(x)dx$. With this in mind, we define the expected error of a prediction $g^{*}$ as

[TABLE]

Note that $g^{*}$ is itself a function of $N$ draws from $\mu$ which make up the training set $D_{N}$. Our quantity of interest, the dataset averaged expected error (DAEE), is Eq. (6) averaged over the ensemble of all possible $N$ sized training sets. We denote this average as $\langle\cdot\rangle_{D_{N}}$, so the DAEE is given by $\langle\left\|g^{*}-g\right\|^{2}\rangle_{D_{N}}$. The learning curve is the dependence of the DAEE on $N$. We see that in order to calculate learning curves, one needs to calculate quantities like $\langle g^{*}\rangle_{D_{N}}$ and $\langle g^{*2}\rangle_{D_{N}}$.

Equation (5) determines the predictions, and therefore the learning-curves, but it is not very convenient for analytic exploration of the expected predictions. This fact is due to the (potentially very) large matrix inversion involved, and the additional averaging over $D_{N}$ required. Nonetheless, there are some approximations for the expected prediction $\langle g^{*}\rangle_{D_{N}}$. The most famous of which is the equivalence kernel (EK) result [12]:

[TABLE]

Where $\lambda_{n}$ and $\phi_{n}(x)$ here are the eigenvalues and eigenfunctions of the kernel w.r.t the input probability measure $\mu$, and $g(x)=\sum_{n}g_{n}\phi_{n}(x)$ is the target function. One notices immediately that this approximation breaks down completely in the noiseless case where Eq. (7) implies perfect estimation of the target with just one data point. To gain some intuition as to why having $\sigma^{2}=0$ hinders predictions of $\langle g^{*}\rangle_{D_{N}}$ one can view it as a hard constraint ($f(x_{n})=g(x_{n})$), and hard constraints are typically less tractable than soft ones. In a related view, finite $\sigma^{2}$ can be seen as a form of averaging which smooths and regulates analytical expressions making them more tractable. Another limitation of the EK result is that (to the best of our knowledge) there is no systematic way to extend it in orders of $1/N$ and get a more detailed picture of generalization in GP regression (GPR).

III.3 From DNNs to GPs through Langevin dynamics

Here we review, for completeness, several recent correspondences between DNNs and GPs. It has long been known [37] that randomly initialized, infinitely wide DNNs with i.i.d weights are equivalent to samples from a GP known as the neural network GP (NNGP). More recently it was shown that training only the last layer of a network with gradient decent is equivalent to posterior sampling of the NNGP [38], and consequently averaging the prediction of many networks trained on the same dataset is equivalent to GPR. Turning to more standard training of the entire DNN, it has been recently established [8] that fully training a network with vanishing learning rate for infinitely long time and MSE loss yields the same predictions as a noiseless GPR with a different kernel, the neural tangent kernel (NTK), along with an additional initialization dependent term. Averaging over many initialization seeds gives an exact correspondence with a GP whose kernel is the NTK.

Recently, another novel correspondence between DNNs and GPs has been introduced [13]. Due to its simplicity we shall re-derive it here. Consider the training of a DNN using gradient descent (full-batch SGD) with weight decay and added white noise, in the limit of vanishing learning rate. For sufficiently small learning rate, and making the reasonable assumption that the gradients of the loss are globally Lifshitz, the SGD equations are ergodic and converge to the same invariant measure (equilibrium distribution) as the following Langevin equation [39, 40]

[TABLE]

where $\xi_{i}(t)$ being a set of Gaussian white noise ($\left\langle\xi_{i}(t)\xi_{j}\left(t^{\prime}\right)\right\rangle=\delta_{ij}\delta\left(t-t^{\prime}\right)$), $T$ accounting for the strength of the noise, $w_{i}$ being the set of networks parameters ($W$), $z_{W}$ is the network output for a given configuration of $W$ and $\mathcal{L}$ is the loss function. The equilibrium-distribution or invariant-measure describing the steady state of the above equation is the Boltzmann distribution [40] $P(W)\propto e^{-\frac{1}{2\sigma^{2}}\mathcal{L}\left[z_{W}\right]-\frac{1}{2\sigma_{w}^{2}}\sum_{i}w_{i}^{2}}$ with $T=2\sigma^{2}$. Notably, various works argue that at low learning rates, discrete SGD dynamics running for long enough times, reaches the above equilibrium, approximately [41, 42].

Next, we adopt the approach of [8] and describe the dynamics in function space ($f$) instead of weight space ($W$). Using the Boltzmann distribution described above, the post-training probability density function for some function $f$ is given by

[TABLE]

where we identify $P_{nd}[f]\propto\int dWe^{-\frac{1}{2\sigma_{w}^{2}}\sum_{i}w_{i}^{2}}\delta[f-z_{W}]$ as the distribution of the output of the network after being trained with no data (or equivalently, a vanishing loss function). $\delta[...]$ is a functional delta function, which can be thought of as the limit of a large product of regular delta-functions on each Fourier component of the argument. As we will discuss, for an infinitely over-parameterized network $P_{nd}$ coincides with the prior of a NNGP with the weights and biases variance determined by training parameters rather than by initialization. However, for finite over-parameterization it becomes a more generic stochastic process determined by the neural network (an NNSP). In [13], the leading finite-width corrections were calculated and shown to result in $f^{4}$ corrections to the prior.

Clearly the practical use of the above result hinges on how quickly the dynamics mixes or reaches ergodicity. While ergodicity in its full sense (for any weight-space observable) seems unrealistic in this non-convex scenario, reaching ergodicity in the mean of the outputs of the DNNs (for low order polynomials in $f(x)$) may be quicker. This milder form of ergodicity was shown numerically for fully-connected DNNs trained on regression problems similar to those studied here as well as CNNs trained on CIFAR-10 [13].

From now on we shall focus on the infinite over-parameterized limit.

IV Field theory formulation of GP learning-curves

IV.1 Rephrasing GPs as a field theory

We begin by phrasing inference with GPs in the language of field theory. To this end we first write a Gaussian distribution over the space of functions that leads to a two point correlation function equal to $K_{xx^{\prime}}$. This is given by

[TABLE]

where $K^{-1}(x,x^{\prime})$ is the inverse kernel function, meaning that $\int d\mu_{x^{\prime}}K(x,x^{\prime})K^{-1}(x^{\prime},x^{\prime\prime})=\delta(x-x^{\prime\prime})/P(x)$ where $d\mu=P(x)dx$. This formalism is sometimes referred to as Information Field Theory (IFT) [43].

A different viewpoint on $P_{0}[f]$ comes from viewing $f(x)$ as the outputs of a wide DNN with weights drawn from an i.i.d Gaussian distribution $P_{0}(W)$. It is well known [44] that correlations between the outputs of random DNNs are Gaussian and governed by some kernel, $K_{xx^{\prime}}$. This kernel is determined, in a tractable manner, by the DNNs architecture. From a field-theory viewpoint this can be stated as

[TABLE]

at infinite width, where $z_{W}(x)$ is the output of a DNN with weights $W$ on an input point $x$. The keen reader may be alarmed by the fact that this definition of $P_{0}[f]$ does not involve the measure $\mu(x)$. However, as shown in [12], the norm $\left\|f\right\|_{K}^{2}$ (called the RKHS norm), and therefore Eq. (10) are in fact the same for any two probability measures with identical support.

Performing Bayesian inference in the context of GPs means conditioning Eq. (10) using Bayes’ theorem and assuming Gaussian noise with amplitude $\sigma^{2}$ on our target function ($g(x)$). This yields the additional factor

[TABLE]

It can be checked that the expectation value of $f(x_{*})$ under the above probability yields Eq. (5).

Notably, by taking into account Eq. (11), the above expression coincides with that obtained via the NNSP correspondence Eq. (9) in the infinite over-parameterization limit where $P_{nd}[f]=P_{0}[f]$ for MSE loss and a suitably chosen $K_{xx^{\prime}}$. In the NNSP context, the data-term came out quadratic when training using MSE loss and more generally it could be replaced with a general loss function $\mathcal{L}\left[f\right]$, so the DNNs predictive distribution becomes

[TABLE]

where $K$ in this context is the kernel of the DNN trained with no data. Though not necessarily Gaussian, this expression can still be treated using mean-field or perturbative approaches. A more detailed treatment of different loss functions, most notably as cross-entropy loss, is left for future work.

Denoting $S[f]=\frac{1}{2}\left\|f\right\|_{K}^{2}+\frac{1}{2\sigma^{2}}\mathcal{L}\left[f\right]$ (the ”Information Hamiltonian”, in IFT terminology), Eq. (12) gives rise to the partition function

[TABLE]

where $\int dx\alpha(x)f(x)$ is a source term used to calculate cumulants of $P[f]$, and specifically the average prediction of the network:

[TABLE]

where $\delta/\delta\alpha$ stands for functional derivative. As shown visually in Fig. 1, this expression leads to a tangible physical picture of how DNNs learn. Their output as a function of the input can be seen as a fluctuating elastic membrane over input space which, in the highly over-parameterized limit, is in its linear elastic regime. The training data appears as isolated points at which this membrane is pinned down to certain value by (loss dependent) springs whose constant is proportional to $1/\sigma^{2}$– the inverse of the noise on the gradients during training. The membrane then interpolates and extrapolate between these pinning points in a way which, on average, minimizes its elastic energy. This elastic energy differs considerably from that of physical membrane and in particular has a non-local dependence on the shape of the membrane. Different DNNs correspond to different elastic energies. Finite networks entail non-linear corrections to the elastic energy which may be beneficial for learning in the case of CNNs [13].

IV.2 Predictions in the grand-canonical ensemble

As mentioned, in order to calculate the learning curve one needs to calculate quantities like $\langle g^{*}\rangle_{D_{N}}$ and $\langle{g^{*}}^{2}\rangle_{D_{N}}$. These averages involve multidimensional integrations over all possible datasets. To facilitate their computation we adopt the approach of [28] and instead consider a related quantity given by the Poisson averaging of the former

[TABLE]

where $...$ can be any quantity, in particular $g^{*}$ and $g^{*2}$. This average can be thought of as a grand-canonical ensemble, though a non-standard one since we average the observables and not the partition function. Taking $\eta=N$ means we are essentially averaging over values of $N$ in an $\sqrt{N}$ vicinity of $N$. This means that as far as the leading asymptotic behavior is concerned, one can safely exchange $N$ and $\eta$ as the differences would be sub-leading. We therefore focus on calculating the grand-canonical DAEE, $\langle\left\|g^{*}-g\right\|^{2}\rangle_{\eta}$. In App. A we compare learning curves as a function of $N$ and $\eta$ and show that they match very well.

By using the grand canonical ensemble, averaging over draws from the dataset can be carried out as follows.

First, using the replica trick:

[TABLE]

where for integer $M$ and assuming that the loss function acts point-wise on the training set $\mathcal{L}[f]=\sum_{i=1}^{n}\mathcal{L}_{f}(x_{i})$, we have

[TABLE]

As shown in App. G.2, a Taylor expansion in $\eta$ of the above r.h.s. yields the $\langle...\rangle_{\eta}$ averaging appearing on the l.h.s.

Second, we notice that the main benefit of Eq. (17) and Eq. (18) over Eq. (5) is that it allows for a controlled expansion in $1/\eta$. At large $\eta$ (or similarly large $N$) we expect the fluctuations in $f_{m}(x)$ to be small and centered around $g(x)$. Indeed, such a behavior is encouraged by the term multiplied by $\eta$ in the exponent. We can therefore systematically Taylor expand the inner exponent

[TABLE]

and each term will yield a higher order of $\langle g^{*}(x_{*})\rangle_{\eta}$ in $1/\eta$.

IV.3 EK as a free theory

Notably, so far the choice of a loss function was largely arbitrary. The advantage of choosing MSE loss, $\mathcal{L}_{f}(x)=(f(x)-g(x))^{2}$, is that $P[f]$ also becomes a GP, or equivalently has a quadratic action. From now on we shall focus on MSE loss.

Aiming for standard pertubative calculations, we wish to perform diagrammatic calculations w.r.t a free quadratic theory. Expanding Eq. (19) to first order and substituting in Eq. (18) we obtain

[TABLE]

where $S_{EK}[f]=\frac{1}{2}\left\|f\right\|_{K}^{2}+\frac{\eta}{2\sigma^{2}}\int d\mu_{x}(f(x)-g(x))^{2}$, which is quadratic in $f$ and therefore induces a Gaussian field.

Substituting Eq. (20) in Eq. (17) we get

[TABLE]

where the second equality is due to the fact that for Gaussian distributions the expectation value coincides with the most probable value, and the third equality is due to [12], with the subtle change that $N$ is being replaced by $\eta$.

Let us denote by $\langle\ldots\rangle_{0}$ the free-theory average, that is an average w.r.t $Z_{EK}$. We therefore get $\langle f\rangle_{0}=g^{*}_{EK,\eta}\neq 0$, meaning that our free theory, though Gaussian, is not centered.

The correlations of the free theory are

[TABLE]

where again, $\lambda_{i}$ and $\phi_{i}$ are the eigenvalues and eigenfunctions of the kernel.

IV.4 Next order corrections

We now wish to perform pertubative calculations w.r.t to the free (Gaussian) EK theory, and obtain a sub-leading (SL) correction for the EK result in the inverse dataset size:

[TABLE]

Expanding Eq. (19) to second order, substituting in Eq. (18) and keeping only $O\left(1/\eta^{2}\right)$ terms, the calculation can be carried using Feynman diagrams w.r.t to the free EK Gaussian theory. Leaving the details to appendix I, the sub-leading correction is

[TABLE]

or explicitly

[TABLE]

As shown App. G.2, similar expressions for $\langle g^{*2}\rangle_{\eta}$ are obtained using two replica indices. Interestingly we find that $\langle g^{*2}\rangle_{\eta}=\langle g^{*}\rangle_{\eta}^{2}+O(1/\eta^{3})$. Hence, up to $O(1/\eta^{3})$ corrections, the averaged MSE error is $(\langle g^{*}(x_{*})\rangle_{\eta}-g(x_{*}))^{2}$ integrated over $x_{*}$. Since the variance of $g^{*}$ came out to be $O(1/\eta^{3})$ one finds that $g^{*}-g$, which is $O(1/\eta)$, is asymptotically much larger than its standard deviation. This implies self averaging at large $\eta$, or equivalently that our dataset-averaged results capture the behavior of a single fixed dataset.

Equations (23), (25) and their application to the calculation of the grand-canonical DAEE are one of our key results. They provide us with closed expressions for the DAEE as a function of $\eta$, namely the fixed-teacher learning curve. They hold without any limitations on the dataset or the kernel and yield a variant of the EK result along with its sub-leading correction. From an analytic perspective, once $\lambda_{i}$ and $\phi_{i}(x)$ are known, the above expressions provide clear insights to how well the GP learns each feature and what cross-talk is generated between features due to the second sub-leading term. Notably for the renormalized NTK introduced below, the number of non-zero $\lambda_{i}$’s is finite, and so the above infinite summations reduce to finite ones. This makes these expressions computationally superior to directly performing the matrix-inversion in Eq. (5) along with an $N-$dimensional integral involved in dataset-averaging. In addition, having the sub-leading correction allows us to estimate the range of validity of our approximation by comparing the sub-leading and leading contributions, as we shall do for the uniform case below.

V Uniform datasets

To make Eq. (25) interpretable, $\phi_{i}(x)$ and $\lambda_{i}$ are required. This can be done most readily for the case of datasets normalized to the hypersphere ($\left\|x_{n}\right\|=1$) with a uniform probability measure and rotation-symmetric kernel functions. By the latter we mean $K_{x,x^{\prime}}=K_{Ox,Ox^{\prime}}$ for any orthogonal matrix $O$ with the same dimension as the inputs. Although beyond the scope of the current work, obvious extensions to consider are datasets which are uniform only in a sub-space of $x$ and/or small perturbations to uniformity.

Importantly, both NNGP and NTK associated with any DNN with a fully connected first layer and weights initialized from a normal distribution, has the above symmetry under rotations (see App. E). It follows that such a kernel can be expanded as $K_{x,x^{\prime}}=\sum_{n}b_{n}(x\cdot x^{\prime})^{n}$. An additional corollary [34] is that its features are hyperspherical harmonics ($Y_{lm}(x)$) as these are the features of all dot product kernels. Hyperspherical harmonics are a complete and orthonormal basis w.r.t a uniform probability measure on the hypersphere. Note that this implies a non-standard normalization for the $Y_{lm}$s in this context, as they are usually normalized w.r.t Lebesgue measure. For each $l$ these can be written as a sum of polynomials in the input coordinates of degree $l$. The extra index $m$ enumerates an orthogonal set of such polynomials (of size $\operatorname{deg}(l)$). For a kernel of the above form the eigenvalues are independent of $m$ and given by [34]

[TABLE]

For ReLU and erf activations, the $b_{n}$’s, can be obtained analytically up to any desirable order [44]. Thus one can semi-analytically obtain the eigenvalues up to any desired accuracy. For the particular case of depth 2 ReLU networks with no biases, we report in App. J closed expression where the above summation can be carried out analytically for the NNGP and NTK kernels. However, as we shall argue soon, it is in fact desirable to trim the NTK in the sense of cutting-off its Taylor expansion at some order $m$, resulting in what we call the renormalized NTK. For such kernels, which would be our main focus next, Eq. (26) can be seen as a closed analytical expression for the eigenvalues.

Interestingly, for any dot-product kernel and uniform data of dimension $d$ on the hypersphere, there is a universal bound given by $\lambda_{l}\leq K_{x,x}/\deg(l)\approx O(d^{-l})$, where $K_{x,x}$ is a constant in $x$. Indeed, $K_{x,x}=\sum_{lm}\lambda_{l}=\sum_{l}\deg(l)\lambda_{l}$. The degeneracy ($\deg(l)$) is fixed from properties of hyper spherical harmonics, and equals $\deg(l)=\frac{2l+d-2}{l+d-2}\binom{l+d-2}{l}$ [45] which goes as $O(d^{l})$ for $l\ll d$. This combined with the positivity of the $\lambda_{l}$’s implies the above bound.

Expressing our target in this feature basis $g(x)=\sum_{l,m}g_{lm}Y_{lm}(x)$, Eq. (25) simplifies to

[TABLE]

where $C_{K,\sigma^{2}/\eta}=\sum_{l}\deg(l)(\lambda_{l}^{-1}+\eta/\sigma^{2})^{-1}$ and notably cross-talk between features has been eliminated at this order since $\sum_{m}Y_{lm}^{2}(x)=\deg(l)$ is independent of $x$, yielding $\sum_{\tilde{m}}\int d\mu_{x}Y_{lm}(x)Y_{l^{\prime}m^{\prime}}(x)Y^{2}_{\tilde{l}\tilde{m}}(x)=\deg(\tilde{l})\delta_{ll^{\prime}}\delta_{mm^{\prime}}$.

By splitting the sum in $C_{K,\sigma^{2}/\eta}$, to cases in which $\lambda_{l}<\sigma^{2}/\eta$ and their complement, one has the bound $C_{K,\sigma^{2}/\eta}<\#F\sigma^{2}/\eta+\sum_{lm|\lambda_{l}<\sigma^{2}/\eta}\lambda_{l}$, where $\#F$ is the number of eigenvalues such that $\lambda_{l}>\sigma^{2}/\eta$. Thus for kernels with a finite number of non-zero $\lambda_{i}$’s (as the renormalized NTK introduced below), and for large enough $\eta$, $\#F$ becomes the number of non-zero eigenvalues and $C_{K,\sigma^{2}/\eta}=\#F\sigma^{2}/\eta$ has a $\eta^{-1}$ asymptotic. This illustrates the fact that the above terms are arranged by their orders in $\eta$.

We can use Eq. (27) to understand the validity of the EK result. We therefore look for sufficient conditions for $g^{*}_{EK,\eta}\gg g^{*}_{SL,\eta}$ to hold. By a term-wise comparison, for some $l$ we obtain $C_{K,\sigma^{2}/\eta}\ll\eta(\lambda_{l}+\sigma^{2}/\eta)$ which holds for $C_{K,\sigma^{2}/\eta}\ll\sigma^{2}$. For trimmed kernels, this yield $\#F\ll\eta$. Notably it means that the original non-trimmed NTK cannot be analyzed perturbatively, since with $\sigma^{2}=0$, $\#F$ becomes infinite. In the next section we tackle this issue.

VI Generalization in the noiseless case and the renormalized NTK

The correspondence between DNNs trained in the NTK regime and GPR using NTK implies noiseless GPR ($\sigma^{2}=0$) for which the perturbative analysis carried in previous sections fails. Here we show that the fluctuations of $f$ associated with small $\lambda_{l}$s can be traded for noise on the fluctuations of $f$ associated with large $\lambda_{l}$s, thereby making our perturbative analysis applicable. As shown in the previous section, for uniform datasets, the smaller $\lambda_{l}$s correspond to higher spherical harmonics (higher $l$) and hence have higher oscillatory components. We argue that these higher oscillatory modes can be marginalized over in a controlled manner to generate both noise and corrections to the large $\lambda_{l}$s. This is very much in spirit of the renormalization group technique, wherein high oscillatory modes are integrated over to generate changes (renormalization) of some parameters in the probability distribution of the low oscillatory modes.

We begin by defining a set of renormalized NTKs. As argued, an NTK of any fully-connected DNN can be expanded as $K_{x,x^{\prime}}=\sum_{q=0}^{\infty}b_{q}(x\cdot x^{\prime})^{q}$. The renormalized NTK at scale $r$ is then simply $K^{(r)}_{x,x^{\prime}}=\sum_{q=0}^{r}b_{q}(x\cdot x^{\prime})^{q}$. Harmoniously with this notation we denote the prediction of GPR with the original kernel as $g_{\infty}^{*}$. Our claim is that GPR with $K$ and a noise of $\sigma^{2}$ can be well approximated by GPR with $K^{(r)}$ and noise $\sigma^{2}+\sigma_{r}^{2}$ (where $\sigma_{r}^{2}=\sum_{q=r+1}^{\infty}b_{q}$), for sufficiently large $r$. Specifically, our claim is that the discrepancy between the original vs. truncated GPR predictions scales as $O(\sqrt{N}d^{-(r+1)/2}/K_{x,x})$, where $d$ is the effective data-input dimension. Importantly, as can be inferred from Eq. (26), the renormalized NTK at scale $r$ has zero eigenvalues for all spherical Harmonics with $l>r$, as well as modified eigenvalues for spherical harmonics with $l\leq r$ (compared to the non-truncated NTK). Thus, as advertised, these high Fourier modes have been removed from the problem in exchange for a renormalized theory with a modified low energy spectrum, and augmented noise. In a related manner, trimming the Taylor expansion after $(x\cdot x^{\prime})^{r}$ effectively reduces our angular resolution and coarse grains the fine angular features captured by these spherical Harmonics with $l>r$.

To justify this approximation we consider the difference matrix $A_{nm}=K_{x_{n},x_{m}}-K^{(r)}_{x_{n},x_{m}}$, given a dataset $\{x_{n}\}_{n=1}^{N}$ drawn from a uniform distribution on a hypersphere of dimension $d$. The terms $b_{q}(x_{n}\cdot x_{m})^{q}$ scale roughly as $d^{-q/2}$ (see App. K for a more accurate expression) due to the tendency of random vectors in high dimensions to be orthogonal. Consequently the above difference diminishes very quickly with $r$. Notably this also applies for the entries of $K_{x_{*},x_{n}}-K^{(r)}_{x_{*},x_{m}}$, provided $x_{*}$ is a test point and not a train point. In contrast, the diagonal part of $A$ is $A_{nn}=\sigma_{r}^{2}$ and may diminish more slowly depending on the coefficients $b_{q>r}$. Upon neglecting $K_{x_{*},x_{n}}-K^{(r)}_{x_{*},x_{m}}$ and the off-diagonal elements of $A$, one finds that Eq. (5) with these two GPRs yields identical predictions. As shown in App. K, these neglected off-diagonal elements yield a discrepancy which scales as $\sqrt{N}d^{-(r+1)/2}$. Consequently, the MSE error between the two GPRs should scale as $N$ times an exponentially small factor ($d^{-r-1}$). This scaling with $N$ should saturate when the accuracy is nearly perfect since then the predictions remain largely constant as $N$ is increased.

Focusing back on the question of how to tackle noiseless GPR, we thus find that as long as the $b_{q}$’s decays slowly enough with $q$, then at any finite $N$ we can choose a large enough $r$ such that two desirable properties are maintained: A. The discrepancy between the GPRs is small and B. $\sigma_{r}^{2}$ is large enough to ensure convergence to our perturbative analysis. The required slow decay of $b_{q}$ is harmonious with the intuition that DNNs should be initialized at the edge of chaos [46] where the output of the network has a fine and multi-scale sensitivity to small changes in the input. As $K_{x,x^{\prime}}$ is the correlation of two outputs with inputs $x$ and $x^{\prime}$, having a power law decaying $b_{q}$ implies such fine and multi-scale sensitivity. Establishing relations between good initialization and effectiveness of our renormalized NTK is left for future work.

We have tested the accuracy of approximating noiseless NTK GPR with renormalized NTK GPR with the appropriate $\sigma_{r}^{2}$, both on artificial datasets (see next section) and on real world dataset such as CIFAR10 (see app. B). In both cases we found an excellent agreement between the two GPRs for $r$’s as small as $3$ and $4$.

VII Generalization in the NTK regime

Collecting the results of all the preceding sections, we can obtain a detailed and clear picture of generalization in fully connected DNNs trained in the NTK-regime on datasets with a uniform distribution normalized to some hypersphere in input space. We begin with a qualitative discussion and consider some renormalized NTK at scale $r$. From Sec. V, we have that the features of this kernel are hyperspherical harmonics and that $\lambda_{l}$ scales as $d^{-l}$. We also recall that $Y_{lm}$ is a polynomial of degree $l$ and that all the hyperspherical harmonics up to degree $l$ span all polynomials on the hypersphere with degree up to $l$. Examining Eq. (27) we find that features with $\lambda_{l}\gg\sigma^{2}/\eta$ are learnable and via the above scaling we find that we learn polynomials of degree $O(\log(\eta/\sigma^{2})/\log(d))$ or less. In particular, a function like parity, which is a polynomial of degree $d$ is very hard to learn whereas a linear function is the easiest to learn. Thus, despite having infinitely more parameters than data-points (due to infinite width) and despite being able to span almost any function (due to the richness of the kernel’s features), the DNN here avoids overfitting by having a strong bias towards low degree polynomials.

To make more quantitative statements we now focus on a specific setting. We consider input data in dimension $d=50$ and a scalar target function $g(x)=\sum_{l=1,2;m}g_{lm}Y_{lm}(x)$ such that the vectors $(g_{l,1},g_{l,2},\ldots,g_{l,\deg(l)})^{T}$ for $l=1,2$ are drawn from a uniform measure on the $\deg(l)$-sphere of radius $1/\sqrt{2}$. We generate several toy datasets $D_{N}$ consisting of $N$ data points ($x_{n}$) uniformly distributed on the hypersphere $S^{d-1}$ and their corresponding targets ($g(x_{n})$). We consider the GP equivalent to training a fully-connected DNN consisting of 4 layer with ReLU activations and width $W$ which we initialize with variance ($\sigma^{2}_{w}=\sigma^{2}_{b}=1/d$) for the input layer and ($\sigma^{2}_{w}=\sigma^{2}_{b}=1/W$) for the hidden layers (see for instance [38] App. C and App. E for how to compute the kernel. Notice there is a factor of $1/W$ between our convention for $\sigma_{w}^{2}$ and [38]). To be in the NTK correspondence regime we consider training such a network at vanishing learning-rate, MSE loss, and with $W\gg N$. One then has that the predictions of the DNN are given by GPR with $\sigma^{2}=0$ and the $K$ given by the NTK kernel [8] (To be more precise, [8] predict correspondence with GPR up to a random initialization factor, so to get exact match with GPR one would also need to average over initialization seeds. Recent research [38] suggests this caveat can be avoided under some conditions).

For each such DNN we obtained the expected MSE loss $\left\lVert g_{\infty}^{*}-g\right\rVert^{2}$ of GPR with the NTK kernel by numerical integration over $x_{*}$. Repeating this process multiple times we obtained the DAEE for $N=1,2,\dots,N_{\max}$ with a relative standard error of less then 5% (this typically required averaging over 10 datasets). For direct comparison with our prediction of the learning curve, we computed the Poisson averaged learning curve $\langle\left\lVert g_{\infty}^{*}-g\right\rVert^{2}\rangle_{\eta}$ in accordance with Eq. (16), neglecting the terms $n>N_{\max}$. We restricted ourselves to $\eta_{\max}\leq N_{\max}-5\sqrt{N_{\max}}$ to make tail effects negligible. Notably the Poisson averaging makes the final statistical error negligible relative to the discrepancies coming from our large $\eta$ approximations (see A). To make it easier to appreciate the power of GPs over simpler regression models we also provide the Poisson averaged DAEE for ordinary least square method (OLS) as a yardstick.

To pick the renormalization scale $r$ we must consider two factors, on the one hand we want the discrepancy between the renormalized and regular NTK to be small, this scales as $O(\sqrt{N}d^{-(r+1)/2}/K_{x,x})$. On the other hand we want the effective noise $\sigma^{2}_{r}$ to be as large as possible to assure the accuracy of the prediction. We found that $r=3$ strikes a good balance for the range of $N$ values used in the experiment, but $r=4,5,6$ also produced adequate predictions since $\sigma^{2}_{r}$ shrinks slowly with $r$ for the architecture used.

Our analytical expressions following Eq. (26) combined with known results [8, 44] about the Taylor coefficients ($b_{n}$) yield $\lambda_{0},...,\lambda_{3}=\{3.19,7.27\cdot 10^{-3},5.98\cdot 10^{-6},1.62\cdot 10^{-7}\}$ and $\sigma_{r}^{2}=0.018$. Since $\lambda_{0},\lambda_{1}\gg\sigma^{2}/\eta\gg\lambda_{2},\lambda_{3}$ for $50<\eta<3500$, $C_{K_{r},\sigma^{2}/\eta}\sigma^{-2}<[\deg(0)+\deg(1)]\sigma^{2}/\eta+O(\deg(2)10^{-6})$, thus $C_{K_{r},\sigma^{2}/\eta}\sigma^{-2}\approx 51/\eta$. Thus we expect perturbation theory to be valid for $\eta\gg 50$. At $\eta=1000$ the $l=1$ features are learned well since $\sigma^{2}/\eta=1.8\cdot 10^{-4}\gg\lambda_{1}$ and the $l=2$ features neglected, at $\eta=1000$ they are learned but suppressed by a factor of about $3$. Had the target contained $l=3$ features, they would have been entirely neglected at these $\eta$ scale. Experimental learning curves along with our leading and sub-leading estimates are shown in Fig. 2 left panel showing an excellent agreement between theory and experiment.

While no actual DNNs were trained in the above experiments, the NTK correspondence means that this would be the exact behavior of a DNN trained in the NTK regime [8, 38, 22]. Furthermore, since our aim was to estimate what the DNNs would predict rather than reach SOTA predictions, we focus on reasonable hyper-parameters but did not perform any hyper-parameter optimization. The complementary case of noisy GPR, which one encounters in the NNSP correspondence, is studied in App. C.

To demonstrate that our results work with more complex functions we repeated the experiment with a different target function which cannot be expressed as a finite order polynomial. We drew a uniformly distributed vector $w$ on the $(d-1)\textrm{-sphere}$ of radius $\sqrt{d}$ and set the target as $g(x)=|w\cdot x|$. Fig. 2 (right) shows good agreement between theory and experiment here as well.

Lastly we argue that the asymptotic behavior of learning-curve we predict is more accurate than the recent PAC based bounds [9, 10, 11]. In App. D we show a log-log plot of the learning-curves contrasted with a $1/\sqrt{\eta}$ which is the most rapidly decaying bound appearing in those works. It can be seen that such an asymptotic cannot be made to fit the experimental learning-curve with any precision close to ours.

VIII Application of results to hyper-parameter optimization

As with most machine learning algorithms, when training a neural network for a particular task one needs to choose a number of hyper-parameters such as the network’s width at each level $W_{l}$, depth $L$, variance of weights at initialization $\sigma_{b}^{2},\sigma_{w}^{2}$, activation function, and optimizer related parameters such as batch size, learning rate etc. There are many considerations for hyper-parameter selection such as training time, memory footprint and optimizer convergence, but here we will focus on the expected loss of the network. While there are some accepted heuristics, there is no a priori way to predict the best preforming architecture other than an expensive process of trial and error. In this section we introduce a scheme for picking theoretically advantageous parameters, given minimal information on the target function and dataset distribution. As recent research suggests [47] $W$ should be increased as much as possible to put the network in the interpolation regime, we will assume that choice was made. We also note that $\sigma_{b}^{2},\sigma_{w}^{2}$ are typically thought to be related more to convergence issues, for example via the exploding/vanishing gradient problem, than to the performance of the network. However as this work as well as [48] suggest, these parameters have an important effect on the network performance by changing the NTK spectrum.

We suggest the following scenario: we have $N=1000$ data points uniformly distributed on $S^{9}$. We are also given the spectral weight of the target in each eigenspace, that is $w_{\ell}^{2}=\|\Pi_{\ell}(g)\|^{2}$ where $\Pi_{\ell}$ is the projection operator on the $\ell$ subspace. For the case of uniform measure on the hyper-sphere the projection operator is simply $\Pi_{\ell}=P_{\ell}(\langle x\cdot y\rangle)\operatorname{deg}(\ell)$ where $P_{\ell}$ is the Legendre polynomial of degree $\ell$ and $\operatorname{deg}(\ell)$ is the dimension of the eigenspace, so finding $w_{\ell}^{2}$ is a much simpler task then finding the $\operatorname{deg}(\ell)$ coefficients of the target (which scale as $d^{\ell}$) and can be accomplished with a few numeric integrals. In this case we focus on a target with $w{{}_{\ell}}^{2}=\frac{1}{3}(\delta_{1,\ell}+\delta_{2,\ell}+\delta_{3,\ell})$. Given this setting, we would like to find a network architecture with minimal expected error. For computational efficiency reasons we decide to focus on ReLU networks with one hidden layer, so we need to choose four hyper-parameters $\bm{\sigma}=\left(\sigma_{w_{1}},\sigma_{w_{2}},\sigma_{b_{1}},\sigma_{b_{2}}\right)$.

We present two typical ways used to choose $\bm{\sigma}$, then propose a better way to do so based our theory. The naive and most prevalent way to choose hyper-parameters is to simply take $\bm{\sigma}_{\textrm{Typical}}=(\sqrt{2},\sqrt{2},0.05,0.05)$ which roughly correspond to He initialization [49], a common heuristic for avoiding gradient propagation issues. A more diligent approach would be to draw some random values in the vicinity of $\bm{\sigma}_{\textrm{Typical}}$, train the network, evaluate the test loss and pick the best preforming hyper-parameters $\bm{\sigma}_{\textrm{Best}}$.

Next we suggest a different approach which utilizes our analytical results. We construct a symbolic expression for the expected loss using the formalism outlined in the paper. By taking $\eta=N$ and applying the renormalization scheme with appropriate $r$ we get an estimator for the expected loss $\hat{L}\left(\bm{\sigma}\right)=\int dx\left\langle\left(f\left(x\right)-g\left(x\right)\right)^{2}\right\rangle$ which we can use to predict the performance of different hyper-parameters without training the network. Moreover, we can use standard numerical optimization algorithms to minimize the predicted loss and obtain $\bm{\sigma}_{\textrm{Optimized}}=\operatornamewithlimits{argmin}_{\bm{\sigma}}\hat{L}\left(\bm{\sigma}\right)$.

To test the scheme experimentally, we drew 21 random hyper-parameter proposals $\{\bm{\sigma}_{i}\}$ uniformly distributed in the hyper-rectangle defined by $\frac{1}{2}\bm{\sigma}_{\textrm{Typical}},\frac{3}{2}\bm{\sigma}_{\textrm{Typical}}$. We constrained the optimization algorithm to this hyper-rectangle as well to avoid convergence issues in the training procedure, as unconstrained optimization leads to very small values of $\sigma_{w}$ ($\ll\sqrt{2}$) which in turn lead to vanishing gradients and long training times. We defined networks corresponding to $\bm{\sigma}_{\textrm{Optimized}},\bm{\sigma}_{\textrm{Typical}},\{\bm{\sigma}_{i}\}$ and trained each network on the same dataset with full-batch gradient descent and learning rate $1.0$ until the train loss was smaller then $0.1\cdot\sigma^{2}_{r}$. A summary of the experiment is outlined in Table 1.

The results clearly demonstrate the effectiveness of our scheme which reduces the test loss by a factor of 4 relative to the typical hyper-parameter choice and a factor of 2 over the best performing random hyper-parameters. In terms of computational complexity, it took approximately 2.5 hours to train each network using Google’s neural tangent package [50] on a 20 core CPU 111While GPUs are generally faster, fully connected DNNs do not gain the full benefit of GPU parallelism and we expect the computation time would only improve by a factor of $O(1)$ with $W=2^{14}$. In comparison, the time it takes to build and optimize $\hat{L}$ is completely negligible at about 30 seconds. The best random hyper-parameters were found on the fifteenth attempt, so had we stopped then we would have wasted 35 computer hours relative to our scheme and gotten inferior test loss. Note also that we focused on $L=2$ in order to speed up training, which scales exponential with depth, but the optimization procedure is not nearly as sensitive to depth and could have been done for any reasonable $L$. Moreover, increasing the depth would have also enlarged the hyper-parameter space, making random search even less effective. For each network we also experimentally obtained the dataset averaged expected loss using GPR with the associated NTK. The fair agreement between the test loss and dataset averaged expected loss (GPR in Table 1) further solidifies previous results and demonstrates our claim of self-averaging.

As expected, the above results required some knowledge of the target function, in particular its spectral weight within each angular momentum space. Alternatively one can capitalize on the fact that our learning curves predictions are quadratic in the target, average them over a target function ensemble, and optimize with respect to this average case. Another option is to consider a min-max optimization scheme in which hyper-parameters are optimized for the worse case target within some domain. The scheme can also be extended to non-uniform datasets and different activation functions as long as some way of computing the eigenvalues is provided.

IX Discussion and Outlook

In this work we laid out a formalism based on field theory tools for predicting learning-curves in the NTK and NNSP correspondence regimes. Despite DNNs’ black-box reputation, well within the validly range of our perturbative analysis, we obtained very low relative mismatch between our best estimate and the experimental curves, with good agreement extending well into regions with low amounts of data compared to that needed to learn the target. A potential use of such learning curves in hyper-parameter optimization was also demonstrated.

Central to our analysis was a renormalization-group transformation leading to effective observation noise on the target and to a simpler renormalized quadratic-action/kernel. Notably this RG transformation implied that wide Fully-Connected networks, even ones working on real-world datasets such as CIFAR10, could be effectively described by very few parameters being the noise level and the $O(1)$-first Taylor expansion parameters of the kernel.

Our analysis provides a lab setting in which deep learning can be understood. In its training phase, DNNs avoid local-minima issues and glassy behavior due to their high over parameterization which makes the optimization problem highly under-constraint [52, 53, 13]. As a result, many different solutions or weights which fit perfectly the training data are possible. While each such solution will behave differently on a test point, this arbitrariness does not entail an erratic behavior. The reason is the implicit bias DNNs have towards simple functions. In the case of the NNSP correspondence, a simple function is, by definition, a function that can be generated, up to some small noise, by a large phase-space of weights.

Simplicity is therefore strongly architecture and dataset dependent. For fully connected DNNs trained in the regime of the NTK or the NNSP correspondences, as well as data uniformly sampled from the hypersphere, simplicity amounts to low order polynomials over that hypersphere. These are the hyper-spherical Harmonics with low $l$, which are the leading eigenfunctions w.r.t. such a uniform measure of a generic kernel associated with a fully connected DNN. As long as the DNN has at least one non-linear layer and biases, depth has only a quantitative effect as it modifies the eigenvalues ($\lambda_{l}$) but does not change their scale. Generally, the eigenvalues and eigenfunctions vary with architecture and data distribution. Convolutional neural networks (CNNs) require further study, however one can argue on a qualitative level that simple functions would be polynomials with certain spatial hierarchy. Moreover, one expects that qualitative details of this hierarchy would depend on depth as it controls the input-fan-in of the hidden activation in the last CNN layer.

It seems unrealistic that a purely analytical approach such as ours would describe well the predictions of state-of-the-art DNNs such as VGG-19 trained on a real-world datasets such as ImageNet. Similarly unrealistic is to expect an analytical computation based on thermodynamics to capture the efficiency of a modern car engine or one based on Naiver-Stoke’s equations to output a better shape of a wing. Still, scientific experience shows that understanding toy-models, especially rich enough ones, has value. Indeed toy-models provide an analytical lab where theories could be refined or refuted, algorithms could be benchmarked and improved, and wider ranging conjectures and intuitions could be formed. Such models are useful whenever domain knowledge possesses some degree of universality or independence from detail. In converse, when all details matter knowledge is nothing more than a log of all experiences. The fact that DNNs work well in variety of different architecture and data-set settings, suggests that some degree of universality worth exploring is present. Further research would thus tell if the tools and methodologies that have enabled us to comprehend our physical world can help us comprehend the artificial world of deep learning.

Many extensions of the current work, aimed at approaching real-world settings, can be considered. First and most, much of the recent excitement about DNNs comes from either CNNs or Long Short Term Memory networks (LSTMs). Considering CNNs, while much of our formalism applies, the spectrum of CNN Kernels is more challenging to obtain as their Kernels are less symmetric compared to Fully-Connected DNNs. From similar reasons the RG approach presented here requires a more elaborate trimming of the CNN kernel since the latter would not consist of only powers of dot-products. Furthermore, CNNs trained with SGD show rather large gaps in performance compared to their NNGP or NTK. The culprit here might very well be the finite-width or finite-number-of-channels corrections to the NNGP or NTK priors. Leading finite-width corrections, considered in Ref. [13], amount to adding quartic terms to $P_{0}[f]$. Those could be dealt with straightforwardly using our perturbation theory formalism. Interestingly, at least for CNNs without pooling, these corrections introduce a qualitative change to the prior, making it reflect the weight-sharing property of CNNs which is lost at the level of the NNGP or NTK [54, 13]. Other viable directions are handling richer datasets distributions, extending EK results to the more common cross-entropy loss, applying RG reasoning on finite-width DNNs, and using the above insights for developing DNN-architecture design principals.

Acknowledgements. Z.R and O.M acknowledge support from ISF grant 2250/19. Both O.M. and O.C. contributed equally to this work.

Appendix A Poisson Averaging Demonstration

Here we demonstrate that Poisson averaging has no substantial effect on the learning curve. To this end the figure below shows the experimental learning curve from the main text pre- and post-averaging. It is evident that other than the unintended consequence of eliminating the experimental noise, the averaged learning curve is equivalent to the original for all practical intents.

Appendix B Comparison of NTK and Renormalized NTK Predictions on Synthetic and Real World Datasets

In section VII we used the renormalized NTK as a proxy for the regular NTK, the purpose of this section is to affirm the validity of this approximation. Moreover, while our lack of knowledge of the NTK eigenvalues and eigenfunctions with respect to a non-uniform measure prevents us from predicting learning curves, we would like to show that the renormalized NTK is a valid approximation in this setting as well.

To this end we used the following procedure. We took the NTK kernel defined in the paper and its associated renormalized kernels at different scales and trained them over the same training set $D_{N}$. In the figure below (top) the training set and target function were the ones defined in the main text. In the figure below (bottom) $D_{N}$ consisted of uniform draws without replacements from the cifar-10 training set, standardized and normalized to unit vectors, and the target function was the one-hot encoding of the labels standardized to have zero mean and $K_{x,x}$ variance.

For each training set we logged the average squared deviation of each renormalized kernel estimation $g^{\star}_{r}$ from the estimation of the non-renormalized kernel $g^{\star}_{\infty}$. This is the quantity $\left\|g^{\star}_{r}-g^{\star}_{\infty}\right\|^{2}$ (where in the cifar-10 case $\|\cdot\|$ implies both the Euclidean norm in $\mathbb{R}^{10}$ and integration over the input measure, which we approximated by averaging over the cifar-10 test set). We averaged this quantity over different draws of training sets to obtain $\left\langle\left\|g^{\star}_{r}-g^{\star}_{\infty}\right\|^{2}\right\rangle_{D_{N}}$. The results show good agreement between $g^{\star}_{\infty}$ and $g^{\star}_{r}$ as $r$ is increased.

Appendix C Learning Curves in the NNSP Protocol

We report here the results of a similar experiment to the one presented in the main text, but with the NTK kernel replaced with the NNGP kernel as appropriate for the NNSP correspondence. In this case we used a kernel simulating a network with a single hidden layer and $\sigma^{2}_{w}=1/W,\sigma^{2}_{b}=0$, and a target function equivalent to the one in the main text. In the NNSP protocol the renormalization group approach is not necessary to introduce noise to the observations, as it comes into play naturally via the temperature dependent fluctuations, so we can choose arbitrary $\sigma^{2}$. Notwithstanding, the renormalization group approach can aid in analyzing low temperature behavior.

Notice, in the figure below, that the sub-leading prediction significantly improves upon the EK prediction. As the inset plot demonstrates, when the dataset size is small the expected error actually increases. Surprisingly, the sub-leading correction manages to capture this behaviour even though the dataset size is small, demonstrating its superiority.

Appendix D Comparison with recent bounds

As mentioned in the main text, various recent bounds, relevant to the NTK regime, have been derived recently. Notwithstanding importance and rigor of these works, their bounds have at best a $1/\sqrt{N}$ asymptotic scaling. The figure below shows that given a functional behavior of the experimental learning curves such a bound cannot be nearly as tight as our predictions.

Appendix E NNGP and NTK are Rotationally Invariant

Let us proof that the NNGP and NTK kernels associated with any network whose first layer is fully-connected, are rotationally invariant. Indeed, let $h_{w}(x)$ be the output vector of the first layer $[h_{w}(x)]_{i}=\phi(\sum_{j}w_{ij}x_{j}+b)$ where $x_{j}$ is the $j$’th component of the input vector $x$. Let $z_{w^{\prime}}(h)$ be the output of the rest of the network given $h$. The covariance function of NNGPs are defined by [44]

[TABLE]

where $P_{0}(w,w^{\prime})$ is a prior over the weights, typically taken to be i.i.d Gaussian for each layer ($P_{0}(w,w^{\prime})=P_{0}(w)P_{0}(w^{\prime})$ and $P_{0}(w)\propto e^{-\sum_{ij}w^{2}_{ij}/(2\sigma^{2})}$). Following this one can show

[TABLE]

where the second equality uses the definition of $h_{w}(x)$, the third results from an orthogonal change of integration variable $w\rightarrow O^{T}w$, and the forth is a property of our prior over $w$. Since the NTK relates to the NNGP kernel in a recursive manner ([8]), it inherits that symmetry as well.

Appendix F Notations for the field theory derivation

For completeness, here we re-state the notations used in the main-text.

$x,x^{\prime},x^{*}$ - Inputs.

$\mu_{x}$ - Measure on input space.

$K\left(x,x^{\prime}\right)$ - Kernel function (covariance) of a Gaussian process. Assumed to be symmetric and positive-semi-definite.

$\phi_{i}\left(x\right)$ - $i$’th eigenfunction of $K\left(x,x^{\prime}\right)$. By the spectral theorem, the set $\left\{\phi_{i}\right\}_{i=1}^{\infty}$ can be assumed to be orthonormal:

[TABLE]

$\lambda_{i}$ - $i$’th eigenvalue of $K\left(x,x^{\prime}\right)$:

[TABLE]

$\|\cdot\|_{K}$ - RKHS norm:

[TABLE]

If $f\left(x\right)=\sum_{i}f_{i}\phi_{i}\left(x\right)$ then $\|f\|_{K}=\sum_{i}\frac{f_{i}^{2}}{\lambda_{i}}$ (where $\phi_{i}$ is an orthonormal set). Note that this norm is independent of $\mu_{x}$ [12].

$g\left(x\right)$ - The target function.

$\sigma^{2}$ - Noise variance.

$N$ - Number of inputs in the data-set.

$D_{N}$ - Data-set of size $N$, $D_{N}=\left\{x_{1},...,x_{N}\right\}$.

$g^{*}$ - The prediction function.

Appendix G Phrasing the Problem as a Field Theory Problem

G.1 Without Data

We start by establishing the exact equivalence between a prior of a centered GP and the corresponding partition function.

For a kernel function $K$, let us define the partition function

[TABLE]

Since the RKHS norm is quadratic in $f$, the distribution over the space of functions induced by $Z$ is Gaussian (a GP). Since a GP is determined by is mean and kernel, it is sufficient to show those equalities.

For the mean we get

[TABLE]

since for Gaussian distributions it holds that the average case is also the most probable case. For the covariance we get

[TABLE]

Indeed, $Z$ is the partition function corresponding to a centered GP with kernel $K$.

G.2 With Data

We continue by establishing the exact equivalence between Bayesian inference on a GP and the corresponding partition function.

From 30 we get that

[TABLE]

For given target function $g$ and a sampled datapoint $\left(x_{1},g(x_{i})\right)$, assuming that $f$ is our prediction it holds that $g(x_{i})\sim\mathcal{N}\left(f\left(x_{i}\right),\sigma^{2}\right)$, since $g$ distributes normally around $f$ with variance $\sigma^{2}$. Therefore, $p\left(g(x_{i})|f\right)\propto\exp\left(-\left(g(x_{i})-f\left(x_{i}\right)\right)^{2}/2\sigma^{2}\right)$, so

[TABLE]

and using Bayes’ theorem we get

[TABLE]

which gives rise to the posterior partition function

[TABLE]

and again, the exponent is quadratic in $f$ leading to a Gaussian distribution over the space of functions. Indeed, for the mean we get

[TABLE]

in agreement with [12].

G.3 Calculating Observables

G.3.1 Averaging $g^{*}$

Applying the replica trick to Eq. 37 and averaging over all the datasets of size $N$ we obtain

[TABLE]

for integer $M$ we get

[TABLE]

and after averaging

[TABLE]

where $\langle\ldots\rangle_{x\sim\mu}=\int\ldots d\mu_{x}$.

Performing the Poissonic averging we get

[TABLE]

so overall

[TABLE]

G.3.2 Averaging ${g^{*}}^{2}$

From 38 we get that

[TABLE]

Therefore

[TABLE]

Appendix H Equivalence Kernel as Free Theory

Expending the nested exponent in Eq. 41 using (first order) Taylor series we get

[TABLE]

where we defined

[TABLE]

under this approximation we get

[TABLE]

Denoting the average w.r.t $Z_{EK}$ as $\langle\dots\rangle_{0}$, The mean of the distribution induced by $Z_{EK}$ is

[TABLE]

where the last equality is due to [12].

The covariance induced by $Z_{EK}$ is

[TABLE]

where in ($*$) the non-centered part of the distribution was deleted, in ($**$) the eigenfunctions of $K$ were chosen as a base for the path integration and in ($***$) we used the fact that $\intop df\cdot f\cdot\exp\left(-\frac{1}{2}\left(\frac{1}{\lambda_{i}}+\frac{\eta}{\sigma^{2}}\right)f^{2}\right)=0$, since it is the mean of a centered (unnormalized) Gaussian distribution.

For a rotationally invariant kernel, the eigenfunctions are $Y_{lm}$ and the eigenvalues are $\lambda_{l}$ (independent of $m$) so Eq. 49 becomes

[TABLE]

which is a constant (independent of $x$ and $y$).

Appendix I Next Order Correction

We now wish to perform the first order correction to the free theory.

Expanding Eq. 41 to the next order (keeping terms up to $O\left(1/\eta^{2}\right)$) we get

[TABLE]

I.1 Calculating $\left\langle g^{*}\right\rangle_{\eta}$

We now wish to calculate the correction to $\left\langle g^{*}\right\rangle_{\eta}$ given by Eq. 51. From Eq. 42 we get

[TABLE]

Simplifying the average of the multiple sums we get

[TABLE]

Since $f$ has a Gaussian distribution ($f\sim EK$), such averages can be calculated using Feynman diagrams.

Let us denote $f(x)-g(x)$ by and $f(x_{*})$ by

. Since our free theory is not centered ($\left\langle f\right\rangle_{0}=g^{*}_{EK,\eta}\neq 0$), we allow edges in the diagrams to be connected at only one side, representing the average of the vertex w.r.t the EK distribution. An edge connected to vertices on both sides represents the covariance. Note that since we divide by $M$ and take the limit $M\to 0$, we do not care about diagrams which are not connected to $f(x_{*})$ since they scale as $M^{2}$.

Calculating the averages we get

[TABLE]

+$$+disconnected diagrams=$$=

[TABLE]

[TABLE]

+$$($$)$$($$)$$+disconnected diagrams=$$=

[TABLE]

[TABLE]

[TABLE]

Taking the limit $M\to 0$ and summing everything together we get

[TABLE]

so finally

[TABLE]

Substituting the expressions for the free variance and the covariance (Eq. 49) we get

[TABLE]

For a rotationally invariant kernel, 60 becomes

[TABLE]

I.2 Calculating $\left\langle{g^{*}}^{2}\right\rangle_{\eta}$

Substituting 39 in 44 we get

[TABLE]

By expanding to the same order we get (all equalities are up to $O\left(1/\eta^{3}\right)$)

[TABLE]

The first integrand was already calculated and is given in 58. For the second integrand we get

[TABLE]

$+$disconnected diagrams=$$=

[TABLE]

so the correction for $\left\langle{g^{*}}^{2}\left(x_{*}\right)\right\rangle_{\eta}$ is

[TABLE]

and for a rotationally invariant kernel we get

[TABLE]

Appendix J Various insights

J.1 Correction means worse generalization

The correction always means worse generalization than what the EK suggests. Indeed, expending equation (61) we get

[TABLE]

J.2 Exact eigenvalues for 2-layer ReLU NNGP and NTK with $\sigma^{2}_{b}=0$

For the NNGP associated with a 2-layer ReLU NTK without bias we were able to fined an exact expression for the eigenvalues for all $l$:

[TABLE]

and for NTK:

[TABLE]

It is interesting to note that for all odd $l>1$ $\lambda_{l}=0$ so the expressive power of the kernel (and hence the neural network) is greatly reduced.

Appendix K Accuracy of the renormalized NTK

Let us consider the random variable $t=x\cdot x^{\prime}$, for two normalized datapoints $x$ and $x^{\prime}$ drawn uniformly from the unit hypersphere $S^{d-1}$. Without loss of generality, $x$ can be assumed to be a unit vector in the direction of the last axis, and therefore $t$ is the last component of $x^{\prime}$. The density at $t\in[-1,1]$ is therefore proportional to the surface area lying at a height between $t$ and $t+dt$ on the unit sphere. That proportion occurs within a belt of height $dt$ and radius $\sqrt{1-t^{2}}$, which is a conical frustum constructed out of a $d-2$ dimensional hypersphere of radius $\sqrt{1-t^{2}}$, of height $dt$, and slope $1/\sqrt{1-t^{2}}$. Hence the probability is proportional to $p(t)dt\propto(1-t^{2})^{(d-3)/2}dt$. Defining $u=(t+1)/2$ it holds that $p(u)du\propto u^{(d-3)/2}(1-u)^{(d-3)/2}$, meaning that $u\sim\mathrm{Beta}\left(\frac{d-1}{2},\frac{d-1}{2}\right)$, and for large $d$ it holds that $\mathrm{Var}[t]=O(d^{-1})$. Since $t$ is bounded to $[-1,1]$, the random variable $t^{r}$ must have a standard deviation which is a decaying function of $r$. Indeed, for $n\ll d$ and large $d$, approximating the integral $\int t^{n}(1-t^{2})^{(d-3)/2}dt$ using saddle point approximation we get that $f(t)=n\ln(t)+\frac{d-3}{2}\ln(1-t^{2})$ is maximal for $t_{0}=[n/(n+d-3)]^{1/2}$, and $f^{\prime\prime}(t_{0})=2[n^{2}-(d-3)^{2}]/(d-3)$ so overall

[TABLE]

This implies that for $r\ll d$, the standard deviation of $t^{r}$ is $O(d^{-r/2})$. Considering next the tail of the Taylor expansion $\sum_{q>r}b_{q}(x\cdot x^{\prime})^{q}$, projected on the dataset ($\sum_{q>r}b_{q}(x_{n}\cdot x_{m})^{q}$). The resulting $N$ by $N$ matrix is $\sum_{q>r}b_{q}$ on the diagonal but $O(d^{-(r+1)/2})$ in all other entries. As we justified in the main text, our renormalization transformation amounts to keeping only the diagonal piece of this matrix and interpreting it as noise.

Consider then (5) for $g^{*}$ in two scenarios: (I) $g^{*}_{\infty}$ with the full NTK ($K(x,x^{\prime})$) and no noise and (II) $g^{*}_{r}$ with the NTK trimmed after the $r$’th power ($K_{r}(x,x^{\prime})$) but with $\sigma_{r}^{2}=\sum_{q>r}b_{q}$. The first $K(x_{*},x_{n})$ piece, for $x_{*}$ drawn from the dataset distribution, obeys $K(x_{*},x_{n})-K_{r}(x_{*},x_{n})=O(d^{-(r+1)/2})$. Next we compare $K_{r}(x_{n},x_{m})+I_{nm}\sigma_{r}^{2}$ and $K(x_{n},x_{n})$. On their diagonal they agree exactly but their off-diagonal terms agree only up to a $O(d^{-(r+1)/2})$ discrepancy. Denoting by $\delta K$ the difference between these two matrices, we may expand $K^{-1}=[K_{r}+\sigma_{m}^{2}{\rm I}+\delta K]^{-1}=[K_{r}+\sigma_{r}^{2}{\rm I}]^{-1}[1-\delta K[K_{r}+\sigma_{r}^{2}{\rm I}]^{-1}+\delta K[K_{r}+\sigma_{r}^{2}{\rm I}]^{-1}\delta K[K_{r}+\sigma_{r}^{2}{\rm I}]^{-1}+...]$.

We next argue that $\delta K[K_{r}+\sigma_{r}^{2}{\rm I}]^{-1}$ multiplied by target vector ($g(x_{n})$) is negligible compared to the identity for large enough $r$ thereby establishing the equivalence of the two scenarios. Indeed consider the eigenvalues of $\delta K[K_{r}+\sigma_{r}^{2}{\rm I}]^{-1}$. As $\delta K_{nm}$ is $O(d^{-(r+1)/2})$ its typical eigenvalues are $O(\sqrt{N}d^{-(r+1)/2})$ and bounded by $O(Nd^{-(r+1)/2})$. The typical eigenvalues of $[K_{r}+\sigma_{m}^{2}{\rm I}]^{-1}$ are of the same order as $K(x_{n},x_{n})=K$ and bounded from below by $\sigma_{r}^{2}$. Thus typical eigenvalues of $\delta K[K_{r}+\sigma_{r}^{2}{\rm I}]^{-1}$ are $O(\sqrt{N}d^{-(r+1)/2}/K)$ and bounded from above by $O(Nd^{-(r+1)/2}/\sigma_{r}^{2})$. The NTK has the desirable property that $\sigma_{r}^{2}$ decays very slowly. Thus certainly in the typical case but even in the worse case scenario we expect good agreement at large $r$. In Fig. 1, right panel, we provide supporting numerical evidence.

We refer to $K_{r}(x,x^{\prime})$ as the renormalized NTKs at the scale $r$. As follows from (26), $\lambda_{l}$’s with $l\geq r$ are zero. Therefore, as advertised, the high-energy-sector has been removed and compensated by noise on the target and a change of the remaining $l<r$ (low-energy) eigenvalues. A proper choice of $r$ involve two considerations. Requiring perturbation theory to hold well ($C_{K_{r},\sigma_{r}^{2}/\eta}<\sigma_{r}^{2}$) which puts an $\eta$-depended upper bound on $r$ and requiring small discrepancy in predictions puts another $\eta$ dependent lower bound on $r$ (typically $\sqrt{N}d^{-(r+1)/2}\ll 1$).

Lastly we comment that our renormalization NTK approach is not limited to uniform datasets. The entire logic relies on having a rapidly decaying ratio of off-diagonal moments ($(x_{n}\cdot x_{m})^{2r}$) and diagonal moments $(x_{n}\cdot x_{n})^{2r}$ as one increases $r$. We expect this to hold in real-world distributions. For instance for a multi-dimension Gaussian data distribution the input dimension ($d$) traded by an effective dimension ($d_{eff}$) defined by the variance of $(x_{m}\cdot x_{n})$.

Appendix L Hyper-parameter optimization experiment results