Year: 2019

Abstract:
A recent line of work studies overparametrized neural networks in the “kernel regime,” i.e. when the network behaves during training as a kernelized linear predictor, and thus training with gradient descent has the effect of finding the minimum RKHS norm solution. This stands in contrast to other studies which demonstrate how gradient descent on overparametrized multilayer networks can induce rich implicit biases that are not RKHS norms. Building on an observation by Chizat and Bach (2018), we show how the scale of the initialization controls the transition between the “kernel” (aka lazy) and “rich” (aka active) regimes and affects generalization properties in multilayer homogeneous models. We provide a complete and detailed analysis for a simple two-layer model that already exhibits an interesting and meaningful transition between the kernel and rich regimes, and we demonstrate the transition for more complex matrix factorization models and multilayer non-linear networks.

heading: 1 INTRODUCTION
text: A string of recent papers study neural networks trained with gradient descent in the “kernel regime.” They observe that, in a certain regime, networks trained with gradient descent behave as kernel methods (Jacot et al., 2018; Daniely et al., 2016; Daniely, 2017). This allows one to prove convergence to zero error solutions in overparametrized settings (Du et al., 2018; 2019; Allen-Zhu et al., 2018), and also implies gradient descent will converge to the minimum norm solution (in the corresponding RKHS) (Chizat and Bach, 2018; Arora et al., 2019b; Mei et al., 2019) and more generally that models inherit the inductive bias and generalization behavior of the RKHS. This suggests that, in a certain regime, deep models can be equivalently replaced by kernel methods with the “right” kernel, and deep learning boils down to a kernel method with a fixed kernel determined by the architecture and initialization, and thus it can only learn problems learnable by some kernel.
This contrasts with other recent results that show how in deep models, including infinitely overparametrized networks, training with gradient descent induces an inductive bias that cannot be represented as an RKHS norm. For example, analytic and/or empirical results suggest that gradient descent on deep linear convolutional networks implicitly biases toward minimizing the Lp bridge penalty, for p = 2/depth ≤ 1, in the frequency domain (Gunasekar et al., 2018b); weight decay on an infinite width single input ReLU implicitly biases towards minimizing the second order total variations ∫ |f ′′(x)| dx of the learned function (Savarese et al., 2019); and gradient descent on a overparametrized matrix factorization, which can be thought of as a two layer linear network, induces nuclear norm minimization of the learned matrix (Gunasekar et al., 2017) and can ensure low rank matrix recovery (Li et al., 2018). All these natural inductive biases (Lp bridge penalty for p < 1, total variation norm, nuclear norm) are not Hilbert norms, and therefore cannot be captured by a kernel. This suggests that training deep models with gradient descent can behave very differently from kernel methods, and have much richer inductive biases.
One might then ask whether the kernel approximation indeed captures the behavior of deep learning in a relevant and interesting regime, or does the success of deep learning come when learning escapes this regime? In order to understand this, we must first carefully understand when each of these regimes hold, and how the transition between the “kernel” regime and the “rich” regime happens.
Some investigations of the kernel regime emphasized the number of parameters (“width”) going to infinity as leading to this regime. However Chizat and Bach (2018) identified the scale of the model
as a quantity controlling entry into the kernel regime. Their results suggest that for any number of parameters (any width), a model can be approximated by a kernel when its scale at initialization goes to infinity (see details in Section 3). Considering models with increasing (or infinite) width, the relevant regime (kernel or rich) is determined by how the scaling at initialization behaves as the width goes to infinity. In this paper we elaborate and expand of this view, carefully studying how the scale of initialization effects the model behaviour for D-homogeneous models.
In Section 4 we provide a complete and detailed study for a simple 2-homogeneous model that can be viewed as linear regression with squared parametrization, or as a “diagonal” linear neural network. For this model we can exactly characterize the implicit bias of training with gradient descent, as a function of the scale α of initialization, and see how this implicit bias becomes the `2 norm in the α → ∞ kernel regime, but the `1 norm in the α → 0 rich regime. We can therefore understand how, e.g. for a high dimensional sparse regression problem, where it is necessary to discover the relevant features, we can get good generalization when the initialization scale α is small, but not when α is large. Deeper networks corresponds to higher orders of homogeneity, and so in Section 5 we extend our study to a D-homogeneous model, studying the effects of D. In Sections 6 and 7, we demonstrate similar transitions experimentally in matrix factorization and non-linear networks.

