Title: ASYMPTOTICS OF WIDE NETWORKS FROM FEYNMAN DIAGRAMS

Abstract: Understanding the asymptotic behavior of wide networks is of considerable interest. In this work, we present a general method for analyzing this large width behavior. The method is an adaptation of Feynman diagrams, a standard tool for computing multivariate Gaussian integrals. We apply our method to study training dynamics, improving existing bounds and deriving new results on wide network evolution during stochastic gradient descent. Going beyond the strict large width limit, we present closed-form expressions for higher-order terms governing wide network training, and test these predictions empirically.

1. INTRODUCTION:

Neural networks achieve remarkable performance on a wide array of machine learning tasks, yet a complete analytic understanding of deep networks remains elusive. One promising approach is to consider the large width limit, in which the number of neurons in one or several layers is taken to be large. In this limit one can use a mean-field approach to better understand the network's properties at initialization (Neal (1996) ; Lee et al. (2018) ), as well as its training dynamics (Daniely (2017) ; Jacot et al. (2018) ). Additional related works are cited below.

Suppose that f (x) is the network function evaluated at an input x. Let us denote the vector of model parameters by θ, whose elements are initially chosen to be i.i.d. Gaussian. In this work we consider a class of functions we call correlation functions, obtained by taking the ensemble averages of f , its products, and its derivatives with respect to the parameters θ, evaluated on arbitrary inputs. Here are a few examples of correlation functions.

E θ [f (x 1 )f (x 2 )] , µ E θ ∂f (x 1 ) ∂θ µ ∂f (x 2 ) ∂θ µ , µ,ν E θ ∂f (x 1 ) ∂θ µ ∂f (x 2 ) ∂θ ν ∂ 2 f (x 3 ) ∂θ µ ∂θ ν f (x 4 ) . (1)

Correlation functions often show up in the study of wide networks. For example, the first correlation function in (1) plays a central role in the Gaussian Process picture of wide networks (Lee et al. (2018) ), and has been used to diagnose signal propagation in wide networks (Pennington et al. (2017) ). The second example in (1) is the ensemble average of the Neural Tangent Kernel (NTK), which controls the evolution of wide networks under gradient flow (Jacot et al. (2018) ), and the third example shows up when computing the time derivative of the NTK with MSE loss.

While correlation functions can be computed analytically in some special cases (Cho & Saul (2009) ), they are not analytically tractable in general. In this work, we present a method for bounding the asymptotic behavior of such functions at large width. In many cases of interest, the bound we obtain is tight. Derivation of the method relies on Feynman diagrams (Feynman (1949) ), a technique for calculating multivariate Gaussian integrals, and specifically on the 't Hooft expansion ('t Hooft (1974) ). However, applying the method is straightforward and does not require any knowledge of Feynman diagrams.

Motivation. The method presented here is a general-purpose tool that can be added to a researcher's toolbox. It allows a researcher to quickly identify ways in which network behavior and training dynamics simplify at large width, for example by deriving which quantities vanish in this limit. As we demonstrate, using this method one is able to dramatically cut down on the amount of effort required to derive several key results related to wide networks, as well as to go beyond the infinite width limit and investigate ordinary neural networks.

Our contribution.

1. We present a general method for bounding the asymptotic behavior of correlation functions.

The method is an adaptation of Feynman diagrams to the case of wide neural networks. The adaptation involves a novel treatment of derivatives of the network function, an element that is not present in the original theoretical physics formulation. 2. We apply the method to the study of wide network evolution under gradient descent. We improve on existing results for gradient flow (Jacot et al. (2018) ) by deriving tighter bounds, and extending the analysis to the case of stochastic gradient descent (SGD). Going beyond the infinite-width limit, we present a formalisn for deriving finite-width corrections to network evolution, and present explicit formulas for the first order correction. To our knowledge, this is the first time this correction has been calculated. 3. As additional applications of our method, in Appendix E.2 we show that in the large width limit the SGD updates are linear in the learning rate, and in Appendix E.3 we discuss finite width corrections to the spectrum of the Hessian.

Limitations of our approach. The main result of this paper is a conjecture. We test our predictions extensively using numerical experiments, and prove the conjecture in some cases, but we do not have a proof that applies to all the cases we tested, including for deep networks with general non-linearities. Furthermore, our method can only be used to derive asymptotic bounds at large width; it does not produce the width-independent coefficient, which is often of interest. et al. (2018) . For additional works discussing the training dynamics of wide networks see Geiger et al. (2019a) ; Arora et al. (2019) . For a previous use of diagrams in this context, see Pennington & Worah (2017) . The Neural Tangent Hierarchy presented in Huang & Yau (2019) , published during the completion of this version, has significant overlap with the recursive differential equations (10) presented below.

The rest of the paper is organized as follows. In Section 2 we present our main conjecture and supporting evidence. In Section 3 we apply the method to gradient descent evolution of wide networks, and in Section 4 we present details on Feynman diagrams, which is the basic technique used in our proofs. We conclude with a Discussion. Proofs, additional applications, and details can be found in the Appendices.

Note: An earlier version of this work appeared in the ICML 2019 workshop, Theoretical Physics for Deep Learning (Dyer & Gur-Ari (2019) ).

2. CORRELATION FUNCTION ASYMPTOTICS:

In this section we present our main result: a method for computing asymptotic bounds on correlation functions of wide networks. We present the result as a conjecture, supported by analytic and empirical evidence.

2.1. NOTATION:

Let f (x) ∈ R be the network output of a deep network with d hidden layers and input x ∈ R Din , defined by

EQUATION

Here U ∈ R n×Din , V ∈ R n , W 1 , . . . , W d-1 are weight matrices of dimension n, and σ : R → R is the non-linearity. 1 We denote the vector of all model parameters by θ. At initialization, the elements of θ are independent Gaussian variables, with each element θ µ ∼ N (0, 1). The corresponding distribution of θ is denoted by P 0 . 2

Let us now define correlation functions, the class of functions that is the focus of this work. These functions involve derivative tensors of the network function. We denote the rank-k derivative tensor by T µ1...µ k (x; f ) := ∂ k f (x)/∂θ µ1 • • • ∂θ µ k . For k = 0 we define T (x; f ) := f (x), and still refer to this as a derivative tensor for consistency. Definition 1. A correlation function is the expectation value of a product of derivative tensors, evaluated at arbitrary inputs, where the tensor indices are summed in pairs over all the model parameters. A general correlation function C takes the form C(x 1 , . . . , x m ) := µ1,...,µ km

∆ (π) µ1...µ km E θ T µ1...µ k 1 (x 1 )T µ k 1 +1 ...µ k 2 (x 2 ) • • • T µ k m-1 +1 ...µ km (x m ) .

(3) 3 m and k m are even, π ∈ S km is a permutation, and ∆

Here, 0 ≤ k 1 ≤ • • • ≤ k m-1 ≤ k m are integers,

(π) µ1...µ km = δ µ π(1) µ π(2) • • • δ µ π(km-1) µ π(km ) . We use δ to denote the Kronecker delta.

If two derivative tensors in a correlation function have matching indices that are summed over, we say that they are contracted. For example, the correlation function µ E θ [∂f (x 1 )/∂θ µ • ∂f (x 2 )/∂θ µ ] has one pair of contracted tensors. See (1) for additional examples of correlation functions.

2.2. ASYMPTOTIC BOUNDS ON WIDE NETWORKS:

We now present our main conjecture, which allows us to place asymptotic bounds on general correlation functions of wide networks. Conjecture 1. Let C(x 1 , . . . , x m ) be a correlation function. The cluster graph G C (V, E) of C is a graph with vertices V = {v 1 , . . . , v m }, where v i = T (x i ) for all i, and edges E = {(v i , v j ) | (T (x i ), T (x j )) contracted in C}. Suppose that the cluster graph G C has n e connected components with an even size (even number of vertices), and n o components of odd size. Then C(x 1 , . . . , x m ) = O(n s C ), where

EQUATION

)

We will refer to the connected components of a cluster graph G C as the clusters of C. Table 1 lists examples of bounds derived using the Conjecture for several correlation functions. The intuition behind Conjecture 1 comes from the following result for deep linear networks. Theorem 1. Conjecture 1 holds for correlation functions of networks with linear activations, including for deep linear networks with biases (affine transformations).

Let us discuss the intuition behind this theorem. Computing correlation functions of deep linear networks amounts to evaluating Gaussian integrals with polynomial integrands in θ. One can evaluate such integrals using Isserlis' theorem, which tells us how to express moments of multivariate Gaussian variables in terms of their second moments. For example, given centered Gaussian variables z 1 , ..., z 4 ,

E z [z 1 z 2 z 3 z 4 ] = E z [z 1 z 2 ] E z [z 3 z 4 ] + E z [z 1 z 3 ] E z [z 2 z 4 ] + E z [z 1 z 4 ] E z [z 2 z 3 ] .

(5) Therefore, correlation functions of deep linear networks can be expressed in terms of the covariances

E θ [U iα U jβ ] = δ ij δ αβ , E θ [V i V j ] = δ ij , and E θ W (l) ij W (l) kl

= δ ik δ jl . For example, for a deep linear network with 2 hidden layers, we have

E θ [f (x 1 )f (x 2 )] = 1 n 2 E θ V T W U x 1 V T W U x 2 = x T 1 x 2 n 2 n i,k δ ik δ ik n j,l

Every correlation function of a deep linear network can be similarly reduced to sums over products of Kronecker delta functions and width-independent functions of the inputs. The asymptotic large width behavior is determined by these sums over delta functions, which are tedious to compute by hand. Feynman diagrams are a graphical tool for computing these sums, allowing us to obtain the asymptotic behavior with minimal effort. This tool, which is described in detail in Section 4, is used to prove Theorem 1.

For networks with non-linear activations we further show the following Theorem 2. Conjecture 1 holds for (1) networks with ReLU activations, where all inputs are set to be equal, and for (2) networks with one hidden layer and smooth activation.

For case (1), the idea behind the proof is to put an asymptotic bound on the ReLU network in terms of a corresponding deep linear network. For case (2), the basic idea is that each network function contains a single sum over the width, and by keeping track of these sums using Feynman diagrams we are able to bound the asymptotic behavior. We refer the reader to Appendix C for details.

2.3. NUMERICAL EXPERIMENTS:

Table 1 lists asymptotic bounds on several correlation functions, derived using Conjecture 1. These are compared against the asymptotic behavior computed using numerical experiments. In addition to the results presented here, we performed experiments using the same correlation functions and experimental setup, but with weights sampled uniformly from {±1} instead of from a Gaussian distribution. The results are shown in Appendix A.1.

In all cases that were tested empirically, we found that Conjecture 1 holds. For networks with smooth activations, we found that Conjecture 1 always gives a tight bound. For networks with linear or ReLU activations, we always find that the Conjecture holds as an upper bound, but that sometimes the bound is not tight. One such case is highlighted in Table 1 . In such cases, a tight bound can be obtained for deep linear networks using the complete Feynman diagram analysis presented below. 4Correlation function C n e , n o s C lin.

