A generalized neural tangent kernel for surrogate gradient learning

Authors: luke  eilers
raoul-martin  memmesheimer
sven  goedeke
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-15 20:22:02

Abstract: State-of-the-art neural network training methods depend on the gradient of the network function. Therefore, they cannot be applied to networks whose activation functions do not have useful derivatives, such as binary and discrete-time spiking neural networks. To overcome this problem, the activation function's derivative is commonly substituted with a surrogate derivative, giving rise to surrogate gradient learning (SGL). This method works well in practice but lacks theoretical foundation.

The neural tangent kernel (NTK) has proven successful in the analysis of gradient descent. Here, we provide a generalization of the NTK, which we call the surrogate gradient NTK, that enables the analysis of SGL. First, we study a naive extension of the NTK to activation functions with jumps, demonstrating that gradient descent for such activation functions is also ill-posed in the infinite-width limit. To address this problem, we generalize the NTK to gradient descent with surrogate derivatives, i.e., SGL. We carefully define this generalization and expand the existing key theorems on the NTK with mathematical rigor. Further, we illustrate our findings with numerical experiments. Finally, we numerically compare SGL in networks with sign activation function and finite width to kernel regression with the surrogate gradient NTK; the results confirm that the surrogate gradient NTK provides a good characterization of SGL.

Main paper: 
 introduction: 
 Artificial neural networks (ANNs) originate from the biologically inspired perceptron [Rosenblatt, 1958] . While the perceptron has a binary output that is faithful to the all-or-none behavior of spiking neurons in the nervous system, most activation functions considered nowadays for ANNs are smooth (like the logistic function) or at least semi-differentiable (like the ReLU function). Differentiable network functions enable the learning of network weights with methods that leverage the gradient of the network function with respect to the network weights like backpropagation [Rumelhart et al., 1986] . These methods are very successful [LeCun et al., 2015] , but cannot be used without welldefined gradients. This is a problem when considering more biologically plausible ANNs, which are typically used in computational neuroscience to understand how the networks of spiking neurons in our nervous system work. These include spiking neural networks (SNNs). Motivated by the energy efficiency of our brain, SNNs and similar networks, such as binary neural networks (BNNs), are considered in the context of neuromorphic computing [Merolla et al., 2014] . Both discrete-time SNNs and BNNs have in common that their activation functions do not have useful derivatives, which renders standard gradient-descent training impossible [Neftci et al., 2019 , Taherkhani et al., 2020 , Tavanaei et al., 2019 , Roy et al., 2019 , Pfeiffer and Pfeil, 2018] . For the scope of this paper, these activation functions can be thought of as step-like functions like the sign function, in which cases the derivative vanishes almost everywhere and is thus no longer informative about the shape of the activation function.Surrogate gradient learning resolves this issue by providing the missing information about the activation function in the form of a surrogate derivative [Hinton, 2012 , Bengio et al., 2013 , Zenke and Ganguli, 2018] . As a result, the gradient-based methods for classical ANNs can be leveraged with great success [Zenke and Vogels, 2021] . However, a theoretical basis underpinning the intuition is missing and it is often unclear which surrogate derivative should be chosen for a particular network. For a review focusing on surrogate gradient learning methods, which we are most interested in, see Neftci et al. [2019] . For a comprehensive review of other learning methods for SNNs, we refer to Taherkhani et al. [2020] , Tavanaei et al. [2019] , and Eshraghian et al. [2021] .The neural tangent kernel (NTK) introduced by Jacot et al. [2018] allows to formulate gradient descent as a kernel method. Just as ANNs with randomly initialized weights converge under certain conditions to Gaussian processes (GPs) in the infinite-width limit [Matthews et al., 2018 , Lee et al., 2018] , the NTK converges at initialization and during training to a deterministic kernel in the same limit. Moreover, the NTK then describes how the network function changes under gradient descent in the infinite-width regime. This led to both a better theoretical understanding of gradient descent and practical applications to neural network training; see Section 1.2 for more details.We study the NTK for networks with sign function as activation function. As the NTK theory is not directly applicable due to an ill-defined derivative, we consider the NTK for a sequence of activation functions that approximates the sign function and then derive a principled way of generalizing the NTK to gain more theoretical insight into surrogate gradient learning. Our contributions are as follows:• We provide a clear definition of the infinite-width limit in Section 2.2, capturing the different notions used in the literature due to the different choices of rates at which the layer widths increase. This definition is used consistently in all mathematical statements and the respective proofs. • In Section 2.3, we demonstrate that the direct extension of the NTK for the sign activation function using an approximating sequence is not well-defined due to the divergence of the kernel's diagonal. This illustrates, from the NTK perspective, that gradient descent is ill-defined for activation functions with jumps and how this problem will be mitigated by surrogate gradient learning. Moreover, we connect this divergence to results by Radhakrishnan et al. [2023] in Theorem 2.3 and show that the direct extension of the NTK can be seen as a well-defined kernel for classification. • We define a generalized version of the NTK in Section 2.4 using so-called quasi-Jacobian matrices and prove the convergence to a deterministic, in general asymmetric, kernel in the infinite width limit at initialization in Theorem 2.4. Using the generalized NTK, we formulate surrogate gradient learning in terms of the generalized NTK for networks with differentiable activation functions. This novel NTK is named the SG-NTK and we prove its convergence to a deterministic kernel during training in Theorem 2.5. • For both the diverging direct extension of the NTK and the SG-NTK with sign activation function and arbitrary surrogate derivative, we derive exact analytical expressions in sections D.1, D.2 and E.2. In particular, we identify the terms that emerge from SGL and that prevent divergence. • In Section 3, we illustrate our findings, in particular Theorem 2.4 and Theorem 2.5, using simulations. Numerical experiments show that the distribution of networks trained with SGL shows agreement with the distribution given by the SG-NTK.Mathematically precise versions of all statements can be found in the appendix, which is selfcontained to ensure rigor and a consistent notation throughout all theorems and proofs.The neural tangent kernel. The convergence of randomly initialized ANNs in the infinite-width limit to Gaussian processes under appropriate scaling of the weights has first been described by Neal [1996] for a single hidden layer and has been extended to multiple hidden layers and other network architectures like CNNs [Matthews et al., 2018 , Garriga-Alonso et al., 2018 , Yang, 2019b , Lee et al., 2018] . The NTK was introduced by Jacot et al. [2018] with first results on its convergence at network initialization and during training. Theoretical results on the implication of the convergence of NTKs on the behaviour of trained wide ANNs were given by Lee et al. [2019] , Arora et al. [2019] , Allen-Zhu et al. [2019] . In Section C.2, we review central theorems on the NTK to enable a clear comparison to our theoretical results.The NTK has been generalized to all kinds of ANN standard architectures Yang [2020] such as CNNs Arora et al. [2019] and attention layers Hron et al. [2020] . In particular, it has been generalized to RNNs Alemohammad et al. [2020] , which are particularly interesting in light of SNNs due to their shared temporal dimension. Bordelon and Pehlevan [2022] derive a generalization of the NTK called effective NTK for different learning rules using an approach similar to ours. Note that all of these extensions require well-defined gradients.The ability of overparameterized neural networks to converge during training and to generalize can be explained using the NTK [Allen-Zhu et al., 2019 , Bordelon et al., 2020 , Bietti and Mairal, 2019 , Du et al., 2019] . The NTK has also been used in more applied areas such as neural architecture search [Chen et al., 2021] and dataset distillation [Nguyen et al., 2021] .Surrogate gradient learning. In the context of computational neuroscience, the idea of replacing the derivative of the output of a spiking neuron with a surrogate derivative was introduced by Bohte [2011]. To deal with the temporal component in SNNs or more generally RNNs, the resulting gradient is usually combined with backpropagation through time (BPTT) . Prominent examples of surrogate gradient approaches include SuperSpike by Zenke and Ganguli [2018] and a number of works with different surrogate derivatives [Wu et al., 2018 , 2019 , Shrestha and Orchard, 2018 , Bellec et al., 2018 , Esser et al., 2016 , Woźniak et al., 2020] . In the general ANN literature, the method is better known as straight-through estimation (STE) and was introduced by Hinton [2012] and by Bengio et al. [2013] in more detail. It was successfully applied by Hubara et al. [2016] and Cai et al. [2017] .Surrogate gradient learning or STE is only heuristically motivated and it is hence desirable to derive a theoretical basis. The influence of different surrogate derivatives on the method has been analyzed through systematic numerical simulations Zenke and Vogels [2021] , revealing that the particular shape has a minor impact compared to the scale. In a more theoretical work, Yin et al. [2019] analyzed the convergence of STE for a Heaviside activation function with three different surrogate derivatives and found that the descent directions of the respective surrogate gradients reduce the population loss when the clipped ReLU function is used as surrogate derivative. Gygax and Zenke [2024] examine how SGL is connected to smoothed probabilistic models [Bengio et al., 2013 , Neftci et al., 2019 , Jang et al., 2019] and to stochastic automatic differentiation [Arya et al., 2022] . In particular, they consider SGL for differentiable activation functions, as we do in our derivation of the SG-NTK.2 Theoretical results 
notation and ntk parametrization: 
 We consider multilayer perceptrons with so-called neural tangent kernel parametrization. For a network with depth L, layer width n l for 0 ≤ l ≤ L, activation function σ, weight matrices W (l) ∈ R n l ×n l-1 , and biases b (l) ∈ R n l , the preactivations with NTK parametrization are given by l+1) for l = 1, . . . , L -1, where σ w > 0, σ b ≥ 0, and we initialize Wh (1) (x) = σ w √ n 0 W (1) x + σ b b (1) , h (l+1) (x) = σ w √ n l W (l+1) σ h (l) (x) + σ b b ((l) ij , b (l) iiid ∼ N (0, 1) for all i, j. The NTK parametrization differs from the standard parametrization by a rescaling factor of 1/ √ n l in layer l + 1. The network function is then given by f ( • ) = h (L) ( • ) : R n0 → R n L and notably the last layer is a preactivation layer. We denote the total number of weights by P . More details can be found in Section C.1.Notation (see Remark C.1 and C.2) By default, we will interpret any vector as a column vector, i.e., we identify R n with R n×1 . This is the case even when writing x = (x 1 , . . . , x n ) ∈ R n for handier notation. Row vectors will be indicated within calculations using the transpose operator, x ⊺ . For a function f : R n1 → R n2 and X = (x 1 , . . . , x d ) ∈ R d•n1 , we define the vector f (X ) := (f (x 1 ),• • • , f (x d )) ∈ R d•n2. For n 2 = 1, we denote the gradient of f by ∇f (x) ∈ R n1 for all x ∈ R n1 . If n 2 > 1, we denote the Jacobian matrix of f by Jf : R n1 → R n2×n1 . A subscript of the form J θ f (x) denotes Jacobian matrices with respect to a subset of variables. We write f (n) ∝ ∼ g(n)to denote that f (n) and g(n) are asymptotically proportional, i.e., f (n) ∼ Kg(n) for some constant K ̸ = 0.We will consider neural networks in the limit of infinitely many hidden layer neurons, i.e., n l → ∞ for all 1 ≤ l ≤ L -1. We will see later, when paraphrasing existing results from the literature, that different ways of taking the number of hidden neurons to infinity can be found. To formalize these notions, we use the definition of a width function from Matthews et al. [2018] with slight modifications: Definition 2.1 (Width function). For every layer l = 0, . . . , L and any m ∈ N, the number of neurons in that layer is given by r l (m), and we call r l : N → N the width function of layer l. We say that a width function r l is strictly increasing if r l (m) < r l (m + 1) for all m ≥ 1. We set R L := (r l ) L-1 l=1 | r l is a strictly increasing width function for all 0 < l < L , the set of collections of strictly increasing width functions for network depth L.Every element of R L provides a way to take the widths of the hidden layers to infinity by setting n l = r l (m) for any 1 ≤ l < L and considering m → ∞. The notions of infinite-width limits used in the literature will now correspond to classes R ⊆ R L for which the respective limiting statements hold. This is captured in the following definition. Definition 2.2 (Types of infinite-width limits). Consider a statement S of the form "Let an ANN have depth L and network layer widths defined by n 0 , n L , and n l := r l (m) for 1 ≤ l < L and some (r l ) L-1 l=1 ∈ R L . Then, for fixed n 0 and any n L , the statement P holds as m → ∞." We also write the statement S as S(r).(i) We say that such a statement S holds strongly, if S(r) holds for any r ∈ R L . This can be interpreted as requiring that the statement holds as min 1≤l<L (n l ) → ∞. We will also write "P holds as n 1 , . . . , n L-1 → ∞ strongly".(ii) We say that such a statement S holds for(n l ) 1≤l≤L-1 ∝ ∼ n, if S holds for all r ∈ R L with r l (m) ∝ ∼ m for all 1 ≤ l < L.This means that S(r) holds for all r ∈ R L such that r p (m)/r q (m) → α p,q ∈ (0, ∞) as m → ∞. We will also write "P holds as (n l ) 1≤l<L ∝ ∼ n".(iii) We say that such a statement S holds weakly, if S holds for at least one r ∈ R L . We will also write "P holds as n 1 , . . . , n L-1 → ∞ weakly". This type of infinite-width limit is tightly connected to the sequential infinite-width limit. Remark 2.1 (Connection between weak and sequential infinite-width limits). In the sequential infinite-width limit, meaning that n 1 → ∞, . . . , n L-1 → ∞ sequentially, the layer widths are not strictly finite. This is opposed to applications, where layer widths may be large but finite. Hence, the infinite-width limits using width functions as explained above are more meaningful in practice. For a sequential limit lim n1→∞ lim n2→∞ f (n 1 , n 2 ) = f * , we can find functions n 1 (m) and n 2 (m) such that lim m→∞ f (n 1 (m), n 2 (m)) = f * . However, the rate at which n 1 (m), n 2 (m) diverge as m → ∞ cannot generally be controlled. Since the weak infinite-width limit allows for arbitrary rates, any sequential infinite-width limit can hence be turned into a weak infinite-width limit as defined in Definition 2.2 (iii).We use Definition 2.2 to paraphrase the convergence of ANNs to GPs in the infinite-width limit: Theorem 2.1 (Theorem 4 from Matthews et al. [2018] ). Any network function f of depth L defined as in Section 2.1 with continuous activation function σ that satisfies the linear envelope property, i.e., there exist c, m ≥ 0 with |σ(u)| ≤ c + m|u| for all u ∈ R, converges in distribution as n 1 , . . . , n L-1 → ∞ strongly to a multidimensional Gaussian process (X j ) n L j=1 for any fixed countable input set L) ) where the covariance function Σ (L) is recursively given byEQUATIONΣ (1) (x, x ′ ) = σ 2 w n 0 ⟨x, x ′ ⟩ + σ 2 b , Σ (L) (x, x ′ ) = σ 2 w E g∼N (0,Σ (L-1) ) [σ(g(x)) σ(g(x ′ ))] + σ 2 b . (1)We also write Σ (L) = Σ σ . The theorem states that the distribution of the network function, which is given by its randomly initialized weights, approaches the distribution of independent GPs as the hidden layer widths increase. Hence, a large-width network will approximately be a realization of the respective GPs. Note that this result can be generalized to non-continuous activation functions without well-defined derivatives that fulfil the linear envelope property, like σ(z) = sign(z). We provide a rigorous proof of this kind of generalization for the case n 1 , . . . , n L-1 → ∞ weakly in Theorem E.3. Σ σ remains well-defined in this case since the expectation in Equation ( 1) does. When approximating the sign function with scaled error function, i.e., σ(z) = erf m (z) := erf(m • z), it holds that lim m→∞ Σ erfm = Σ sign due to the dominated convergence theorem.We consider a datasetD = (X , Y) with X = (x 1 , . . . , x d ) ∈ R d•n0 and Y = (y 1 , . . . , y d ) ∈ R d•n L .To solve the regression problem, i.e., to find weights θ ∈ R P such that f (x i ; θ) = y i for all i = 1, . . . , d, we apply gradient descent in continuous time, also know as gradient flow, with learning rate η and loss function L :R d•n L × R d•n L → R.This means that, using the chain rule, the learning rule can then be written asEQUATIONTo derive the NTK, we observe that the network function then evolves according tod dt f (x; θ t ) = J θ f (x; θ t ) d dt θ t (2) = -η J θ f (x; θ t )J θ f (X ; θ t ) ⊺ ∇ f (X ;θt) L(f (X ; θ t ); Y) (3) =: -η Θt (x, X )∇ f (X ;θt) L(f (X ; θ t ); Y),where we implicitly defined the empirical neural tangent kernel as Θ(x, x ′ ) := J θ f (x; θ)J θ f (x ′ ; θ) ⊺ , which depends on the current weights θ = θ t . This means that the learning dynamics of the network functions during training are given by a kernel whose entries are the scalar products between the gradients of the network's output neuron activity, Θi j (x,x ′ ) = ⟨∇ θ f i (x; θ), ∇ θ f j (x ′ ; θ)⟩.The key result on the NTK is that the empirical NTK converges in the infinite-width limit to a constant kernel, the analytic NTK, at initialization, θ = θ 0 , and during training, θ = θ t : Theorem 2.2 (Theorem 1 from Jacot et al. [2018] for general σ w > 0). For any network function of depth L defined as in Section 2.1 with Lipschitz continuous activation function σ, Θ =: Θ(L) converges in probability to a constant kernelΘ (L) ⊗ I n L as n 1 , . . . , n L-1 → ∞ weakly. This means that for all x, x ′ ∈ R n0 and 1 ≤ i, j ≤ n L , it holds Θ(L) i j (x, x ′ ) - → δ ij Θ (L) (x, x ′) in probability, where δ ij denotes the Kronecker delta. We call Θ (L) the analytic neural tangent kernel of the network, which is recursively given byΘ (1) (x, x ′ ) = Σ (1) (x, x ′ ), Θ (L) (x, x ′ ) = Σ (L) (x, x ′ ) + Θ (L-1) (x, x ′ ) • Σ(L) (x, x ′ ),where Σ (l) are defined as in Theorem 2.1 and we defineEQUATIONWe also write Θ (L) = Θ σ . If θ t are the weights during gradient flow learning at time t ≥ 0 as before, Theorem 2.2 shows that Θ(L) t → Θ (L) ⊗ I n L for t = 0 in the infinite-width limit. This reveals that the gradients of the output neurons, ∇ θ f i (x; θ), are mutually orthogonal. Under additional assumptions this convergence also holds for the whole training duration, t > 0, see Theorem 2 of Jacot et al. [2018] for the case n 1 , . . . , n L-1 → ∞ weakly, and Chapter G of Lee et al. [2019] for (n l ) 1≤l<L ∝ ∼ n. Hence, the kernel that describes the learning dynamics stays constant in the infinite-width limit. This implies that the distribution of the network function during training also converges to GPs [Lee et al., 2019, Theorem 2.2] . Then, any output neuron after training has mean m(x) = Θ (L) (x, X )Θ (L) (X , X ) -1 Y under the assumption of a mean squared error (MSE) loss, see Section C.2.1. The expression for the mean is equivalent to kernel regression with the NTK.The above results show that gradient flow for networks with randomly initialized weights is characterized by the analytic NTK Θ (L) in the infinite-width limit. Since we are interested in gradient flow for networks with activation functions inadequate for gradient descent training, we want to know whether the analytic NTK can be extended to such activation functions. We see that the derivative of the activation function does not have to be defined point-wise for Equation (4). In particular, activation functions with distributional derivatives, e.g., d dz sign(z) = 2 δ(z) can be taken into consideration, where δ denotes the Dirac delta distribution. If we again approximate the sign function with scaled error functions erf m , m ∈ N, it holdsEQUATIONin case the right-hand side exists. A rigorous analysis of this limit can be found in Section D.1. Heuristically, a simple observation suffices: if x = y, we have a one-dimensional integral over two delta distributions, which yields infinity. On the other hand, if x ̸ = y and Σ sign is non-degenerate, a two-dimensional integral over two delta distributions yields a finite value. We derive analytic expressions for lim m→∞ Θ erfm =: Θ sign in Lemma D.3. We call this kernel singular, becauseΘ sign (x, x) = ∞ and Θ sign (x, y) ∈ R if x ̸ = y.By the same reasoning, this divergence occurs for any activation function with jumps. Note that for the mean given by kernel regression, this implies m(x i ) = y i and m(x) = 0 if x is not a training point. Intuitively, this is because the network is initialized with zero mean and different input points are uncorrelated under the singular kernel, so only the training points are learned. A limit kernel with this property also arises if the activation function is fixed but the depth of the network goes to infinity as was shown by Radhakrishnan et al. [2023] . We adopt the ideas of their proof to show that the sign of m is still useful for classification:Theorem 2.3 (Inspired by Lemma 5 of Radhakrishnan et al. [2023] ; see Theorem D.4). Let σ 2 b > 0 or let all x i ∈ R n0 be pairwise non-parallel. Let L ≥ 2 and x i ∈ S n0-1 R for all i = 1, . . . , d, whereS n0-1 R ⊆ R n0 is the sphere of radius R. Assuming that Θ (L) ∞ (x, X )Y ̸ = 0 for almost all x ∈ S n0-1 R , it holds lim m→∞ sign Θ (L) m (x, X )Θ (L) m (X , X ) -1 Y = sign Θ (L) ∞ (x, X )Y a.e. on S n0-1 R . The estimator Θ (L) ∞ (x, X )Y = n i=1 Θ (L)∞ (x, x i ) y i has the form of a so-called Nadaraya-Watson estimator and is well-defined for singular kernels such as Θ (L) . If we assume a classification problem in the sense thatY ⊆ {-1, 1} n , it holds sign Θ (L) ∞ (x i , X )Y = Θ (L) ∞ (x i , X )Y at training point x i .The above singularity of the limit kernel can be avoided by considering surrogate gradient learning instead of gradient descent. First, we introduce a generalization of the NTK that later leads to the surrogate gradient NTK.By definition and originally due to the chain rule, the empirical NTK consists of two Jacobian matrices of the network function with respect to the weights. The Jacobian matrix can be thought of as a recursive formula, J θ f (x; θ) = G(x, θ; σ, σ), which is given by the input x and the architecture of the network, including the activation function σ and its derivative σ. This formula can be modified to define a quasi-Jacobian matrix, J σ1,σ1 (x; θ) := G(x, θ; σ 1 , σ1 ), where σ1 does not have to be the derivative of σ 1 . Analogous to the definition of the empirical NTK we define the empirical generalized NTK to be Î(x, x ′ ) := J σ1,σ1 (x; θ) J σ2,σ2 (x ′ ; θ) ⊺ , where σ 1 , σ1 , σ 2 , σ2 are specified in the following theorem.Theorem 2.4 (Generalized version of Theorem 1 by Jacot et al. [2018] ; see Theorem E.4). For activation functions σ 1 , σ 2 and so-called surrogate derivatives σ1 , σ2 such that σ 1 , σ 2 , σ1 , and σ2 satisfy the linear envelope property and are continuous except for finitely many jump points, denote the empirical generalized neural tangent kernelÎ(L) (x, x ′ ) = J (L),σ1,σ1 (x; θ) J (L),σ2,σ2 (x ′ ; θ) ⊺ for x, x ′ ∈ R n0 , as before. Then, for any x, x ′ ∈ R n0 and 1 ≤ i, j ≤ n L , it holds Î(L) ij (x, x ′ ) P -→ δ ij I (L) (x, x ′), as n 1 , . . . , n L-1 → ∞ weakly. We call I (L) the analytic generalized neural tangent kernel, which is recursively given byI (1) (x, x ′ ) = Σ (1) 1,2 (x, x ′ ), I (L) (x, x ′ ) = Σ (L) 1,2 (x, x ′ ) + I (L-1) (x, x ′ ) • Σ(L) 1,2 (x, x ′ ) for L ≥ 2, with Σ (L) 1,2 (x, x ′ ) = σ 2 w E[σ 1 (Z 1 ) σ 2 (Z 2 )] + σ 2 b for L ≥ 2 and Σ 1,2 (x, x ′ ) = σ 2 w n 0 ⟨x, x ′ ⟩ + σ 2 b ,whereΣ(L) 1,2 (x, x ′ ) = σ 2 w E[σ 1 (Z 1 ) σ2 (Z 2 )] and (Z 1 , Z 2 ) ∼ N 0, Σ (L-1) 1 (x,x) Σ (L-1) 1,2 (x,x ′ ) Σ (L-1) 1,2 (x,x ′ ) Σ (L-1) 2 (x ′ ,x ′ ).Σ 1 and Σ 2 denote the covariance functions of the Gaussian processes that arise from network functions f 1 , f 2 with activation functions σ 1 , σ 2 , respectively, in the infinite-width limit. The covariance between these two Gaussian processes is denoted by Σ 1,2 . This covariance function is asymmetric in the sense thatCov[f 1 (x 1 ), f 2 (x 2 )] ̸ = Cov[f 1 (x 2 ), f 2 (x 1 )] in general.We show in Theorem 2.4 that the generalized empirical NTK tends to a generalized analytic NTK at initialization in an infinite-width limit. Now, the SGL rule can be written using the generalized NTK, similar to Equation (3):EQUATIONHere, we chose σ 1 = σ 2 = σ, σ1 = σ and σ2 = σ in the previous definition of the generalized NTK. This we call the surrogate gradient NTK (SG-NTK) with activation function σ and surrogate derivative σ. Compared to the classical NTK, one of the true Jacobian matrices is replaced by the quasi-Jacobian matrix, since the learning rule is SGL instead of gradient flow. Note that we assume that the derivative of the activation function, σ, exists. Theorem 2.4 shows convergence at time t = 0.We extend this convergence to t > 0 in the following theorem: Theorem 2.5 (Based on Theorem G.2 from Lee et al. [2019] ; see Theorem E.5). Let σ, σ, σ as before, all Lipschitz continuous, and σ, σ bounded. Let f t be a network function with depth L initialized as in Section 2.1 and trained with MSE loss and SGL with surrogate derivative σ. Assume that the generalized NTK converges in probability to the analytic generalized NTK of Theorem 2.4, Î(L) → I (L) ⊗ I n L , as (n l ) L-1 l=1 ∝ ∼ n. Furthermore, assume that the smallest and largest eigenvalue of the symmetrization of I (L) (X , X ), S (L) := I (L) (X , X ) + I (L) (X , X ) ⊺ /2, are given by 0 < λ min ≤ λ max < ∞ and that the learning rate is given by η > 0. Then, for any δ > 0 there exist R > 0, N ∈ N and K > 1 such that for every n ≥ N , the following holds with probability at least 1 -δ over random initialization:sup t∈[0,∞) Î(L) t (X , X ) -I (L) (X , X ) F ≤ 6K 3 R λ min n -1 2, where ∥ • ∥ F denotes the Frobenius norm.We also write I (L) = I σ,σ or simply I σ . The explicit analytic expression of the SG-NTK is derived in Section E.2 for activation function σ = erf m , m ∈ N and surrogate derivative σ = erf as well as general surrogate derivatives. I erfm,σ converges as m → ∞, because compared to Equation (5) we obtain E[ ė rf m (g(x)) σ(g(y))] → E[2δ(g(x)) σ(g(y))] and the delta distribution yields a finite value. In Section E.2, we show this rigorously and derive the analytic expressions. Hence, we define I sign,σ := lim m→∞ I erfm,σ .For any m ∈ N we can consider the SGL dynamics given by Equation ( 7) in the infinite-width limit to obtaind dt f (x ; θ t ) = -η I erfm (x, X ) ∇ f (X ;θt) L(f (X ; θ t ); Y).With MSE error, this equation is solved by a GP with mean m(x) = I erfm (x, X )I erfm (X , X ) -1 Y for t → ∞ and it is natural to assume that the networks trained with SGL converge in distribution to this GP in the infinite-width limit, analogous to the standard NTK case [Lee et al., 2019, Theorem 2.2 ]. However, the empirical counterpart to I sign does not exist as we cannot formulate an empirical SG-NTK for the sign activation function due to the missing Jacobian matrix, compare Equation (6). We suggest that even in this case the network trained with SGL will converge in distribution to the GP given by I sign , since SGL with activation function erf m will approach SGL with activation function sign as m → ∞ and the GP given by I erfm will approach the GP given by I sign as m → ∞. Numerical experiments in Section 3 indicate that this is indeed the case. We conclude that, remarkably, the SG-NTK can be defined for network functions without well-defined gradients and is informative about their learning dynamics under SGL. 