heading: 2 SETUP AND PRELIMINARIES
text: We consider models f : Rp × X → R which map parameters w ∈ Rp and examples x ∈ X to predictions f(w,x) ∈ R. We denote the predictor implemented by the parameters w as hw = F (w) such that hw(x) = f(w,x). Much of our focus will on models, such a linear networks, which are linear in x (but not on the parameters w!), in which case F (w) ∈ X ∗ is a linear predictor and can be represented as a vector βw with f(w,x) = 〈βw, x〉. Such models are essentially alternate parametrizations of linear models, but as we shall see that change of parametrization is crucial.
In this paper, we consider models that are D-positive homogeneous in the parameters w, for some integerD ≥ 1, meaning that for any c ∈ R+, F (c·w) = cDF (w) and f(c·w,x) = cDf(w,x). We refer to such models simply as D-homogeneous. Many interesting model classes have this property, including multi-layer ReLU networks with fully connected and convolutional layers, layered linear neural networks, and matrix factorization where D corresponds to the depth of the network.
Consider a training set {(x(n), y(n))}Nn=1 consisting of N examples of input label pairs. For a given loss function ` : R × R → R, the loss of the model parametrized by w is L(w) = L(F (w)) = ∑N n=1 `(f(w,x
(n)), y(n)). We will focus mostly on the squared loss `sq(ŷ, y) = (ŷ − y)2. We slightly abuse notation and use f(w, X) ∈ RN to denote the vector of predictions [f(w,x(1)), . . . , f(w,x(N))] and so for the squared loss we can write L(w) = ‖f(w, X)− y‖22, where y ∈ RN is the vector of target labels. Minimizing the loss L(w) using gradient descent amounts to iteratively updating the parameters
w(k + 1) = w(k)− η∇L(w(k)). (1)
We consider gradient descent with infinitesimally small stepsize η, i.e. gradient flow dynamics
ẇ(t) = −∇L(w(t)). (2)
We are particularly interested in the scale of initialization and we capture it through a scalar parameter α ∈ R+. For scale α, we will denote by wα(t) the gradient flow path (2) with the initial condition wα(0) = αw0 for some fixed w0. We use hα(t) = F (wα(t)), and for linear predictors βα(t) = βwα(t), to denote the dynamics on the predictor F (w) induced by the gradient flow on w.
In many cases, we expect the dynamics to converge to a minimizer of L(w), though proving this happens will not be our main focus. Rather, we are interested in the underdetermined case, N p, where there are generally many minimizers of L(w), all with f(w, X) = y and L(w) = 0. Our main focus is which of the many minimizers does gradient flow converge to. That is, we want to characterize wα(∞) = limt→∞wα(t) or, more importantly, the predictor hα(∞) = F (wα(∞)) or βα(∞) = βwα(∞) we converge to, and how these depend on the scale α. In underdetermined problems, where there are many zero error solutions, simply fitting the data using the model does not provide enough inductive bias to ensure generalization. But in many cases, the specific solution
reached by gradient flow (or some other optimization procedure) has special structure, or minimizes some implicit regularizer, and this structure or regularizer provides the needed inductive bias (Gunasekar et al., 2018b;a; Soudry et al., 2018; Ji and Telgarsky, 2018).

heading: 3 THE KERNEL REGIME
text: Gradient descent/flow considers only the first-order approximation of the model w.r.t. w:
f(w,x) = f(w(t), x) + 〈w −w(t), ∇wf(w(t),x)〉+O(‖w −w(t)‖2). (3) That is, locally around any w(t), gradient flow operates on the model as if it were an affine model f(w,x) ≈ f0(x) + 〈 w, φw(t)(x) 〉 with feature map φw(t)(x) = ∇wf(w(t),x), corresponding to the tangent kernel Kw(t)(x, x′) = 〈∇wf(w(t),x), ∇wf(w(t),x)〉 (Jacot et al., 2018; Zou et al., 2018; Yang, 2019; Lee et al., 2019). Of particular interest is the tangent kernel at initialization, Kwα(0) = α 2(D−1)K0 where we denote K0 = Kw0 .
The “kernel regime” refers to a situation in which the tangent kernel Kw(t) does not change over the course of optimization, and less formally to the regime where it does not change significantly, i.e. where ∀tKw(t) ≈ Kw(0). In this regime, training the model is exactly equivalent to training an affine model f̃(w,x) = αDf(w0,x) + 〈 w, αD−1φ0(x) 〉 with kernelized gradient descent/flow with the kernel α2(D−1)K0 and a “bias term” of αDf(w0,x). To avoid handling this bias term, and in particular its scaling, Chizat and Bach (2018) suggest using “unbiased” initializations such that F (w0) = 0, so that the bias term vanishes. This can often be achieved by replicating units or components with opposite signs at initialization, which is the approach we use here (see Sections 4–6 for examples and details).
For underdetermined problem with multiple solutions f(w, X) = y, unbiased kernel gradient flow (or gradient descent) converges to the minimum norm solution ĥK = arg minh(X)=y ‖h‖K , where ‖h‖K is the RKHS norm corresponding to the kernel. And so, in the kernel regime, we will have that h(∞) = ĥK0 , and the implicit bias of training is precisely given by the kernel. When does the “kernel regime” happen? Chizat and Bach (2018) showed that for any homogeneous1 model satisfying some technical conditions, the kernel regime is reached as α→∞. That is, as we increase the scale of initialization, the dynamics converge to the kernel gradient flow dynamics with the kernel K0, and we have limα→∞ hα(∞) = ĥK . In Sections 4 and 5 we prove this limit directly for our specific models, and we also demonstrate it empirically for matrix factorization and deep networks in Sections 6 and 7.
In contrast, and as we shall see in later sections, the α→ 0 small initialization limit often leads to a very different and rich inductive bias, e.g. inducing sparsity or low-rank structure (Gunasekar et al., 2017; Li et al., 2018; Gunasekar et al., 2018b), that allows for generalization in many settings where kernel methods would not. We refer to this limit reached as α→ 0 as the “rich regime.” This regime is also referred to as the “active” or “adaptive” regime (Chizat and Bach, 2018) since the tangent kernel Kw(t) changes over the course of training, in a sense adapting to the data. We argue that this regime is the one that truly allows us to exploit the power of depth, and thus is the more relevant regime for understanding the success of deep learning.

heading: 4 DETAILED STUDY OF A SIMPLE DEPTH-2 MODEL
text: We study in detail a simple 2-homogeneous model. Consider the class of linear functions over X = Rd, with squared parameterization as follows:
f(w,x) = d∑ i=1 (w2+,i−w2−,i)xi = 〈βw, x〉 , where w = [ w+ w− ] ∈ R2d and βw = w2+−w2− (4)
1Chizat and Bach did not consider only homogenous models, and instead of studying the scale of initialization they studied scaling the output of the model. For homogeneous models, the dynamics obtained by scaling the initialization are equivalent to those obtained by scaling the output, and so here we focus on homogenous models and on scaling the initialization.
where we use the notation z2 for z ∈ Rd to denote elementwise squaring. We consider initializing all weights equally with w0 = 1.
This is nothing but a linear regression model, except with an unconventional parametrization. The model can also be thought of as a “diagonal” linear neural network (i.e. where the weight matrices have diagonal structure) with 2d units. A standard diagonal linear network would have d units, with each unit connected to just a single input unit with weights ui and the output with weight vi, thus implementing the model f((u,v),x) = ∑ i uivixi. But if at initialization |ui| = |vi|, their magnitude will remain equal and their signs will not flip throughout training, and so we can equivalently replace both with a single weight wi, yielding the model f(w,x) = 〈 w2, x 〉 .
The reason for using both w+ and w− (or 2d units) is two-fold. First, it ensures that the image of F (w) is all (signed) linear functions, and thus the model is truly equivalent to standard linear regression. Second, it allows for initialization at F (αw0) = 0 without this being a saddle point from which gradient flow will never escape.2
The model (4) is perhaps the simplest non-trivial D-homogeneous model for D > 1, and we chose it for this reason, as it already exhibits distinct and interesting kernel and rich regimes. Furthermore, we can completely understand both the implicit regularization driving this model and the transition between the regimes analytically.
Consider the behavior of the limit of gradient flow (2) as a function of the initialization, in the underdeterminedN d case where there are many possible solutionsXβ = y. The tangent kernel at initialization is K0(x,x′) = 8α2 〈x, x′〉, i.e. a scaling of the standard inner product kernel, so ‖β‖K0 ∝ ‖β‖2. Thus, in the kernel regime, gradient flow leads to the minimum `2 norm solution, β∗L2 . = arg minXβ=y ‖β‖2. Following Chizat and Bach (2018) and the discussion in Section 3, we thus expect that limα→∞ βα(∞) = β∗L2, and we also show this below. In contrast, Gunasekar et al. (2017) shows that as α → 0, gradient flow leads instead to the minimum `1 norm solution limα→0 βα(∞) = β∗L1 . = arg minXβ=y ‖β‖1. This is the “rich regime.” Comparing this with the kernel regime, we already see two very distinct behaviors and, in high dimensions, two very different inductive biases. In particular, the rich regime’s bias is not an RKHS norm for any choice of kernel. Can we charactarize and understand the transition between the two regimes as α transitions from very small to very large? The following theorem does just that. Theorem 1. For any 0 < α <∞,
βα(∞) = β̂α . = arg min
β Qα (β) s.t. Xβ = y, (5)
where Qα (β) = ∑d i=1 q ( βi α2 ) and q(z) = ∫ z 0 arcsinh ( u 2 ) du = 2− √ 4 + z2 + z arcsinh ( z 2 ) Proof sketch The proof in Appendix A proceeds by showing the gradient flow dynamics on w lead to a solution of the form
βα(∞) = α2 ( exp ( −4X> ∫ ∞ 0 rα(t)dt ) − exp ( 4X> ∫ ∞ 0 rα(t)dt )) (6)
where rα(t) = Xβα(t)− y. While evaluating the integral would be very difficult, the fact that βα(∞) ∈ { α2 ( exp ( −X>r̄ ) − exp ( X>r̄ ) : r̄ ∈ RN )} (7)
already provides a dual certificate for the KKT conditions for minβ Qα (β) s.t. Xβ = y.
In light of Theorem 1, the function Qα (referred to elsewhere as the “hypentropy” function (Ghai et al., 2019)) can be understood as an implicit regularizer which biases the gradient flow solution towards one particular zero-error solution out of the many possibilities. As α ranges from 0 to ∞, the Qα regularizer interpolates between the `1 and `2 norms, as illustrated (labelled D = 2) in Figure 2a, which shows the coordinate function q. As α → ∞ we have that βi/α2 → 0, and so the behaviour of Qα(β) is controlled by the behaviour of q(z) around z = 0. In this regime
2Our results can be generalized to non-uniform initialization, “biased initiliation” (i.e. where w− 6= w+ at initialization), or the asymmetric parameterization f((u, v),x) = ∑ i uivixi, however this complicates the presentation without adding much insight.
q(z) = z 2
4 + O(z 4) is quadratic, and so Qα(β)∝̃
∑
i β 2 i = ‖β‖ 2 2. On the other hand when α → 0,
βi/α 2 →∞ is governed by the asymptotic behaviour q(z) = Θ(z log z) as z →∞. In this regime Qα(β)∝̃ ∑ i βi α2 log βi α2 ∝ ‖β‖1 + o(1). For any initialization scale α, the function Qα describes exactly how training will interpolate between the kernel and rich regimes. The following Theorems, proven in Appendix B, provide a quantitative statement of how the `1 and `2 norms are approached as α→ 0 and α→∞ respectively: Theorem 2. For any 0 < < d
α ≤ min {
(2(1 + ) ‖β∗L1‖1) − 2+ 2 , exp (−d/( ‖β∗L1‖1)) } =⇒ ‖β̂α‖1 ≤ (1 + ) ‖β∗L1‖1
α ≥ √
2(1 + )(1 + 2/ ) ‖β∗L2‖2 =⇒ ‖β̂α‖ 2 2 ≤ (1 + ) ‖β∗L2‖ 2 2
Theorem 2 indicates a certain asymmetry between reaching the rich and kernel regimes: polynomially large α suffices to approximate β∗L2 to a very high degree of accuracy. On the other hand, exponentially small α is sufficient to approximate β∗L1, and Lemma 2 in Appendix B proves that α ≤ d−Ω(1/ ) is necessary in order for Qα to be proportional to the `1 norm, which indicates that α must be exceedingly small to approximate β∗L1 for certain problems.
This suggestive of an explanation for the difficulty of demonstrating rich regime behavior empirically in matrix factorization problems (Gunasekar et al., 2017; Arora et al., 2019a). If the initialization really needs to be exponentially small, then conducting experiments in this regime may be infeasible for practical reasons.
In order to understand the effects of initialization on generalization, consider a simple sparse regression problem, where x(1), . . . ,x(N) ∼ N(0, I) and y(n) = 〈 β∗, x(n) 〉 + N(0, 0.01) where β∗ is r∗-sparse and its non-zero entries are 1/ √ r∗. When N ≤ d, gradient flow will reach a zero training error solution, however, not all of these solutions will generalize the same. With N = Θ(r∗ log d) samples, the rich regime, i.e. the minimum `1 norm solution will generalize well. However, even though we can fit the training data perfectly well, we should not expect any generalization in the kernel regime with this sample size (N = Ω(d) samples would be needed that regime), see Figure 1c. In this case, to generalize well may require using very small initialization, and generalization will improve as we decrease α. From an optimization perspective this is unfortunate because w = 0 is a saddle point, so taking α→ 0 drastically increases the time needed to escape the saddle point. Thus, there is a tension here between generalization and optimization: a smaller α might improve generalization, but it makes optimization trickier. This suggests that in practice we would want to compromise, and operate just at the edge of the rich regime, using the largest α that still allows for generalization. This is borne out in our neural network experiments in Section 7, where standard initialization schemes correspond to being right on the edge of entering the kernel regime, where we expect models to both generalize well and avoid serious optimization difficulties.
The tension between optimization and generalization can also be seen through a tradeoff between the sample size and the largest α we can use and still generalize. In Figure 1c, for each sample size N , we plot the largest α for which the gradient flow solution β̂α achieves population risk below
some threshold. As N approaches the number of samples needed for β∗L1 to generalize (the vertical dashed line), α must become extremely small. However, generalization is much easier when the number of samples is only slightly larger, and we can use much more moderate initialization.

heading: 5 HIGHER ORDER MODELS
text: In the previous Section, we considered a 2-homogeneous model, corresponding to a simple depth-2 “diagonal” network. Deeper models correspond to higher order homogeneity (a depth-D ReLU or linear network is D-homogeneous), motivating us to understand the effect of the order of homogeneity on the transition between the regimes. We therefore generalize our model and consider:
FD(w) = βw,D = w D + −wD− and fD(w,x) = 〈 wD+ −wD− , x 〉 (8)
We again consider initializing all weights equally so w0 = 1. As before, this is just a linear regression model with an unconventional parametrization. It is equivalent to a depth-D matrix factorization model with commutative measurement matrices, as studied by Arora et al. (2019a), and can be thought of as a depth-D diagonal linear network.
We can again study the effect of the scale of the initialization α on the implicit bias. Let βα,D(∞) denote the limit of gradient flow on w when initialized at α1 for the D-homogeneous model. Using the same approach as in Section 4, in Appendix C we show: Theorem 3. For any α and D ≥ 3, if gradient flow reaches a solution Xβα,D(∞) = y, then
βα,D(∞) = arg minβ QDα (β) s.t. Xβ = y
where QDα (β) = ∑d i=1 qD(βi/α D) and qD = ∫ h−1D is the antiderivative of the unique inverse of hD(z) = (1 − z)− D D−2 − (1 + z)− D D−2 on [−1, 1]. Furthermore, limα→0 βα,D(∞) = β∗L1 and limα→∞ βα,D(∞) = β∗L2.
In the two extremes we see that we again get the minimum `2 solution in the kernel regime, and more interestingly, for any depth D ≥ 2, we get the same minimum `1 norm solution in the rich regime, as has also been observed by Arora et al. (2019a). The fact that the rich regime solution does not change with depth is perhaps surprising, and does not agree from what is obtained with explicit regularization (regularizing ‖w‖2 is equivalent to ‖β‖2/D regularization), nor with implicit regularization on logistic-type loss Gunasekar et al. (2017).
Although the two extremes do not change as we go beyondD = 2, what does change is the intermediate regime (the regularizer QDα is unique and cannot be obtained with any other order D
′ 6= D), as well as the sharpness of the transition into the extreme regimes, as illustrated in Figures 2a-2c. The most striking difference is that for order D > 2 the scale of α needed to approximate the `1 is polynomial rather then exponential, yielding a much quicker transition to the “rich regime”, and
allowing near-optimal sparse regression with reasonable initialization scales. Increasing D further hastens the transition. This might also help explain some of the empirical observations about the benefit of depth in deep matrix factorization Arora et al. (2019a).