ReLU tanh

E θ [f (x 1 )f (x 2 )] 0,2 0 -0.02 0.003 -0.02 E θ [f (x 1 )f (x 2 )f (x 3 )f (x 4 )] 0,4 0 -0.01 0.03 -0.03 µ E θ [∂ µ f (x 1 )∂ µ f (x 2 )]

1,0 0 0.00 0.00 0.00

µ,ν E θ [∂ µ f (x 1 )∂ ν f (x 2 )∂ µ,ν f (x 3 )f (x 4 )] 0,2 -1 -0.98 -1.03 -1.01 µ,ν,ρ E θ [∂ µ f (x 1 )∂ ν f (x 2 )∂ ρ f (x 3 )∂ µ,ν,ρ f (x 4 )] 1,0 -1 -2.01 * -2.01 * -0.98 µ,ν,ρ,σ E θ ∂µf (x 1 )∂ν f (x 2 )∂µ,ν f (x 3 )∂ρf (x 4 )∂σ f (x 5 )∂ρ,σ f (x 6 ) 0,2 -2 -2.05 -2.01 -1.99

Table 1 : Examples of bounds on correlation functions obtained from Conjecture 1. The 3 rightmost columns list numerical results for fully-connected networks with 3 hidden layers and with linear, ReLU, and tanh activations. The values listed in these columns are the fitted exponents for networks with the corresponding activations. The entries marked with an asterisk are those for which the predicted bound is not tight. The numerical results are obtained by computing the correlation functions for networks with widths 2 7 , 2 8 , . . . , 2 13 , each averaged over 1,000 initializations, and fitting the exponent. Inputs are chosen to be random vectors of dimension 4.

3. APPLICATIONS TO TRAINING DYNAMICS:

In this section we apply Conjecture 1 to study the evolution of wide networks under gradient flow and gradient descent. We begin by briefly reviewing existing results. Let D tr be a training set of size M , and let L = (x,y)∈Dtr (x, y) be the MSE loss, with single sample loss (x, y) = 1 2 (f (x)y) 2 . The gradient flow equation is dθ dt = -∇ θ L. The evolution of the network function under gradient flow is given by

EQUATION

Here, Θ is the Neural Tangent Kernel (NTK), defined by Θ(x 1 , x 2 ) := ∇ θ f T (x 1 )∇ θ f (x 2 ). The authors of Jacot et al. (2018) showed that the kernel is constant during training up to O(n -1/2 ) corrections. This leads to a dramatic simplification in training dynamics (Jacot et al. (2018) ; Lee et al. (2019) ). In particular, for MSE loss the network map evaluated on the training data evolves as f (t) = y + e -tΘ (0) (f (0)y). 5 We will use our technology to derive a tighter bound on finite-width corrections to the kernel during training and present explicit formulas for the leading correction.

The following result is useful in analyzing the behavior of correlation functions under gradient flow.

Lemma 1. Let C( x) = E θ [F ( x)

] be a correlation function, where x = (x 1 , . . . , x m ), and suppose

that C = O(n s C ) for s C as defined in Conjecture 1. Then E θ d k F ( x) dt k = O(n s C ) for all k.

Here we prove the statement for k = 1. Appendix D.2 contains a proof for the general case.

proof (k = 1). Let C have n e even clusters and With this result, it is easy to understand the constancy of the NTK at large width. The first derivative of the NTK is given by

n o odd clusters. Consider E θ dF ( x) dt = -µ x ∈Dtr E θ ∂F ( x) ∂θ µ ∂f (x ) ∂θ µ f (x ) .

E θ dΘ(x 1 , x 2 ) dt = - M x ∈Dtr µ,ν E θ ∂ 2 f (x 1 ) ∂θ µ ∂θ ν ∂f (x 2 ) ∂θ µ ∂f (x ) ∂θ ν f (x ) + (x 1 ↔ x 2 ) . (8)

Here, (x 1 ↔ x 2 ) means "add the same term, but exchange x 1 and x 2 ". This correlation function has n e = 0, n o = 2, and m = 4, and is therefore O(n -1 ) by Conjecture 1. By Lemma 1, all higher-order time derivatives of the NTK are O(n -1 ) as well. If we now assume that the time-evolved kernel Θ(t) is analytic in training time t, and that we are free to exchange the Taylor expansion in time with the large width limit, then we find that

E θ [Θ(t) -Θ(0)] = ∞ k=1 t k k! E θ d k Θ(0) dt k = O(n -1

) for any fixed t. This bound, which is tighter than that found in Jacot et al. (2018) , was noticed empirically in Lee et al. (2019) as well as in our own experiments, see Figure 1a .

This analysis can be extended to show that E θ [Θ(t) -Θ(0)] = O(n -1 ) for SGD as well. The technique is similar, and again relies on Conjecture 1. We refer the reader to Appendix E.1 for details. These results improve on existing ones in several ways. Our method applies to the case of SGD as well as to networks where all layer widths are increased simultaneously -a setup that has proven to be difficult to analyze. In addition, the O(n -1 ) bound we derive on kernel corrections during SGD is empirically tight, and improves on the existing bound of O(n -1/2 ) which was derived for gradient flow (Jacot et al. (2018) ).

3.1. FINITE WIDTH CORRECTIONS:

Next, we will compute the explicit time dependence of Θ and f at order O(n -1 ) under gradient flow. This is the leading correction to the infinite width result. We define the functions O 1 (x) := f (x) and

EQUATION

)

Notice that O 2 = Θ is the kernel. It is easy to check that

EQUATION

)

We can use equations ( 9) and ( 10) to solve for the evolution of the kernel and network map. Notice that each function O s has s derivative tensors and a single cluster. As a result, correlation functions involving operators with larger s are increasingly suppressed in width. In particular,

E θ [dO 4 /dt] = -x∈Dtr E θ [O 5 (x)f (x)] = O(n -2

) at all times, using Conjecture 1 and Lemma 1. Thus to solve for f and Θ at O(n -1 ) we can set O s = 0 for all s ≥ 5.

Let us denote the time-evolved kernel by Θ(t) = Θ (0) + Θ 1 (t) + O(n -2 ), where Θ (0) is the kernel at initialization, and Θ 1 (t) is the O(n -1 ) correction we are seeking. Integrating equations (10) starting with s = 4, we find

EQUATION

Here we have introduced the notation ∆f (x; t) = e -tΘ0 (f (0)y). A detailed derivation can be found in Appendix E.4. There we also evaluate the integrals in (11) in terms of the NTK spectrum.

To obtain the O(n -1 ) correction to the network map (evaluated for simplicity on the training data), we further integrate (10) for s = 1 and find

EQUATION

Here we have denote the infinite width evolution by f 0 (t) = y + e -tΘ (0) (f (0)y). Figures 1b and 1c compare these predictions against empirical results. 

4. FEYNMAN DIAGRAMS FOR DEEP LINEAR NETWORKS:

In this section we present the Feynman diagram technique, and show how it allows us to compute the asymptotic behavior of correlation functions. We end this section with a proof of Theorem 1 for the case of networks with a single hidden layer and linear activations.

Given a correlation function C, we map it to a family of graphs called Feynman diagrams. The graphs are independent of the inputs, and are defined as follows.

Definition 2. Let C(x 1 , . . . , x m ) be a correlation function for a network with d hidden layers. The family Γ(C) is the set of all graphs that have the following properties.

1. There are m vertices v 1 , . . . , v m , each of degree d + 1.

2. Each edge has a type t ∈ {U, W 1 , . . . , W d-1 , V }. Every vertex has one edge of each type.

3. If two derivative tensors T µ1,...,µ (x i ), T ν1,...,ν (x j ) are contracted k times in C, the graph must have at least k edges (of any type) connecting the vertices v i , v j .

The graphs in Γ(C) are called the Feynman diagrams of C.

Single hidden layer. For the rest of this section we focus on networks with a single hidden layer. We refer the reader to Appendix B for a full treatment of deep linear networks. For networks with one hidden layer and linear activation, the network output is f

EQUATION

)

As the factors of x 1 and x 2 are independent of n, we see that C(x 1 , x 2 ) = O(n 0 ). Notice that there are two relevant contributions to this answer: each factor of the network function in the integrand contributes n -1/2 , and the summed-over product of Kronecker deltas contributes n. Other details, such as the input dependence, are irrelevant. Feynman diagrams allow us to encode only those details that affect the n scaling, ignoring the rest.

The set Γ(C) for the correlation function (13) consists of a single Feynman diagram, shown in Figure 2a . The asymptotic bound on a correlation function is obtained by the following result, which is due to 't Hooft (1974) . Theorem 3. Let C(x 1 , . . . , x m ) be a correlation function with one hidden layer and linear activation.

Then C = O(n s ) where s = max γ∈Γ(C) l γ -m 2 , and l γ is the number of loops in γ.7 

Let us give some intuition for Theorem 3. Each Feynman diagram γ encodes a subset of the terms contributing to the correlation function. To get the asymptotic bound on the correlation function, we sum over the contributions of individual diagrams. We can compute the asymptotic behavior of a single diagram γ using the following Feynman rules: (1) each vertex contributes a factor of n -1/2 , and (2) each loop contributes a factor of n. Therefore, if a diagram has l γ loops, its contribution to the correlation function scales as n lγ -m/2 . Rule (1) is due to the explicit n -1/2 factor in the network definition. Rule (2) follows from applying Isserlis' theorem, as follows. Each covariance factor (such as the factor E [V i V j ] = δ ij in eq. ( 13)) corresponds to an edge in a Feynman diagram. A loop in the diagram corresponds to a sum over a product of Kronecker deltas, yielding a factor of n.

Returning to our example, the graph in Figure 2a has two vertices and one closed loop, so we recover the asymptotic behavior O(n 0 ) from Theorem 3. As another example, consider the correlation function

C(x 1 , x 2 , x 3 , x 4 ) = E θ [f (x 1 )f (x 2 )f (x 3 )f (x 4 )].

It has two Feynman diagrams, shown in Figures 2b and 2c . Using the Feynman rules, we find that the disconnected graph represents terms that scale as O(n 0 ), while the connected graph represents terms that scale as O(n -1 ). Therefore,

C(x 1 , x 2 , x 3 , x 4 ) = O(n 0

). This is an example of a more general phenomenon, that connected graphs vanish faster at large n compared with disconnected graphs.

Correlation functions with derivatives. We now extend the Feynman diagram technique to correlation functions that include derivatives of f .8 As a concrete example, consider the correlation function

C(x, x ) = E θ [Θ(x, x )],

where Θ is the kernel defined in Section 3. For a single hidden layer, we have

EQUATION

)

The two derivative tensors in this correlation function are contracted: their indices are set to be equal and summed over. Therefore, according to Definition 2, Γ(C) includes all diagrams in which the corresponding vertices share at least one edge. The resulting diagrams are shown in Figure 3 . The edges forced by the contraction are explicitly marked by dashed lines for clarity, but mathematically they are ordinary edges. We see that in fact there is only one diagram contributing to this correlation function -the same one shown in Figure 2a . Following the Feynman rules, we find that