numerical experiments: 
 We numerically illustrate the divergence of the analytic NTK, Θ erfm , shown in Section 2.3 and the convergence of the SG-NTK in the infinite-width limit, Î(L) → I [Radhakrishnan et al., 2023] and thus the angle between its two arguments, ∆α = ∢(x, x ′ ). Hence, we plot the NTKs as functions of this angle, where ∆α = 0 corresponds to x = x ′ .In Figure 1 , the empirical and analytic NTKs for the networks described above and trained with gradient descent are plotted for m ∈ {2, 5, 20} and hidden layer widths n ∈ {10, 100, 500, 1000}.In addition, the analytic NTK for m → ∞ is plotted. Note that the steep slope of erf m for m = 20 results in erf m being very close to the sign function. For any m, we observe that the empirical NTKs converge to the analytic NTK both at initialization and after training as NTK theory states. For Figure 2 , we use the same setup as before, but train the networks using SGL with surrogate derivative σ = ė rf, and compare the empirical and analytic SG-NTKs instead of NTKs. We observe that the empirical SG-NTKs converge to the analytic SG-NTK as n → ∞ both at initialization and after training in accordance with Theorem 2.4 and Theorem 2.5. Figure B.4 illustrates this further by displaying the mean squared errors between empirical SG-NTKs and respective analytic SG-NTK.Moreover, we observe that the analytic SG-NTKs indeed converge to a finite limit as m → ∞, as shown in Section 2.4. Finally, we consider SGL for networks with the same architecture and training objective as before, but with sign activation function, which can be seen as the case m → ∞ of the setups above. We examine whether the distribution of network functions trained with SGL agrees with the distribution of the GP given by the limit kernel I sign . Specifically, we compare 500 networks trained with SGL for t = 30000 time steps, which represent the distribution of the network function after training, to the mean and confidence band of the GP. The mean of the GP is given by kernel regression using the SG-NTK, m(x) = I sign (x, X )I sign (X , X ) -1 Y, and the confidence band is given by m(x) ± 2σ GP (x), where σ GP (x) is the standard deviation at x. We observe in Figure 3a that the mean of the trained networks is close to the GP's mean for network width n = 500 and that most networks lie within the confidence band. The mean of the networks differs from the kernel regression using the kernel Σ sign .Figure 3b shows that this agreement between SGL and the SG-NTK already exists for a network width of n = 100, demonstrating that the SG-NTK predicts the SGL dynamics of networks with moderate width. Note that the variance in the networks' output and the confidence band can be reduced (see Arora et al. [2019] and Section A).  
conclusion: 
 Gradient descent training is not applicable to networks with sign activation function. In the present study, we have first shown that this is even true for the infinite-width limit in the sense that the NTK diverges to a singular kernel. We found that this singular kernel still has interesting properties and allows for classification, but it is unusable for regression.We then studied SGL, which is applicable to networks with sign activation function. We defined a generalized version of the NTK that can be applied to SGL and derived a novel SG-NTK. We proved that the convergence of the NTK in the infinite-width limit extends to the SG-NTK, both at initialization and during training. Strikingly, we were able to derive an SG-NTK for the sign activation function, I sign , by approximating the sign function with error functions. We suggest that this SG-NTK predicts the learning dynamics of SGL, and support this claim with heuristic arguments and numerical simulations.A limitation of our work is that due to the considered NTK framework, our results are naturally only applicable to sufficiently wide networks with random initialization. Further, we only prove the convergence of the SG-NTK during training for activation functions with well-defined derivatives. More rigorous analysis should be carried out on how the connection between SGL and the SG-NTK carries over to activation functions with jumps, as shown by our simulations.Our derivation of the SG-NTK opens a novel path towards addressing the many unanswered questions regarding the training of binary networks, in particular regarding the class of functions that SGL learns for wide networks and how that class differs for different activation functions and surrogate derivatives.  
experiments compute resources: 
 Question: For each experiment, does the paper provide sufficient information on the computer resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?Answer: [Yes] Justification: The compute resources used for the numerical simulation are described in Section A.Guidelines:• The answer NA means that the paper does not include experiments.• The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage. • The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute. • The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make it into the paper). 
code of ethics: 
 Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https://neurips.cc/public/EthicsGuidelines?Answer: [Yes] Justification: We reviewed the NeurIPS Code of Ethics and the paper conforms with it.Guidelines:• The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.• If the authors answer No, they should explain the special circumstances that require a deviation from the Code of Ethics. • The authors should make sure to preserve anonymity (e.g., if there is a special consideration due to laws or regulations in their jurisdiction). 
broader impacts: 
 Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?• The authors should state which version of the asset is used and, if possible, include a URL. • The name of the license (e.g., CC-BY 4.0) should be included for each asset.• For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided. • If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswithcode.com/datasets has curated licenses for some datasets. Their licensing guide can help determine the license of a dataset. • For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided. • If this information is not available online, the authors are encouraged to reach out to the asset's creators. 
new assets: 
 Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?Answer: [NA] Justification:Guidelines:• The answer NA means that the paper does not release new assets.• Researchers should communicate the details of the dataset/code/model as part of their submissions via structured templates. This includes details about training, license, limitations, etc. • The paper should discuss whether and how consent was obtained from people whose asset is used. • At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file. 
crowdsourcing and research with human subjects: 
 Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well as details about compensation (if any)?Answer: [NA] Justification:Guidelines:• The answer NA means that the paper does not involve crowdsourcing nor research with human subjects. • Including this information in the supplemental material is fine, but if the main contribution of the paper involves human subjects, then as much detail as possible should be included in the main paper. • According to the NeurIPS Code of Ethics, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector. 