heading: 6 DEMONSTRATION IN MATRIX COMPLETION
text: We now turn to a more complex depth two model, namely a matrix factorization model, and demonstrate similar transitions empirically. Specifically, we consider the model over matrix inputs X ∈ Rd×d defined by f((U, V ), X) = 〈 UV >, X 〉 , where U, V ∈ Rd×k. This corresponds to linear predictors over matrix arguments specified by F (U, V ) = UV >. In overparameterized regime with k ≥ d, the parameterization itself does not introduce any explicit rank constraints. We consider here a random low rank matrix completion problem where Xn = eine > jn
represents an uniform random observation of entry (in, jn) of a planted low rank matrix M∗ of rank r∗ d: yn = 〈M∗, Xn〉 = M∗ij . For underdetermined problems where N d2, there are many trivial global minimizers, most of which are not low rank and hence will not guarantee recovery. As was demonstrated empirically by Gunasekar et al. (2017) and also proven rigorously for Gaussian measurements by Li et al. (2018), as α→ 0 gradient flow implicitly regularizes the nuclear norm, which for random measurements leads to recovery of ground truth (Candès and Recht, 2009; Recht et al., 2010): these are very different and rich implicit biases that are not RKHS norms.
Crucially, the reconstruction results in Gunasekar et al. (2017) and Li et al. (2018) are dependent on initialization with scale α → 0. Here we further explore the role of initialization. Similar to Section 4, in order to get unbiased 0 initialization, we consider k ≥ 2d and initialization of the form U(0) = α [U0,−U0] and V (0) = α [V0, V0], where U0, V0 ∈ Rd×k/2. We will study implicit bias of gradient flow over the factorized parameterization with above initialization.
For matrix completion problems with X = eiXe > jX , the tangent kernel at initialization is given by K0(X,X ′) ∝ 〈U0[iX , :], U0[iX′ , :]〉1(jX = jX′) + 〈V0[jX , :], V0[jX′ , :]〉1(iX = iX′). This defaults to the trivial delta kernel K(X,X ′) = 1(iX = iX′) ·1(jX = jX′) for the two special cases (a) U0, V0 have orthogonal columns (e.g. U0 = V0 = I), or (b) U0, V0 have independent Gaussian entries and k → ∞. In these cases, minimizing the RKHS norm of the tangent kernel corresponds to returning a zero imputed matrix (minimum Frobenius norm solution). Figure 3 demonstrates the behaviour of gradient flow updates in the “rich" regime (where for α→ 0 recovers the ground truth) and in the “kernel" regime (where for large α, there are no updates to the unobserved entries).