E θ [Θ(x, x )] = O(n 0 ).

The fact that contracted derivatives should be mapped to forced edges in the Feynman diagrams is proved in Appendix B. The basic reason behind this rule is the relation

k ∂Vi ∂V k ∂Vj ∂V k = δ ij = E [V i V j ].

Namely, when derivatives act in pairs they yield a Kronecker delta factor (δ ij ), which is equal to the factor obtained from a covariance (E [V i V j ]). While Isserlis' theorem instructs us to sum over all possible covariance configurations (and therefore over all possible edge configurations), a pair of summed derivative leads to a particular covariance factor. Therefore, we should only consider graphs that include the edge corresponding to this covariance factor.

We are now ready to prove Theorem 1 for the case of single hidden layer with linear activations. A proof for the general case can be found in Appendix B.

Proof (Theorem 1, one hidden layer). Let C(x 1 , . . . , x m ) be a correlation function for a network with a single hidden layer and linear activation. Let G C (V, E) be the cluster graph of C, and let γ(V, E ) ∈ Γ(C) be a Feynman diagram. Notice that E ⊂ E . Indeed, E contains an edge corresponding to each pair of contracted derivative tensors in C, and

E = {(v i , v j ) | (T (x i ), T (x j )) contracted in C}.

In addition, notice that γ only contains connected components (i.e. loops) with an even number of vertices, because every vertex has exactly one edge of each type. Therefore, n γ ≤ n e + no 2 , where n γ is the number of loops in γ, and n e (n o ) is the number of even (odd) components in G C . The bound is saturated when every even component of G C belongs to a different component of γ, and pairs of odd components in G C belong to different components of γ. From Theorem 3, we have that

C = O(n s ) where s = max γ∈Γ(C) n γ -m 2 ≤ n e + no 2 -m 2 .

5. DISCUSSION:

Ensemble averages of the network function and its derivatives are an important class of functions that often show up in the study of wide neural networks. Examples include the ensemble average of the train and test losses, the covariance of the network function, and the Neural Tangent Kernel (Jacot et al. (2018) ). In this work we presented Conjecture 1, which allows one to derive the asymptotic behavior of such functions at large width.

For the case of deep linear networks, we presented a complete analytic understanding of the Conjecture based on Feynman diagrams. In addition, we presented empirical and analytic evidence showing that the Conjecture also holds for deep networks with non-linear activations, as well as for networks with non-Gaussian initialization. We found that the Conjecture holds in all cases we tested.

The basic tools presented in this work can be applied to many aspects of wide network research, greatly simplifying theoretical calculations. We presented several applications of our method to the asymptotic behavior of wide networks during stochastic gradient descent, and additional applications are presented in Appendix E. We were able to improve upon known results by tightening existing bounds, and by applying the technique to SGD as well as to gradient flow. In addition, we took a step beyond the infinite width limit, deriving closed-form expressions for the first finite-width correction to the network evolution. These novel results open the door to studying finite-width networks by systematically expanding around the infinite width limit.

A central question in the study of wide networks is whether the infinite width limit is a good model for describing the behavior of realistic deep networks (Chizat et al. (2018) ; Ghorbani et al. (2019) ; Geiger et al. (2019b) ). In this work we take a step toward answering this question, by working out the next order in a perturbative expansion around the infinite width limit, potentially bringing us closer to an analytic description of finite-width networks. We hope that the techniques presented here provide a basis to systematically answering these and other questions about the behavior of wide networks.

A EXPERIMENTAL DETAILS AND ADDITIONAL RESULTS:



A.1 NON-GAUSSIAN INITIALIZATION:

In Table 1 we listed asymptotic bounds on several correlation functions, where the model parameters were initialized from a Gaussian distribution. Table 2 shows additional results using the same experimental setup, but with weights sampled uniformly from {±1}. We again find good agreement with the predictions of Conjecture 1.

Correlation function C n e , n o s C lin.

ReLU tanh

E θ [f (x 1 )f (x 2 )] 0,2 0 0.00 0.00 0.02 E θ [f (x 1 )f (x 2 )f (x 3 )f (x 4 )] 0,4 0 -0.07 0.06 -0.02 µ E θ [∂ µ f (x 1 )∂ µ f (x 2 )] 1,0 0 0.00 0.00 0.00 µ,ν E θ [∂ µ f (x 1 )∂ ν f (x 2 )∂ µ,ν f (x 3 )f (x 4 )] 0,2 -1 -1.02 -1.01 -0.97 µ,ν,ρ E θ [∂ µ f (x 1 )∂ ν f (x 2 )∂ ρ f (x 3 )∂ µ,ν,ρ f (x 4 )] 1,0 -1 -2.00 -1.99 -2.02 µ,ν,ρ,σ E θ ∂µf (x 1 )∂ν f (x 2 )∂µ,ν f (x 3 )∂ρf (x 4 )∂σ f (x 5 )∂ρ,σ f (x 6 ) 0,2 -2 -2.05 -2.01 -1.99

Table 2 : Examples of bounds on correlation functions obtained from Conjecture 1. The experimental setup is the same as in Table 1 , but the model parameters are sampled uniformly from {±1} instead of from a Gaussian distribution. We find good agreement with the theoretical predictoins, and in many cases the bound is tight.

A.2 EXPERIMENTAL DETAILS:

The experiments in Figure 1 

EQUATION

B FEYNMAN DIAGRAMS FOR DEEP LINEAR NETWORKS:

Feynman diagrams can be used to derive asymptotic upper bounds on deep linear networks in the large width limit. In this section we describe the method in detail, and use it to prove Theorem 1.

B.1 FEYNMAN DIAGRAMS AND DOUBLE-LINE DIAGRAMS:

In this section we build on the results of Section 4 and consider correlation functions of deep linear networks with d hidden layers. The network function was defined in (2), and here we set the activation σ to be the identity. Definition 2 describes how to map a correlation function C to Γ(C), a family of graphs called Feynman diagrams. The Feynman diagram method relies on Isserlis' theorem, which allows us to express arbitrary moments of multivariate Gaussian variables in terms of their covariance.

Theorem 4 (Isserlis). Let z = (z 1 , . . . , z l ) be a centered multivariate Gaussian variable. For any positive integer k,

EQUATION

)

In particular, if the covariance matrix of z is the identity then

EQUATION

)

Using this theorem, a correlation function C can be expressed as a sum over permutations as in ( 16). Each term in this sum maps to a Feynman diagram in Γ(C). As an example, consider the correlation function

C(x, x ) = E θ [f (x)f (x )]

for a network with 2 hidden layers. An explicit calculation gives

EQUATION

)

To get from ( 19) to (20), we applied Isserlis' theorem for every choice of indices i, j, k, l, α, β. We find that there is at most one permutation π of the network parameters such that the covariances do not vanish. This is because the covariance of parameters across different layers vanishes identically (for example E θ V i W 1 jk = 0). Correspondingly, this correlation function has a single Feynman diagram, shown in Figure 4a .

For networks with one hidden layer, we saw in Section 4 that the asymptotic behavior is determined by the number of loops in a graph. This observation does not immediately generalize to networks with general depth, because it is not obvious how to count loops in diagrams such as the one in Figure 4a . The problem can be traced back to the fact that weight matrices have covariances of the form E [W ij W kl ] = δ ik δ jl involving two Kronecker deltas, but the procedure described so far assumes that each covariance (and each edge in the graph) corresponds to a single Kronecker delta.

A similar question appeared in the context of theoretical physics, and the correct generalization is due to 't Hooft (1974) . The idea is to treat each Feynman diagram as the triangulation of a Riemann surface, and to define the number of loops in a graph to be the number of faces of the triangulation. In practice, this involves mapping each Feynman diagram to a new double-line diagram: A graph in which each edge corresponds to a single Kronecker delta factor, and loops correspond to triangulation faces of the original diagram. Definition 3. Let γ ∈ Γ(C) be a Feynman diagram for a correlation function C involving k derivative tensors for a network of depth d. Its double-line graph, DL(γ) is a graph with kd vertices of degree 2, defined by the following blow-up procedure.

• Each vertex v (i) in γ is mapped to d vertices v (i) 1 , . . . , v (i) d in DL(γ). • Each edge (v (i) , v (j) ) in γ of type U is mapped to a single edge (v (i) 1 , v (j) 1 ). • Each edge (v (i) , v (j) ) in γ of type W l is mapped to two edges (v (i) l , v (j) l ), (v (i) l+1 , v (j) l+1 ). • Each edge (v (i) , v (j) ) in γ of type V is mapped to a single edge (v (i) d , v (j) d ).

The number of faces in γ is given by the number of loops in the double-line graph DL(γ). Generally, every Feynman diagram of a deep linear network can be interpreted as the triangulation of a 2-dimensional manifold with at least one boundary. Intuitively, the presence of the boundary is due to the fact that the U, V weights at both ends of network have only one dimension that scales linearly with the width. As a result, in the double-line diagrams the U, V edges become single lines rather than double lines. These 'missing' lines translate into fewer faces in the triangulation, and the 'missing' faces can be geometrically interpreted as boundaries in the corresponding surface.9 

The discussion above is summarized by the following result, due to 't Hooft (1974) , that describes how the asymptotic behavior of a general correlation function can be computed using the Feynman rules for deep linear networks. The intuition for the formula s γ = l γ -dm 2 is similar to the single hidden layer case, Theorem 3. The term l γ counts the number of factors of the form i1,...,i k δ i1i2 δ i2i3 • • • δ i k i1 = n that appear in the Correlation function after applying Isserlis' theorem. The term -dm 2 is due to the explicit n -d/2 normalization of the network function.

B.2 ASYMPTOTICS OF DEEP LINEAR NETWORKS:

We now prove Theorem 1. The theorem follows from the following lemma, again due to 't Hooft (1974) , that relates the asymptotic behavior to the number of connected components in a Feynman diagram. Lemma 2. Let C(x 1 , . . . , x m ) be a correlation function for a deep linear network. Let c γ be the number of connected components of a graph γ ∈ Γ(C). Then C = O(n s ), where

EQUATION

)

Proof. We prove the result for 1D inputs (D in = 1), and it is easy to generalize to arbitrary input dimension. Let γ ∈ Γ(C) be a Feynman diagram and let γ be a connected component of γ with v γ vertices. Notice that we can apply the Feynman rules of Theorem 5 separately to each component γ and find a bound O(n s γ ). Then, γ = O(n sγ ) where s γ = γ s γ , and the sum runs over the connected components of s γ . We will show below that s γ ≤ 1 -

v γ

2 , and therefore

s γ ≤ c γ -m 2 , which is sufficient to prove (22).