institutional review board (irb) approvals or equivalent for research with human subjects: 
 Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained?Answer: [NA] Justification:Guidelines:• The answer NA means that the paper does not involve crowdsourcing nor research with human subjects. • Depending on the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper. • We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution. • For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review. 
 
 Appendix: Weight initialization and implementation. All networks are initialized with σ w = 1, σ b = 0.1 For the implementation of the NTK and SG-NTK we use the JAX package [Bradbury et al., 2018] and Neural Tangents package [Novak et al., 2020 , 2022 , Han et al., 2022 , Sohl-Dickstein et al., 2020 , Hron et al., 2020] with modifications. Computations were done using an Intel Core i7-1355U CPU and 16 GB RAM. The simulations for Figure 1 and Figure 2 took two hours each. The simulations for Figure 3 took 12 hours.Variance reduction trick. The variance in the outputs of the networks trained with gradient descent or SGL and the confidence band given the NTK or SG-NTK respectively can be reduced by multiplying the weights of the last layer with a constant κ < 1 at initialization [Arora et al., 2019] . This is explained in detail at the end of Section C.2.1. We can see from Figure B .2 that the agreement between the distribution of the trained networks and the distribution given by the SG-NTK still holds; however, network width and training time have to be increased. We begin by defining an artificial neural network with a parameterization suitable for considering the limit of infinitely many hidden neurons. This parameterization is called neural tangent kernel parameterization and differs from standard parameterization of multilayer perceptrons by a rescaling factor of 1/ √ n l in layer l + 1, where n l is the width of layer l. We follow the slightly more general definition in [Lee et al., 2019, Equation (1) ]. A discussion of this kind of parameterization can be found in [Jacot et al., 2018, Remark 1] .Definition C.1 (Artificial neural network with NTK parameterization). Let L be the depth of the network, n k for k = 0, . . . , L the widths of the layers, and σ : R → R an activation function. We draw network weight matricesW (l) ∈ R n l ×n l-1 and biases b (l) ∈ R n l for l = 1, . . . , L from a probability distribution such that W (l) ij , b (l) i iid ∼ N (0, 1).For the parameters, we denote by θ (l) = (W (l) , b (l) ) the parameters of layer l and by θ (1 : l) = (θ (1) , θ (2) , . . . , θ (l) ) the parameters of layers 1 up to and including l. Given some σ w > 0 and σ b ≥ 0 we then define for all x ∈ R n0h (1) x; θ (1) = σ w √ n 0 W (1) x + σ b b (1) , h (l+1) x; θ (1 : l+1) = σ w √ n l W (l+1) σ h (l) x; θ (1 : l) + σ b b (l+1) for l = 1, . . . , L -1. Therefore, h (l) ( • ; θ (1 : l)) is a map from R n0 to R n l and we use the short-hand notation h (l) (x; θ) = h (l) (x; θ (1 : l) ). Finally, we define our network functionf ( • ; θ) = h (L) ( • ; θ) : R n0 → R n L .With this definition, the total number of parameters, P = |θ|, isP = L l=1 θ (l) = L l=1 n l (n l-1 + 1).Given such a network function with NTK parameterization, we consider a datasetD = (X , Y) with X = (x 1 , . . . , x d ) ∈ R d•n0 and Y = (y 1 , . . . , y d ) ∈ R d•n L .First, we want to solve the regression problem, i.e., find parameters θ ′ such that f (x i ; θ ′ ) = y i for all i = 1, . . . , d. Later, we will also consider the classification problem, i.e., we assume y i ∈ {-1, 1} and want to solve sign(f (x i ; θ ′ )) = y i for all i = 1, . . . , d. Tackling both cases from the regression perspective, we define the so-called loss functional and loss function. The following two definitions are similarly formulated by Jacot et al. [2018] .Definition C.2 (Loss functional and loss function). Let F = {g | g : R n0 → R n L } and L : F → R a so-called loss functional. In addition, let L be convex, i.e., for all λ ∈ [0, 1] and g 1 , g 2 ∈ F it holdsL(λg 1 + (1 -λ)g 2 ) ≤ λL(g 1 ) + (1 -λ)L(g 2 ).We will assume that L can be written asL(f ) = 1 d d i=1 ℓ(f (x i ; θ); y i ) =: L(f (X ); Y), so that the so-called loss function L( • ; Y) : R d•n L → R is differentiable.Remark C.1 (Function evaluation at sets and vector notation). We want to detail the notation used in Definition C.2. For a function f : R n0 → R n L and X = (x 1 , . . . , x d ) ∈ R d•n0 we define the vectorf (X ) := (f (x 1 ), • • • , f (x d )) ∈ R d•n L .By default, we will interpret any vector as a column vector, i.e., we identify R n with R n×1 . This is the case even when writing x = (x 1 , . . . , x n ) ∈ R n for handier notation. Row vectors will be indicated within calculations using the transpose operator, x ⊺ .Let us first consider regular gradient descent learning in continuous time, also known as gradient flow.For this, we assume that our network function is differentiable with respect to its parameters.Remark C.2 (Gradient and Jacobian matrix notation). Let f : R n → R m . For m = 1, we denote the gradient by ∇f (x) ∈ R n for all x ∈ R n . If m > 1 we denote the Jacobian matrix of f by Jf : R n → R m×n . Therefore, Jf (x) is a m × n matrix for all x ∈ R n . We do not always want to consider the gradient or Jacobian matrix with respect to all variables. We indicate this with subscripts as follows. Let f :R n × R P → R m and g x (θ) = f (x; θ) for fixed x ∈ R n . Then, we write ∇ θ f (x; θ) := ∇g x (θ) and J θ f (x; θ) := Jg x (θ).In particular, for a map f : R P → R, the gradient with respect to the j-th variable, 1 ≤ j ≤ P , is a scalar and denoted by ∂ j f (θ) := ∇ θj f (θ). This is called the partial derivative of f with respect to the j-th variable.With this notation, we consider the gradient flow method with learning rate η > 0 and recall the derivation of the neural tangent kernel. We move the weights in the opposite direction of the gradient of the loss function with respect to the parameters of the network evaluated at the training points:d dt θ t = -η ∇ θ L(f (X ; θ t ); Y)) = -η J θ f (X ; θ t ) ⊺ ∇ f (X ;θt) L(f (X ; θ t ); Y) (S8) = -η 1 d d i=1 J θ f (x i ; θ t ) ⊺ ∇ f (xi;θt) ℓ(f (x i ; θ t ); y i ),using the chain rule for the second equality and with A ⊺ denoting the transpose of a matrix A. Again using the chain rule, this then implies for anyx ∈ R n0 d dt f (x; θ t ) = J θ f (x; θ t ) d dt θ t (S8) = -η J θ f (x; θ t )J θ f (X ; θ t ) ⊺ ∇ f (X ;θt) L(f (X ; θ t ); Y) (S9) = -η 1 d d i=1 J θ f (x; θ t )J θ f (x i ; θ t ) ⊺ ∇ f (xi;θt) ℓ(f (x i ; θ t ); y i ). (S10)We therefore define the neural tangent kernel as follows: Definition C.3 (Neural tangent kernel). Let f be a network function of depth L as in Definition C.1 with parameters θ, not necessarily drawn randomly. Then, we define the neural tangent kernel (NTK) as:Θ(L) : R n0 × R n0 → R n L ×n L (x, y) → J θ f (x; θ)J θ f (y; θ) ⊺ .Therefore, it holds for all x, y ∈ R n0 and1 ≤ i, j ≤ n L Θ(L) i j (x, y) = P p=1 ∂ θp f i (x; θ) ∂ θp f j (y; θ).Notice that the NTK depends on the parameters of the networks. It is therefore initialized randomly and varies over the course of the training. With notation f t (x) = f (x; θ t ) and Θ(L) t = Θ(L) for parameters θ t at training time t we can now rewrite Equations ( S9) and (S10) as follows:d dt f t (x) = -η Θ(L) t (x, X )∇ ft(X ) L(f t (X ); Y) (S11) = -η 1 d d i=1 Θ(L) t (x, x i ) ∇ ft(xi) ℓ(f t (x i ); y i ).We are hence able to express the change of the network function during training in a kernel fashion. Later, we will consider this change of the network function in the infinite-width limit, i.e., n 1 , . . . , n L-1 → ∞.Before doing so, we will generalize the NTK definition in order to apply the NTK to surrogate gradient learning later. In particular, we will break the symmetry of the above definition and generalize the Jacobian matrices to quasi-Jacobian matrices by replacing the derivatives of the activation function by surrogate derivatives. Let us write the recursive formula of the Jacobian matrix of the network function given by the chain rule as J θ f (x; θ) = G(σ; σ; x; θ), where σ is the derivative of the activation function σ. Then, surrogate gradient learning replaces J θ f (x; θ) with G(σ; σ; x; θ) for the surrogate derivative σ of the activation function σ. We call this the quasi-Jacobian matrix: Definition C.4 (Quasi-Jacobian matrices for neural networks). Let L be the depth of the network, n k for k = 0, . . . , L the width of the layers, σ : R → R the activation function, and σ : R → R the so-called surrogate derivative of the activation function. Let f be the network function, h (l) , l = 1, . . . , L -1, the intermediate layers as in Definition C.1 and θ the network parameters. We then define the quasi-Jacobian matrix J (L) of f at point x recursively as follows:J (1) x; θ (1) ∈ R n1×|θ (1) | with J (1) k θp x; θ (1) = δ ki σw √ n0 x j if θ p = W (1) ij δ ki σ b if θ p = b (1) i (S12) J (l) x; θ (1 : l) ∈ R n l ×|θ (1 : l) | for 2 ≤ l ≤ L with J (l) k θp x; θ (1 : l) =        δ ki σw √ n l-1 σ h (l-1) j (x; θ) if θ p = W (l) ij δ ki σ b if θ p = b (l) i σw √ n l-1 n l-1 j=1 W (l) ij σ h (l-1) j (x; θ) J (l-1) j θp x; θ (l-1) if θ p ∈ θ (l-1) . (S13)Remark C.3 (Notations for the quasi-Jacobian). With the above definition of the quasi-Jacobian matrix of the network function f with activation function σ and surrogate derivative σ we write ) x; θ (L) .J (L),σ,σ (x; θ) := J (LIt then holds J (L),σ, σ (x; θ) = J θ f (x; θ). For a data set X instead of a single point x, we concatenate the matrices row-wise as before, namelyJ (L) (X ; θ) ∈ R d•n L ×|θ| .Definition C.5 (The generalized neural tangent kernel). Let σ 1 , σ 2 be activation functions and σ1 , σ2 the surrogate derivatives respectively. Given a network depth L and parameters θ we define the generalized neural tangent kernel as:Î(L) : R n0 × R n0 → R n L ×n L (S14) (x, y) → J (L),σ1,σ1 (x; θ) J (L),σ2,σ2 (y; θ) ⊺ .Remark C.4. The generalized neural tangent kernel agrees with the neural tangent kernel in the case where σ = σ 1 = σ 2 and σ = σ1 = σ2 .C.2 Notation for the infinite-width limit and review of key theorems for the NTKIn this section we will formulate all important theorems on the NTK that we will need for our later analysis using the introduced notation. Furthermore, we will discuss and remark their proofs, in particular in view of the generalizations that will be proved in Section E.1.Convergence of networks to Gaussian processes in the infinite-width limit. We will consider neural networks in the limit of infinitely many hidden layer neurons. The fact that such networks converge to Gaussian processes was first mentioned by Neal [1996] . We follow and present the formulations of Matthews et al. [2018] for the general mathematical statement. First, we formalize the limit of infinitely many hidden neurons. Definition C.6 (Width function, as in Definition 3 of Matthews et al. [2018] with modifications). For every layer l = 0, . . . , L and any m ∈ N, the number of neurons at that layer is given by r l (m), and we call r l : N → N the width function of layer l. We say that a width function r l is strictly increasing if r l (m) < r l (m + 1) for all m ≥ 1. We set R L := (r l ) L-1 l=1 | r l is a strictly increasing width function for all 1 ≤ l < L , the set of collections of strictly increasing width functions for network depth L.Every element of R L provides a way to take the widths of the hidden layers to infinity by setting n 0 and n L to some constant, setting n l = r l (m) for any 1 ≤ l < L and considering m → ∞. To formally define for which ways of taking the widths of hidden layers to infinity a statement holds, we can now state the set R ⊆ R L such that the statement holds for widths given by any r ∈ R as m → ∞. Clearly, the claim "The statement holds for all r ∈ R 1 ." is stronger than "The statement holds for all r ∈ R 2 ." if R 2 ⊆ R 1 . On the basis of these considerations, we define three types of infinite-width limits using the previous definition in order to structure the different types of limits in this thesis as well as in the literature. Definition C.7 (Types of infinite-width limits). Consider a statement S of the form "Let an ANN have depth L and network layer widths defined by n 0 , n L and n l := r l (m) for 1 ≤ l < L and some (r l ) L-1 l=1 ∈ R L . Then, for fixed n 0 and any n L , the statement P holds as m → ∞." We also write the statement S as S(r).(i) We say that such a statement S holds strongly, if S(r) holds for any r ∈ R L . This can be interpreted as requiring that the statement holds as min 1≤l<L (n l ) → ∞. We will also write "P holds as n 1 , . . . , n L-1 → ∞ strongly".(ii) We say that such a statement S holds for (n l ) 1≤l≤L-1 ∝ ∼ n, if S holds for all r ∈ R L with r l (m) ∝ ∼ m for all 1 ≤ l < L. In other words S(r) holds for all r ∈ R L such that r p (m)/r q (m) → α p,q ∈ (0, ∞) as m → ∞. We will also write "P holds as (n l ) 1≤l<L ∝ ∼ n".(iii) We say that such a statement S holds weakly, if S holds for at least one r ∈ R L . This can be read as requiring that the statement holds as n 1 → ∞, . . . , n L-1 → ∞ sequentially. We will also write "P holds as n 1 , . . . , n L-1 → ∞ weakly". Theorem C.1 (Theorem 4 from Matthews et al. [2018] ). Any network function f of depth L defined as in Definition C.1 with continuous activation function σ that satisfies the linear envelope property, i.e., there exist c, m ≥ 0 with |σ(u)| ≤ c + m|u| ∀u ∈ R, converges in distribution as n 1 , . . . , n L-1 → ∞ strongly to a multidimensional Gaussian process (X j ) n L j=1 for any fixed countable input set(x i ) ∞ i=1 . It holds X j iid ∼ N (0, Σ (L)) where the covariance function Σ (L) is recursively given byΣ (1) (x, x ′ ) = σ 2 w n 0 ⟨x, x ′ ⟩ + σ 2 b , Σ (L) (x, x ′ ) = σ 2 w E g∼N (0,Σ (L-1) ) [σ(g(x)) σ(g(x ′ ))] + σ 2 b . Remark C.5.First, a proof of the above theorem can be found in the paper of Matthews et al. [2018] . While it takes a lot of effort to show that the statement holds strongly in the sense of Definition C.7, the weak version of the statement can be proved via induction. This has been done by Jacot et al.[2018] and we will later adapt their proof to show a generalized version, Theorem E.3.Second, in the context of analyzing the network behavior, we are interested in the finite-dimensional distributions first of all, since neural networks are trained and tested on a finite number of data points. From the convergence of the marginal distributions, we can infer the convergence to an stochastic process via the Kolmogorov extension theorem. However, this assumes the product σ-algebra, which is why Theorem C.1 assumes a fixed countable input set. Matthews et al. [2018] have discussed these formal restrictions in more detail (Chapter 2.2). If one does not want to be restricted to a countable index set, one could, for example, consider the condition by Prokhorov [1956, Theorem 2.1]. A similar approach was taken by Bracale et al. [2021] , which applied the Kolmogorov-Chentsov criterion [Kallenberg, 2021, Theorem 4.23] .Finally, note that the theorem assumes continuity of the activation function. In the proof of Matthews et al. [2018] this is only used in order to apply the continuous mapping theorem. However, it is sufficient for the limiting process to attain possible points of discontinuity with probability zero for the continuous mapping theorem to be applicable. The theorem is thus also valid for activation functions that are continuous except at finitely many jump points, such as step-like activation functions.Convergence of the NTK at initialization in the infinite-width limit. Jacot et al. [2018] showed that the previously defined empirical NTK converges to a deterministic limit, which we will call the analytic NTK. Theorem C.2 (Theorem 1 from Jacot et al. [2018] , slightly generalized). For any network function of depth L defined as in Definition C.1 with Lipschitz continuous activation function σ, the empirical neural tangent kernel Θ(L) converges in probability to a constant kernelΘ (L) ⊗ I n L as n 1 , . . . , n L-1 → ∞ weakly. For all x, x ′ ∈ R n0 and 1 ≤ i, j ≤ n L , it holds Θ(L) i j (x, x ′ ) P --→ δ ij Θ (L) (x, x ′ ),which we also write asΘ(L) P --→ Θ (L) ⊗ I n L .We call Θ (L) the analytic neural tangent kernel of the network, which is recursively given byΘ (1) (x, x ′ ) = Σ (1) (x, x ′ ) Θ (L) (x, x ′ ) = Σ (L) (x, x ′ ) + Θ (L-1) (x, x ′ ) • Σ(L) (x, x ′ ),where Σ (l) are defined as in Theorem C.1 and we defineΣ(L) (x, x ′ ) = σ 2 w E g∼N (0,Σ (L-1) ) [ σ(g(x)) σ(g(x ′ ))] .Compared to Theorem 1 of Jacot et al. [2018] , the statement is slightly generalized in the sense that it allows for arbitrary σ w > 0. The arguments in the proof work the same way.Remark C.6 (Versions of Theorem C.2 in the literature). A proof of this theorem for (n l ) 1≤l<L ∝ ∼ n is given by Yang [2019a] and his proof is also referenced by Lee et al. [2019] . However, the proof is given in terms of so-called tensor programs and therefore harder to follow. For the ReLU activation function, a proof for n 1 , . . . , n L-1 → ∞ strongly is provided by Arora et al. [2019, Theorem 3.1 ].We will later prove a version of this theorem for the generalized NTK.Convergence of the NTK during training in the infinite-width limit. Not only does the NTK converge to a constant kernel in the infinite-width limit, even the kernel during training, Θ(L) t , converges to this constant kernel. This was also discovered by Jacot et al. [2018] . Theorem C.3 (Theorem 2 by Jacot et al. [2018] ). Assume any network function of depth L defined as in Definition C.1 with Lipschitz continuous activation function σ, twice differentiable with bounded second derivative, and trained with gradient flow as in Equation S10. Let T > 0 such thatT 0 ∥∇ℓ(f t ( • ); f * ( • ))∥ pemp dt = T 0 √ d ∇ ft(X ) L(f t (X ); f * (X )) 2 dt ∈ O p (1), ( * )where X ∈ O p (1) denotes that X is stochastically bounded. Then, as n 1 , . . . , n L-1 → ∞ weakly, the empirical NTK Θ(L) t converges in probability to the analytic NTK Θ (L) ⊗ I n L in probability uniformly for t ∈ [0, T ]. We therefore writeΘ(L) t P --→ Θ (L) ⊗ I n L .Remark C.7 (Versions of Theorem C.3 in the literature). The proof of Jacot et al. [2018] relies heavily on a function space perspective. Since this formulation tends to lack mathematical rigor, we will rely on the proof of the theorem for the case (n l ) 1≤l<L ∝ ∼ n given by Lee et al. [2019, Chapter G] . In particular, the first inequality of (S51) in Theorem G.2 of Lee et al. [2019] implies the condition ( * ). Furthermore, a different approach to proving the above statement for the case (n l ) 1≤l<L ∝ ∼ n using the Hessian matrix of the network function was taken by Liu et al. [2020, Proposition 2.3, Theorem 3.2] . A partial proof of a version of this theorem for the generalized NTK will be given later. Only an auxiliary lemma remains to be proved.Given the results of the previous section, we can formulate an infinite width version of Equation (S10) by replacing the empirical with the analytic NTK. This allows us to analyze the learning dynamics of networks in the infinite-width limit, which yields connections to kernel methods and reproducing kernel Hilbert spaces. We then discuss how far the resulting functions and solutions in the infinite-width limit deviate from the finite width networks. This is essential to evaluate to what extend the results in the infinite-width limit can inform us about the behavior of gradient flow in the finite-width networks. First, we state the infinite-width version of Equation (S11) using Theorem C.3,d dt f t (x) = -η Θ (L) ⊗ I n L (x, X ) ∇ ft(X ) L(f t (X ); Y) = -η 1 d d i=1 Θ (L) (x, x i ) • I n L ∇ℓ(f t (x i ); y i ) = -η d i=1 Θ (L) (x, x i ) 1 d ∇ℓ(f t (x i ); y i ) (S15) = -η Θ (L) (x, X ) 1 d [∇ℓ(f t (x 1 ); y 1 ), . . . , ∇ℓ(f t (x d ); y d )] ⊺ =: -η Θ (L) (x, X ) ∇ ft(X ) L(f t (X ); Y), (S16)where in the last line we interpret∇ ft(X ) L(f t (X ); Y) ∈ R d×n L as a matrix of size d × n L with entries ∇L ft(X ) (f t (X ); Y) i j = 1/d • ∂ j ℓ(f t (x i ); y i ). Recall that ∂ j ℓ(f t (x i ); y i )is the partial derivative of ℓ( • ; y i ) with respect to its j-th entry, i.e., with respect to f t,j (x i ). Note that the last line is a row vector, which we can identify as a column vector. The fact that the NTK is now time-independent and non-random has two interesting implications:• Equation ( S16) is now an differential equation that can be solved explicitly or numerically for certain loss functions. • According to Equation (S15), the time derivative of f t (x) can now be expressed elementwise as a linear combination of functions of the type Θ (L) ( • , x) : R n0 → R. For an arbitrary symmetric and positive definite kernel k( • , • ), the completion of the linear span of functions of this type is called the reproducing kernel Hilbert space (RKHS) of k. Assuming that the solution of Equation ( S15) is an element of the RKHS of Θ (L) , one can ask what the space looks like.The ODE of Equation (S16) has already been considered by Jacot et al. [2018, Chapter 5] and Lee et al. [2019, Chapter 2 .2], and we will follow the observations made there. To do this, we will assume the mean squared error (MSE) loss, L( Ỹ; Y) = 1 2 ∥ Ỹ -Y∥ 2 2 , implying ∇ ft(X ) L(f t (X ); Y) = f t (X ) -Y,(x) = -η Θ (L) (x, X ) (f t (X ) -Y) .Now, for simplicity, we denote Θ(x, y) := Θ (L) (x, y) and Θ := Θ(X , X ). Furthermore, we consider an arbitrary set of test points X T . The solution of the ODE is then given byf t (X T ) = µ t (X T ) + γ t (X T ) for µ t (X T ) = Θ(X T , X )Θ -1 I d -e -ηΘt Y and γ t (X T ) = f 0 (X T ) -Θ(X T , X )Θ -1 I d -e -ηΘt f 0 (X ).Recall that by Theorem C.1, the components f 0,j are independent and identically distributed Gaussian processes with mean zero and covariance function Σ := Σ (L) . Hence, γ t has mean zero and the mean of f t is given by µ t . By looking at the components of f t ,f t,j (X T ) = f 0,j (X T ) -Θ(X T , X )Θ -1 I d -e -ηΘt (f 0,j (X T ) -Y j ) ,we can conclude that they are independent and identically distributed as well. One can show that the components are indeed Gaussian processes again with mean µ t . Using γ t we can also compute the covariance matrix for our arbitrary set of test points X T ,Γ t (X T , X T ) := E [γ t,j (X T ) γ t,j (X T ) ⊺ ] = E [f 0,j (X T ) f 0,j (X T ) ⊺ ] -E f 0,j (X T )f 0,j (X ) ⊺ I d -e -ηΘt Θ -1 Θ(X , X T ) -E Θ(X T , X )Θ -1 I d -e -ηΘt f 0,j (X )f 0,j (X T ) ⊺ + E Θ(X T , X )Θ -1 I d -e -ηΘt f 0,j (X )f 0,j (X ) ⊺ I d -e -ηΘt Θ -1 Θ(X , X T ) = Σ(X T , X T ) -Σ(X T , X ) I d -e -ηΘt Θ -1 Θ(X , X T ) -Θ(X T , X )Θ -1 I d -e -ηΘt Σ(X , X T ) + Θ(X T , X )Θ -1 I d -e -ηΘt Σ(X , X ) I d -e -ηΘt Θ -1 Θ(X , X T ).Assuming that Θ is positive definite immediately leads to pointwise convergence of the mean and covariance functions. This implies that the gradient flow solution for networks in the infinite-width limit converges to a Gaussian process as t → ∞ with mean function µ ∞ and covariance function Γ ∞ given below. This follows from the weak convergence of the finite-dimensional marginal distributions by Lévy's convergence theorem [Williams, 1991, Section 18.1] . Again, the discussions of Remark C.5 applies. We have Lee et al. [2019] state that a network trained with gradient flow will indeed converge in distribution to this Gaussian process as the width goes to infinity:µ ∞ (X T ) = Θ(X T , X )Θ -1 Y and (S17) Γ ∞ (X T , X T ) = Σ(X T , X T ) -Σ(X T , X )Θ -1 Θ(X , X T ) -Θ(X T , X )Θ -1 Σ(X , X T ) + Θ(X T , X )Θ -1 Σ(X , X )Θ -1 Θ(X , X T ). (S18)Theorem C.4 (Theorem 2.2 from Lee et al. [2019] ). Let the learning rate η < η critical forη critical := 2(λ min (Θ (L) (X , X )) + λ max (Θ (L) (X , X )))with a network function f t as in Theorem C.3 with hidden layer widthsn 1 = • • • = n L-1 = n and restricted to x ∈ R n0 with ∥x∥ 2 ≤ 1. If λ min (Θ (L) (X , X )) > 0, then the components of f t converge in distribution to independent, identically distributed Gaussian processes N (µ t , Γ t ) as n → ∞ for all t ∈ [0, ∞) ∪ {∞}.Hence, the result of training a finite-width network with gradient flow for an infinite amount of time will be arbitrarily close in distribution to a Gaussian process with mean function µ ∞ and covariance function Γ ∞ , if the width is sufficiently large. Note that by Equations ( S17) and (S18) the variance at the training points X is zero and the mean at the training points is exactly Y.Since we will focus on the mean, we will first sketch a trick introduced in Chapter 3 of Arora et al. [2019] to make the variance term arbitrarily small. If f 0 had mean variance, this would consequently also be the case for all f t and for the solution f ∞ := lim t→∞ f t . This can be achieved by multiplying f 0 by a small constant κ > 0 and considering the network function g 0 = κf 0 instead. It then holdsΘg (x, y) = J θ g 0 (x; θ)J θ g 0 (y; θ) ⊺ = J θ κf 0 (x; θ)J θ κf 0 (y; θ) ⊺ = κ 2 Θf (x, y),and thus we have Θ g = κ 2 Θ f . In the infinite-width limit, the derivative of g t is then given byd dt g t (x) = -η Θ g (x, X ) (g t (X ) -Y) = -η κ 2 Θ f (x, X ) (κf t (X ) -Y) ,which implies as beforeg t (x) = g 0 (x) -Θ g (x, X )Θ g (X , X ) -1 I d -e -η Θg(X ,X )t (g 0 (X ) -Y) = Θ f (x, X )Θ f (X , X ) -1 I d -e -ηκ 2 Θ f (X ,X )t Y + κ f 0 (x) -Θ f (x, X )Θ f (X , X ) -1 I d -e -ηκ 2 Θ f (X ,X )t f 0 (X ) .Note that the term in the second last line corresponds to the non-random mean of f trained with learning rate ηκ 2 , and that the term in the last line is random, but can be made arbitrarily small using κ. We can think of this as a trade-off between learning rate and variance. This justifies why we can focus on the mean in the next section.To sum up, we are interested in network functions in the infinite-width limit that are trained over time according toEQUATIONwhere we change from a row vector to a column vector in the last equation. The mean of such network functions after infinite training time is given byf NTK (x) := Θ(x, X )Θ(X , X ) -1 Y = µ ∞ (x). (S20)The first observation to make in our attempt to apply the neural tangent kernel to networks with the sign function as activation function is that the sign function has a zero derivative almost everywhere. Thus, the derivative of the network function with respect to the network weights is zero for all weights that are not part of the last layer. The case where the weights θ (1:L-1) are frozen after initialization and only θ (L) is trained has already been discussed by Lee et al. [2019, Chapter 2.3.1 and Chapter D] . For a network in the infinite-width limit, this approach is equivalent to applying Gaussian process regression, i.e., knowing that f ∼ N 0, Σ (L) for infinite width, one considers f | f (X ) = Y. This can be seen by realizing that Θ (L) = Σ (L) if σ = 0 almost everywhere and applying Theorem C.4.While this is an interesting observation, and the strategy of optimizing only the last layer can also be transferred to finite width networks, we would prefer to train the whole network and not identify the derivative of the sign function with zero, since this discards all information about the jump discontinuities in our networks. An obvious alternative would be to use the distributional derivative of the sign function, which is given by 2 δ 0 , where δ 0 denotes the delta distribution. We will see that Σ(L) still exists when the distributional derivative is substituted into its formula. Alternatively, we can obtain the same expression by approximating the sign function with scaled error functions,erf m (z) = erf(m • z) = 2 √ π m•z 0 e -t 2 dt,and considering the limit m → ∞.Due to the previous considerations, we begin by deriving the analytic NTK for the error function.Following the notation of Lee et al. [2019] , we need to find analytic expressions for the termsT m (Σ) := E (X,Y )∼N (0,Σ) [erf m (X) erf m (Y )] and Ṫm (Σ) := E (X,Y )∼N (0,Σ) [ ė rf m (X) ė rf m (Y )].Note that by a change of variables we can alternatively consider the termsT (m 2 • Σ) := E (X,Y )∼N (0,m 2 •Σ) [erf(X) erf(Y )] = T m (Σ) and Ṫ (m 2 • Σ) := E (X,Y )∼N (0,m 2 •Σ) [ ė rf(X) ė rf(Y )] = 1 m 2 Ṫm (Σ). For Σ ′ = ( x•x x•y x•y y•y ), T (Σ ′) and Ṫ (Σ ′ ) are given in Chapter C of the supplementary material of Lee et al. [2019] . However, we cannot assume that Σ ′ always has this form. While Ṫ can be easily calculated, T is harder to deal with and a reference to Williams [1996, Chapter 3 .1] is used. There, the main idea of the proof, how to evaluate a more general expression, is given without further details. We will derive analytic expressions for both terms explicitly.We start by evaluating Ṫ . Note thatd dz erf(z) = 2 √ π e -z 2 and d dz erf m (z) = 2m √ π e -m 2 z 2 . Lemma D.1. Given U ∼ N (0, Σ) with invertible covariance matrix Σ ∈ R d×d and x, y ∈ R d , it holds E[ ė rf(U ⊺ x) ė rf(U ⊺ y)] = 4 π (1 + 2x ⊺ Σx)(1 + 2y ⊺ Σy) -(2x ⊺ Σy) 2 -1/2 . (S21)In particular, given (X, Y ) ∼ N (0, Σ) with invertible covariance matrix Σ ∈ R 2×2 or with X = Y and singular covariance matrix Σ ∈ R 2×2 , it holdsEQUATION)where |A| denotes the determinant of a matrix A.Proof. It holds for U ∼ N (0, Σ) with covariance matrix Σ ∈ R d×d and x, y ∈ R d :E[ ė rf(U ⊺ x) ė rf(U ⊺ y)] = R d 1 (2π) d/2 |Σ| 1/2 2 √ π e -(u ⊺ x) 2 2 √ π e -(u ⊺ y) 2 e -1 2 u ⊺ Σ -1 u du (⋆) = 4 π R d 1 (2π) d/2 exp - 1 2 v ⊺ (I d + 2Σ 1/2 xx ⊺ Σ 1/2 + 2Σ 1/2 yy ⊺ Σ 1/2 )v dv = 4 π I d + 2Σ 1/2 xx ⊺ Σ 1/2 + 2Σ 1/2 yy ⊺ Σ 1/2 -1 1/2 = 4 π I d + 2Σ 1/2 xx ⊺ Σ 1/2 + 2Σ 1/2 yy ⊺ Σ 1/2 -1/2 ,using a change of variable Σ 1/2 u = v for Equation (⋆) and using basic properties of the determinant. We can evaluate the determinant in the last line by applying the Sylvester's determinant theorem [Pozrikidis, 2014, (B.1.16 )], i.e.,|I n + AB ⊺ | = |I m + B ⊺ A| for any matrices A, B ∈ R n×m . We then define A = B = ( √ 2Σ 1/2 x, √ 2Σ 1/2 y) ∈ R d×2 , which yields I d + 2Σ 1/2 xx ⊺ Σ 1/2 + 2Σ 1/2 yy ⊺ Σ 1/2 = |I d + AB ⊺ | = |I 2 + B ⊺ A| = 1 + 2x ⊺ Σx 2x ⊺ Σy 2x ⊺ Σy 1 + 2y ⊺ Σy .This directly implies Equation (S21). Furthermore, Equation (S22) follows with Σ ∈ R 2×2 and x = ( 1 0 ), y = ( 0 1 ) or with x = y = ( 1 0 ) and an arbitrary invertible covariance matrix Σ ′ such thatΣ ′ 11 = Σ 11 . Corollary D.1. Given (X, Y ) ∼ N (0, Σ) with invertible covariance matrix Σ, it holds Ṫm (Σ) = E[ ė rf m (X) ė rf m (Y )] = 2 π Σ + I 2 /(2m 2 ) -1/2 m→∞ ----→ 2 π |Σ| -1/2 . (S23) If (X, Y ) ∼ N (0, Σ) with X = Y and singular covariance matrix Σ ∈ R 2×2 , it holds Ṫm (Σ) m→∞ ----→ ∞.Proof.Ṫm (Σ) = m 2 Ṫ (m 2 • Σ) Lemma D.1 = 4m 2 π |I 2 + 2m 2 Σ| -1/2 = 2m 2 π 4m 4 |I 2 /(2m 2 ) + Σ| -1/2 = 2 π Σ + I 2 /(2m 2 ) -1/2 m→∞ ----→ 2 π |Σ| -1/2 if Σ is invertible, ∞ if Σ is singular.Remark D.1. As mentioned at the beginning of this chapter, we can get the same result by considering the distributional derivative of the sign function, 2δ 0 . It holds for (X, Y ) ∼ N (0, Σ) with invertible covariance matrix Σ as beforeE[δ 0 (X) δ 0 (Y )] = R 2 1 2π|Σ| 1/2 2δ 0 (z 1 ) 2δ 0 (z 2 ) e -1 2 z ⊺ Σ -1 z dz = 2 π |Σ| -1/2 .In the case of X = Y , the integral is no longer well-defined.Next, we consider T by solving a more general problem, which was formulated in slightly less general form by Williams [1996, Chapter 3.1] . Lemma D.2. Given U ∼ N (0, Σ) with invertible covariance matrix Σ ∈ R d×d and x, y ∈ R d , it holdsV := E[erf(U ⊺ x) erf(U ⊺ y)] = 2 π arcsin 2 x ⊺ Σy √ 1 + 2 x ⊺ Σx √ 1 + 2 y ⊺ Σy .In particular, given (X, Y ) ∼ N (0, Σ) with invertible covariance matrix Σ = Σ1 Σ3 Σ3 Σ2 ∈ R 2×2 or with X = Y and singular covariance matrixΣ = Σ1 Σ3 Σ3 Σ2 ∈ R 2×2 , Σ 1 = Σ 2 = Σ 3 , it holds E[erf(X) erf(Y )] = 2 π arcsin 2 Σ 3 √ 1 + 2 Σ 1 √ 1 + 2 Σ 2 .Proof. We follow the proof idea given by Williams [1996, Chapter 3.1] , that is, we define V (λ), differentiate the expectation, and integrate by parts. We can then see that d dλ V (λ) = (1 -γ 2 ) -1/2 dγ dλ , which gives the desired arcsin. So, we defineV (λ) = E[erf(λ • U ⊺ x) erf(U ⊺ y)] = R d 1 (2π) d/2 |Σ| 1/2 erf(λ • u ⊺ x) erf(u ⊺ y) e -1 2 u ⊺ Σ -1 u du,and then differentiate with respect to λ on both sidesEQUATION)using a change of variables u = Σ 1/2 v in Equation (⋆) and using partial integration and Gauss' divergence theorem in the last equation. In addition, we used that∇ v e -1 2 v ⊺ (Id+2λ 2 Σ 1/2 xx ⊺ Σ 1/2 )v = -(I d + 2λ 2 Σ 1/2 xx ⊺ Σ 1/2 )v e -1 2 v ⊺ (Id+2λ 2 Σ 1/2 xx ⊺ Σ 1/2 )vand the partial integration rule for scalar functions. To be precise, for differentiable scalar functions f and g it holds∇(f • g) = f • ∇g + g • ∇f,which then implies the partial integration rule. Furthermore, the left-hand side vanishes in our case due to Gauss' divergence theorem:R d ∇ v erf(v ⊺ Σ 1/2 y) • e -1 2 v ⊺ (Id+2λ 2 Σ 1/2 xx ⊺ Σ 1/2 )v dv = lim R→∞ S d-1 (R) erf(v ⊺ Σ 1/2 y) • e -1 2 v ⊺ (Id+2λ 2 Σ 1/2 xx ⊺ Σ 1/2 )v dv ≤ lim R→∞ S d-1 (R) e -1 2 v ⊺ (Id+2λ 2 Σ 1/2 xx ⊺ Σ 1/2 )v dv ≤ lim R→∞ S d-1 (R) • e -1 2 R 2 λmin(I d +2λ 2 Σ 1/2 xx ⊺ Σ 1/2 ) = 0,where S d-1 (R) is the surface area of the sphere in R d with radius R. Continuing with our previous calculations, we see that(S24) = 2x ⊺ Σ 1/2 √ π(2π) d/2 (I d + 2λ 2 Σ 1/2 xx ⊺ Σ 1/2 ) -1 × R d e -1 2 v ⊺ (Id+2λ 2 Σ 1/2 xx ⊺ Σ 1/2 )v • Σ 1/2 y 2 √ π e -(v ⊺ Σ 1/2 y) 2 dv = 4 π x ⊺ Σ 1/2 (I d + 2λ 2 Σ 1/2 xx ⊺ Σ 1/2 ) -1 Σ 1/2 y (2π) d/2 R d e -1 2 v ⊺ (Id+2λ 2 Σ 1/2 xx ⊺ Σ 1/2 +2Σ 1/2 yy ⊺ Σ 1/2 )v dv. (S25)We evaluate the expression outside of the integral in the last line by applying the Sherman-Morrison-Woodbury formula. For A ∈ R d×d and w 1 , w 2 ∈ R d it holds [Golub and Van Loan, 1996, (2.4.1) ](C + w 1 w ⊺ 2 ) -1 = C -1 - C -1 w 1 w ⊺ 2 C -1 1 + w ⊺ 2 C -1 w 1 . For C = I d and w = w 1 = w 2 = √ 2λΣ 1/2 x this yields (I d + 2λ 2 Σ 1/2 xx ⊺ Σ 1/2 ) -1 = (I d + ww ⊺ ) -1 = I d - ww ⊺ 1 + w ⊺ w = I d - 2λ 2 Σ 1/2 xx ⊺ Σ 1/2 1 + 2λ 2 x ⊺ Σx .With this we see thatx ⊺ Σ 1/2 (I d + 2λ 2 Σ 1/2 xx ⊺ Σ 1/2 ) -1 Σ 1/2 y = x ⊺ Σ 1/2 I d - 2λ 2 Σ 1/2 xx ⊺ Σ 1/2 1 + 2λ 2 x ⊺ Σx Σ 1/2 y =x ⊺ Σy - 2λ 2 (x ⊺ Σx)(x ⊺ Σy) 1 + 2λ 2 x ⊺ Σx = x ⊺ Σy 1 - 2λ 2 x ⊺ Σx 1 + 2λ 2 x ⊺ Σx = x ⊺ Σy 1 + 2λ 2 x ⊺ ΣxInserting this into Equation (S25), we obtain(S25) = 2 π 2x ⊺ Σy 1 + 2λ 2 x ⊺ Σx R d 1 (2π) d/2 e -1 2 v ⊺ (Id+2λ 2 Σ 1/2 xx ⊺ Σ 1/2 +2Σ 1/2 yy ⊺ Σ 1/2 )v dv = 2 π 2x ⊺ Σy 1 + 2λ 2 x ⊺ Σx I d + 2λ 2 Σ 1/2 xx ⊺ Σ 1/2 + 2Σ 1/2 yy ⊺ Σ 1/2 -1 1/2 = 2 π 2x ⊺ Σy 1 + 2λ 2 x ⊺ Σx I d + 2λ 2 Σ 1/2 xx ⊺ Σ 1/2 + 2Σ 1/2 yy ⊺ Σ 1/2 -1/2 .As in the proof of Lemma D.1, we evaluate this using Sylvester's determinant theorem. With the same notation as before, we can defineA = B = ( √ 2λΣ 1/2 x, √ 2Σ 1/2 y) ∈ R d×2. This yieldsI d + 2λ 2 Σ 1/2 xx ⊺ Σ 1/2 + 2Σ 1/2 yy ⊺ Σ 1/2 = |I d + AB ⊺ | = |I 2 + B ⊺ A| = 1 + 2λ 2 x ⊺ Σx 2λx ⊺ Σy 2λx ⊺ Σy 1 + 2y ⊺ Σy = (1 + 2λ 2 x ⊺ Σx)(1 + 2y ⊺ Σy) -4λ 2 (x ⊺ Σy) 2 .If we insert this this again, we have so far shownEQUATION)We now defineγ(λ) := 2λx ⊺ Σy (1 + 2λ 2 x ⊺ Σx)(1 + 2y ⊺ Σy) ,and claim that2 π 1 -γ(λ) 2 -1/2 d dλ γ(λ) = d dλ V (λ). (S27)We can find a solution to the claimed equation by finding a function Ṽ that satisfiesd dγ(λ) Ṽ (γ(λ)) = 2 π 1 -γ(λ) 2 -1/2 ,and by setting V (λ) := Ṽ (γ(λ)). This follows from the chain rule. Ṽ is thus simply given by Ṽ (γ(λ)) = 2 π arcsin(γ(λ)). In particular, this yieldsV = V (1) = Ṽ (γ(1)) = 2 π arcsin 2x ⊺ Σy (1 + 2x ⊺ Σx)(1 + 2y ⊺ Σy) .It is now left to show Equation (S27) using Equation (S26). First, see that1 -γ(λ) 2 -1/2 = 1 - (2λx ⊺ Σy) 2 (1 + 2λx ⊺ Σx)(1 + 2y ⊺ Σy) -1/2 = (1 + 2λx ⊺ Σx)(1 + 2y ⊺ Σy) (1 + 2λx ⊺ Σx)(1 + 2y ⊺ Σy) -4λ 2 (x ⊺ Σy) 2 1/2 . Second, it holds d dλ λ(1 + 2λ 2 x ⊺ Σx) -1/2 = d dλ (λ -2 + 2x ⊺ Σx) -1/2 = λ -3 (λ -2 + 2x ⊺ Σx) -3/2 = (1 + 2λ 2 x ⊺ Σx) -3/2 = 1 √ 1 + 2λ 2 x ⊺ Σx 3 .With the results of both calculations, we get2 π 1 -γ(λ) 2 -1/2 d dλ γ(λ) = 2 π 1 -γ(λ) 2 -1/2 2x ⊺ Σy √ 1 + 2y ⊺ Σy d dλ λ(1 + 2λ 2 x ⊺ Σx) -1/2 = 2 π √ 1 + 2λx ⊺ Σx √ 1 + 2y ⊺ Σy (1 + 2λx ⊺ Σx)(1 + 2y ⊺ Σy) -4λ 2 (x ⊺ Σy) 2 2x ⊺ Σy √ 1 + 2y ⊺ Σy 1 √ 1 + 2λ 2 x ⊺ Σx 3 = 2 π 1 (1 + 2λx ⊺ Σx)(1 + 2y ⊺ Σy) -4λ 2 (x ⊺ Σy) 2 2x ⊺ Σy 1 + 2λ 2 x ⊺ Σx (S26) = d dλ V (λ),which yields Equation (S27) and concludes the proof. The special case Σ ∈ R 2×2 with invertible covariance matrix or X = Y and singular covariance matrix follows as in the proof of Lemma D.1.Corollary D.2. If (X, Y ) ∼ N (0, Σ) with invertible covariance matrix Σ = Σ1 Σ3 Σ3 Σ2 ∈ R 2×2 or with X = Y and singular covariance matrix Σ = Σ1 Σ3 Σ3 Σ2 ∈ R 2×2 , Σ 1 = Σ 2 = Σ 3 , it holds T m (Σ) = 2 π arcsin   Σ 3 1 2m 2 + Σ 1 1 2m 2 + Σ 2   m→∞ ----→ 2 π arcsin Σ 3 √ Σ 1 √ Σ 2 . (S28)Proof. It holdsT m (Σ) = T (m 2 • Σ) Lemma D.2 = 2 π arcsin 2m 2 Σ 3 √ 1 + 2m 2 Σ 1 √ 1 + 2m 2 Σ 2 = 2 π arcsin   Σ 3 1 2m 2 + Σ 1 1 2m 2 + Σ 2   m→∞ ----→ 2 π arcsin Σ 3 √ Σ 1 √ Σ 2 .We can now use Corollary D.1 and Corollary D.2 to evaluate Σ (L) , Σ(L) and Θ (L) for activation function erf m , which we denote byEQUATIONΘ (L)m respectively. We then are interested in the limit m → ∞. First recall thatEQUATIONm is independent of m, and it holds for any x, y ∈ R n0Σ (1) m (x, y) = σ 2 w n 0 ⟨x, y⟩ + σ 2 b =: Σ (1) ∞ (x, y).Assuming that the limit lim m→∞ Σ(L) m (x, y) =: Σ (L) ∞ (x, y) exists and that Σ (L) ∞ (x, x) ̸ = 0, Σ (L) ∞ (y, y) ̸ = 0, we can define Σ (L+1) m (x, y) = σ 2 w E (X,Y )∼N (0,Σ (L) m;x,y ) [erf m (X) erf m (Y )] + σ 2 b = σ 2 w T m Σ (L) m;x,y + σ 2 b Cor. D.2 = 2σ 2 w π arcsin   Σ (L) m (x, y) 1 2m 2 + Σ (L) m (x, x) 1 2m 2 + Σ (L) m (y, y)   + σ 2 b m→∞ ----→ 2σ 2 w π arcsin   Σ (L) ∞ (x, y) Σ (L) ∞ (x, x) Σ (L) ∞ (y, y)   + σ 2 b =: Σ (L+1) ∞ (x, y).(S29) Hence, it follows via induction and from the continuity of the arcsin function thatlim m→∞ Σ (L) m (x, y) =: Σ (L)∞ (x, y) is well defined. We discuss the resulting kernelΣ (L) ∞ in Re- mark D.3. For Σ(L+1) m , L ≥ 1, it holds Σ(L+1) m (x, y) = σ 2 w E (X,Y )∼N (0,Σ (L) m;x,y ) [ ė rf m (X) ė rf m (Y )] = σ 2 w Ṫm Σ (L) m;x,y Cor. D.1 = 2σ 2 w π Σ (L) m;x,y + 1 2m 2 I 2 -1/2 = 2σ 2 w π Σ (L) m (x, x) + 1 2m 2 Σ (L) m (y, y) + 1 2m 2 -Σ (L) m (x, y) 2 -1/2.(S30) As we will see later, the limitlim m→∞ Σ(L) m (x, y) =: Σ(L) ∞ (x, y) (S31)exists for x ̸ = y apart from a few exceptions. In the case x = y, we obtainEQUATIONwhere ∼ denotes asymptotic equality. Therefore, for x ̸ = y, the limit lim m→∞ Θ (L) m (x, y) =: Θ (L) ∞ (x, y) exists. However, due to Equation (S32), the NTK diverges for x = y as m → ∞. We will call a kernel with this property a singular kernel. We also say that a kernel with this property is singular along the diagonal. Remark D.2 (Distributional neural tangent kernel). An alternative conceivable approach to obtain the same kernel Θ (L) ∞ would have been to consider the distributional Jacobian matrix of the network function with step-like activation function. Whether or not a distributional NTK can be formulated is a question for further research. Here, we only want to point out that the corresponding formulas for the recursive definition of the analytical distributional NTK would then naturally read as follows,Θ (1) ∞ (x, x ′ ) = Σ (1) ∞ (x, x ′ ) Θ (L) ∞ (x, x ′ ) = Σ (L) ∞ (x, x ′ ) + Θ (L-1) ∞ (x, x ′ ) • Σ(L) ∞ (x, x ′ ) for L ≥ 2, where Σ (L) ∞ (x, x ′ ) = σ 2 w E g∼N (0,Σ (L-1) ∞ ) [sign(g(x)) sign(g(y))] + σ 2 b for L ≥ 2, (S33) Σ(L) ∞ (x, x ′ ) = σ 2 w E g∼N (0,Σ (L-1) ∞ ) [2δ 0 (g(x)) 2δ 0 (g(y))] for L ≥ 2. (S34)Note that Equation (S34) can be derived from Remark D.1, and that Equation (S33) is also easy to show.We now want to explore the implications of our findings in this section for f NTK , which we introduced at the end of Section C.2. Equation (S20) yieldslim m→∞ f m NTK (x) = lim m→∞ Θ m (x, X )Θ m (X , X ) -1 Y =: f ∞ NTK (x). Note that Θ m (X , X ) := Θ (L)m (X , X ) is invertible for sufficiently large m, as the matrix will be dominated by the diagonal for large m. Using this and the fact Θ m (x, X ) m→∞ ----→ Θ ∞ (x, X ) if x ̸ = x i for all i = 1, . . . , d, we obtain for such x:f m NTK (x) = Θ m (x, X )Θ m (X , X ) -1 Y m→∞ ----→ 0However, if x = x i ∈ X and denoting the i-th basic vector by e i ∈ R d we get:f m NTK (x i ) = Θ m (x i , X )Θ m (X , X ) -1 Y = e ⊺ i Y = y i Therefore, f mNTK converges pointwise to a function that is zero almost everywhere, but interpolates exactly at the data points. The singular kernel Θ (L) ∞ , which we can see as the analytic NTK for the sign function, is for that reason not suitable for regression. We will deal with this problem in the following section. D.2 From singular kernel to Nadaraya-Watson estimator Radhakrishnan et al. [2023] considered the limit of infinite depth, lim L→∞ Θ (L) , and a singular kernel emerged similar to the previous section. It was shown that the resulting estimator using this singular kernel behaves like a Nadaraya-Watson estimator when classification tasks are considered instead of regression tasks. We will adapt the ideas of Radhakrishnan et al. [2023] to show similar results. To consider a classification task, we assume that Y ∈ {-1, 1} d . The network is then trained with data points (X , Y) as before, but we apply the sign function at the end to obtain the final classifier. As Radhakrishnan et al. [2023] , we assume that we are operating in the infinite-width limit after infinite training. So, we are interested in limm→∞ sign Θ m (x, X )Θ m (X , X ) -1 Y . (S35)We will see that the resulting classifier is equal tosign (Θ ∞ (x, X )Y) = sign d i=1 Θ(x, x i )y i .The form of this classifier is similar to a Nadaraya-Watson estimator. To be precise, for a singular kernel k and data points (X , Y), the Nadaraya-Watson estimator is defined asc(x) := d i=1 k(x, x i )y i d i=1 k(x, x i ) for all x ̸ = x i , 1 ≤ i ≤ d.Since k(x, x i ) → ∞ as x → x i , it holds that c(x) → y i as x → x i . Thus, the continuous extension of c to R n0 interpolates the data points.We now begin to further evaluate Σ (L)∞ to analyze Equation (S35):EQUATIONΣ (L+1) ∞ (x, y) (S29) = 2σ 2 w π arcsin   Σ (L) ∞ (x, y) Σ (L) ∞ (x, x) Σ (L) ∞ (y, y)   + σ 2 b (S37) = 2σ 2 w π arcsin Σ (L) ∞ (x, y) σ 2 w + σ 2 b + σ 2 b for all L ≥ 2. (S39)Regarding the assumptions made for Equation (S29), note that that Σ(L) ∞ (x, x) ̸ = 0 if L ≥ 2 and Σ (1) ∞ (x, x) ̸ = 0 if x ̸ = 0 or σ 2 b > 0.After Equation (S31), we mentioned exceptions to the existence of Σ(L) ∞ (x, y) for x ̸ = y that were postponed at that time. For Σ(2) m with x ̸ = y we can see that Σ(2) m (x, y)EQUATIONassuming that either σ 2 b > 0 or that x, y are not parallel, i.e., |⟨x, y⟩| ̸ = ∥x∥∥y∥. The term 1 2m 2 A m is given by1 2m 2 A m = 1 2m 2 σ 2 w n 0 ∥x∥ 2 + σ 2 b + σ 2 w n 0 ∥y∥ 2 + σ 2 b + 1 2m 2 m→∞ ----→ 0.For L ≥ 3 the expression simplifies:EQUATIONwith the term 1 2m 2 B m given by1 2m 2 B m = 1 2m 2 Σ (L-1) m (x, x) + Σ (L-1) m (y, y) + 1 2m 2 m→∞ ----→ 0.From Equation (S39) and Equation (S41) we see thatd d Σ (L-1) ∞ (x, y) Σ (L) ∞ (x, y) = Σ(L) ∞ (x, y).Although we did not use dual activation functions to denote the NTK as [Jacot et al., 2018, Section A.4 ], we still find that the property (σ) ′ = (σ ′ ) applies, i.e., the derivative of the dual activation function is the dual of the derivative of the activation function. Here σ denotes the dual activation function of σ.In the case x = y we getΣ(2) m (x, x) (S32) ∼ 2σ 2 w π m σ 2 w n 0 ∥x∥ 2 + σ 2 b -1/2, and using (S37) (S42)Σ(L) m (x, x) (S32) ∼ 2σ 2 w π m σ 2 w + σ 2 b -1/2 for L ≥ 3. (S43)We summarize the above calculations in the following lemma:Lemma D.3. For m, L ∈ N let Σ (L)m and Σ(L) m be as in Theorem C.2 for activation function erf m . It then holds for any x ̸ = y with x, y ̸ = 0:Σ (1) ∞ (x, y) = lim m→∞ Σ (1) m (x, y) (S36) = σ 2 w n 0 ⟨x, y⟩ + σ 2 b , Σ (2) ∞ (x, y) = lim m→∞ Σ (2) m (x, y) (S38) = 2σ 2 w π arcsin   σ 2 w n0 ⟨x, y⟩ + σ 2 b σ 2 w n0 ∥x∥ 2 + σ 2 b σ 2 w n0 ∥y∥ 2 + σ 2 b   + σ 2 b , Σ (L) ∞ (x, y) = lim m→∞ Σ (L) m (x, y) (S39) = 2σ 2 w π arcsin Σ (L-1) ∞ (x, y) σ 2 w + σ 2 b + σ 2 b for all L ≥ 3, Σ (L) ∞ (x, x) = lim m→∞ Σ (L) m (x, x) (S37) = σ 2 w + σ 2 b for all L ≥ 2,and, assuming that x, y are not parallel or thatσ 2 b > 0, Σ(2) ∞ (x, y) = lim m→∞ Σ(2) m (x, y) (S40) = 2σ 2 w π σ 4 w n 2 0 ∥x∥ 2 ∥y∥ 2 -⟨x, y⟩ 2 + σ 2 w σ 2 b n 0 ∥x -y∥ 2 -1 2 = 2σ 2 w π Σ (1) ∞ (x, x) Σ (1) ∞ (y, y)) -Σ (1) ∞ (x, y) 2 -1 2 , Σ(L) ∞ (x, y) = lim m→∞ Σ(L) m (x, y) (S41) = 2σ 2 w π (σ 2 w + σ 2 b ) 2 -Σ (L-1) ∞ (x, y) 2 -1 2 for all L ≥ 3, Σ(2) m (x, x) (S42) ∼ 2σ 2 w π m σ 2 w n 0 ∥x∥ 2 + σ 2 b -1 2 , Σ(L) m (x, x) (S43) ∼ 2σ 2 w π m σ 2 w + σ 2 b -1 2 for L ≥ 3.Remark D.3 (Addendum to Remark C.5). In Remark C.5 we discussed the topology of the space of the Gaussian processes to which ANNs with continuous activation functions converge in the infinite-width limit. The product σ-algebra restricts us to a countable input set, so it is not possible to check for properties such as continuity or even differentiability. While Theorem C.4 is stated only for continuous activation functions with linear envelope property, we will see in Theorem E.3 that the convergence also holds in the (weak) infinite-width limit even for step-like activation functions.For the sign function as activation function, the covariance function of this Gaussian process is then given by Σ (L)∞ . Since we know explicitly what this covariance function looks like, we can examine the sample-continuity of the process. Note that to get the full picture, one would still have to show functional convergence of the network to this process, as Bracale et al. [2021] have done.∞ is isotropic when restricted to a sphere, as will be discussed in the next section. In the case of isotropic covariance functions, k(x, y) = k(∥x -y∥), the simplest way to show sample-continuity using the Kolmogorov-Chentsov criterion is to show Lipschitz continuity of k at zero. This can be seen from the proof of Lemma 4.3 by Lang and Schwab [2015] . In our case, the covariance function is basically given by a composition of arcsin functions. Lipschitz continuity of k at zero is therefore equivalent to Lipschitz continuity of the arcsin function at 1, which does not hold. In conclusion, it is not possible to show sample-continuity of the Gaussian process given by Σ (L) ∞ in the established way using the Kolmogorov-Chentsov criterion.Clearly, the kernel can be rewritten in terms of the arccos function. Cho and Saul [2009] analyzed arc-cosine kernels in the context of deep learning in detail. b > 0, then the limitΘ (L) ∞ (x, y) = lim m→∞ Θ (L) m (x, y)exists. Furthermore, it holds, asymptotically as m → ∞,Θ (L) m (x, x) ∼ 2σ 2 w π σ 2 w n 0 ∥x∥ 2 + σ 2 b 1 2 2σ 2 w π σ 2 w + σ 2 b -1 2 L-2 m L-1 for L ≥ 2.Proof. The first statement directly follows from Lemma D.3 and the definition ofΘ (L)m . The recursive definition can be resolved to the following formula:Θ (L) m (x, y) = L k=1 Σ (k) m (x, y) • L-1 l=k Σ(l+1) m (x, y). With K(z) := 2σ 2 w π σ 2 w n 0 z + σ 2 b -1 2 , we get from Lemma D.3 that Σ(2) m (x, x) ∼ K ∥x∥ 2 •m and Σ(L) m (x, x) ∼ K(n 0 )•m for L ≥ 3. This implies L-1 l=1 Σ(l+1) m (x, x) ∼ K(∥x∥ 2 )K(n 0 ) L-2 • m L-1 and L-1 l=k Σ(l+1) m (x, x) ∼ K(n 0 ) L-k • m L-k for any k ≥ 2. For L ≥ 2, we get that Θ (L) m (x, x) ∼ Σ (1) ∞ (x, y)K ∥x∥ 2 K(n 0 ) L-2 • m L-1 + L k=2 Σ k ∞ (x, x) • K(n 0 ) L-k • m L-k ∼ σ 2 w n 0 ∥x∥ 2 + σ 2 b 2σ 2 w π σ 2 w n 0 ∥x∥ 2 + σ 2 b -1 2 K(n 0 ) L-2 m L-1 = 2σ 2 w π σ 2 w n 0 ∥x∥ 2 + σ 2 b 1 2 2σ 2 w π σ 2 w + σ 2 b -1 2 L-1 m L-2 .Theorem D.4 (Inspired by Lemma 5 of Radhakrishnan et al. [2023] ). Let σ 2 b > 0 or let all x i ∈ R n0 be pairwise non-parallel. Let L ≥ 2 and x i ∈ S n0-1 R for all i = 1, . . . , d, where S n0-1 R ⊆ R n0 is the sphere of radius R. Then, withc (L) (x) := lim m→∞ sign Θ (L) m (x, X )Θ (L) m (X , X ) -1 Y ,and assuming that Θ(L) ∞ (x, X )Y ̸ = 0 for almost all x ∈ S n0-1 R , it holds c (L) (x) = sign Θ (L) ∞ (x, X )Y a.e. on S n0-1 R .Proof. First note that almost all x ∈ S n0-1 R are not parallel to any x i ∈ X . We denote Θ (L) m (X , X ) =: Θ m , m ∈ N ∪ {∞}, for convenience. Let a m > 0 be a positive constant that we will choose later. It then holds for almost allx ∈ S n0-1 R sign Θ (L) m (x, X )Θ -1 m Y = sign a m Θ (L) m (x, X )Θ -1 m Y =sign a m Θ (L) m (x, X )Θ -1 m Y -Θ (L) m (x, X )Y =:Am + Θ (L) m (x, X )Y -Θ (L) ∞ (x, X )Y =:Bm + Θ (L) ∞ (x, X )Y .We now show that A m and B m go to zero as m → ∞ for a suitable choice of a m . First, note that|B m | ≤ Θ (L) m (x, X ) -Θ (L) ∞ (x, X ) 2 ∥Y∥ 2 m→∞ ----→ 0, since ∥Y∥ 2 < ∞ and by Corollary D.3 it holds that Θ (L) m (x, x i ) → Θ (L) ∞ (x, x i ) for all x i ∈ X as m → ∞. Second, we have |A m | ≤ a m Θ (L) m (x, X )Θ -1 m -Θ (L) m (x, X ) 2 ∥Y∥ 2 ≤ Θ (L) m (x, X ) 2 a m Θ -1 m -I d 2 ∥Y∥ 2 . (S44)Again by Corollary D.3 it holds that ∥Θ(L) m (x, X )∥ 2 → ∥Θ (L) ∞ (x, X )∥ 2 as m → ∞. Furthermore, for all 1 ≤ i ≤ d we have Θ (L) m (x i , x i ) ∼ C(R) • m L-1 for some constant C(R) which depends on R = ∥x i ∥. Since it holds Θ (L) m (x j , x i ) → Θ (L) ∞ (x j , x i ) as m → ∞ for any i ̸ = j, we choose a m = C(R) • m L-1 and conclude a -1 m Θ (L) m (x i , x j ) m→∞ ----→ 1 if i = j 0 if i ̸ = j. This implies a -1 m Θ m → I d as m → ∞. In particular, {I d } ∪ {a -1 m Θ m | m ∈ N} is a compact set with respect to ∥ • ∥ 2 . Since D → D -1 is a continuous function, {I d } ∪ {a m Θ -1 m | m ∈ N} is bounded. We get a m Θ -1 m -I d 2 ≤ a m Θ -1 2 I d -a -1 m Θ m 2 m→∞ ----→ 0, and thus |A m | (S44) ≤ Θ (L) m (x, X ) 2 a m Θ -1 m -I d 2 ∥Y∥ 2 m→∞ ----→ 0.Recall that Θ (L)∞ (x, X )Y ̸ = 0, which now concludes the proof:lim m→∞ sign Θ (L) m (x, X )Θ -1 m Y = lim m→∞ sign A m + B m + Θ (L) ∞ (x, X )Y = sign Θ (L) ∞ (x, X )Y .Remark D.4. One can generalize the above theorem by dropping the restriction toS n0-1 R . By Lemma D.3 it holds Θ (L) m (x i , x i ) ∼ σ 2 w n0 ∥x i ∥ 2 + σ 2 b • C • m L-1 . For some constant C independent of L and ∥x i ∥. If we now define a m := C • m L-1 , we get a -1 m Θ (L) m (X , X ) m→∞ ----→ diag σ 2 w ∥x i ∥ 2 /n 0 + σ 2 b 1≤i≤d =: D,where diag{v} of a vector v is a square matrix with diagonal v. As in the proof above, and using the continuity of the inverse map D → D -1 , this implieslim m→∞ sign Θ (L) m (x, X )Θ -1 m Y = sign Θ (L) ∞ (x, X )D -1 Y = sign Θ(L) ∞ (x, X )Y , if we define Θ(L) ∞ (x, y) := σ 2 w ∥x∥ 2 /n 0 + σ 2 b -1/2 σ 2 w ∥y∥ 2 /n 0 + σ 2 b -1/2 Θ (L) ∞ (x, y).In conclusion, dropping the restriction leads to a similar form for our classifier, but with a modified kernel.lim m→∞ sign Θ (L) m (x, X )Θ -1 m Y = lim m→∞ sign (f m NTK (x)) has an easily interpretable form that is close to the Nadaraya-Watson estimator with singular kernel Θ (L) ∞ . This is despite the fact that the pointwise limit of f m NTK is trivial, i.e., regression is no longer possibleIn addition to the ideas used in the above section, Radhakrishnan et al. [2023] proved that the resulting estimator is Bayes optimal under certain conditions. This property is also referred to as consistency. Let the probability distribution P data on our data space be given by a random variable(X, Y ) ∈ R n0 × {-1, 1} and let c(x) = arg max ỹ∈{-1,1} P(Y = ỹ | X = x)be the Bayes classifier with respect to this distribution. Denote X d = (x 1 , . . . , x d ) and Y d = (y 1 , . . . , y d ) for (x i , y i ) drawn independently from P data . We then define Bayes optimality as follows:Definition D.1 (Bayes optimality). Let c d ( • ) = c d ( • ; X d , Y) be classifiers for d ∈ N and estimators of c( • ). We then say that (c d ) d∈N is Bayes optimal, if it is a consistent estimator of c, i.e., for all ε > 0 and X-almost all x limd→∞ P |c d (x) -c(x)| > ε = 0.First, we summarize the results of Radhakrishnan et al. [2023] and consider Θ (L) . If the singular limiting kernel for L → ∞ behaves like a singular kernel of the form k(x, y) = ∥x -y∥ -α , then the classifier for L → ∞ will be of the formsign d i=1 y i /∥x -x i ∥ α d i=1 1/∥x -x i ∥ α , assuming d i=1 1/∥x -x i ∥ α > 0.This classifier satisfies Bayes optimality for α = n 0 by Devroye et al. [1998] . Radhakrishnan et al. [2023] generalized the results of Devroye et al. [1998] , expressed α in terms of the activation function and its derivative, and chose them in such a way as to achieve α = n 0 . Going back to our setup, we want to see if we can writelim m→∞ Θ (L) m (x, y) = R(∥x -y∥) ∥x -y∥ α ,for some constant α and for some function R : R + → R bounded away from zero as ∥x -y∥→ 0.We start by proving that Θ (L) m (x, y) = G(∥x -y∥) for some function G. Recall that we have to restrict ourselves to a sphere by Theorem D.4. For simplicity, we restrict ourselves to the unit sphere S n0-1 . Then, it holds ⟨x, x⟩ = 1 for all x ∈ S n0-1 , which gives us∥x -y∥ 2 = ⟨x -y, x -y⟩ = ⟨x, x⟩ -2⟨x, y⟩ + ⟨y, y⟩ = 2(1 -⟨x, y⟩).In the following, we will substitute z = ⟨x, y⟩. Taking ∥x -y∥→ 0 is then equivalent to taking z → 1. We can conclude from Lemma D.3 that we can indeed writeΣ (L) ∞ (x, y) = Σ (L) ∞ (z), Σ(L) ∞ (x, y) = Σ(L) ∞ (z), and hence Θ (L) ∞ (x, y) = Θ (L) ∞ (z) for all L ≥ 1. Note that Σ (L) ∞ (1) = σ 2 w +σ 2 b for L ≥ 2 and Σ (L) ∞ (1) = σ 2 w /n 0 + σ 2 b . Next, we consider Σ(2) ∞ (z) for z → 1 using Lemma D.3 Σ(2) ∞ (z) = 2σ 2 w π σ 4 w n 2 0 (1 -z 2 ) + σ 2 w σ 2 b n 0 • 2(1 -z) -1 2 = 2σ 2 w π σ 4 w n 2 0 (1 -z)(1 + z) + 2σ 2 w σ 2 b n 0 (1 -z) -1 2 z→1 ∼ 2σ 2 w π 2σ 4 w n 2 0 (1 -z) + 2σ 2 w σ 2 b n 0 (1 -z) -1 2 = 2σ 2 w π 2σ 2 w n 0 σ 2 w n 0 + σ 2 b -1 2 (1 -z) -1 2 = n 0 π 2σ 2 w n 0 1 2 σ 2 w n 0 + σ 2 b -1 2 (1 -z) -1 2 . (S45)For L ≥ 3, we can observe thatEQUATION)We will now prove a lemma that allows us to analyze the behavior of Σ(L) ∞ (z) as z → 1. Lemma D.5. It holds:lim z→1 (1 -z) 1/2 1 -2 π arcsin (z) = π 2 √ 2 .Proof. Since the numerator and the denominator both go to zero as z → 1, we apply l'Hôpital's rule and differentiate both. This yieldsd dz (1 -z) 1 2 = - 1 2 (1 -z) -1 2 and d dz 1 - 2 π arcsin (z) = - 2 π 1 -z 2 -1 2 .Thus,lim z→1 (1 -z) 1/2 1 -2 π arcsin (z) = lim z→1 1 2 √ 1 -z √ 1 + z 2 π √ 1 -z = π 4 √ 2 = π 2 √ 2 ,concluding the proof.For L = 3 it now holdsΣ(3) ∞ (z) (S46) ∼ 2σ 2 w π 2 σ 2 w + σ 2 b -1 2 σ 2 w + σ 2 b -Σ (2) ∞ (z) -1 2 = 2σ 2 w π 2 σ 2 w + σ 2 b -1 2 σ -1 w 1 - 2 π arcsin σ 2 w z/n 0 + σ 2 b σ 2 w /n 0 + σ 2 b -1 2 (⋆) ∼ 2σ 2 w π 2 σ 2 w + σ 2 b -1 2 σ -1 w π 2 √ 2 1 - σ 2 w z/n 0 + σ 2 b σ 2 w /n 0 + σ 2 b -1 4 = 2σ 2 w π 2 σ 2 w + σ 2 b -1 2 σ -1 w π 2 √ 2 σ 2 w /n 0 • (1 -z) σ 2 w /n 0 + σ 2 b -1 4 = 2σ 2 w π 2 √ 2 σ 2 w + σ 2 b -1 2 σ 2 w /n 0 σ 2 w /n 0 + σ 2 b -1 4 (1 -z) -1 4 ,where we used Lemma D.5 at (⋆). For L ≥ 4 it holdsΣ(L) ∞ (z) (S46) ∼ 2σ 2 w π 2 σ 2 w + σ 2 b -1 2 σ 2 w + σ 2 b -Σ (L-1) ∞ (z) -1 2 L-1≥3 = 2σ 2 w π 2 σ 2 w + σ 2 b -1 2 σ -1 w 1 - 2 π arcsin Σ (L-2) ∞ (z) σ 2 w + σ 2 b -1 2 (⋆) ∼ 2σ 2 w π 2 σ 2 w + σ 2 b -1 2 σ -1 w π 2 √ 2 1 - Σ (L-2) ∞ (z) σ 2 w + σ 2 b -1 4 = 2σ 2 w π 2 σ 2 w + σ 2 b -1 2 σ -1 w π 2 √ 2 σ 2 w + σ 2 b 1 4 σ 2 w + σ 2 b -Σ (L-2) ∞ (z) -1 4 (S46) ∼ 2σ 2 w π 2 σ 2 w + σ 2 b -1 4 σ -1 w π 2 √ 2 σ 2 w + σ 2 b 1 4 Σ(L-1) ∞ (z) = Σ(L-1) ∞ (z)/2,again using Lemma D.5 at (⋆). This implies for arbitraryL ≥ 2 Σ(L) ∞ (z) ∼ K(L) • (1 -z) -1/2 L-1 ,for some constant K(L) depending on L. Recall the formula for the NTK, that was used in the proof of Corollary D.3Θ (L) ∞ (z) = L k=1 Σ (k) ∞ (z) • L-1 l=k Σ(l+1) ∞ (z) ∼ L k=1 Σ (k) ∞ (z) • L-1 l=k K(l + 1) • (1 -z) -1/2 l ∼ Σ (1) ∞ (1) L-1 l=1 K(l + 1) • (1 -z) -1/2 l = K(L) • (1 -z) -(1-1/2 L-1 ) = K ′ (L) • ∥x -y∥ -(2-1/2 L-2 ) =: K ′ (L) • ∥x -y∥ -α(L) ,for some constants K(L), K ′ (L). It is therefore to possible find a function R that is bounded away from zero near ∥x -y∥ → 0, such that L) .Θ (L) ∞ (x, y) = K ′ (L) • R(∥x -y∥) ∥x -y∥ α(However, we have that α(1) = 1 and α(L) ↑ 2 as L → ∞. So it is only possible to choose L such that α(L) = n 0 , if n 0 = 1. But this is a trivial case, since we have restricted ourselves to the unit sphere. In conclusion, we cannot prove that c (L) is a Bayes optimal classifier for any choice of L. Chapter 4 of Devroye et al. [1998] suggests that, in fact, the estimator will not be universally consistent.In this chapter we explore surrogate gradient learning, introduced in Section 1, by connecting it to the NTK. Recall that the standard gradient flow dynamics of the parameters are given byd dt θ t = -η ∇ θ L(f (X ; θ t ); Y) = -η J θ f (X ; θ t ) ⊺ ∇ f (X ;θt) L(f (X ; θ t ); Y). (S8)As mentioned several times before, the Jacobian matrix J θ will vanish for all parameters except those in the last layer if we consider the sign activation function due to its zero derivative. The idea of surrogate gradient learning is to circumvent the zero derivative of the sign function by replacing the derivative with a surrogate derivative [Neftci et al., 2019] . We can replace the derivative in two main ways:• The activation function in the full network can be replaced by a differentiable surrogate activation function. In particular, we obtain a non-vanishing surrogate derivative of the activation function and thus non-vanishing network gradients. Let g denote the network with surrogate activation function. Therefore, we can train the weights according tod dt θ t = -η J θ g(X ; θ t ) ⊺ ∇ g(X ;θt) L(g(X ; θ t ); Y),or consider the loss with respect to f and train according toEQUATION• Instead of replacing the activation function, we can only replace the derivative of the activation function σ with a surrogate derivative σ in Equation (S8). Let J σ,σ be the quasi-Jacobian matrix as in Definition C.4 with activation function σ and surrogate derivative σ.Then the training is given byEQUATIONwith Î(L) as in Definition C.5 for σ 1 = σ, σ1 = σ, σ 2 = σ and σ2 = σ. For Equations (S48) and (S49) we assume that σ exists and is non-vanishing, but we deliberately train with a surrogate gradient. For example, we can again consider erf m as the activation function.Its derivative explodes at zero and vanishes everywhere else as m → ∞. The hope is thatlim n1,...n L-1 →∞ Î(L) m = I (L)m exists andlim m→∞ I (L)m is not a singular kernel as before.In this chapter we will deal with the second approach, because J σ,σ (x; θ t ) =: G(σ; σ; x; θ t ) is closer to J θ f (x; θ t ) = G(σ; σ; x; θ t ) than J θ g(x; θ t ) = G(η; σ; x; θ t ) as a formula if J θ f (x; θ t ) exists and is non-vanishing. Here, η denotes the surrogate activation function with derivative σ. To do this, we provide asymmetric generalizations of Theorem C.1, Theorem C.2, and Theorem C.3. These generalizations would, in principle, even allow us to compare the two approaches.For this section we adopt the so-called linear envelope property from [Matthews et al., 2018, Definition 1] to ensure that all expectations exist in the following theorems:∃ m, c ≥ 0 ∀ u ∈ R : |σ(u)| ≤ c + m|u| (S50)First, we consider networks under the weak infinite-width limit, and are thus interested in taking the number of hidden neurons to infinity sequentially. This is done inductively while using the central limit theorem and the weak law of large numbers. In order to do this rigorously, we state and prove two lemmata.The first lemma is stated in terms of Gaussian measures on Hilbert spaces. An introduction to Gaussian measures on Hilbert spaces can be found in Chapter 1 of Da Prato [2006] . A rigorous derivation and definition of convergence in distribution on arbitrary metric spaces is given by Heyer [2009, Chapter 1.2] . This includes a definition of weak convergence in terms of continuous and bounded functions (Remark 1.2.5 (b)), a version of the Portemanteau theorem (Theorem 1.2.7), and the fact that convergence in probability implies convergence in distribution (Application 1.2.15). Lemma E.1. Let H m i and Z i , i, m ∈ N, be random variables with values in a separable Hilbert space H. Furthermore, let Z i be independent and identically distributed with finite mean and covariance operatorEQUATIONfor an H-valued Gaussian random variable Z with mean and covariance operator like Z 1 . Similarly, there existsk ′ : N → N such that k ′ (m) → ∞ monotonically as m → ∞ and 1 k ′ (m) k ′ (m) i=1 H m i D ----→ m→∞ E[Z]. (S52)Proof. Since H is separable and complete, convergence in distribution and convergence with respect to the Prokhorov metric d (also known as the Lévy-Prokhorov metric) are equivalent [Billingsley, 1999, Theorem 6.8 ]. By the central limit theorem for separable Hilbert spaces [Zalesskiȋ et al., 1991] , this implieslim k→∞ d 1 √ k k i=1 Z i , Z = 0. (S53)In addition, the assumption together with the continuous mapping theorem gives thatlim m→∞ d 1 √ k k i=1 H m i , 1 √ k k i=1 Z i = 0 for all k ∈ N.In particular, for any k ∈ N, there exists some m k ∈ N such thatd 1 √ k k i=1 H m i , 1 √ k k i=1 Z i ≤ 1 k for all m ≥ m k . (S54)We now want to choose k(m) as large as possible for any m, but small enough to ensure Inequality (S54), i.e., m ≥ m k(m) . So we definek(m) := sup{k | m ≥ m k }. First note that {k | m ≥ m k } ̸ = ∅, if m ≥ m 1 .The map k : N → N is therefore well-defined, as we consider m → ∞. Similarly, we can find a m ≥ m k for any given k. This yieldslim m→∞ k(m) = lim m→∞ sup{k | m ≥ m k } = ∞. By definition of k(m), it holds m ≥ m k(m) for all m ∈ N and thus d   1 k(m) k(m) i=1 H m i , 1 k(m) k(m) i=1 Z i   ≤ 1 k(m) for all m ∈ N. (S55)Together with Equation (S53) this yields the claim, (S51):d   1 k(m) k(m) i=1 H m i , Z   ≤ d   1 k(m) k(m) i=1 H m i , 1 k(m) k(m) i=1 Z i   + d   1 k(m) k(m) i=1 Z i , Z   (S55) ≤ 1 k(m) + d   1 k(m) k(m) i=1 Z i , Z   m→∞ ----→ 0.For the second claim, (S52), one can follow the same procedure but use the law of large numbers for Banach spaces instead of the central limit theorem. Suitable results are given by Ledoux and Talagrand [1991, Corollary 7 .10] and Hoffmann-Jørgensen and Pisier [1976, Theorem 2.1] . Note that even the strong law of large numbers holds, but the weak law is sufficient.In the second lemma, some properties about convergence in distribution and convergence in probability are stated. Lemma E.2 (Theorem 2.7 from van der Vaart [1998] , modified and (iv) added). Let X n , X and Y n be random vectors. Then(i) X n P -→ X implies X n D -→ X; (ii) X n P -→ c for a constant c if and only if X n D -→ c; (iii) if X n D -→ X and Y n P -→ c for a constant c, then (X n , Y n ) D -→ (X, c); (iv) if X n D -→ X and W is a random vector independent of (X n ) n∈N , then (X n , W ) D -→ (X, W ).Proof. (i) -(iii). The proofs are given by van der Vaart [1998] .(iv). Let f be a bounded and continuous function. Then,lim n→∞ E[f (X n , W )] = lim n→∞ E [f (X n , W ) | W ] (x) dP(w) (⋆) = lim n→∞ E [f (X n , W (w))] dP(w) = lim n→∞ E [f (X n , W (w))] dP(w) = E [f (X, W (w))] dP(w) = E [f (X, W ) | W ] (w) dP(w) = E[f (X, W )],where we used the independence of W and (X n ) n∈N in Equation (⋆) and the boundedness of f for the interchange of limit and integration. This proves convergence in distribution.Remark E.1. The above theorem can be generalized to metric spaces. One can easily check that the proofs in [van der Vaart, 1998, Theorem 2.7 ] also work for metric spaces using the Portemanteau theorem provided by Heyer [2009, Theorem 1.2.7 ]. However, it is necessary to derive some more equivalent characterizations of convergence in distribution, which are given and used by van der Vaart [1998] but are missing in the work of Heyer [2009] . Theorem E.3 (Generalized version of Proposition 1 by Jacot et al. [2018] ). For activation functions σ 1 and σ 2 with property (S50), which are continuous except for finitely many jump points, let f 1 ( • ; θ) and f 2 ( • ; θ) be network functions with hidden layers hEQUATION2 ( • ; θ), for 1 ≤ l ≤ L, respectively as in Definition C.1 and with shared weights θ.Then (f 1 ( • ; θ), f 2 ( • ; θ)) converges in distribution to a multidimensional Gaussian process (X (L) j , Y (L) j ) j=1,...,n L as (n l ) 1≤l≤L-1 → ∞ weakly for any fixed countable input set (z i ) ∞ i=1 .The Gaussian process is defined byX (L) j iid ∼ N 0, Σ (L) 1 , Y (L) j iid ∼ N 0, Σ (L) 2, where we have forEQUATION)where(Z 1 , Z 2 ) ∼ N 0, Σ (L-1) 1 (x,x) Σ (L-1) 1,2 (x,x ′ ) Σ (L-1) 1,2 (x,x ′ ) Σ (L-1) 2 (x ′ ,x ′ ).Proof. We write h (l) (x) = h (l) (x; θ). We prove the theorem by induction, as in the proof of Proposition 1 of Jacot et al. [2018] , but expand on the technical details.L = 1. By definition, we haveh (1) 1 (x) = h (1) 2 (x) = σ w √ n 0 W (1) x + σ b b (1) .This implies that h(1) k2,j (x ′ ; θ), the i-th and j-th component of h(1) k1 (x; θ) and h(1)EQUATION2 ( • ; θ)) is a Gaussian process, it is thus sufficient to show that vectors of the form ( h (1)1,i (x1) ... h (1) 1,i (xn) h (1) 2,i (x1) ... h (1) 2,i (xn) )⊺ have a multivariate Gaussian distribution for anyEQUATION1,i (x 1 ) . . . hEQUATION, with 1 d a column vector of ones with length d. Now sinceW (1) i • ⊺ b (1) i ∼ N (0, I n0+1 ),it holds( h (1) 1,i (x1) ... h (1) 1,i (xn) h (1) 2,i (x1) ... h (1) 2,i (xn) ) ⊺ ∼ N 0, σwX ⊺ /n0 σ b 1 d σwX ⊺ /n0 σ b 1 d σwX ⊺ /n0 σ b 1 d σwX ⊺ /n0 σ b 1 d ⊺ .Therefore, the vector has a multivariate Gaussian distribution with the covariances required for Equation (S56) and the first case of Equations (S58).L → L + 1. We assume that the convergence holds for depth L. This means that there exists some r ∈ R L such that, for given constant width n 0 , any width n L , and widthsn l = r l (m) , 1 ≤ l < L, it holds h (L) 1 (•), h (L) 2 (•) D ----→ m→∞ X (L) j , X (L) j j=1,...,n L .To be precise, there is no such r in the case L = 1 → L + 1. However, this only makes the proof simpler and one can still follow the same steps as for L ≥ 2.By the continuous mapping theorem [Billingsley, 1999, Theorem 2.7 ] it holdsσ 1 h (L) 1 (•) , σ 2 h (L) 2 (•) D ----→ m→∞ σ 1 X (L) j , σ 2 Y (L) j j=1,...,n L . (S59)The theorem is applicable despite the finitely many jump points, since (X (L) , Y (L) ) assumes the values of the jump points with zero probability.We now need to find an increasing width function r L : N → N such that, if we additionally set n L = r L (m), it holds for any fixed n L+1h (L+1) 1 (•), h (L+1) 2 (•) D ----→ m→∞ X (L) j , X (L) j j=1,...,n L+1.Note that by Remark C.5 we consider the product σ-algebra. Therefore, to show convergence in distribution for the whole process, it is sufficient to show convergence in distribution for the marginal distributions. By definition, we have for k ∈ {1, 2} thath (L+1) k (x) = σ w √ n L W (L+1) σ k h (L) k (x) + σ b b (L+1) = σ w √ n L n L i=1 σ k h (L) k,i (x) W (L+1) • i + σ b b (L) .The marginal vector for points x 1 , . . . , x n as before can thus be written as           h (L+1) 1 (x 1 ) . . . h (L+1) 1 (x n ) h (L+1) 2 (x 1 ) . . . h (L+1) 2 (x n )            = σ w √ n L n L i=1            σ 1 h (L) 1,i (x 1 ) I n L+1 . . . σ 1 h (L) 1,i (x n ) I n L+1 σ 2 h (L) 2,i (x 1 ) I n L+1 . . . σ 2 h (L) 2,i (x n ) I n L+1            W (L+1) • i + σ b           b (L+1) . . . b (L+1) b (L+1) . . . b (L+1)          .With the same arguments as before and using the continuous mapping theorem in combination with Lemma E.2 (iv) and the independence ofW (L+1) , it holds                       σ 1 h (L) 1,i (x 1 ) I n L+1 . . . σ 1 h (L) 1,i (x n ) I n L+1 σ 2 h (L) 2,i (x 1 ) I n L+1 . . . σ 2 h (L) 2,i (x n ) I n L+1            W (L+1) • i            n L i=1 D ----→ m→∞                       σ 1 X (L) i (x 1 ) I n L+1 . . . σ 1 X (L) i (x n ) I n L+1 σ 2 Y (L) i (x 1 ) I n L+1 . . . σ 2 Y (L) i (x n ) I n L+1            W (L+1) • i            n L i=1 .Since (X i , Y i ) and (X j , Y j ) are independent for i ̸ = j, the conditions of Lemma E.1 are satisfied. Now, there exists k : N → N such that k(m) → ∞ monotonically as m → ∞ and, when settingn L := r L (m) := k(m), it holds            h (L+1) 1 (x 1 ) . . . h (L+1) 1 (x n ) h (L+1) 2 (x 1 ) . . . h (L+1) 2 (x n )            D ----→ m→∞ σ w            R (L+1) x1 . . . R (L+1) xn S (L+1) x1 . . . S (L+1) xn            + σ b           b (L+1) . . . b (L+1) b (L+1) . . . b (L+1)           (⋆) =           X (L+1) (x 1 ) . . . X (L+1) (x n ) Y (L+1) (x 1 ) . . . Y (L+1) (x n )           , for a Gaussian random variable R (L+1) x1 . . . , R (L+1) xn , S (L+1) x1 , . . . , S (L+1) xn ∈ R 2•n L+1 •d with            R (L+1) x1 . . . R (L+1) xn S (L+1) x1 . . . S (L+1) xn            ∼            σ 1 X (L) 1 (x 1 ) I n L+1 . . . σ 1 X (L) 1 (x n ) I n L+1 σ 2 Y (L) 1 (x 1 ) I n L+1 . . . σ 2 Y (L) 1 (x n ) I n L+1            W (L+1) • 1.Before considering the covariances, we want to comment on r L . First, note that this sequence may not initially be strictly increasing, but this can be circumvented by considering a strictly increasing subsequence. Second, r L could theoretically depend on X . However, this can be resolved by not evaluating the pair (h(L) 1 , h (L)2 ) at certain data points, but doing the same calculation as above for (h(L) 1 , h (L)2 ). The starting point for this is Equation (S59). To apply Lemma E.1, note additionally that (h(L) 1 , h (L) 2 ) ∈ ∞ i=1 R 2•n L ,which is a separable Hilbert space because we are considering a countable input set. Above, we worked with the marginal distribution because this makes the following calculation of the covariances easier. Finally, the choice of r L should be independent of n L+1 . This follows from the independence of W (L+1) ji and W (L+1) j ′ i for j ̸ = j ′ .To verify Equation (⋆), we check the covariances of the random vector. They are given byCov R (L+1) xi , R (L+1) xj = E σ 1 X (L) 1 (x i ) W (L+1) • 1 W (L+1) • 1 ⊺ σ 1 X (L) 1 (x j ) =E σ 1 X (L) 1 (x i ) E W (L+1) • 1 W (L+1) • 1 ⊺ X (L) 1 (x i ), X (L) 1 (x j ) σ 1 X (L) 1 (x j ) =E σ 1 (X (L) 1 (x i )) I n L+1 σ 1 X (L) 1 (x j ) = E σ 1 X (L) 1 (x i ) σ 1 X (L) 1 (x j ) I n L+1 .Similarly, we get thatCov S (L+1) xi , S (L+1) xj = E σ 2 Y (L) 1 (x i ) σ 2 Y (L) 1 (x j ) I n L+1 ,and together this implies using the independence of biases b (L+1) and weight matrices W (L+1)Cov σ w R (L+1) xi,k + σ b b (L+1) k , σ w R (L+1) xj ,l + σ b b (L+1) l =δ kl σ 2 w E σ 1 X (L) 1 (x i ) σ 1 X (L) 1 (x j ) + σ 2 b = Cov X (L+1) k (x i ), X (L+1) l (x j ) , Cov σ w S (L+1) xi,k + σ b b (L+1) k , σ w S (L+1) xj ,l + σ b b (L+1) l =δ kl σ 2 w E σ 2 Y (L) 1 (x i ) σ 2 Y (L) 1 (x j ) + σ 2 b = Cov Y (L+1) k (x i ), Y (L+1) l (x j ) .We therefore proved Equation (S57). For the second case of (S58), we see thatCov R (L+1) xi , S (L+1) xj = E σ 1 X (L) 1 (x i ) W (L+1) • 1 W (L+1) • 1 ⊺ σ 2 Y (L) 1 (x j ) =E σ 1 X (L) 1 (x i ) E W (L+1) • 1 W (L+1) • 1 ⊺ X (L) 1 (x i ), Y (L) 1 (x j ) σ 2 Y (L) 1 (x j ) =E σ 1 X (L) 1 (x i ) σ 2 Y (L) 1 (x j ) I n L+1 ,with, by induction hypothesis,X (L) 1 (x i ), Y (L) 1 (x j ) ∼ N 0, Σ (L) 1 (x i , x i ) Σ (L) 1,2 (x i , x j ) Σ (L) 1,2 (x i , x j ) Σ (L) 2 (x j , x j ). This finished the proof, since it now holdsCov σ w R (L+1) xi,k + σ b b (L+1) k , σ w S (L+1) xj ,l + σ b b (L+1) l =δ kl σ 2 w E σ 1 X (L) 1 (x i ) σ 2 Y (L) 1 (x j ) + σ 2 b = Cov X (L+1) k (x i ), Y (L+1) l (x j ) ,Remark E.2. In the preceding proof, we checked the marginal distributions of arbitrary size in order to give a complete proof of the convergence to a Gaussian process. However, the covariances can be derived by considering only a pair of data points (x 1 , x 2 ), which drastically simplifies the notation. Also, we only need the distributions of pairs for the next theorems.Theorem E.4 (Generalized version of Theorem 1 by Jacot et al. [2018] ). For activation functions σ 1 , σ 2 and so-called surrogate derivatives σ1 , σ2 such that σ 1 , σ 2 , σ1 , and σ2 are continuous except for finitely many jump points with property (S50), let f 1 ( • ; θ) and f 2 ( • ; θ) be network functions with hidden layers hEQUATION2 ( • ; θ), 1 ≤ l ≤ L, respectively as in Definition C.1 with shared weights θ. Denote the empirical generalized neural tangent kernel Î(L) (x, x ′ ) = J (L),σ1,σ1 (x; θ) J (L),σ2,σ2 (x ′ ; θ) ⊺ for x, x ′ ∈ R n0 , as in Definition C.5. Then, for any x, x ′ ∈ R n0 and 1 ≤ i, j ≤ n L , it holdsÎ(L) ij (x, x ′ ) P -→ δ ij I (L) (x, x ′ ),as n 1 , . . . , n L-1 → ∞ weakly. We call I (L) the analytic generalized neural tangent kernel, which is recursively given byI (1) (x, x ′ ) = Σ (1)1,2 (x, x ′ ) (S60)I (L) (x, x ′ ) = Σ (L) 1,2 (x, x ′ ) + I (L-1) (x, x ′ ) • Σ(L) 1,2 (x, x ′ ) for L ≥ 2, (S61) with Σ (L)1,2 as in Theorem E.3 andΣ(L) 1,2 (x, x ′ ) = σ 2 w E[σ 1 (Z 1 ) σ2 (Z 2 )] for L ≥ 2, where (Z 1 , Z 2 ) ∼ N 0, Σ (L-1) 1 (x,x) Σ (L-1) 1,2 (x,x ′ ) Σ (L-1) 1,2 (x,x ′ ) Σ (L-1) 2 (x ′ ,x ′ ).Proof. We prove the theorem by induction over L, as in the proof of Theorem 1 by Jacot et al. [2018] . We denoteJ (L),k (z) = J (L),σ k ,σ k (z; θ) for k ∈ {1, 2}, z ∈ R n0 . L = 1. For 1 ≤ i, j ≤ n 1 it holds Î(1) ij (x, x ′ ) = J (1),1 i • (x) J (1),2 j • (x ′ ) ⊺ = θ ′ ∈θ (1)(1),1 i θ ′ (x) J(1),2j θ ′ (x ′ ) (S12) = 1≤k≤n1 1≤l≤n0 δ ki σ w √ n 0 x l δ kj σ w √ n 0 x ′ l + 1≤k≤n1 δ ki σ b δ kj σ b = σ 2 w n 0 δ ij 1≤k≤n0 x l x ′ l + σ 2 b δ ij = δ ij σ 2 w n 0 ⟨x, x ′ ⟩ + σ 2 b = δ ij I (1) (x, x ′ ).This proves Equation (S60).L → L + 1. Now we assume that the statement is true for L and need to prove it for L + 1. Instead of considering an explicit r ∈ R L as in the proof of Theorem E.3 and taking m → ∞, we will just write "n 1 , . . . , n L-1 → ∞ weakly".In the induction step, we would like to use Theorem E.3. However, it is not obvious that this is possible. In the setting of Definition C.7, let S 1 and S 2 be two statements that hold weakly. In our setting, these are the induction hypothesis and Theorem E.3 for depth L. Then there exist s, t ∈ R L such that S 1 (s) and S 2 (t) are true. It is not clear that there exists some r ∈ R L such that S 1 (r) and S 2 (r) are true. It would be natural to define r by r l (m) := max{s l (m), t l (m)} for all 1 ≤ l < L and m ∈ N, but it is still unclear that this implies S 1 (r) and S 2 (r). Instead, one can consider the combined statement "S 1 and S 2 " as a new statement S. In our case, for any depth L, we would need to find a r ∈ R L such that the statements in Theorem E.3 and Theorem E.4 are both true, which can be done. Since we used the first part of Lemma E.1 to prove Theorem E.3 and will use the second part of Lemma E.1 for this proof, we would have to define r L (m) := min{k(m), k ′ (m)}. For simplicity, we will assume that the r ∈ R L we get by the induction hypothesis also gives us the convergence statement of Theorem E.3.Using the definition of the generalized NTK and the quasi-Jacobian matrices, we obtain for1 ≤ i, j ≤ n L+1 I (L+1) ij (x, x ′ ) = J (L+1),1 i • (x) J (L+1),2 j • (x ′ ) ⊺ = θ ′ ∈θ (1 : L+1) J (1),1 i θ ′ (x) J (1),2 j θ ′ (x ′ ) = 1≤k≤n L+1 1≤l≤n L δ ki σ w √ n L σ 1 h (L) 1,l (x) δ kj σ w √ n L σ 2 h (L) 2,l (x ′ ) + 1≤k≤n L+1 δ ki σ b δ kj σ b + θ ′ ∈θ (1 : L) σ 2 w n L n L m=1 W (L+1) i,m σ1 h (L) 1,m (x) J (L),1 m θ ′ (x) n L r=1 W (L+1) i,r σ2 h (L) 2,r (x ′ ) J (L),2 r θ ′ (x ′ ) =δ ij σ 2 w n L n L l=1 σ 1 h (L) 1,l (x) σ 2 h (L) 2,l (x ′ ) + σ 2 b (S62) + σ 2 w n L n L m,r=1 W (L+1) i,m W (L+1) j,r σ1 h (L) 1,m (x) σ2 h (L) 2,r (x ′ ) • θ ′ ∈θ (1 : L) J (L),1 m θ ′ (x)J (L),2 r θ ′ (x ′ ) (S63)We want to apply the second part of Lemma E.1. We will consider the terms (S62) and (S63) separately. First note thatEQUATIONas n 1 , . . . , n L-1 → ∞ weakly by Theorem E.3 and the continuous mapping theorem withX (L) l (x), Y (L) l (x ′ ) iid ∼ N 0, Σ (L) 1,1 (x, x) Σ (L) 1,2 (x, x ′ ) Σ (L) 1,2 (x, x ′ ) Σ (L) 2,2 (x ′ , x ′ ). Again, we used that the values of the jump points are assumed with zero probability. Thus, again by the continuous mapping theorem and by the second part of Lemma E.1, it holdsσ 2 w n L n L l=1 σ 1 h (L) 1,l (x) σ 2 h (L) 2,l (x ′ ) + σ 2 b D --------→ n1,...,n L →∞ weakly E σ 1 X (L) 1 (x) σ 2 Y (L) 1 (x ′ ) + σ 2 b .Here, as in the proof of Theorem E.3, n L → ∞ is given by n L := r L (m) := k ′ (m), which is in turn is given by Lemma E.1. Note also that the limit is a constant, which implies convergence in probability according to Lemma E.2 (ii). In conclusion, we obtainδ ij σ 2 w n L n L l=1 σ 1 h (L) 1,l (x) σ 2 h (L) 2,l (x ′ ) + σ 2 b P --------→ n1,...,n L →∞ weakly δ ij Σ (L+1) 1,2 (x, x ′ ). (S65)For term (S63), we can apply the induction hypothesis to obtainθ ′ ∈θ (1 : L) J (L),1 θ ′ ,m (x)J (L),2 θ ′ ,r (y) = Î(L) m,r (x, y) P ----------→ n1,...,n L-1 →∞ weakly δ mr I (L) (x, y). (S66)Also, we can again use the convergence given by (S64), but with σ 1 , σ 2 replaced by σ1 , σ2 respectively. This gives using Lemma E.2 (iv) L+1) .EQUATIONTogether with (S66) and Lemma E.2 (iii) this implies L+1) , δ mr I (L) (x, x ′ ) .EQUATIONFor the summands in (S63) we therefore have by the continuous mapping theoremn L r=1 W (L+1) i,m W (L+1) j,r σ1 h (L) 1,m (x) σ2 h (L) 2,r (x ′ ) • θ ′ ∈θ (1 : L) J (L),1 θ ′ ,m (x)J (L),2 θ ′ ,r (x ′ ) D -------→ n1,...,n L-1 weakly n L r=1 W (L+1) i,m W (L+1) j,r σ1 X (L) m (x) σ2 Y (L) r (x ′ ) • δ mr I (L) (x, x ′ ) =W (L+1) i,m W (L+1) j,m σ1 X (L) m (x) σ2 Y (L) m (x ′ ) • I (L) (x, x ′ )These terms are independent for different 1 ≤ m ≤ n L . Therefore, the second part of Lemma E.1 can be applied as before. This yieldsσ 2 w n L n L r=1 W (L+1) i,m W (L+1) j,r σ1 h (L) 1,m (x) σ2 h (L) 2,r (x ′ ) • θ ′ ∈θ (1 : L) J (L),1 θ ′ ,m (x)J (L),2 θ ′ ,r (x ′ ) D --------→ n1,...,n L →∞ weakly I (L) (x, x ′ ) • σ 2 w E W (L+1) i,1 W (L+1) j,1 σ1 X L) 1 (x) σ2 Y (L) 1 (x ′ ) =I (L) (x, x ′ ) • δ ij σ 2 w E σ1 X L) 1 (x) σ2 Y (L) 1 (x ′ ) = δ ij I (L) (x, x ′ ) • Σ(L) 1,2 (x, x ′). Together with (S65) this yields Equation (S61) and concludes the proof.The two theorems above are complemented by the convergence of the generalized NTK in the infinitewidth limit during training, which we will prove below. The theorem is a generalization of Theorem G.2 from Lee et al. [2019] . There, the convergence of the NTK at initialization as (n l ) L-1 l=1 ∝ ∼ n is assumed. More precisely, they refer to a result of Yang [2019a] . In consequence, we have to assume the same for the next theorem. The proof of Yang [2019a] could also be generalizable to our case. Alternatively, one could also try to generalize the statement of Jacot et al. [2018] on stability during training in the weak infinite-width limit. Theorem E.5 (Based on Theorem G.2 from Lee et al. [2019] ). Let σ be a Lipschitz continuous and differentiable activation function. Let the derivative of the activation function σ and a socalled surrogate derivative σ be Lipschitz continuous and bounded. Let f t ( • ; θ) be the corresponding network function with depth L initialized as in Definition C.1 and trained with MSE loss and surrogate gradient learning, i.e., according to Equation (S49) with surrogate derivative σ. The hidden layers are denoted by h (n l ) L-1 l=1 ∝ ∼ n, J (L),σ, σ J (L),σ,σ ⊺ = Î(L) P ----→ n→∞ I (L) ⊗ I n L .Furthermore, assume that the smallest and largest eigenvalue of the symmetrization of I (L) (X , X ),S (L) := 1 2 I (L) (X , X ) + I (L) (X , X ) ⊺ ,are given by 0 < λ min ≤ λ max < ∞ and that the learning rate is given by η > 0. Then, for any δ > 0 there exist R > 0, N ∈ N and K > 1 such that for every n ≥ N , the following holds with probability at least 1 -δ over random initializationsup t∈[0,∞) Î(L) (X , X ) -I (L) (X , X ) F ≤ 6K 3 R λ min n -1 2 ,where ∥ • ∥ F denotes the Frobenius norm.Proof. We follow the proofs in Chapter G of Lee et al. [2019] . The analogous statement is given by Theorem G.2 of Lee et al. [2019] . Since we are considering the infinite-width limit (n l ) L-1 l=1 ∝ ∼ n, we will assume that the width of the hidden layers are given by n ∈ N. The more general case can easily be proved using the same steps. Recall thatEQUATIONwith probability at least 1 -δ 1 using Lemma E.6 as before. With the same calculations as for Inequality (S71), this implies for the symmetrizations thatS (L) 0 -S (L) t 2 ≤ Î(L) 0 (X , X ) - Î(L) t (X , X ) 2 ≤ 6(K ′ ) 2 KR 0 √ n ≤ 1 3 λ min , (S73) for all n ≥ N := max m 0 , m 1 , λmin 18(K ′ ) 2 KR0 2 . Now Inequalities (S70), (S71) and (S73) give us z, S (L) t z = z, S (L) + S (L) 0 -S (L) + S (L) t -S (L) 0 z = z, S (L) z + z, S (L) 0 -S (L) z + z, S (L) t -S (L) 0 z (S70) ≥ λ min ∥z∥ 2 -z, S (L) 0 -S (L) z -z, S (L) t -S (L) 0 z (S71)+(S73) ≥ λ min ∥z∥ 2 - 1 3 λ min ∥z∥ 2 - 1 3 λ min ∥z∥ 2 = 1 3 z, 1 3 λ min I z ,and thus imply Inequality (S70). To summarize, it holds with probability at least 1 -δ 2 -δ 3 for any n ≥ N and θ t ∈ B (θ 0 , C) that∥g(θ 0 )∥ 2 2 (S68) ≤ R 2 0 and d dt ∥g(θ t )∥ 2 2 (S69) ≤ - 2 3 ηλ min ∥g(θ t )∥ 2 2 .Grönwall's inequality now implies∥g(θ t )∥ 2 2 ≤ e -2 3 ηλmint ∥g(θ 0 )∥ 2 2 ≤ e -2 3 ηλmint R 2 0 .We can now return to Inequality (S67) to obtain with probability at least 1-δ 1 -δ 2 -δ 3 d dt ∥θ t -θ 0 ∥ 2 (S67) ≤ ηK∥g(θ t )∥ 2 ≤ ηKR 0 e -1 3 ηλmint .Integrating the inequality on both sides yields for all θ t ∈ B(θ 0 , C)∥θ t -θ 0 ∥ 2 ≤ 3KR 0 λ min 1 -e -1 3 ηλmin t < 3KR 0 λ min = C. (S74) Let t 1 := inf {t : ∥θ t -θ 0 ∥ 2 < C}. Now, if t 1 < ∞, it holds C = lim t↑t1 ∥θ t -θ 0 ∥ 2 (S74) < C.This is a contradiction, so we conclude that t 1 = ∞. In particular, Inequality (S74) holds for all t > 0. We can now repeat the calculations for Inequality (S72) with the Frobenius norm to finally obtainÎ(L) 0 (X , X ) - Î(L) t (X , X ) F ≤ 2K 2 √ n ∥θ t -θ 0 ∥ 2 ≤ 6K 3 R 0 λ min n -1 2 .Therefore, if we choose δ 1 = δ 2 = δ 3 = 1 3 δ, the desired inequality holds for any n ≥ N with probability at least 1 -δ.Lemma E.6 (Based on Lemma 1 and Lemma 2 from Lee et al. [2019] ). In the setting of Theorem E.5, for any δ > 0 there exists a K > 0 such that for any C > 0, with probability at least 1 -δ it holds J (L),σ,σ (X ; θ) -J (L),σ,σ (X ; θ)F ≤ K √ n θ -θ 2 and J (L),σ,σ (X ; θ) F ≤ K,for all θ, θ ∈ B (θ 0 , C) and σ ∈ { σ, σ}. Due to the equivalence of matrix norms, the same inequalities hold for the norm ∥ • ∥ 2 and some constant K ′ > 0.Proof. For σ = σ and a different but equivalent parameterization of the network parameters [Lee et al., 2019, Chapter F] , the proof can be found in Section G.2 of Lee et al. [2019] . The surrogate derivative σ is bounded and Lipschitz continuous. Also, J (L),σ,σ (X ; θ) and J (L),σ, σ (X ; θ) = J θ f (X ; θ) share the same recursive formula. Therefore, the proof, which builds on the so-called Gaussian conditioning technique [Yang, 2019a, Section E.1], should also hold for the surrogate derivative.Remark E.3 (The analytic generalized NTK for surrogate gradient learning). Since in Theorem E.5 we consider only one activation function σ instead of two different ones σ 1 , σ 2 as in Theorem E.3 and E.4, it holdsΣ (L) 1,2 = Σ (L) 1 = Σ (L) 2 = Σ (L) .The analytic generalized NTK in the setting of Theorem E.4 is thus given byI (1) = Σ (1) and I (L+1) = Σ (L+1) + I (L) • Σ(L+1) 1,2.We therefore still obtain an asymmetric kernel due to the contribution of Σ(L) 1,2 . Remark E.4 (Positive definiteness of the generalized NTK). In Theorem E.5 we require that the matrix I (L) (X , X ) be positive definite. Equivalently, its symmetrization,S (L) = 1 2 I (L) (X , X ) + I (L) (X , X ) ⊺ ,should be positive definite. For applications this can be checked numerically. More generally, it would be interesting to know whether the symmetric kernel given byS (L) (x, x ′ ) = 1 2 I (L) (x, x ′ ) + I (L) (x ′ , x)is positive definite. One could try to approach this question inductively. However, this leads to problems in the induction step. We want to present a different ansatz, which reduces the question to a question about Σ(L) 1,2 . We use the closed form of the analytic NTK and the symmetry of Σ (l) for all l ∈ N to see thatS (L) (x, x ′ ) = 1 2 I (L) (x, x ′ ) + I (L) (x ′ , x) = 1 2 L k=1 Σ (k) (x, x ′ ) • L-1 l=k Σ(l+1) 1,2 (x, x ′ ) + L k=1 Σ (k) (x ′ , x) • L-1 l=k Σ(l+1) 1,2 (x ′ , x) = L k=1 Σ (k) (x, x ′ ) • 1 2 L-1 l=k Σ(l+1) 1,2 (x, x ′ ) + L-1 l=k Σ(l+1) 1,2 (x ′ , x) .This defines a symmetric positive definite kernel if Σ (L) is positive definite (see Jacot et al. [2018, Section A.4 ] for comparison) and the symmetrized kernels,1 2 L-1 l=k Σ(l+1) 1,2 (x, x ′ ) + L-1 l=k Σ(l+1) 1,2 (x ′ , x) ,are positive semi-definite for all k = 1, . . . , L -1.We denote e 1 = ( 1 0 ), e 2 = ( 0 1 ) and UEQUATIONusing Lemma D.1 in Equation (⋆). For the penultimate equality we used Equation (S36) and Equation (S37). Compared to Equation (S41),Σ(L) m (x, y) m→∞ ----→ 2σ 2 w π Σ (L-1) ∞ (x, x) • Σ (L-1) ∞ (y, y) -Σ (L-1) ∞ (x, y) 2 = 2σ 2 w π Σ (L-1) ∞;x,y -1 2 , (S41)which in fact holds for both L = 2 and L ≥ 3, an additional term appeared in Equation (S76). It always holds σ 2 w + σ 2 b > 0 and it holds σ 2 w ∥x∥ 2 /n 0 + σ 2 b > 0 if σ 2 b > 0 or x ̸ = 0. As discussed earlier, we always assume that x ̸ = 0 is satisfied. It follows that this asymmetric NTK is not singular, sinceΣ(L) 1,2;∞ (x, x) = 2σ 2 w π Σ (L-1) ∞;x,x + 1 2 Σ (L-1) ∞ (x, x) -1 2 = 2σ 2 w π 0 + 1 2 Σ (L-1) ∞ (x, x) -1 2 = √ 2 2σ 2 w π Σ (L-1) ∞ (x, x) -1 2 ∈ R. (S77)Note that this is reminiscent of the constant factor depending on x in the asymptotics of Σ(L) m (x, x) as m → ∞, given by (S42) and (S43),Σ(L) m (x, x) ∼ 2σ 2 w π m Σ (L-1) ∞ (x, x) -1 2 .According to Remark E.3 it holdsI (L) m (x, y) = Σ (L) m (x, y) + I (L-1) m (x, y) • Σ1,2;m (x, y). Since Σ (L)m (x, y) and Σ(L) 1,2;m are continuous functions of the entries of the matrix Σ (L-1) m;x,y , it follows by induction using Equations (S75) and (S76) that the limit of the analytic NTK is well defined for m → ∞. We can write I (L) ∞ (x, y) = lim m→∞ I (L) m (x, y). Thus, by approximating the sign function with error functions, we found that the analytic NTK for surrogate gradient learning with sign function and error function as surrogate derivative is welldefined and non-singular as a kernel on R n0 × R n0 . Furthermore, comparing this NTK with the NTK we derived in Section D, the term Σ (L) ∞ does not change. This is a direct consequence of not replacing the activation function with a surrogate activation function. In this sense, we inherit more properties with this approach than by replacing not only the derivative but the entire activation function including its derivative with a surrogate. Comparing Equation (S41) with Equation (S76), we see that Σ(L) 1,2;∞ (x, y) can be obtained from Σ(L) ∞ (x, y) by adding a regularizing term depending on x, Σ (L-1) ∞ (x, x)/2.Now we turn to the general case with a general surrogate derivative σ. Similar to before, for m ∈ N ∪ {∞} we denote(Z m 1 , Z m 2 ) ∼ N 0, Σ (L-1) 1;m (x,x) Σ (L-1)1,2;m (x,y) Σ (L-1) 1,2;m (x,y) Σ EQUATIONwhere Equation (⋆) seems natural, but requires further reasoning. We will prove the above rigorously in Lemma E.7. For now, we assume that the equality holds. Since σ is bounded and continuous, this yields This agrees with the fact thatΣ(L) 1,2;∞ (x, x) = lim m→∞ σ 2 w E ė rf m (Z m 1 ) σ(Z m 1 ) = σ 2 w E [2δ 0 (Z ∞ 1 ) σ(Z ∞ 1 )] (S79) =2σ 2 w 1 √ 2π Σ (L-1) ∞ (x, x) -1 2 σ(0) = σ 2 w 2 π Σ (L-1) ∞ (x, x) -1 2 σ(0),where we again assumed an equality very similar to Equation (⋆). It is easy to check that Equations (S76) and (S77) can be recovered by inserting σ = ė rf into the derived formulas.We now prove the missing part of Equation (S78).Lemma E.7. Let σ be a bounded and continuous function and ((Z m 1 , Z m 2 )) m∈N a sequence of random variables, (Z m 1 , Z m 2 ) ∼ N (0, Σ m ). If the covariance matrices are invertible and converge to an invertible matrix, Σ m → Σ ∞ ∈ R 2×2 as m → ∞, then it holdsEQUATION)If Z m 1 = Z m 2 for all m ∈ N so that the covariance matrices are given by Σm = Σ m 1 Σ m 1 Σ m 1 Σ m 1 , and if Σ m 1 → Σ ∞ 1 ̸ = 0 as m → ∞, then it holds lim m→∞ E ė rf m (Z m 1 ) σ(Z m 2 ) = 2 π (Σ ∞ 1 )-1 2 σ(0). (S81)Proof. We begin with the case of invertible covariance matrices. Again denoting u = ( z1 z2 ), it holds by assumptionEQUATION)forB m := 1 |Σ m | Σ m 2,2 /m 2 -Σ m 1,2 /m -Σ m 1,2 /m Σ m 1,1 + 2 0 0 m→∞ ----→ 2 0 0 Σ ∞ 1,1 /|Σ ∞ | =: B ∞ . (S83)The determinant of B m is given by|B m | = Σ m 2,2 m 2 |Σ m | + 2 Σ m 1,1 |Σ m | - -Σ m 1,2 2 m 2 |Σ m | 2 = 1 m 2 |Σ m | + 2 Σ m 1,1 |Σ m | .This now yields(S82) = 2 √ π |Σ m | -1 2 (B m ) -1 1 2 R 2 1 2π (B m ) -1 -1 2 • σ(z 2 ) • e -1 2 u ⊺ B m u du = 2 √ π |Σ m | -1 2 |B m | -1 2 E (Y m 1 ,Y m 2 )∼N 0,(B m ) -1 [σ(Y m 2 )].The continuity of matrix inversion implies B m -1 → B ∞ -1 as m → ∞. The characteristic functions of finite-dimensional Gaussian random variables are fully defined by their means and covariance matrices. Thus, the convergence of the covariance matrices implies convergence of the characteristic functions, which in turn implies convergence in distribution. Since σ is continuous and bounded and the determinant is continuous, we obtain2 √ π |Σ m | -1 2 |B m | -1 2 E (Y m 1 ,Y m 2 )∼N 0,(B m ) -1 [σ(Y m 2 )] m→∞ ----→ 2 √ π |Σ ∞ | -1 2 |B ∞ | -1 2 E (Y ∞ 1 ,Y ∞ 2 )∼N 0,(B ∞ ) -1 [σ(Y ∞ 2 )] (S83) = 2 √ π |Σ ∞ | -1 2 2Σ ∞ 1,1 |Σ ∞ | -1 2 E Y ∼N 0,|Σ ∞ |/Σ ∞ 1,1 [σ(Y )] = 2 π Σ ∞ 1,1 -1 2 E Y ∼N 0,|Σ ∞ |/Σ ∞ 1,1 [σ(Y )].This proves Equation (S80). In the case Z m 1 = Z m 2 , we haveE ė rf m (Z m 1 ) σ(Z m 2 ) = E Y ∼N (0,Σ m 1 ) [ ė rf m (Y ) σ(Y )] = R 1 √ 2π 1 Σ m 1 2 √ π m • e -m 2 y 2 • σ(y) • e -1 2 y 2 Σ m 1 dy = R 1 √ 2π 1 Σ m 1 2 √ π m • σ(y) • e -1 2 y 2 (1/Σ m 1 +2m 2 ) dy = 2 √ π 1 Σ m 1 m 1/Σ m 1 + 2m 2 R 1 √ 2π (1/Σ m 1 + 2m 2 ) 1 2 • σ(y) • e -1 2 y 2 (1/Σ m 1 +2m 2 ) dy = 2 √ π 1 1/m 2 + 2Σ m 1 E Y ∼N (0,(1/Σ m 1 +2m 2 ) -1 ) [σ(Y )] m→∞ ----→ 2 √ π 1 2Σ ∞ 1 σ(0) = 2 π (Σ ∞ 1 ) -1 2 σ(0),again using the boundedness and continuity of σ in the last line. This proves Equation (S81).With Σ m = Σ In conclusion, we found that the analytic NTK for surrogate gradient learning with sign function is well-defined and non-singular for any bounded and Lipschitz continuous surrogate derivative. As in the special case of the error function, the term Σ(L) 1,2;∞ (x, y) can be expressed in terms of Σ Guidelines:• The answer NA means that the paper does not include experiments.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions). • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors).• It should be clear whether the error bar is the standard deviation or the standard error of the mean. • It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have a 96% CI, if the hypothesis of Normality of errors is not verified. • For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative error rates). • If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text. 