heading: 7 NEURAL NETWORK EXPERIMENTS
text: In the preceding sections, we intentionally focused on the simplest possible models in which a kernel-to-rich transition can be observed, in order to isolate this phenomena and understand it in
detail. In those simple models, we were able to obtain a complete analytic description of the transition. Obtaining such a precise description in more complex models is somewhat optimistic at this point, as we do not yet have a satisfying description of even just the rich regime. Instead, we now provide empirical evidence suggesting that also for non-linear and realistic networks, the scale of initialization induces a transition into and out of a “kernel” regime, and that to reach good generalization we must operate outside of the “kernel” regime. To track whether we are in the “kernel” regime, we track how much the gradient ∇wf(w(t),x) changes throughout training. In particular, we define the gradient distance to be the cosine distance between the initial tangent kernel feature map∇wf(w(0),x) and the final tangent kernel feature map∇wf(w(T ),x). In Figures 4a and 4b, we see that also for a non-linear ReLU network, we remain in the kernel regime when the initialization is large, and that exiting from the kernel regime is necessary in order to achieve small test error on the synthetic data. Interestingly, when αD ≈ 1, the models achieve good test error but have smaller gradient distance which, not coincidentally, corresponds to using the out-of-the-box Uniform He initialization. This lies on the boundary between the rich and kernel regimes, which is desirable due to the learning vs. optimization tradeoffs discussed in Section 4. On MNIST data, Figure 4e shows that previously published successes with training overly wide depth-2 ReLU networks without explicit regularization (e.g. Neyshabur et al., 2014) relies on the initialization being small, i.e. being outside of the “kernel regime”. In fact, the 2.4% test error reached for large initialization is no better than what can be achieved with a linear model over a random feature map. Turning to a more realistic network, 4f shows similar behavior when training a VGG11-like network on CIFAR10.
So far, we attempted to use the best fixed stepsize for each initialization (i.e. achieving the best test error). But as demonstrated in Figures 4c and 4d, the stepsize choice can also have a significant effect, with larger stepsizes allowing one to exit the kernel regime even at an initialization scale where a smaller stepsize would remain trapped in the kernel regime. Further analytic and empirical studies are necessary in order to understand the joint behavior of the stepsize and initialization scale.

heading: 8 DISCUSSION
text: The main point of this paper is to emphasize the distinction between the “kernel” regime in training overparametrized multi-layered networks, and the “rich” (active, adaptive) regime, show how the scaling of the initialization can transition between them, and understand this transition in detail. We argue that rich inductive bias that enables generalization may arise in the rich regime, but that focusing on the kernel regime restricts us to only what can be done with an RKHS. By studying the transition we also see a tension between generalization and optimization, which suggests we would tend to operate just on the edge of the rich regime, and so understanding this transition, rather then just the extremes, is important. Furthermore, we see that at the edge of the rich regime, the implicit bias of gradient descent differs substantively from that of explicit regularization. Although in our theoretical study we focused on a simple model so that we can carry out a complete and exact analysis analytically, our experiments show that this is representative of the behaviour also in other homogeneous models, and serve as a basis of a more general understanding.
Effect of Width Our treatment focused on the effect of scale on the transition between the regimes, and we saw that, as pointed out by Chizat and Bach, we can observe a very meaningful transition between a kernel and rich regime even for finite width parametric models. The transition becomes even more interesting if the width of the model (the number of units per layer, and so also the number of parameters) increases towards infinity. In this case, we must be careful as to how the initialization of each individual unit scales when the total number of units increase, and which regime we fall in to is controlled by the relative scaling of the width and the scale of individual units at initialization. This is demonstrated, for example, in Figure 5, which shows the regime change in matrix factorization problems, from minimum Frobenius norm recovery (the kernel regime) to minimum nuclear norm recovery (the rich regime), as a function of both the number of factors k and the scale of initialization of each factor α . As is expected, the scale α at which we see the transition decreases as the model becomes wider, but further study is necessary to obtain a complete understanding of this scaling.
A particularly interesting aspect of infinite width networks is that, unlike for fixed-width networks, it may be possible to scale α relative to the width k such that at the infinite-width limit we would
have an (asymptotically) unbiased predictor at initialization limk→∞ Fk(w(0)) = 0, or at least a non-exploding initialization lim supk→∞ ‖Fk(w(0))‖ = O(1), even with random initialization (without a doubling trick leading to artificially unbiased initialization), while still being in the kernel regime. For two-layer networks with ReLU activation, Arora et al. (2019b) showed that with width k ≥ poly(1/max‖x‖2≤1,N ‖f0(x)‖) the gradient dynamics stay in the kernel regime forever.