Let us prove the remaining statement about s γ . The Euler character of the graph γ with v = v γ vertices, e edges and f faces is χ = ve + f . The degree of each vertex in the graph is d + 1, and therefore e = (d+1)v 2 . Using Theorem 5 the graph is O(n s γ ) where s γ = f -dv 2 . We therefore find that χ = v 2 + s γ . The graph γ is a triangulation of some connected surface with at least one boundary. The Euler character for such a surface is bounded by χ ≤ 1, and therefore

s γ ≤ 1-v 2 .

Let us now prove Theorem 1.

Proof. Let C(x 1 , . . . , x m ) be a correlation function for a deep linear network. Suppose that the cluster graph G C has n e even size components and n o odd size components. Let γ ∈ Γ(C) be a Feynman diagram with c γ connected components. We will show that c γ ≤ n e + no 2 . It then follows immediately from Lemma 2 that C = O(n s ) where s = n e + no 2 -m 2 , concluding the proof. Let us derive the bound on c γ . First, all vertices that belong to a given cluster (a component of G C ) will also belong to the same connected component in γ. This is because every edge in G C is also an edge in γ (note that G C and γ have the same set of vertices). Therefore, c γ ≤ n e + n o . Second, note that every connected component of the graph γ has an even number of vertices. Indeed, each edge has a type t, and each vertex has exactly one edge of each type. Therefore, a connected component with v vertices has v 2 edges of each type, and so v must be even. It follows that the vertices of even clusters can form their own connected components in a Feynman diagrams, while odd clusters must be connected in sets of 2 or more to form connected components. The bound on c γ then follows.

For deep linear networks with bias the proof is the same, except we use Lemma 3 below instead of Lemma 2.

B.3 DEEP LINEAR NETWORKS WITH BIAS:

In this section we argue that Conjecture 1 holds for deep linear networks with bias. Let us add bias to the definition (2) of a deep linear network (setting σ to be the identity). We choose x = 1 without loss of generality. The activations z l at layer l and the network function f are given by

EQUATION

EQUATION

EQUATION

)

Here b 1 , . . . , b d-1 ∈ R n are biases whose elements are chosen independently from N (0, 1) at initialization. The network function can be written as 1. There are m vertices v 1 , . . . , v m , where each v i has degree k i ∈ {2, . . . , d + 1}.

EQUATION

2. Each edge has a type t. Every vertex of degree k has one edge of each type t ∈ {b d+1-k , W d+2-k , . . . , W d-1 , V }.

3. If two derivative tensors T µ1,...,µ (x i ), T ν1,...,ν (x j ) are contracted s times in C, the graph must have at least s edges (of any type) connecting the vertices v i , v j .

The intuition behind this definition is that term k in the sum in ( 26) is equal to the network function of a deep linear network with depth k -1. Correlation functions of such a function can be represented diagramatically by including vertices of degree k. To represent the entire sum we must consider all diagrams with vertices of varying degree.

We now extend Definition 3 of double-line graphs to accommodate vertices with varying degree. Definition 5. Let γ ∈ Γ(C) be a Feynman diagram for a correlation function C involving m derivative tensors for a network of depth d with biases. Its double-line graph, DL(γ) is a graph with vertices of degree 2, defined by the following blow-up procedure.

EQUATION

ki-1 in DL(γ).

• Each edge (v (i) , v (j) ) in γ of type W l is mapped to two edges, (v

(i) l+ki-d-1 , v (j) l+kj -d-1 ) and (v (i) l+ki-d , v (j) l+kj -d ). • Each edge (v (i) , v (j) ) in γ of type b d+1-ki is mapped to a single edge (v (i) 1 , v (j) 1 ) (notice that k i = k j ). • Each edge (v (i) , v (j) ) in γ of type V is mapped to a single edge (v (i) ki-1 , v (j) kj -1 ).

As in Definition 3, the number of faces in γ is given by the number of loops in the double-line graph DL(γ). 

E θ [f (x 1 )f (x 2 )f (x 3 )f (x 4 )]

for a two-hidden-layer network with biases.

Let us now extend Theorem 5 (the Feynman rules) to accomodate bias. The only difference is in the explicit normalization. Previously, each vertex had degree d + 1 and contributed a factor of n -d/2 due to the explicit normalization of the function map. With biases, as a result of (26), a vertex with degree k contributes a factor of n (k-1)/2 . Theorem 6. Let C(x 1 , . . . , x m ) be a correlation function of a deep linear network with d hidden layers and with biases. Let γ ∈ Γ(C) be a Feynman diagram. The diagram represents a subset of terms that contribute to C, and its asymptotic behavior is determined by the Feynman rules: the subset is O(n sγ ) where s γ = l γ -m i=1 ki-1 2 , with k i the degree of the i-th vertex, and l γ the number of loops in the double-line diagram DL(γ). Furthremore, the correlation function is C = O(n s ), where s = max γ∈Γ(C) s γ .

We now introduce the following result, which generalizes Lemma 2 to the case of networks with bias. This Lemma is used above in the proof of Theorem 1. Lemma 3. Let C(x 1 , . . . , x m ) be a correlation function for a deep linear network with biases. Let c γ be the number of connected components of a graph γ ∈ Γ(C). Then C = O(n s ), where

EQUATION

)

Proof. It is enough to show that each connected component γ in γ is bounded as O(n s γ ) where

s γ ≤ 1 - v γ 2 .

The graph γ is a triangulation of a Riemann surface with f faces e edges, and v γ vertices of degrees k 1 , . . . , k v γ . The Feynman rules give

EQUATION

)

Using the relation e = 

EQUATION

)

The diagram γ is a triangulation of a Riemann surface with at least one boundary, thus χ ≤ 1 and

s γ ≤ 1 - v γ 2 .

C NON-LINEARITIES:

In previous sections we presented the Feynman diagram method for computing the large width asymptotics of correlation functions. In this section we show that the method applies as-is for deep networks with ReLU non-linearities and all-equal inputs, as well as to networks with a single hidden layer, a broader class of non-linearities, and arbitrary inputs. Theorem 2 follows immediately from the results presented in this section.

C.1 DEEP NETWORKS WITH RELU NON-LINEARITIES:

The following result guarantees that the presence of ReLU non-linearities does not change the asymptotic upper bound compared with linear activations, when all inputs are the same. Theorem 7. Let f nl be a network function of the form (2) with ReLU activation. Let f be a network with the same architecture but with linear activation. Let C(x 1 , . . . , x m ; f nl ) be a correlation function of the deep network f nl with width n, and let f be a deep linear network with the same width and depth. Suppose that all inputs are the same,

EQUATION

)

Intuitively, we will rely on the fact that for ReLU networks we can, in some cases, treat the binary neuron activations as being statistically independent of the weights. This result is due to (Hanin & Nica (2018) ). Given this result, we can bound the contribution of the binary activations, and the remaining Gaussian integral is equivalent to that found in a deep linear network. The proof does not work for correlation functions with non-equal inputs, because in that case the independence result of Hanin & Nica (2018) no longer holds.

Proof. We can write the network function as

EQUATION

EQUATION

Here, H is the Heaviside step function acting elementwise on its vector argument. We now introduce the construction from Hanin & Nica (2018) . Let ξ j , η j , j = 1, . . . , d be diagonal matrices of dimension n, whose diagonal elements are ±1-Bernoulli(p) variables with p = 1 2 . We define the new variables Û := ξ 1 U , V := η d V , Ŵ j := ξ j+1 W j η j , j = 1, . . . , d .

Let fnl be a network function with the same architecture as f nl but using the re-defined weights. We define fnl

EQUATION

ρ j := η j ξ j , j = 1, . . . , d .

(37) The following was shown in Hanin & Nica (2018) .

1. f nl and fnl are equal in distribution.

2. { Dj (x), ρ j , j = 1, . . . , d} are independent of {U, V, W 1 , . . . , W d-1 }.

3. { Dj (x), j = 1, . . . , d} are independent of each other for fixed x. The diagonal entries of each diagonal matrix Dj (x) are independent, and take the values {0, 1} with probability 1 2 .

Now, the correlation function can be written as

EQUATION

The second equality follows from the fact that f nl d = fnl . Let us assume for now that the correlation function has no contracted derivatives. The we have

C(x 1 , . . . , x m ; f nl ) = E θ,η,ξ fnl (x 1 ) • • • fnl (x m ) (39) = 1 n dm/2 n i Din α E m l=1 V i l,d d-1 l =1 E m l =1 W l i l ,l +1 ,i l ,l E m l=1 U i l,1 ,α l E i, α . (40)

Here, i = {i 1,1 , . . . , i m,d } and α = {α 1 , . . . , α m }, and

EQUATION

In writing the equation ( 40) we used the facts that D, ρ are independent of the parameters, and that parameters in different layers are independent. We now use Theorem 4 (Isserlis), which says that each of the expectation values over products of V , W j , and U elements is equal to a sum over permutations, where each term is a product over Kronecker delta functions -the covariance matrices of the parameters.

EQUATION

Here, ∆ σ i, α is a product of Kronecker delta functions. The precise form of this object will not be important for our purpose, though it can be easily derived using Theorem 4.

We can now bound the correlation function as follows.

EQUATION

EQUATION

EQUATION

Here, E max = max i, α E i, α , and v ∈ R Din is the all-ones vector. The diagonal elements of ρ l , Dl (x) are identical variables. Therefore, E i, α only takes an O(1) number of values in the large width limit. Note that each of these values are independent of n, and therefore E max = O(1).

We managed to bound the correlation function for f nl in terms of the correlation function for the corresponding linear network f with fixed inputs. The asymptotics of the linear-network correlation function do not depend on the inputs, and therefore this concludes the proof.

C.2 SINGLE HIDDEN LAYER NETWORKS:

For networks with a single hidden layer, defined by

EQUATION

we can extend our asymptotic analysis to smooth non-linearities. We will show in Theorem 8 that for any correlation function C, we have

EQUATION

where f is a deep linear network of equal width and sufficient depth. Therefore, computing the asymptotics using Feynman diagrams for deep linear networks yields a bound on networks with a single hidden layer and smooth non-linearities.

Before delving into the proof of this claim, we consider a few simple examples. Let us begin with the correlation function:

EQUATION

Here we used two facts. First, the weights V i are unaffected by the non-linearity, so we can carry out the V integral. Second, the summand in the last equation, E σ(U T i x)σ(U T i x ) , is independent of both i and n because U i are i.i.d. variables.

Next, consider the following correlation function. For simplicity, here we set the input dimension to be D in = 1, and we set all inputs to 1; this does not change the asymptotics.

E θ dΘ dt = n j,k E θ f nl ∂ 2 f nl ∂U j ∂V k ∂f nl ∂U j ∂f nl ∂V k = 1 n 2 n i1,i2 E θ [V i1 V i2 σ(U i1 )σ (U i2 )σ (U i2 )σ(U i2 )] = 1 n 2 n i1 E θ σ(U i1 ) 2 σ (U i1 ) 2 = O(n -1

) .

(50) In the last line, we again used the fact that the summands are equal and independent of n.