Decision: The paper does a thorough analysis of the NTK regime in in the presence of surrogate gradients. The technique seem general enough that they might be relevant to other models. Most reviewers are excited about the result and consider the paper to be worth of a spotlight.

Reviews: review: paper_summary: paper_summary: The paper adapts the neural tangent kernel framework to surrogate gradient learning (and so to learning in spiking neural networks).
strength_weakness: strength_weakness: The main weakness is the NTK approach itself: NTK is known to be a poor approximation to learning in deep networks in many cases, which limits what kind of conclusions we can make with NTK-driven theoretical analysis. However, this weakness is not specific to this paper.
questions: questions: - I don’t see how definition 2.2 iii can work for strictly increasing width functions. Sequential limits mean some of the widths stay fixed, right? 
- Fig.1 should have the black line at the back I think, otherwise it’s not clear if the peaks always overlap
limitations: limitations: The authors have addressed the limitations of the work.
paper_summary: paper_summary: This paper explore the neural tangent kernel (NTK) with regards to surrogate gradient learning for non-differentiable activation functions. The authors show that the standard neural tangent kernel is not equipped to deal with such activation functions and causes the kernel function to become singular. They provide a derivation of a surrogate derivative based NTK that can handle such activations thus further generalizing the NTK for a larger subset of networks.
strength_weakness: strength_weakness: While the paper is well developed overall, there are a few spots that can benefit from additional clarity regarding notation:

Line 144 (Definition 2.1) - $r_l(m)$ refers to the number of neurons for a particular layer however $m$ is not directly described. My understanding is that it is referring to the output size dimension but this is not immediately clear.

Line 215 - I am assuming that $\delta(z)$ refers to the Dirac delta function (delta distribution). 

Line 244 - Here we have $\delta_{ij}$ which I assume is *not* the Diract delta and instead a constant for a given kernel matrix entry.

Line 289 - Since we are talking about divergence, notationally it would be better to use $\Theta_{\text{erf}_m} \not\rightarrow \Theta\_{\text{sign}}$ or specify the formal definition of divergence *or* avoid it all together.

In regard to the experiments line 303 mentions that the NTK diverge as $m \rightarrow \infty$, however this does not seem to be clear from the plots especially since $m=20$ does not seem to be sufficient to show the divergence occurring. In addition, since we are discussing divergence, I believe the paper could benefit by including the average error values between the analytic and empirical kernels since graphs tend to do a poor job of illustrating this.
questions: questions: Are $\delta(z)$ and $\delta_{ij}$ the same?
limitations: limitations: Everything is adequately addressed.
paper_summary: paper_summary: The paper considers an extension of neural tangent kernel methods for the analysis of training of neural networks with non-differentiable activations with surrogate gradient learning.  The basic approach is the define a generalized NTK (the SG-NTK) based on the quasi-Jacobian matrix (that is, the Jacobian constructed using the surrogate gradients) rather than the (ill-defined) Jacobian based on the actual gradients of the activations.  This construct is shown to be deterministic (distributionally determined) in the infinite width limit.  For the sign activation (with erfm acting as surrogate gradient) both the NTK and SG-NTK are derived.  Finally some experimental results are given that show that the distribution of networks trained with SG-NTK matches reality.
strength_weakness: strength_weakness: - it would perhaps be more useful to present figure 1 with a log y scale as the linear scale used, combined with the divergence at $\Delta \alpha = 0$, tends to flatten all features of the NTK, masking the differences between the NTKs away from the point of divergence.
- also regarding figure 1 is it feasible to plot the limiting case $m \to \infty$ with the point of divergence elided?
questions: questions: See previous sections.
limitations: limitations: The authors have addressed limitations adequately and there appear to be no obvious negative societal impacts.
paper_summary: paper_summary: The paper addresses the challenge of applying gradient-based training methods to neural networks with non-differentiable activation functions, such as binary and spiking neural networks. These networks use surrogate derivatives to enable gradient descent, but this approach lacks theoretical foundation. The authors propose a generalization of the NTK, called the surrogate gradient NTK (SG-NTK), to provide a rigorous theoretical basis for SGL.
strength_weakness: strength_weakness: 1. The complexity of the theoretical framework may limit accessibility for practitioners
2. Additional experiments on real datasets (MNIST, CIFAR10, ImageNet) could further validate the approach.
questions: questions: What is the application of the results? The activation function, the ReLU function is commonly used and I believe it can be solved in previous works. Are there any commonly used activation functions that can solved in your paper but not in previous works?
limitations: limitations: The results in this paper are only applicable to sufficiently wide networks with random initialization.

Comments: comment: Thank you for taking the time to consider the other reviews and rebuttals. We are glad that we were able to answer your questions and we appreciate the score improvement.
Thank you for your quick response, we are glad that we were able to answer your questions.
Thank you for your quick response. We agree that the continuity of $f$ does not guarantee the exchangeability of limits. A small note to avoid misunderstandings: We are interested in finding a parameterization that corresponds to a particular sequential limit, not in changing the order of that sequential limit.
Given the enthusiasm for this paper from the reviewers, the question
of accepting this paper as a spotlight/oral comes up. To make this
decision, it is particularly important to consider the novelty and
general applicability of this paper. It would be very helpful if the
reviewers could comment whether they believe this paper is
novel/important enough to get a spotlight/oral, and why.
Thank you for your rebuttal as well as your clarifications and additions to your manuscript. I read through the other reviews and your rebuttals as well. Given the additional changes and the clarity of your presentation, I am comfortable to bump up my score to an 8.
Thank you for the response, I'm happy to keep my score of 8.
Thank you for the answer! This response addressed my questions, so I'm keeping the score of 8. (The only note is that $f(n_1, n_2)$ would need to have additional conditions to guarantee exchangeability of limits, but that's a small point.)