heading: A PROOF OF THEOREM 1
text: It is straightforward, given the expression for Qα, to prove that β̂α is the minimum Qα solution to Xβ = y. In other words, if we had been able to guess from the beginning that the implicit bias would be governed by Qα, it would have been easy to prove this fact. However, a key contribution of this work is in developing a method for determining what the implicit bias is when we do not already have a good guess.
We will first describe a general approach for deriving the implied regularizer reached at the limit of some gradient flow dynamics which minimize the loss. In particular, our approach will apply whenever solutions of the gradient flow dynamics can be written as βα(t) = fα(X>ν(t)) for some ν(t).
The first step is to analyze the gradient flow dynamics of βα and show that
Xβα(∞) = y βα(∞) = fα(X>ν(∞))
(9)
for some function fα and some vector ν(∞). It is not important for our approach to know exactly what ν(∞) is. This is useful because calculating this ν(∞) will often be difficult, even for the simple examples we consider.
The next step is to suppose that there is some function Qα such that
βα(∞) = arg min β Qα(β) s.t. Xβ = y (10)
and to write down the KKT optimality conditions for (10):
Xβ∗ = y
∃ν ∇Qα (β∗) = X>ν (11)
Finally, we connect (9) with (11). Specifically, βα(∞) satisfies the first KKT condition, and, if βα(∞) = β∗ then by taking ν = ν(∞) we have
∇Qα (βα(∞)) = ∇Qα ( fα(X >ν) ) = X>ν (12)
It follows that ∇Qα(β) = f−1α (β). Therefore, to derive Qα we simply invert fα and integrate its inverse. We prove Theorems 1 and 3 using this method. Theorem 1. For any 0 < α <∞,
βα(∞) = β̂α . = arg min
β Qα (β) s.t. Xβ = y, (5)
where Qα (β) = ∑d i=1 q ( βi α2 ) and q(z) = ∫ z 0 arcsinh ( u 2 ) du = 2− √ 4 + z2 + z arcsinh ( z 2 ) Proof. We begin by calculating the gradient flow dynamics on w, since the linear predictor βα(∞) is given by F applied to the limit of the gradient flow dynamics on w. Recalling that X̃ = [X −X],
ẇα(t) = −∇L(wα(t)) = −∇ (∥∥∥X̃wα(t)2 − y∥∥∥2
2
) = −2X̃>rα(t) ◦wα(t) (13)
where the residual rα(t) , X̃wα(t)2 − y, and a ◦ b denotes the element-wise product of a and b. It is easily confirmed that these dynamics have a solution:
wα(t) = wα(0) ◦ exp ( −2X̃> ∫ t 0 rα(s)ds ) (14)
This immediately gives an expression for βα(t):
βα(t) = wα,+(t) 2 −wα,−(t)2 (15) = α2 ( exp ( −4X> ∫ t 0 rα(s)ds ) − exp ( 4X> ∫ t 0 rα(s)ds )) (16)
= 2α2 sinh ( −4X> ∫ t 0 rα(s)ds ) (17)
In addition, this problem satisfies the strict saddle property (Ge et al., 2015) (Zhao et al., 2019, Lemma 2.1), therefore gradient flow will converge to a zero-error solution, i.e. Xβα(∞) = y. Thus, we conclude that βα(∞) is a global minimum with zero error, i.e. Xβα(∞) = y. Thus, we have established that
Xβα(∞) = y βα(∞) = fα(X>ν(∞))
(18)
for fα(z) = 2α2 sinh(z) and ν(∞) = −4 ∫ t
0 rα(s)ds. This corresponds to (9) from our general
approach detailed above. Continuing from (12) we have that
∇Qα(β) = f−1α (β) = arcsinh ( β/2α2 ) (19)
Integrating this expression completes the proof.