In general, a correlation function C(x 1 , . . . , x m ; f nl ) can be reduced to the form

EQUATION

EQUATION

)

This form is obtained by carrying out the V i integrals, as well as all the sums that can be trivially carried out due to the presence of Kronecker deltas. Here, the α sum represents a sum over all K terms that appear from performing the V i integrals using Isserlis' theorem; σ ( ) denotes the -th derivative of the non-linearity; k α s is the number of σ (•) (U T is x) factors sharing the i s index; σα ∈ S m is a permutation; and F (α) is a function whose exact form is not important for us. We will denote by r max the maximum number of index sums appearing in (51), namely r max := max α r α .

Our approach to establishing the asymptotic scaling will be to first bound the maximum number of index sums appearing in any term in our correlation function, written in the form (51), and then to argue that the summands are bounded by an n-independent constant.

Let us introduce a family of diagrams, Γ (C), which are different in general than the Feynman diagrams. A given diagram g ∈ Γ (C) is constructed as follows. Definition 6. Let C(x 1 , . . . , x m ; f nl ) be a correlation function, where f nl is the output of a network with one hidden layer, defined in (46). The family Γ (C) is the set of all graphs that have the following properties.

• Each derivative tensor in C is mapped to a vertex in the graph.

• Each edge has a type that corresponds to one of the weight vectors U or V .

• Each vertex has exactly one edge of V type.

• If two derivative tensors are contracted in C, the graph must have at least one edge (of any type) connecting the corresponding vertices for each contraction.

The reason for introducing this graphical structure is that it allows us to make two important statements.

If we define cg the number of connected components in a graph g ∈ Γ (C) and cmax := max g∈Γ (C) cg then the following holds.

1. The maximal number of sums appearing in a correlation function of the form (51) is cmax .

2. cmax = c max , where c max is the maximal number of connected components in the family of Feynman diagrams corresponding to a correlation function C(x 1 , . . . , x m ; f d ), where f d is a deep linear network with d hidden layers, and d is large enough that none of the derivative tensors vanish.10 

These two results, combined with a bound on the summands occurring in C(x 1 , . . . , x m ; f nl ) will establish the bound (47). Lemma 4. A correlation function C has r max = cmax .

Proof. Each vertex corresponds to a derivative tensor, T µ1...µ k , which contains a single sum over paired U i and V i indices. If two vertices are connected by one or more edges, there is a Kronecker delta factor that sets the corresponding indices to be equal (due either to a derivative contraction, or to a V covariance), reducing the number of sums by one. The result is that any connected component of a graph g ∈ Γ (C) corresponds to a single sum, and the total numer of index sums corresponding to a particular graph is cg , the number of connected components. As a result, the maximum number of sums corresponding to the collection of graphs Γ (C) is cmax . Proof. Let n e (n o ) be the number of even(odd) clusters in the cluster graph G C of C(x 1 , . . . , x m ; f d ).

The cluster graph, G C is a subgraph of any graph g ∈ Γ (C). We can thus think about the embedding of even and odd clusters into g. In any graph g ∈ Γ (C), an even cluster may belong to its own connected component, while for odd clusters there must be an even number of them in any connected component. This is because an even (odd) cluster contains an even (odd) number of factors of V , which must be paired up in any connected component. We find that

EQUATION

)

The last inequality was used below Lemma 2 in the proof of Theorem 1.

Theorem 8. Let C(x 1 , . . . , x m ; f nl ) be a correlation function for a one hidden layer network. Let c max be the maximal number of connected components over the collection of graphs Γ(C) corresponding to C(x 1 , . . . , x m ; f d ), where f d (x) is a deep linear network map, with with depth d greater than or equal to the maximum number of derivatives appearing in any single derivative tensor in C.

EQUATION

)

Proof. The correlation function in ( 51) can be bound as

EQUATION

)

For fixed inputs, S (α) i1,...,ir α can take at most O(1) different values as a function of its indices, and the values are independent of n. This is because the variables U i are identical. We define s max as the maximum value of |S (α) i1,...,ir α | as a function of α and the indices. Combining this with the above lemmas we can then write.

EQUATION

)

The result of the theorem follows from Lemma 2.

D CORRELATION FUNCTION ASYMPTOTICS:

In this section we prove several general results about correlation function asymptotics in the large width limit. Throughout this section, we assume that Conjecture 1 holds.

D.1 VARIANCE ASYMPTOTICS:

Conjecture 1 can be used to bound the variance of the integrands that appear inside correlation functions. 

Lemma 6. Let C(x 1 , . . . , x m ) = E θ [F θ (x 1 , . . . , x m )]

C = n e + no 2 -m 2 . Then Var θ [F θ (x 1 , . . . , x m )] = O(n 2s C ).

Proof. To bound the variance, it is enough to bound the correlation function As a corrolary, notice that if C = O(n s ) according to Conjecture 1, then typical realizations of the integrand will also be O(n s ). In other words, the asymptotic bound of Conjecture 1 holds for typical initializations, not just in expectation.

EQUATION

Table 3 shows empirical results for the variance of several functions. In all cases we find that the bound of Lemma 6 holds. For µ Var θ [∂ µ f (x 1 )∂ µ f (x 2 )] we prove a tight bound below in Appendix E.1 for deep linear networks.

Function n e , n o 2s C lin.

ReLU tanh:

Var θ [f (x 1 )f (x 2 )] 0,2 0 -0.08 -0.00 -0.02

Var θ [f (x 1 )f (x 2 )f (x 3 )f (x 4 )] 0,4 0 -0.03 0.02 -0.05 µ Var θ [∂ µ f (x 1 )∂ µ f (x 2 )] 1,0 0 -1.01 -1.02 -0.99 µ,ν Var θ [∂ µ f (x 1 )∂ ν f (x 2 )∂ µ,ν f (x 3 )f (x 4 )] 0,2 -2 -2.1 -2.13 -2.14 µ,ν,ρ Var θ [∂ µ f (x 1 )∂ ν f (x 2 )∂ ρ f (x 3 )∂ µ,ν,ρ f (x 4 )] 1,0 -2 -4.02 -4.1 -3.05 µ,ν,ρ,σ Var θ ∂µf (x 1 )∂ν f (x 2 )∂µ,ν f (x 3 )∂ρf (x 4 )∂σ f (x 5 )∂ρ,σ f (x 6 )

0,2 -4 -4.09 -4.14 -4.01

Table 3 : Bounds on variances obtained from Lemma 6, where the predicted exponent is 2s C , compared with empirical results. The predicted exponent is 2s C , and the 3 right-most columns list the empirical exponents. The experimental setup is the same as that of Table 1 .

D.2 GRADIENT FLOW:

The following results are used in the gradient flow calculations of Section 3. We now prove Lemma 1 giving the scaling of time derivatives of correlation functions at initialization. We prove the result for polynomial loss functions. Extension to more general loss functions requires interchanging the large width limit and the Taylor expansion of the loss, which we do not discuss.

proof (Lemma 1). Notice that

EQUATION

)

On the right-hand side we have a linear combination A α A C A of correlation functions C A , where the coefficients depend on the training set labels. Here we used the polynomial loss assumption. For each correlation function C A , its integrand can be obtained from F by finitely many operations of the form (1) multiply the integrand by f (x) for some input x, and (2) act with a pair of contracted derivatives on two of the derivative tensors. In the cluster graph representation, these two operations correspond to (1) adding a vertex, and (2) adding an edge. Therefore, the cluster graph

G C of C = E θ [F ]

is a subgraph of the cluster graph G C A of each one of the correlation functions C A . By Lemma 7 we have that C A = O(n s C ) for all A, and therefore the bound applies to E θ d k F/dt k as well.

D.3 STOCHASTIC GRADIENT DESCENT:

In this section we show that the asymptotics of a correlation function do not change under stochastic gradient descent (SGD) updates. Let C(x 1 , . . . , x m ) = E θ∼P0 [F θ (x 1 , . . . , x m )] be a correlation function, where F θ is the integrand which explicitly depends on the parameters θ. Under a single SGD step, the parameters are updated as θ t+1 = θ t -η∇L t (θ t ), where L is the mini-batch loss at time t. Let P t denote the distribution of parameters at step t, where P 0 is the initial distribution. We define the evolved correlation function at step t to be

EQUATION

We have the following Theorem 9. Let C be a correlation function for a network with linear activations, and assume that Conjecture 1 holds, namely C = O(n s C ) where s C is defined in the Conjecture. If C t is the evolved correlation function after t SGD steps, then

C t = O(n s C ).

Proof. Notice that

EQUATION

The integrand F θ can be written as a product of derivative tensors of the form T µ1...µ k (x; θ), with contracted derivatives. Suppose that the network has d hidden layers. Then under an SGD step, we have

EQUATION

Here D B is the mini-batch, and x a are mini-batch samples. The k sum is truncated because higherorder derivatives of f vanish.

We can now see how taking a gradient descent step affects the cluster graph. After taking a step, each derivative tensor in C t is replaced by a sum (over k, x 1 , . . . , x k ). Each term in the combination of these sums is a correlation function, whose cluster graph is a subgraph of C. Therefore, by Lemma 7,

C t = O(n s C ).

We note that for general activation functions the sum in (62) may be infinite. In this case, to complete the proof we would need to show that the infinite sum obeys the same bound as each individual term in the sum. We leave this for future work.

E APPLICATIONS:

Here we present several applications of our Feynman diagram method for computing large width asymptotics. We assume throughout this section that Conjecture 1 holds.

E.1 NEURAL TANGENT KERNEL:

In this section we prove two results regarding the NTK at large width. We show that the kernel converges in probability, and that during gradient descent it is constant up to O(n -1 ) corrections.

Theorem 10. The Neural Tangent Kernel Θ of a deep linear network converges in probability in the large width limit, and its variance is O(n -1 ).

Conjecture 1 is not sufficient for proving this theorem, as we need to use a more detailed Feynman diagram argument. Therefore, we only prove the theorem for the case of deep linear networks.

Proof. First, we will show that Var[Θ] = O(n -1 ). Given a model function f , the variance is given by Var

EQUATION

EQUATION

EQUATION

Here

f (1) = • • • = f (4)

= f ; we introduced the superscripts so we can easily refer to individual factors. The crux of the argument is that the set of Feynman diagrams representing the expression (63) includes only connected graphs. Assuming this is the case, according to the bound ( 22) these diagrams will all scale as O(n (2-4)/2 ) = O(n -1 ) and so will the variance, completing the proof. To see why only connected graphs contribute, notice that

• A(x, x ) includes all diagrams in which the vertices corresponding to f (1) , f (2) share an edge, and also the vertices f (3) , f (4) share an edge (due to the explicit derivatives);

• B(x, x ) includes all diagrams in which all edges are either between f (1) , f (2) or between f (3) , f (4) .

Therefore, in the full expression (63), the only remaining diagrams (i.e. the diagrams that do not cancel between the two terms) are those that include

• an edge connecting f (1) , f (2) ,

• an edge connecting f (3) , f (4) , and

• an edge connecting one of f (1) , f (2) to one of f (3) , f (4) .

These diagrams are all connected graphs, and this completes the proof that Var[Θ] = O(n -1 ).

Here and above we have established that

Var[Θ] = O(n -1

) and E [Θ] = O(n 0 ). More generally, let us consider a random variable O equal to the product of f and its derivatives, where the derivatives indices are fully summed in pairs. As we will now show, if  and  (2) O converges in probability to this limit. In particular, the NTK is an example of such a random variable, and so this will conclude the proof that the NTK converges in probability.

E [O] = O(n 0 ) and Var[O] = O(n -1 ) then (1) the limit lim n→∞ E [O] exists,

First, consider the mean. There is a finite number of diagrams contributing to E [O], and each has a well-defined n scaling. Therefore, we can write

EQUATION

for some values of k max ≥ 0 (integer) and c 0 , c 1 , . . . ∈ R. We see that the mean has a well-defined large n limit,

EQUATION

Next, let us show that O -E [O] converges in probability using Chebyshev's inequality. Indeed, by the variance assumption there exists c > 0 such that

EQUATION

Combining the facts that O -E [O] → p 0 and E [O] → c 0 , we find that O converges in probability to c 0 .

Next, we show that the large width NTK is constant during training, and compute the asymptotics of the higher-order terms. The following argument is phrased for deep linear networks. More generally, the same argument holds for deep networks with smooth non-linear activations under the additional assumption that the large width limit and Taylor series can be exchanged (note that for such networks, the network function is analytic in the weights).

Theorem 11. Let f (x) be the network output of a deep linear network with MSE loss L. Let Θ t (x, x ) be the Neural Tangent Kernel at SGD step t, for some inputs x, x . Then in the large width limit, the kernel is constant in t in expectation, and

E θ∼P0 [Θ t (x, x ) -Θ 0 (x, x )] = O(n -1 ) , (69) Var θ∼P0 [Θ t (x, x ) -Θ 0 (x, x )] = O(n -2 ) .

(70)

Proof. It is enough to show that E θ∼P0 [Θ 1 (x, x ) -Θ 0 (x, x )] = O(n -1 ) . It then follows from Theorem 9 that E θ∼P0 [Θ t+1 (x, x ) -Θ t (x, x )] = O(n -1

) for all t, and therefore

EQUATION

concluding the proof. We have

EQUATION

EQUATION

The fact that the k, l sums are truncated at d + 1 follows from using linear activations, as higher-order derivatives of the network function vanish in this case. All terms in the sum over k, l include a tensor product of the form µ, µ, ν T µµ1...

µ k (x)T µν1...ν l (x )(• • • ) with either k ≥ 1 or l ≥ 1, where (• • • )

stands for additional derivative tensor factors. Therefore, all terms in the k, l sum are correlation functions that have a cluster of size at least 3, including T (x), T (x ), and at least one other tensor contracted through the µ 1 or ν 1 index. It follows from the Conjecture that each term in the sum is O(n -1 ). Lemma 6 then implies that the variance of these updates is O(n -2 ).

E.2 WIDE NETWORK EVOLUTION IS LINEAR IN THE LEARNING RATE:

In this section we prove that, at large width, the NTK determines the evolution of the network function not just for continuous-time gradient descent but also for discrete-time gradient descent. A similar result holds for stochastic gradient descent, using a stochastic kernel.11 Again we prove the deep linear case explicitly, but the result holds for deep networks with smooth non-linear activations under the additional assumption that the large width limit and Taylor series can be exchanged.

Theorem 12. Let f (x) be the network output of a deep linear network, and let f t (x) be the evolved function after t gradient descent steps, defined by θ t+1 = θ t -η∇L(θ t ). In the large width limit, each gradient descent step update of f t is linear in the learning rate η. Furthermore,

E θ [f t+1 (x) -f t (x)] = -η (x ,y )∈Dtr E θ Θ 0 (x, x ) ∂ t (x , y ) ∂f + O(n -1 ) . (73)

Here, Θ 0 is the NTK at initialization and t is the single-sample loss at time t.

Proof. For a deep linear network, under a single gradient descent step we have

EQUATION

)

First, consider the sum over k in (73). Each term in the sum is a correlation function for which the cluster graph contains a connected subgraph of size at least 3, and is therefore O(n -1 ) by Lemma 7 and Theorem 1. In the remaining O(η) term, by the same argument as Theorem 11 we can replace

Θ t = Θ 0 + O(n -1

) in the correlation function. The result is equation ( 73).

As mentioned above, we note that the proof goes through when using stochastic gradient descent updates, with the difference that in (73) we should sum over mini-batch samples instead of over the entire training set.

E.3 SPECTRAL PROPERTIES OF THE NTK AND THE HESSIAN:

With an eye towards understanding the structure of the loss landscape at large width and as another example use case of our approach, we investigate the relation between the spectra of the Hessian, and the NTK.

Among other observations, we present an argument that for a network with d hidden layers and MSE loss, the top D in ordered eigenvalues of the Hessian, λ (H) i

, are related to those of the kernel λ

EQUATION

)

The remaining Hessian eigenvalues vanish at large width as O(n -1/2 ) for one hidden layer networks and O(n -1 ) for deeper networks. 12 This is experimentally corroborated in Figure 9 .

The Hessian of a general loss takes the form,

EQUATION

)

The moments of H are a useful way to understand its spectrum. As in other examples that we have seen, traces of powers of B lead to connected diagrams with higher and higher numbers of vertices, and so scale as increasingly negative powers of n. More explicitly, traces involving powers of B contain multiple contracted derivative tensors. As a result, taking the average of these traces over the 12 Empirically we actually find the even stronger bound E θ λ

(H) i -λ (Θ) i = O(n -1

) for the top Din eigenvalue differences and O(n -1/2 ) for the remaining eigenvalues in the one hidden layer case. We can gain insight into this improved scaling through the perspective of degenerate eigenvalue perturbation theory (Trefethen & Bau (1997) ), but this is outside the scope of the current presentation. 

EQUATION

What's more, we can relate moments of A to those of the kernel, as both are built out of two logit derivatives. Explicitly,

EQUATION

and so the moments of the Hessian are also related to those of the kernel.14 

EQUATION

Here, we have defined,

EQUATION

For the case of MSE loss, we can go even further. In that case, M (x a , x b ) = δ ab and B decays to zero during training. We thus have,

EQUATION

These results indicate that the only difference between the spectra of the Hessian and the NTK come from B and that B must have eigenvalues which scale as 1/ dim(B). As dim(B) = O(n) for one hidden layer networks and O(n 2 ) for deep networks we are left with the relation (76) between the eigenvalues of Θ and H. As the network trains, the difference between these eigenvalues gets even smaller. These results are confirmed experimentally in Figure 9 and Figure 10 .

Multi-class. This story can be generalized to multi-class neural network maps. In this case the logit is a map f A : R Din → R Nclass , and the NTK has two class indices.

EQUATION

Expression (80) relating the moments of the Hessian to the NTK still holds in this context provided we take Tr ((M Θ) m ) → Tr (M AB Θ AB ) m , with the more general matrix, Eigenvalues

EQUATION

2 -4 2 -3 2 -2

Top eigenvalues At large width, the NTK approaches the identity matrix in class space, Θ AB → N -1 class δ AB Tr (Θ) (Jacot et al. (2018) ). This implies that the NTK spectrum consists of N class repeated copies. This has consequences for the spectrum of the Hessian at large width. For instance, for the case of MSE error it also implies N class repeated copies of the Hessian spectrum (See Figure 9b ). It is intriguing to think this could serve as a path towards understanding the emergence of the N class (or N class -1) large eigenvalues and corresponding subspace observed in the Hessian spectrum (Sagun et al. (2016; 2017) ; Gur-Ari et al. (2018) ; Ghorbani et al. (2019) ; Papyan (2019)).

E.4 HIGHER-ORDER NETWORK EVOLUTION:

In this section we will explain how to compute higher-order corrections to training dynamics. In principle, this prescription allows one to compute model dynamics as an expansion in 1/n to arbitrary order. We apply this explicitly to compute the O(1/n) training dynamics. In Figure 11 , we present experimental confirmation of our predictions for the evolution of the NTK. In the main text, we presented these results for the special case of gradient flow with MSE loss. So far we have mostly been using diagrammatic techniques to understanding the leading order scaling of correlation functions. It is interesting to try and understand finite width networks by asking how training dynamics are modified beyond the leading order asymptotics. There are three clear sources of corrections to the leading behavior.

1. The initial kernel, Θ 0 , receives finite width corrections.

2. The linear (in learning rate) update for the network function, equation ( 73), is modified away from infinite width.

3. The kernel is not constant at finite width.

The first source of corrections is automatically taken into account in typical empirical settings, as the finite-width Θ 0 can be computed explicitly. The other sources of corrections are non-trivial, and will be the focus of this section. We begin by explaining how to take into account the non-constancy of the kernel order by order in 1/n, while maintaining the continuous time approximation. Next, we back away from the continuous time limit and write down the full discrete evolution. We introduce a method to compute arbitrary order corrections, and explicitly work out the corrections at order 1/n for MSE loss.

E.4.1 CONTINUOUS TIME.:

Continuous time evolution of the model function in neural networks is governed by the differential equation

EQUATION

As we have discussed at length, in the large width limit, this equation simplifies and Θ(t) is asymptotically constant (Jacot et al. (2018) ). Our goal is to move beyond this leading behavior at large width and solve (85) order-by-order in a 1/n expansion. In Section 3 we described how to compute the O(n -1 ) corrections for the case of MSE loss. Here we explain how to handle corrections more generally as well as giving a more detailed discussion of the MSE case.

The network map f and kernel, Θ are members of the following family of operators.

EQUATION

with O 1 := f . Here, as above, O 2 = Θ is the kernel. With a general loss function, these operators satisfy.

EQUATION

Equations ( 86) and (87). Give an infinite tower of first order ODEs, the solution of which gives the time evolution of the network map and the kernel.

EQUATION

. . . Solving this infinite tower is not feasible. If we wish to work to a given order in an expansion in 1/n, however, there is a dramatic simplification, which makes a solution possible. Firstly, we can truncate these equations. To see this note that the operators O s contain s contracted derivative tensors. As a result, by Lemma 7 and Conjecture 1, correlation functions involving the operators, O s satisfy

EQUATION

where F ( x; t) is arbitrary additional contribution to the integrand.

Thus, if we wish to work to solve for the time evolution up to corrections which scale as O(n -r ) we can truncate the tower at s = 2r and set O 2r (t) to be equal to its initial value. Note that the leading order solution, (85), is the result of this procedure with r = 1 and the results presented in the main text for the O(n -1 ) evolution correspond to r = 2.

The truncation provides a dramatic simplification, however it is not immediately clear how to solve even the truncated differential equations in (87). We now describe how to organize the perturbative expansion of the operators O s (t) (including f and Θ) in such a way that the differential equations become tractable.