heading: B PROOF OF THEOREM 2
text: Lemma 1. For any β ∈ Rd,
α ≤ α1 ( , ‖β‖1 , d) := min { 1, √ ‖β‖1, (2 ‖β‖1) − 12 , exp ( − d
2 ‖β‖1 )} guarantees that
‖β‖1 (1− ) ≤ α2
ln(1/α2) Q(β/α2) ≤ ‖β‖1 (1 + )
Proof. First, we show that Q(β/α2) = Qα(|β| /α2). Observe that g(x) = x arcsin(x/2) is even because x and arcsin(x/2) are odd. Therefore,
Q(β/α2) = d∑ i=1 2− √ 4 + β2i α4 + βi α2 arcsinh ( βi 2α2 ) (20)
= d∑ i=1 2− √ 4 + β2i α4 + g ( βi α2 ) (21)
= d∑ i=1 2−
√ 4 + |βi|2
α4 + g (∣∣∣∣ βiα2 ∣∣∣∣) (22)
= Qα(|β| /α2) (23)
Therefore, we can rewrite α2
ln(1/α2) Q(β/α2) =
α2
ln(1/α2) Q(|β| /α2) (24)
= d∑ i=1 2α2 ln(1/α2) − √ 4α4 + β2i ln(1/α2) + |βi| ln(1/α2) arcsinh ( |βi| 2α2 ) (25)
= d∑ i=1 2α2 ln(1/α2) − √ 4α4 + β2i ln(1/α2) + |βi| ln(1/α2) ln ( |βi| 2α2 + √ 1 + β2i 4α4 ) (26)
= d∑ i=1 2α2 ln(1/α2) − √ 4α4 + β2i ln(1/α2) + |βi| 1 + ln ( |βi| 2 + √ α4 + β2i 4 ) ln(1/α2)  (27)
Using the fact that |a| ≤ √ a2 + b2 ≤ |a|+ |b| (28)
we can bound for α < 1
α2
ln(1/α2) Q(β/α2) ≤ d∑ i=1 2α2 ln(1/α2) − 2α 2 ln(1/α2) + |βi|
1 + ln ( |βi| 2 + α 2 + |βi|2 ) ln(1/α2)  (29) =
d∑ i=1 |βi|
( 1 + ln ( |βi|+ α2 ) ln(1/α2) ) (30)
≤ ‖β‖1
( 1 + max
i∈[d]
ln ( |βi|+ α2 ) ln(1/α2) ) (31)
≤ ‖β‖1
( 1 + ln ( ‖β‖1 + α2 ) ln(1/α2) ) (32)
(33) So, for any α ≤ min { 1, √ ‖β‖1, (2 ‖β‖1) − 12 } , then
α2
ln(1/α2) Q(β/α2) ≤ ‖β‖1
( 1 + ln ( ‖β‖1 + α2 ) ln(1/α2) ) (34)
≤ ‖β‖1
( 1 +
ln (2 ‖β‖1) ln(1/α2)
) (35)
≤ ‖β‖1 (1 + ) (36)
On the other hand, using (27) and (28) again,
α2
ln(1/α2) Q(β/α2) ≥ d∑ i=1 2α2 ln(1/α2) − |βi|+ 2α 2 ln(1/α2) + |βi| ( 1 + ln (|βi|) ln(1/α2) ) (37)
= d∑ i=1 |βi| ( 1 + ln (|βi|)− 1 ln(1/α2) ) (38)
Using the inequality ln(x) ≥ 1− 1x , this can be further lower bounded by
α2
ln(1/α2) Q(β/α2) ≥ d∑ i=1 |βi| − 1 ln(1/α2) (39)
= ‖β‖1 − d
ln(1/α2) (40)
Therefore, for any α ≤ exp ( − d2 ‖β‖1 ) then
α2
ln(1/α2) Q(β/α2) ≥ ‖β‖1 (1− ) (41)
We conclude that for α ≤ min { 1, √ ‖β‖1, (2 ‖β‖1) − 12 , exp ( − d2 ‖β‖1 )} that
‖β‖1 (1− ) ≤ α2
ln(1/α2) Q(β/α2) ≤ ‖β‖1 (1 + ) (42)
Lemma 2. Fix any > 0 and d ≥ max { e, 124 } . Then for any α ≥ d− 14− 18 , Q(β/α2) 6∝ ‖β‖1 in the sense that there exist vectors v, w such that
Q ( 1 α2 v )
‖v‖1 ≥ (1 + )
Q (
1 α2w ) ‖w‖1
Proof. First, recall that
q
( 1
cα2
) = 2− √ 4 + 1
c2α4 +
1
cα2 arcsinh
( 1
2cα2
) (43)
= 1
cα2
( 2cα2 − √ 4c2α4 + 1 + ln ( 1
2cα2 +
√ 1 + 1
4c2α4
)) (44)
= 1
cα2
( 2cα2 − √ 4c2α4 + 1 + ln ( 1
cα2
) + ln ( 1
2 +
√ c2α4 + 1
2
)) (45)
Thus,
− 1 + ln ( 1
cα2
) ≤ cα2q ( 1
cα2
) ≤ 3cα2 − 1 + ln ( 1
cα2
) (46)
Now, consider the ratio
Q (
1 α2 e1 ) Q (
1 α2 1d ‖1d‖2
) = q ( 1α2 ) dq (
1 α2 √ d ) (47) =
1√ d
α2q (
1 α2 ) α2 √ dq (
1 α2 √ d ) (48) Using (46), we conclude
√ d
Q (
1 α2 e1 ) Q (
1 α2 1d ‖1d‖2
) ≥ −1 + ln ( 1α2 ) 3 √ dα2 − 1 + ln ( 1
α2 √ d ) (49) = −1 + ln ( 1 α2 ) 3 √ dα2 − 1 + ln ( 1 α2 ) − 12 ln(d) (50)
= 1 + ln(d)− 6
√ dα2
6 √ dα2 − 2 + 2 ln
( 1
α2 √ d
) (51)
Fix any > 0 and d ≥ max { e, 124 } , and set α = d− 1 4− 1 8 . Then
2d 1 4 ≥ 6 and d
1 4 2 ln d ≥ 6 (52)
=⇒ 2d 14 − 6 + d 1 4
2 ln d− 6 ≥ 0 (53)
=⇒ ( 1 − 1
2
) ln d− 6 d− 1 4 ≥ 6d− 14 − 2 (54)
=⇒ 1 ln d− 6 d− 1 4 ≥ 6d− 14 − 2 + 1
2 ln d (55)
=⇒ ln (d)− 6α2 √ d ≥ ( 6α2 √ d− 2 + 2 ln ( 1
α2 √ d
)) (56)
This implies that the second term of (51) is at least . We conclude that for any > 0 and d ≥ max { e, 124 } , α = d− 1 4− 1 8 implies that
Q (
1 α2 e1 ) Q (
1 α2 1d ‖1d‖2 ) ≥ (1 + ) ‖e1‖1∥∥∥ 1d‖1d‖2 ∥∥∥1 (57) Consequently, for at least one of these two vectors, Q is not proportional to the `1 norm up to accuracy O( ) for this value of α.
It is straightforward to confirm that
d
dα
Q (
1 α2 e1 ) Q (
1 α2 1d ‖1d‖2 ) ≥ 0 (58) which concludes the proof.
Lemma 3. For any β ∈ Rd, α ≥ α2( , ‖w‖2) := √ ‖β‖2 ( 1 + − 1 4 ) guarantees that
(1− ) ‖β‖22 ≤ 4α 4Q(β/α2) ≤ (1 + ) ‖β‖22
Proof. The regularizer Q can be written
Q(β/α2) = d∑ i=1 ∫ βi/α2 0 arcsinh ( t 2 ) dt (59)
Let φ(z) = ∫ z/α2
0 arcsinh ( t 2 ) dt, then
φ(0) = 0 (60)
φ′(0) = 1
α2 arcsinh ( z 2α2 )∣∣∣∣ z=0 = 0 (61)
φ′′(0) = 1 α4 √ 4 + z 2
α4 ∣∣∣∣∣∣ z=0 = 1 2α4 (62)
φ′′′(0) = −z α8 ( 4 + z 2
α4 )3/2 ∣∣∣∣∣∣ z=0 = 0 (63)
φ′′′′(z) = 3z2 α12 ( 4 + z 2
α4
)5/2 − 1 α8 ( 4 + z 2
α4
)3/2 (64)
Also, note that
|φ′′′′(z)| = ∣∣2z2 − 4α4∣∣
α12 ( 4 + z 2
α4 )5/2 (65) ≤ z 2 + 2α4
16α12 (66)
Therefore, by Taylor’s theorem, for some ξ with |ξ| ≤ |z|∣∣∣∣φ(z)− z24α4 ∣∣∣∣ = φ′′′′(ξ)4! z4 (67)
=⇒ ∣∣∣∣φ(z)− z24α4 ∣∣∣∣ ≤ sup |ξ|≤|z| φ′′′′(ξ) 4! z4 ≤ z 6 + 2α4z4 384α12 = z2 4α4 z4 + 2α4z2 96α8 (68)
Therefore, for any β ∈ Rd
∣∣∣4α4Qα(β)− ‖β‖22∣∣∣ = 4α4 ∣∣∣∣∣ d∑ i=1 φ(βi)− β2i 4α4 ∣∣∣∣∣ (69) ≤ 4α4
d∑ i=1 ∣∣∣∣φ(βi)− β2i4α4 ∣∣∣∣ (70)
≤ d∑ i=1 β2i · β4i + 2α 4β2i 96α8
(71)
≤ ‖β‖22 maxi β4i + 2α 4β2i 96α8
(72)
Therefore, α ≥ √ ‖β‖2 ( 1 + − 1 4 ) ensures
(1− ) ‖β‖22 ≤ 4α 4Q(β/α2) ≤ (1 + ) ‖β‖22 (73)
Theorem 2. For any 0 < < d
α ≤ min {
(2(1 + ) ‖β∗L1‖1) − 2+ 2 , exp (−d/( ‖β∗L1‖1)) } =⇒ ‖β̂α‖1 ≤ (1 + ) ‖β∗L1‖1
α ≥ √
2(1 + )(1 + 2/ ) ‖β∗L2‖2 =⇒ ‖β̂α‖ 2 2 ≤ (1 + ) ‖β∗L2‖ 2 2
Proof. We prove the `1 and `2 statements separately. `1 approximation First, we will prove that ∥∥∥β̂α∥∥∥
1 < (1 + 2 ) ‖β∗L1‖1. By Lemma 1, since α ≤
α1
(
2+ , (1 + 2 ) ‖β ∗ L1‖1 , d ) , for all β with ‖β‖1 ≤ (1 + 2 ) ‖β∗L1‖1 we have
‖β‖1
( 1−
2 +
) ≤ α 2
ln(1/α2) Q(β/α2) ≤ ‖β‖1
( 1 +
2 +
) (74)
Let β be such that Xβ = y and ‖β‖1 = (1 + 2 ) ‖β∗L1‖1. Then α2
ln(1/α2) Q(β/α2) ≥
( 1−
2 +
) ‖β‖1 (75)
= ( 1−
2 +
) (1 + 2 ) ‖β∗L1‖1 (76)
≥
( 1− 2+ ) (
1 + 2+
) (1 + 2 ) α2 ln(1/α2) Q(β∗L1/α 2) (77)
= 1 + 2
1 +
α2
ln(1/α2) Q(β∗L1/α 2) (78)
> α2
ln(1/α2) Q(β∗L1/α 2) (79)
≥ α 2
ln(1/α2) Q(β̂α/α
2) (80)
Therefore, β 6= β̂α. Furthermore, let β be any solution Xβ = y with ‖β‖1 > (1 + 2 ) ‖β∗L1‖1. It is easily confirmed that there exists c ∈ (0, 1) such that the point β′ = (1 − c)β + cβ∗L1 is satisfies both Xβ′ = y and ‖β′‖1 = (1 + 2 ) ‖β∗L1‖1. By the convexity of Q, this implies Q(β/α2) ≥ Q(β′/α2) > Qα(β̂α/α2). Thus a β with a large `1 norm cannot be a solution, even if
α2
ln(1/α2)Q(β/α 2) 6≈ ‖β‖1. Since ∥∥∥β̂α∥∥∥
1 < (1 + 2 ) ‖β∗L1‖1, we conclude∥∥∥β̂α∥∥∥
1 ≤ 1 1− 2+ α2 ln(1/α2) Q(β̂α/α 2) (81)
≤ 1 1− 2+
α2
ln(1/α2) Q(β∗L1/α 2) (82)
≤ 1 + 2+ 1− 2+ ‖β∗L1‖1 (83) = (1 + ) ‖β∗L1‖1 (84)
First, we will prove that ∥∥∥β̂α∥∥∥
2 < (1 + 2 ) ‖β∗L2‖2. By Lemma 3, since α ≥
α2
(
2+ , (1 + 2 ) ‖β ∗ L2‖2 ) , for all β with ‖β‖2 ≤ (1 + 2 ) ‖β∗L2‖2 we have
‖β‖22
( 1−
2 +
) ≤ 4α4Q(β/α2) ≤ ‖β‖22 ( 1 +
2 +
) (85)
Let β be such that Xβ = y and ‖β‖2 = (1 + 2 ) ‖β∗L2‖2. Then 4α4Q(β/α2) ≥ (
1− 2 +
) ‖β‖22 (86)
= ( 1−
2 +
) (1 + 2 ) ‖β∗L2‖ 2 2 (87)
≥
( 1− 2+ ) (
1 + 2+ ) (1 + 2 )4α4Q(β∗L2/α2) (88) = 1 + 2
1 + 4α4Q(β∗L2/α 2) (89) > 4α4Q(β∗L2/α 2) (90)
≥ 4α4Q(β̂α/α2) (91)
Therefore, β 6= β̂α. Furthermore, let β be any solution Xβ = y with ‖β‖2 > (1 + 2 ) ‖β∗L2‖2. It is easily confirmed that there exists c ∈ (0, 1) such that the point β′ = (1−c)β+cβ∗L2 satisfiesXβ′ = y and ‖β′‖2 = (1 + 2 ) ‖β∗L2‖2. By the convexity of Q, this implies Q(β/α
2) ≥ Q(β′/α2) > Q(β∗L2/α 2). Thus a β with a large `2 norm cannot be a solution, even if 4α4Q(β/α2) 6≈ ‖β‖22.
Since ∥∥∥β̂α∥∥∥
2 < (1 + 2 ) ‖β∗2‖2, we conclude∥∥∥β̂α∥∥∥2
2 ≤ 1 1− 2+ 4α4Q(β̂α/α 2) (92)
≤ 1 1− 2+ 4α4Q(β∗L2/α 2) (93)
≤ 1 + 2+ 1− 2+ ‖β∗L2‖ 2 2 (94) = (1 + ) ‖β∗L2‖ 2 2 (95)