The central idea is to write each operator, O s (t), as an expansion.

EQUATION

where each order O (r) s order by order. For example, in Section 3, we used the leading order solution for f and Θ to solve for Θ (1) .

EQUATION

In principle this procedure can be extended to arbitrary order in 1/n. Before going onto explain the finite step corrections to this procedure, we reproduce the results of section 3 in more detail. Here, we have introduced the eigenvectors, êi , which are vectors over the training dataset and eigenvalues λ i of the initial kernel, Θ 0 . The vector x = (x 1 , x 2 ) as Θ depends on two inputs. O 3 ( x; 0) is a vector over the dataset while O 4 ( x; 0) is a square matrix and ∆f 0 := f 0y is a vector over the training data.

Finally, we can plug in the expression (100) into (93) and give the sub-leading behavior for f (t) df (x; t) dt = - Note, one consequence of ( 101) is that the O(1/n) corrections to Θ go to a constant at late times.

EQUATION

)

And the asymptotic evolution of ∆f (t) is again constant kernel evolution, but with a corrected kernel.

f (t) → y + e -t(Θ0+Θ (1) ∞ ) (f 0y) .

(105)

It is worth noting that these predictions are quite detailed, giving the full time dependence, rather than just the overall scaling with width. In deriving these results and in the discrete time analysis below, we implicitly make an assumption that we can control the discrete time corrections or in the continuous time case that f (t), Θ(t) and all the O s (t) are differentiable. This relies on being able to differentiate through the activation function. For non-smooth activations such as ReLU, this assumption is suspect, and indeed the evolution is more subtle in the ReLU case and a full analysis is left to future work. For this reason we present experimental results for networks with smooth activation functions.

E.4.2 DISCRETE TIME.

With the 1/n corrections to continuous time evolution under our belts, it is possible to also keep track of corrections coming from the discrete update step. The essential point, that it is possible to solve iteratively, order by order in 1/n is not changed. To see this we can look at how the update equations are modified. Beginning with the network output update equation, we have,

EQUATION

In Section E.2 we argued that the O(η 2 ) terms vanish at large width. In more detail, each factor of the gradient, g µ , contains a contracted derivative tensor. Thus, the cluster graph for correlation functions involving the η 2 term will always contain three contracted vertices and thus the correlation functions will scale as O(n -1 ). The η 3 term will give four vertices, and thus also scale as O(1/n), the η 4 contribution will be O(1/n 2 ) and so on.

Just as in the continuous time case this equation is still solvable order by order in 1/n. What we mean by solvable is the following. All terms on the RHS appearing with f t are already higher order as a result of their derivative structure. This means, just as in the continuous case, we can use the lower order solution of f t in the gradient terms on the RHS to solve for f t on the LHS.

Similarly, the update equation for Θ receives discrete time modifications.

Θ t+1 (x 1 , x 2 ) -Θ t (x 1 , x 2 ) = -η a O 3;t (x 1 , x 2 , x) ∂ t (x, y) ∂f + η 2 2 µν ∂ 2 Θ t (x 1 , x 2 ) ∂θ µ ∂θ ν g µ t g ν t + • • • .

(107)

Here we have adopted a notation where, x = (x 1 , . . . , x s ). Just as for f , increased powers in η are increasingly suppressed in n and we can iteratively solve this equation using solutions at leading orders in n to solve for sub-leading behavior.

More generally the tower of differential equations defined recursively in (87) becomes As in the continuous time case, at a given order in n, we can truncate this tower and use lower order solutions to solve for the time evolution up to the desired order. To ground this discussion in a concrete use case, we again walk through the 1/n corrections for MSE loss, now in discrete time.

MSE discrete time example.:

For discrete time evolution, the leading order behavior of the model map is (112)

EQUATION

To order 1/n we can drop the order η 2 and higher terms, as they contain expressions with 5 or more contracted derivative tensors. Thus at this order, we are left with the discrete analogue of (96). 

Explicitly, this takes the form

EQUATION

Here we have adopted notation similar to above, where O 4;0 ( x) is a vector over the training dataset.

So far, this procedure has been a discrete analogue of what we have done in the continuous time case, however as we move on to compute Θ and f we will have to keep track of the novel corrections, which vanish in continuous time. Explicitly, at order 1/n the discrete update for Θ t is given by,

Θ:

(1) t+1 ( x) -Θ 

This can be summed to give Θ

(1) t . To move up to the neural network map itself there is one additional complication. The O(η 2 /n) term, (120)

Here (x a , y a ) are elements of the training set D tr and summed over. This equation can be solved to give f

(1)

t .

f

(1) t

= (1 -ηΘ 0 ) t t t =1

(1 -ηΘ 0 ) -(t +1) -ηΘ NLO, t (1 -ηΘ 0 ) t ∆f 0 + Disc t .

Where Disc t contains the discrete derivative updates at O(1/n). (122)

These expressions may look fairly intimidating. The key point is that all terms in the summand in (116) and thus (121) are known functions of time and initial data, just as in the continuous time setting.

:

For deep linear networks, the proof of Theorem 1 relies on mapping diagrams to triangulations of 2D surfaces with at least 1 boundary. There are cases, such as the one highlighted in Table1, in which all diagrams one can write down have at least 2 boundaries. In such cases, the full diagrammatic calculation gives a tight bound, while Conjecture 1 provides a valid upper bound that is not tight.

Here we are using condensed notation: Θ (0) and f (0) are values at initialization, and f , f (0) , y are treated as vectors in training set space. The kernel Θ (0) is a square matrix in the same space.

Here we are extending the use of the term cluster to refer to derivative tensors in the integrand itself.

For networks with a single hidden layer, the number of loops in a diagram is equal to the number of connected components. This is not true for a general deep linear network.

We are not aware of a physics application in which such derivatives are included in a Feynman diagram description of correlation functions. Therefore, to our knowledge our treatment of these derivatives is novel both in machine learning and in physics.

One can consider the cyclic model f(x) = n -(d+1)/2 Tr Ṽ W d-1 • • • W 1 Ũ x with 1D input x, in whichall the weight tensors Ũ , Ṽ , W 1 , . . . , W d-1 are n × n matrices. The Feynman diagram construction for this model is similar to the deep linear case, except that all edges in the Feynman diagram map to double edges in the double-line diagrams. Such diagrams can be interpreted as triangulations of surfaces with no boundaries. Unlike the deep network diagrams, they have no 'missing' loops.

Note that any derivative tensor of f that has rank greater than d vanishes.

The perspective presented here helps understand the results of(Lee et al. (2019)) where it was observed empirically that linearized evolution is a good description of wide networks even for relatively large learning rates.

For the first moment in linear or ReLU networks, Tr B = 0.

Here we have argued for relations relating the mean of moments. It is not too difficult to see that these relations will also hold for typical realizations. This follows from Lemma 6.

Note that in Section 3 we used the alternate notation Θ1 for Θ(1) .

Figures and Tables:

Related work. For additional works on wide networks, including relating them to Gaussian processes, see de G. Matthews et al. (2018); Novak et al. (2019); Garriga-Alonso et al. (2019); Yang (2019); Pennington & Bahri (2017); Pennington & Worah (2018); Schoenholz et al. (2017); Xiao et al. (2018); Chen et al. (2018); Daniely et al. (2016); Lee et al. (2019); Du et al. (2018b;a); Allen-Zhu

Denote by n e (n o ) the number of even (odd) clusters in this correlation function, which has m = m + 2 derivative tensors. One can check that either (n e , n o ) = (n e + 1, n o ) or (n e , n o ) = (n e -1, n o + 2), depending on whether the ∂ µ F derivative is acting on an odd or even cluster in F . 6 Therefore, n e + n o 2 -m 2 ≤ s C .

Figure 1: Empirical verification of predicted asymptotics. (a) The mean deviation of the NTK from its initial value for a variety of widths and activation functions. The fit (dashed) matches well with the predicted O(n -1 ) asymptotics. (b-c) Comparison between the empirical evolution (solid) and the O(n -1 ) predicted evolution (dashed) for the network function and the kernel. All experiments were performed on two-class MNIST, computing a single randomly-chosen component of Θ or f . See Appendix A for additional experimental details.

Figure 2: Feynman diagram examples. (a) The Feynman diagram of E θ [f (x)f (x )]. (b)-(c) The Feynman diagrams of E θ [f (x 1 )f (x 2 )f (x 3 )f (x 4 )]; additional, equivalent diagrams are not shown.

Figure 3: Feynman diagrams for E θ [Θ(x, x )] with one hidden layer. The dashed vertical line represents vertices forced by contracted derivatives.

were performed on two-class MNIST, computing a single randomlychosen component of the kernel Θ. Sub-figure (a) uses networks trained for 1024 steps with learning rate 1.0 and 1000 samples per class, averaged over 100 initializations. Each curve in figure (b) represents a single instance of the network map evaluated on a random image over the corse of training. The models were trained with 10 samples per class and learning rate 0.1. The input to the network is normalized by the square root of the input dimension as in(Jacot et al. (2018))

Figure 4: Feynman diagrams for E θ [f (x)f (x )] with 2 hidden layers. Notice that the U, V edges in the Feynman diagram map to single edges in the double-line diagram, while the W edge maps to a double edge.

Figure 4 shows the Feynman diagram and corresponding double-line diagram for E θ [f (x)f (x )] with 2 hidden layers. We can interpret this Feynman diagram as a triangulation of the disc: a 2-dimensional

Figure 6: Double-line graphs describing the expectation value of the NTK, E θ [Θ], for a deep linear network with two hidden layers. Crossed edges mark the edges that are forced by contracted derivatives. The derivative can act on either U , V , or W 1 , and therefore there are three diagrams. Each diagram has 2 vertices and 2 faces, and therefore the correlation function is O(n 0 ) according to the Feynman rules.

Theorem 5. Let C(x 1 , . . . , x m ) be a correlation function of a deep linear network with d hidden layers, and let γ ∈ Γ(C) be a Feynman diagram. The diagram represents a subset of terms that contribute to C, and its asymptotic behavior is determined by the Feynman rules: the subset is O(n sγ ) where s γ = l γ -dm 2 , and l γ is the number of loops in the double-line diagram DL(γ). Furthremore, the correlation function is C = O(n s ), where s = max γ∈Γ(C) s γ .

Figure 7: New Feynman diagrams for E θ [f (x 1 )f (x 2 )f (x 3 )f (x 4 )] with 2 hidden layers and bias.

Figure 7 is an example of the new single line and double line diagrams that appear in the correlation functionE θ [f (x 1 )f (x 2 )f (x 3 )f (x 4 )]for a two-hidden-layer network with biases.

definition of the Euler character, χ = v γe + f we have

The maximal number of connected components, cmax , over the collection of graphs Γ (C) corresponding to C(x 1 , . . . , x m ; f nl ) is bounded by the maximal number of components, c max , over the collection Γ(C) corresponding to C(x 1 , . . . , x m ; f d ) of Feynman diagrams, where f d is a deep linear network of sufficient depth d.