heading: C PROOF OF THEOREM 3
text: Lemma 4. For the D-homogeneous model (8), ∀t ∥∥∥∥X> ∫ t
0
r(τ)dτ ∥∥∥∥ ∞ ≤ α 2−D D(D − 2)
Proof. For the order-D unbiased model β(t) = wD+ −wD− , the gradient flow dynamics are
ẇ+(t) = − dL
dw+ = −DX>r(t) ◦wD−1+ (t), w+(0) = α1 (96)
=⇒ w+(t) = ( α2−D1 +D(D − 2)X> ∫ t 0 r(τ)dτ )− 1D−2 (97)
Where ◦ denotes elementwise multiplication, r(t) = Xβ(t) − y, and where all exponentiation is elementwise. Similarly,
ẇ−(t) = − dL
dw− = DX>r(t) ◦wD−1− (t), w−(0) = α1 (98)
=⇒ w−(t) = ( α2−D1−D(D − 2)X> ∫ t 0 r(τ)dτ )− 1D−2 (99)
First, we observe that ∀t∀i w+(t)i ≥ 0 and ∀t∀i w−(t)i ≥ 0. This is because at time 0, w+(0)i = w−(0)i = α > 0; the gradient flow dynamics are continuous; and w+(t)i = 0 =⇒ ẇ+(t)i = 0 and w−(t)i = 0 =⇒ ẇ−(t)i = 0. Consequently,
0 ≤ w+(t)2−Di = α 2−D +D(D − 2)
[ X> ∫ t 0 r(τ)dτ ] i
(100)
0 ≤ w−(t)2−Di = α 2−D −D(D − 2)
[ X> ∫ t 0 r(τ)dτ ] i
(101)
=⇒ −α2−D ≤ D(D − 2) [ X> ∫ t 0 r(τ)dτ ] i ≤ α2−D (102)
which concludes the proof.
Theorem 3. For any α and D ≥ 3, if gradient flow reaches a solution Xβα,D(∞) = y, then
βα,D(∞) = arg minβ QDα (β) s.t. Xβ = y
where QDα (β) = ∑d i=1 qD(βi/α D) and qD = ∫ h−1D is the antiderivative of the unique inverse of hD(z) = (1 − z)− D D−2 − (1 + z)− D D−2 on [−1, 1]. Furthermore, limα→0 βα,D(∞) = β∗L1 and limα→∞ βα,D(∞) = β∗L2.
Proof. For the order-D unbiased model β(t) = wD+ −wD− , the gradient flow dynamics are
ẇ(t) = dL
dw = −DX̃>r(t) ◦wD−1, w(0) = α1
=⇒ w(t) = ( α2−D +D(D − 2)X̃> ∫ t 0 r(τ)dτ )− 1D−2 =⇒ β(t) = αD ( 1 + αD−2D(D − 2)X>
∫ t 0 r(τ)dτ )− DD−2 −αD ( 1− αD−2D(D − 2)X>
∫ t 0 r(τ)dτ )− DD−2 where X̃ = [X −X] and r(t) = Xβ(t)− y. Supposing β(t) converges to a zero-error solution,
Xβ(∞) = yβ(∞) = αDhD(X>ν(∞))
where ν(∞) = −αD−2D(D − 2) ∫∞
0 r(τ)dτ and the function hD is applied elementwise and is
defined hD(z) = (1− z)− D D−2 − (1 + z)− D D−2
By Lemma 4, ∥∥X>ν∥∥∞ ≤ 1, so the domain of hD is the interval [−1, 1], upon which it is monotonically increasing from hD(−1) = −∞ to hD(1) = ∞. Therefore, there exists an inverse mapping h−1D (t) with domain [−∞,∞] and range [−1, 1]. This inverse mapping unfortunately does not have a simple closed form. Nevertheless, it is the root of a rational equation. Using the approach outlined in Appendix A (12), we conclude:
QDα (β) = ∑ i ∫ βi/αD 0 h−1D (t)dt
Rich regime Next, we show that if gradient flow reaches a solution Xβα,D(∞) = y, then limα→0 βα,D(∞) = β∗L1 for any D. This is implied by the work of Arora et al. (2019a), but we include it here for an alternative, simpler proof for our special case, and for completeness’s sake.
The KKT conditions for β = β∗L1 are Xβ = y and ∃ν sign(β) = X>ν (where sign(0) = [−1, 1]). The first condition is satisfied by assumption. Define ν as above. We will demonstrate that the second condition holds too in the limit as α→ 0. First, by Lemma 4,
∥∥X>ν∥∥∞ ≤ 1 for all α and D. Thus, for any coordinates i such that limα→0[βα,D(∞)]i = 0, the second KKT condition holds. Consider now i for which limα→0[βα,D(∞)]i > 0. As shown above,
lim α→0 [βα,D(∞)]i = lim α→0
αD ( 1− [X>ν]i )− DD−2 − αD (1 + [X>ν]i)− DD−2 > 0 (103) =⇒ lim α→0 αD ( 1− [X>ν]i
)− DD−2 > 0 (104) This and [X>ν]i ≤ 1 implies limα→0[X>ν]i = 1, and thus the positive coordinates satisfy the second KKT condition. An identical argument can be made for the negative coordinates.
Kernel Regime Finally, we show that if gradient flow reaches a solution Xβα,D(∞) = y, then limα→∞ βα,D(∞) = β∗L2 for any D. First, since X and y are finite, there exists a solution β∗ whose entries are all finite, and thus all the entries of βα,D(∞), which is the QDα -minimizing solution, will be finite.
The KKT conditions for β = β∗L2 are Xβ = y and ∃µ β = X>µ. The first condition is satisfied by assumption. Defining ν as above, we have
lim α→∞ [βα,D(∞)]i = lim α→0
αD ( 1− [X>ν]i )− DD−2 − αD (1 + [X>ν]i)− DD−2 <∞ (105) =⇒ lim
α→∞ [X>ν]i = 0 (106)
Consequently, defining µ = 2Dα D
D−2 ν, and observing that for small z,
(1− z)− D D−2 − (1 + z)− D D−2 =
2D
D − 2 z +O(z3) (107)
we conclude
lim α→∞ [βα,D(∞)]i [X>µ]i = lim α→0
αD ( 1− [X>ν]i )− DD−2 − αD (1 + [X>ν]i)− DD−2 [X>µ]i
(108)
= lim α→0
αD (
2D D−2 [X >ν]i +O([X >ν]3i ) ) 2DαD
D−2 [X >ν]i
(109)
= 1 + lim α→0
O([X>ν]2i ) (110)
= 1 (111) Thus, the KKT conditions are satisfied for limα→∞ βα,D(∞) = β∗L2.

heading: D MATRIX EXPERIMENTS
text: In the appendix, we provide additional results similar to those in Section 6. First, in Figure 5, we plot the implicit regularization behavior of gradient flow limits with identity initialization as in Figure 3(a): in “rich" regime (small α) we recover the minimum nuclear norm solution M∗NN = arg minPΩ(M)=y ‖M‖?, while “kernel" regime recovers the minimum Frobenius norm solution M∗L@ = arg minPΩ(M)=y ‖M‖2

heading: E NEURAL NETWORK EXPERIMENTS
text: Synthetic Experiments We construct a synthetic training set withN = 10 points drawn uniformly from the unit circle in R2 and labelled by a teacher model with 1 hidden layer of 3 units. We train fully connected ReLU networks with depths 2, 3, and 5 with 30 units per layer to minimize the square loss using full gradient descent with constant stepsize 0.01 until the training loss is below 10−9. We use Uniform He initialization for the weights and then multiply them by α.
Here, we describe the details of the neural network implementations for the MNIST and CIFAR10 experiments.
MNIST Since our theoretical results hold for the squared loss and gradient flow dynamics, here we empirically assess whether different regimes can be observed when training neural networks following standard practices.
We train a fully-connected neural network with a single hidden layer composed of 5000 units on the MNIST dataset, where weights are initialized as αw0, w0 ∼ N ( 0, √
2 nin
) , nin denoting the
number of units in the previous layer, as suggested by He et al. (2015). SGD with a batch size of 256 is used to minimize the cross-entropy loss over the 60000 training points, and error over the 10000 test samples are used as measure of generalization. For each value of α, we search over learning rates (0.5, 0.01, 0.05, . . . ) and use the one which resulted in best generalization.
There is a visible phase transition in Figure 4e in terms of generalization (≈ 1.4% error for α ≤ 2, and ≈ 2.4% error for α ≥ 50), even though every network reached 100% training accuracy and less than 10−5 cross-entropy loss. The black line indicates the test error (2.7%) when training only the output layer of the network, as a proxy for the performance of a linear predictor with features given by a fixed, randomly-initialized hidden layer.
CIFAR10 We trained a VGG11-like architecture, which is as follows: 64-M-128-M-256-256-M512-512-M-512-512-M-FC (numbers represent the number of channels in a convolution layers with no bias, M is a maxpooling layer, and FC is a fully connected layer). Weights were initialized using Uniform He initialization multiplied by α. No data augmentation was used, and training done using SGD with batch size of 128 and learning rate of 0.0001. All experiments ran for 2000 epochs, and reached 100% train accuracy except when training only the last layer, which reached 50.38% train accuracy with LR = 0.001 (chosen after hyperparameter tuning).
In addition, to approximate the test error in the kernel regime, we experimented with freezing the bottom layers and only training the output layer for both datasets (the solid lines in Figures 4e and 4f).
Figure 6 illustrates some of the optimization difficulties that arise from using smaller α as discussed in Section 4.

Review summaries:
This paper analyzes an inductive bias of the gradient flow for diagonal two-or higher-homogeneous models and characterizes a limit point depending on the initialization scale of parameters. Concretely, the paper shows that the gradient flow converges to an interpolator attaining minimum L1- (or L2-norm) when the scale is small (or large). In addition, these analyses are well verified empirically on MNIST and CIFAR-10 datasets.

I really appreciated this paper. It discusses a very complex question ("Are we learning in a kernel regime, or in a rich regime where features are identified") by looking at perhaps the simplest model the authors could think of, and then study in detail the model. And how simple it turns out to be: just a linear regression with a twist. All in all, the paper is indeed is a clear demonstration that the differences between"Kernel" regime and one where some actual Learning is done can be demonstrated on simple examples. It is also the simplest model where one can observe  a non-trivial inductive bias and 'implicit regularisation'