be a correlation function, and let G C be the cluster graph of C with n e even components and n o odd components. Assume Conjecture 1 holds such that C = O(n s C ), where s

57) because Var θ [F θ (x 1 , . . . , x m )] ≤ C(x 1 , . . . , x m , x 1 , . . . , x m ). The correlation function C has 2n e even clusters and 2n o odd clusters. It follows from Conjecture 1 that C = O(n 2s C ).

Lemma 7. Let G be a graph with m vertices, n e even components, and n o odd components. Let G be a subgraph of G with m vertices, n e even components, and n o odd components. Then s(n e , n o , m) ≥ s(n e , n o , m ) where s(a, b, c) := a + b-c 2 . Proof. It is enough to show that s(n e , n o , m) does not increase if we (1) add a vertex to G, or (2) add an edge to G, because G can be obtained from G by performing such operations finitely many times. Adding a vertex to G changes n e → n e , n o → n o + 1, and m → m + 1, leaving s(n e , n o , m) unchanged. Next, if we add an edge to G then m does not change, and there are 4 possibilities for how n e and n o change. 1. The edge connects two even components. Then n e → n e -1, n o → n o , and s(n e , n o , m) decreases by 1. 2. The edge connects two odd components. Then n e → n e + 1, n o → n o -2, and s(n e , n o , m) does not change. 3. The edge connects an even component and an odd component. Then n e → n e -1, n o → n o , and s(n e , n o , m) decreases by 1. 4. The edge connects two vertices that belong to the same component. In this case n e , n o , and s(n e , n o , m) do not change.

Figure 8: A prototypical diagram corresponding to E θ Tr B 2

Figure 9: (a) The mean difference between the top eigenvalues of the NTK and the Hessian at initializations for a two hidden layer ReLU network of varying width match well with the predicted O(n -1 ) behavior. The mean is taken over 100 instances per width on two-class MNIST with 10 images per class. (b) The top 10 eigenvalues of both the Hessian and the NTK for a two-layer ReLU network of width 2048 trained on a three-class subset of MNIST. The top eigenvalues match well and show aspects of the repeated N class structure predicted at large width. The eigenvalues of H are computed using the lanczos algorithm.

Figure 10: Relation between the spectrum of the Hessian and NTK. (a) The first 4 moments of H and Θ for two hidden layer ReLU networks. At initialization the first third and fourth moments are numerically similar, while the relative difference of the second moment is O(1). After training all moments are numerically close. (b) The root mean square eigenvalue of the B matrix shows good agreement with the predicted 1/ √ n fall off for single hidden layer ReLU networks. Experiments were performed on ReLU networks trained on two-class MNIST with 100 images per class. Moments of the NTK were computed exactly, while moments of the Hessian and B were computed by randomly sampling 1000 vectors and using Hutchinson's trick.

n -r ) and, by Lemma 6, that typical realizations of the operators scale as O(n -r ).15   Once the operators O s are organized in this way, solving the differential equations (85) is tractable. As the differential equations describing the evolution of O (r) s for r > 0 only depend on the time dependent solutions of O (r-1) s , we can iteratively solve for the O (r)

Figure 11: Single instance evolution for a randomly selected element of the NTK and network map. Theoretical predictions at O(n -1 ) are represented by dashed lines, while the true network evolution is shown in solid. (a-b) Wide networks show excellent agreement with predictions. (c-d) Predictions match reasonably well even down to modest widths. All plots correspond to single training runs on 2-class MNIST with 10 samples per class. The learning rate is 0.1 for all runs. All runs reach training accuracy 1.

x ; 0) + Θ (1) (x, x ; t) (f (x ; t)y) + O(n -2 ) (102) f (t) = y + e -Θ0t 1 -t 0 dt e t Θ0 Θ (1) (t )e -t Θ0 (f 0y) . (103)This completes the full O(n -1 ) time dependance.This prediction for the O(1/n) time dependence is confirmed experimentally in Figure11.

s;t+1 ( x) -O s;t ( x) = -(x ,y )∈Dtr O s+1;t ( x, x ) ∂ (x , y ; ,µ2,...,µ k ∂ k O s;t+1 ( x) ∂θ µ1 ∂θ µ2 . . . ∂θ µ k g µ1 g µ2 • • • g µ k , s ≥ 1 . (108)

t+1 ( x) -O 3;t ( x) = -(x ,y )∈Dtr O 4;t ( x, x ) (f t (x )y ) . (113)which we can sum to get the discrete version of (98),O (1) 3;t ( x) = O 3;0 ( x)η t t =1 (x ,y )∈DtrO 4;0 ( x, x ) (f 0;t (x )y ) .

∂θ ν , in (106) has non trivial time dependence. We can deal with this just as we have been doing, by taking an extra time derivative and integrating. Defining,O 1,1;t (x 1 , x 2 , x 2 ) := µν ∂ 2 f t (x 1 ) ∂θ µ ∂θ ν ∂f t (x 2 ) t ( x) = O 1,1;0 ( x)µ ∂O 1,1;0 ( x) ∂θ µ ∂f 0 ∂θ µ • Θ -1 0 (1 -ηΘ 0 )(1 -(1 -ηΘ 0 ) t )∆f 0 a . (119)This can in turn be plugged into the update equation for f t .f t+1 (x)f t (x) = -η a Θ 0 (x, x a ) + Θ (1) t (x, x a ) (f t (x a )y a ) 0 (x) ∂θ µ ∂θ ν ∂θ ρ ∂f 0 (x a ) ∂θ µ ∂f 0 (x b ) ∂θ ν ∂f 0 (x c ) ∂θ ρ (f t (x a )y a ) (f t (x b )y b ) (f t (x c )y c ) .

t (x a , x b ) (f t (x a )y a ) (f t (x b )y b ) 0 ∂θ µ ∂θ ν ∂θ ρ ∂f 0 (x a ) ∂θ µ ∂f 0 (x b ) ∂θ ν ∂f 0 (x c ) ∂θ ρ (f t (x a )y a ) (f t (x b )y b ) (f t (x c )y c ) .

)

Table content: <table><tr><td/><td/><td colspan="3">Top eigenvalue differences</td><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/></tr><tr><td>2 5</td><td/><td/><td/><td colspan="3">top e.v. diff. fit (slope = -0.98)</td><td/><td/><td/><td/><td/><td/><td/><td/><td/></tr><tr><td>2 7</td><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/></tr><tr><td>2 9</td><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/></tr><tr><td>2 5</td><td>2 6</td><td>2 7</td><td>2 8 n</td><td>2 9</td><td>2 10</td><td>2 11</td><td>10th 9th</td><td>8th</td><td>7th</td><td>6th</td><td>5th</td><td>4th</td><td>3rd</td><td>2nd</td><td>1st</td></tr><tr><td colspan="7">(a) Hessian and NTK e.v. difference.</td><td/><td/><td/><td/><td/><td/><td/><td/><td/></tr></table>

At next order we can proceed by solving the truncated set of equations in (108). The first equation, for O 4;t still gives a constant solution, Here, both the discrete and continuous time terms that would appear on the RHS are O(1/n 2 ).Next, we must solve for O 3;t . The discrete time update is, O 3;t+1 ( x) -O 3;t ( x) = -η

Table content: <table><tr><td/><td/><td>O 4;t = O 4;0 . (1)</td><td/><td>(111)</td></tr><tr><td colspan="2">(x ,y )∈Dtr</td><td>O 4;t ( x, x ) (f t (x ) -y ) +</td><td>η 2 2 µν</td><td>∂ 2 O 3;t ( x) ∂θ µ ∂θ ν g µ t g ν t + • • • .</td></tr><tr><td>Θ</td><td colspan="2">(0) t = Θ 0 .</td><td/><td>(110)</td></tr></table>

Peer review (iclr2020):
This is a positive review. Feel free to skip to the feedback. SUMMARY OF PAPER This paper explains how to use cluster graphs to easily compute the asymptotic behaviour of any given correlation function (Definition 1) for deep linear networks. By *asymptotic behaviour* I mean that it upper-bounds correlation functions by c·n^s, where s is a nonnegative integer given by the particular cluster graph, and c a *constant* (which I think depends on the particular input, x, to correlation function, among other things). The authors then conjecture (Conjecture 1) that these bounds transfer to deep nonlinear networks, and that they are tight: - Appendix C proves that these upper bounds also hold for deep ReLU networks, and 1-hidden-layer networks with smooth nonlinearity. This is also mentioned in page 3. - Section 2.3 empirically shows that these bounds are pretty tight. (in terms of the exponent, none of the theory here gives a value for the constant c) The tool provided above is the main result. The authors then use it to provide some results about wide networks - They give a different proof that for large width, the Neural Tangent Kernel stays constant during training. This is because its derivative wrt. time as a function of width n is O(n^-1), and thus 0 for n->infinity. - Using the ease of calculation from the tool, they approximate the change in the NTK over training time for any network. They do this using its value at initialization + a term that depends on n^-1. These results are numerically verified in Figure 1. - They present numerical evidence for the accuracy of this approximation to the change in NTK over time. The authors spend the last 2 pages explaining how cluster graphs derive from Feynman diagrams (FDs), and why these help compute asymptotics. WHY I AM ACCEPTING THIS PAPER The paper adapts FDs and cluster graphs, which is a potentially very useful tool for other wide-network researchers, and could accelerate research in this whole sub-field. It also shows their power by providing a surprisingly large amount of novel theoretical results. FEEDBACK At a very high level, there is only one thing that I think isn*t made quite clear by the presentation, and it should be. If I understand correctly, Feynman diagrams (or cluster graphs) are only used here to calculate correlation functions for deep *linear* networks. Then, other results establish that the width-dependent asymptotic behaviour for linear networks holds as-is for nonlinear networks, and these results with FDs constitute Conjecture 1. There are proofs for ReLU and 1-layer smooth networks, mentioned in pg. 3; and the experiments in the paper support it for common nonlinear deep networks as well. I think that asymptotics for linear networks transfer to nonlinear ones is an interesting result, which doesn*t depend on FDs. What follows are details. It is unclear to me whether cluster graphs are as *powerful* as FDs, i.e. whether the bound at the end of the Proof in page 8 is always saturated. Are there some cases in which you need to use FDs to get a tighter upper bound? In Table 1 you should say that the values under *lin. ReLU tanh* are the fitted exponent s_C. This is not explained. Perhaps you can mark the only 2 cases (in the 5th row) where the bound is not tight. It would be nice to know how much error remains between the fitted c·n^(s_C) and the empirical values. Please explain x1 <-> x2 in eq. 8 In figure 1b, consider adding the finite-width limit prediction for the training dynamics. You have already done so for the prediction of the NTK during training in figure 5c, you could indicate it in the same way in 5b. Typos: Figure 2 caption: feynman -> Feynman pg 8. anlytic evicence -> analytic