Implicit Acceleration of Gradient Flow in Overparameterized Linear Models
Authors: Salma Tarmoun
Guilherme França
Benjamin David Haeffele
Rene Vidal
Venue:  (2021)

Abstract: We study the implicit acceleration of gradient flow in over-parameterized two-layer linear models. We show that implicit acceleration emerges from a conservation law that constrains the dynamics to follow certain trajectories. More precisely, gradient flow preserves the difference of the Gramian~matrices of the input and output weights and we show that the amount of acceleration depends on both the magnitude of that difference (which is fixed at initialization) and the spectrum of the data. In addition, and generalizing prior work, we prove our results without assuming small, balanced or spectral initialization for the weights, and establish interesting connections between the matrix factorization problem and Riccati type differential equations.

Under review as a conference paper at ICLR 2021
IMPLICIT ACCELERATION OF GRADIENT FLOW IN OVERPARAMETERIZED LINEAR MODELS
A Anonymous authors
Paper under double-blind review
A ABSTRACT
We study the implicit acceleration of gradient flow in wve-ppaaamcterieed two-layer linear models. We show that implicit acceleration emerges from a conservation law that constrains the dynamics to follow certain trajectories. More precisely, gradient flow preserves the difference of the Gramian ma- trices of the input and output weights and show that the amount of acceleration depends on both the magnitude of that difference (which is fixed at initialization) and the spectrum of the data. In addition, and generalizing prior work, we prove our results without assuming small, balanced or spectral initialization for the weights, and establish interesting connections between the matrix factorization problem and Riccati lype differential equations.
INTRODUCTION
Understanding peerpprrameteriaaton in deep learning a puzzling question. Contrary to the common belief that over parammterrtiiinoonnn may hurt generalization and optimization, recent work suggests that over-parameterization may actually bias the optimization algorithm towards solutions that generalize well, phenomenon known implicit regularization or implicit bias, and even accelerate convergence, a phenomenon known as implicit acceleration.
Recent work on the implicil bias in the overeparaeeeerrized regime (e.g. Gunasekar et al. (2018a b); Chizal & Bach (2020); Ji & Telgarsky (2019b)) shows that gradient descent on unregularized problems finds minomum norm solutions. For instance, [Soudry et al. (2018);
Ji & Telgarsky (2019a) analyze linear networks trained for binary classification on linearly separable data and show that the predictor converges to a max-margin solution. Similar idcas have been developed for matrix factorization, yiclding solutions with minimum muclear norm (Gunasckar et al., 2017 Li et al., 2018) or low-rank (Arora al. (2019a). It has also been shown that optimization methods which introduce multiplicative stochastic noise, such as dropout and dropblock, induce nuclear norm regularization (Cavazza et al. 2018 and spectral k-support norm regularization PPll et al. 2020), respectively.
Recent work on the implicit acceleration of gradient descent for matrix factorization and deep linear networks (Arora et al. 2018) shows that when the initialization is sufficiently small and balanced (sce Definition 2), oxer-parameterization acts as a pre-conditioning of the gradient that can be interpreted rate. They claim that acceleration for f--regression is possible only if > 2, though there is no theory supporting such claim. Saxe al. (2014) focused on fz-regression with balanced spectral inilializationss (see Definition 3 and similarly concluded that depth may act ually slow down the convergence. For two-layer linear networks, [aae al. (2019): Gidel et al. (2019) analyzed the dynamics of gradient flow and obtained explicit solutions under the assumption of vanishing spectral initialization, highlighting the sequential learning of the hicrarchical components
a combination of momentum and an adaptive learning
phenomenon that could improve generalization. Several recent papers (e.g Arora et al. (2019b;; Du & Hu (2019): Du et al. (2018b)) have also analyzed the convergence behaviour of gradient descent in the wee-pparamttrized setting, particularly for very wide networks and have concluded linear convergence when the initialization is gaussian balanced. While precise study of the connections between gradient descent and gradient flow dynamics for non-convex problems remains elusive, recent work (Franca et al. 2020) shows that. discrete-time convergence rates can be derived from their continuous-time

---

Under review as a conference paper at ICLR 2021
Table 1: Relationships between work and the state of the art
Small and Balanced
Imbalanced
Spectral
Non spectral
Gradient flow
Our work
Our work
Our work
Our work
Saxe et al.
(2014)
Saxe ct al.
2014
Saxe ct al.
(2019
Saxe ct al.
2019
Gradient descent
Arora et al.
2018
None
Gidel et al.
[2019 Arora et al. 20118)
Gidel et al.
2019
counterparts via symplectic integrators. Therefore, our work focuses on the analysis of gradient flow as a stepping stone for future analysis of gradient descent.
In this paper, we present new analysis of the implicit acceleration of gradient flow for over- parametrized two-layer neural networks that applies not only in the case of small, balanced, spectral initialization but also extends to imbalanced and nou-spectral initializations. We show that the key reason for the implicit acceleration of gradient flow is the existence a conservation law that constrains the dynamics to follow porticular path [ More prociscly, the quantity that is preserved by gradient flow is the difference of the Gramians of the input and output weight matrices, which in turn implies that the difference of the square of the norm of the weight matrices is preserved. The particular case where this difference is zero corresponds to the case of balanced weights, but the more general case of imbalanced weights also emerges as a conserved quantity and plays an important role. In particular, we show that acceleration can occur even in the case of Ez-regression as result of imbalanced initialization.
The reason this phenomenon was not previously observed in (Saxe et al. 2014 2019 Gidel ct al. 2019) is preciscly due to the assumption of balancod initialization, which follows as a particular case of our analysis. Our work also establishes interesting connections with Riccati type differential equations. Indeed, some of our results have a similar flavor to those in (Fukumizu 1998 while others are more general and provide explicit characterization of
the contintous time convergence rate. In short, our work makes the following contributions.
1. In Section 2 we analyze the implicit acceleration propertics of gradient flow for symmetric matrix factorization, providing closed form solution and a convergence rate that depends on the eigenvalues of the data without the assumptions of spectral and small initialization.
2. In Section D we analyze the implicit acceleration properties gradient flow for asymemetric matrix factorization with spectral initialization. We show that implicit acceleration emerges as a consequence of conservation laws that only appear in pverpparameterized settings due to an underlying rotational symmetry.
3. In Section ! we analyze the implicit acceleration propertics gradient flow for asymmetric matrix factorization with an arbitrary initialization. We make connections with Riecati differentia ial equations, obtaining a more general characterization of the convergence rale and establish an interesting link with explicit regularization.
2 GRADIENT FLOW DYNAMICS FOR SYMMETRIC MATRIX FACTORIZATION
In this section, we analyze and compare the dynamics of gradient flow
(t) -xx((((t))
(1)
when applied to two problems.
The first one is learning symmetric one-layer linear model -
min {((X) XERmXm
(2)
where Y E mm is given data matrix that one wishes to approximate by X E m The second one is learning its over-paree wver-parameterized symmetric matrix factorization counterpart
min UERmx
(E(U)= = 2 -U-
(3)
'A quantity Q(x(1)) is said to be conserved by the flow a(1) = s(((t) if it remains constant through dynamical evolution, ie.. di Q(x(t)) 0. For example, in mechanics the of potential and kinetic energies remains constant for conservative system. A conservation law usually a consequence of an underlying symmetry (Noether's theorem). In optimization, this can be seen as a constraint Q(x)= = Qo that is automatically satisfied without having to explicitly enforce it. "Gradient descent, Fnt1 = Zn V)(,,,) is simply an explicit Fuler discretization of
2

---

Under review as a conference paper at ICLR 2021
We show that the dynamics of the linear model converge at rate O(e"),, while the over- parameterized model has rate of ((e=4t))))) where Ø; is the ith eigenvalue of the data matrix Y. Therefore, different spectral components are learned at different rates, which can be faster slower than the non-overparamettrized model depending
the eigenvalues of Y.
Linear model. Let us start with the trivial problem of learning a zon-overparrmettrried
linear model (2 3 Applying the gradient flow ① to problem (2) yields (() + x(1) Y
with (0) Xo. This is a linear differential equation whose unique solution is given by
X(t) Y + (Xo Y)e
(4)
Thus, |xx(t) - Yl1r = - 1I Y and lim ((t) Y at an exponential rate of ((e ').
For completeness, it will be interesting to consider the particular case in which X is constrained to be positive semidefinite (PSD), i.. X = 0. In this casc, notice that if Xo ł 0 and Y ž 0, then (t) ; 0 for all t > O, hence the same dynamics and convergence rate still apply without having to enforce the PSD constraint. Otherwise, if Y is not PSD, gradient flow is not directly applicable.
Symmetric matrix factorization model. Consider the more interesting case of learning a two-layer linear model with tied weights U E R formulated as the symmetric matrix factorization problem in (3). In classical low-rank matrix factorization, one assumes & < m. Here, we consider an over- paramcterized formulation where k > m plays the role of the number of hidden units (width). The gradient flow
(1) on problem B for U, now yiclds Ú = 2(Y - UU with (0)) = Ug. Letting X(1) = (()U((1)) 2 and X()) = UpU = one can easily verify that
== + UÙT 2 X)X + 2X(Y - x) 2Y + 2XY -X
(5)
Equation is known to be rank preserving, i.c. if rank = m and X* lim (t) exists then rank(X")) S. It is thus impossible to recover solution of rank higher than the rank of the initialization.
We also note that (5) is (matrix) differential equation of the Ricoali type. Such equations often characterize dynamical systems behind least squares problems and have been extensively studied in the context of optimal control. Using results from this literature, we obtain (see
Appendix A for the proof):
Proposition 1. For any Xo E Rm the solution to 5 exists and is given by
X(1) = 20 Xo (I
2
provided that Y and the malrir inside the parenthesis above are inverlible.
This solution is derived for any Xo, while the over-parameterizee model requires Xo UgUg > 0. Thus in using ( as an analysis lool, it is important keep in mind the set of allowable initializations. In what follows, we consider the spectral initialization (Saxe et al. 2019: Gidel et al.) 2019). and show that the eigenspace of the data is preserved throughout the entire evolution of the learning dynamics.
Definition (Symmetric Spectral initialization). Let Y $M be the eigendecompo- sition of the data. A spectral initialization is defined as Uo= =I and Up=U TT where Eg 0 is a diagonal malria.
From the explicit solution 6 we can readily obtain a convergence rate for a spectral initialization (see Appendix for short proof).
Corollary 1. 1fY = x Dm xi a'p.o?" is invertible and Xo 22. En oooofoof
is spectral inilialization, the solution to (5) is X(1) +Σ(1)T" Em (D)ro where
o. ( -;)
(t σ; + op,i(e
= oi
(7)
σ; (e -1)
1)'
IIn a linear neural network, Y plays the role of the input-output data correlation matrix and X plays the role of the model's input-output map. In this trivial model, the input correlation matrix is assumed to be the identity is the case when the data is whitened.
3

---

Under review as a conference paper at ICLR 2021
provided the denominator is nonzero. Moreover, if E max(o,; 0) = +ΣΦ" is the projection of Y onto the PSD cone, then for all initializations Xo = 0 such that rank( ( rank(E), (t) converges to Y, at rate Ole (e) where minn min; |o is the smallest eigenvalue of Y in magnilude.
Note from 7 that the ith eigencomponent of converges at rale of O(e so that components with 4/0,1 > 1 are accelerated compared to the non-overpparameleerieel case, components with 4 < are slowed down, and minin( > all components are accelerated. This result about different components of the network being learned at different rates by gradient flow is related in spirit to the result in (Saxe et al. 2019: Gidel et al., 2019) about sequential learning in the asymmetric case with spectral balanced initialization. Here the balancedness is enforced by construction.
Next, we derive the same convergence rate with a more general (non-spectral) initialization. The proof is in Appendix C and makes use of several interesting relations for Riccali differential equations.
Proposition (Convergence rate). Consider the eigenvalue decomposition Y E 0.010f Let Y = > max(o,, 0) be the projection of Y onto the PSD cone and Y = = la.l0:0? For any initialization Xo ž 0, assume that I Y -Y) and Y are nonsingular. Then the solution (t) of converyes exponentially to Y at a rate ( S Ce Ce- (ee-tommn(Y))
(8)
where Omin (Y) the smallest eigenwalue of Y in absolute value and C > constant.
It follows from Proposition 2 that the implicit acceleration for symmetric matrix factorization with spectral initialization can be extended to any positive semidefinite initialization Xo provided I y-1 (Xo - Y) is invertible, which is an exlension the previous assumption on Xg, namely that rank(E, = rank(E). The primary difference is that in the spectral initialization case we can derive the convergence rate for cach cigenvalue of X(t), while in general can only obtain a global convergence rate of the solution.
3 GRADIENT FLOW DYNAMICS FOR ASYMMETRIC MATRIX FACTORIZATION WITH SPECTRAL INITIALIZATION
In this section, analyze the dynamics of gradient flow for the more general asymmetric matrix factorization problem. We show that the implicit acceleration phenomenon is still present and provide explanation for it based the Gramians of the factors. We transform the dynamics to a canonical form and show that the solutions under the spectral initialization are diagonal and can be computed in closed
conservation law for the difference of
form. The closed form solution reveals convergence rate of O(e -Voo7+ where is the ith singular value Y and Ao defines the level of imbalance in the initialization for the ith componeut. As in the symmetric case, data matrices with large singular values induce implicit acceleration. However, in the asymmetric formulation, additional acceleration can be gained by choosing an imbalanced initialization.
Asymmetric matrix factorization model. Consider the asymmetric factorization prob- lem:
min (((U = - - -VVT U.V
(9)
where U E Rmnk V E Raxk and k 2 n 2 m. The gradient flow for this problem takes the form
(10)
In what follows, We will make use of conservation law for the difference of the Gramian matrices U and VVV which has been previously identified in (Du et al. 2018a Arora et al. 2018). Previous works have used this conservation law to ensure balancedness under vanishingl small initialization. In contrast, our analysis is the first to highlight the role of imbalance and the resulting acceleration of gradient flow.
4

---

Under review as a conference paper at ICLR 2021
Conservation law. A straightforward calculation shows that
10) admits invariant:
dQ 1ー
Q=U"U~
2(1 = 200). (11)
dt
The origin behind this conserved quantity Q is a global rotational symmetry of the system in 10,, i.e., the system is invariant under the orthogonal group O(k). To see this, consider the singular value decompasition == = and, following (Saxe et al. 2019 define matrices U and V through
U == "-GG
V=wVG"
(12)
where G is an arbitrary element of O(k). These transformed variables obey
--E-UV-Y' Note that these equations have eractly the same form as 11 up to a gauge freedom on the choice of G. Since Q is real and symmetric, it is diagonalizable by an orthogonal matrix.
Ú=
(13)
Therefore, we can choose G to be the matrix that diagonalizes Q this is a gauge choice. Hence, from 11p we have
- ÚTU = Q(TO( = A Ae. - V.TV
(14)
where Ag, is the (constant) diagonal matrir containing the k eigenvalues of Q(0) = Qo (or 2) which is completely specilied by the initial conditions Un and Vú alone. Note that the number of conserved quantities in Aq. depends k which is equal to the degree of over- over- prappaarameteiiitionn. Though we do not assume balanced initialization in this paper, for further reference and comparison with prior work (Arora et al. 2018 Saxe et al. 2014 2019 let us state its precise meaning since it relates to the conservation law.
Definition (Balanced initialization). (Up, Vo) is said to be balanced if | [|ooll F < c for sufficiently small e > 0, i.e. the conserved quantity in or equiva- lently almost vanishes.
Under the above transformation, the matrix problem with spectral initialization can be reduced to solving one-dimensional systems (one for each component). Proposition provides a closed form solution and explicitly characterizes the evolution of each component along with its implicit acceleration (See Appendix D for the proof).
Definition 3 (Asymmetric Spectral initialization). Let Y DΣST be the SVD of the data. The spectral initialization. is defined as Un = DUnG, Vn = VVG and xo= UVV. where Un and Vo are rectangular diagonal matrices and G is any orthogonal matrix.
Proposition 3 (Implicit acceleration of asymmetric factorization with spectral initialization). Let Y DΣ Σ o'ONp?" be the SVD of the data. The solution to I0 with spectral imitialization Xo D'N Vr -mm To,ibaa yields X(1)
- UI)
+ <>(() Em= o.(t)oru? where
(t) ofe 4407+22 - 20,1 - /4n - -0, C?
(15)
24 107 + 80, Sie /10 - Cx2.C²
A = diag(UfU - VoTV) and C. = =C(as, Tox) is a constant. Moreover, the ith
eigencomponent of x() converges to the ith eigencomponent Y at a rate 01177
Implicit acceleration. Recall from that the convergence rate for the non-overparrmeterried problem in (2) is O(e-t), which does not depend on
Table 2: Regimes where implicit acceler- ation is achieved (V) not achieved (x) depending on the data singular values and the level of imbalance An of the initial
the data or the initialization. It follows from (15) that the asymptotic behavior of the singular values of the over-parameterized solution is:
conditions (see (16)).
Balanced Au = 0
Imbalanced
Ày >
(t)(t) a 2A, + 2e -V (16) small a
X
which depends on both o; (singular values of the data) and A (level of initialization imbalance). parameterization leads to implicit acceleration when Ao X3,3 > 1. As illustrated in
large σ
It immediate see that over-
5

---

Under review as a conference paper at ICLR 2021
Table 2 and previously observed in the literature (aax al. 2014 acceleration is possible under a balanced initialization i.., λn, s when the singular values are sufficiently large. This is similar to the symmetric case, which is not surprising since a symmetric factorization is by construction balanced. A new phenomenon that emerges from our analysis is that acceleration also be achieved by increasing the level of imbalance, and that acceleration can always be achieved regardless of the data when Do,s > 1. Moreover, faster convergence can be achieved by using more imbalanced initialization.
4 DYNAMICS OF GRADIENT FLOW FOR ASYMMETRIC MATRIX FACTORIZATION WITHOUT SPECTRAL INITIALIZATION
We now relax the assumption of spectral initialization (Definition
B Defining the quantities
(0))
Σ
5= 1
Lm
0
R(t)= =
S
(17)
[
yr
0
2
In
one can immediately obtain from 13 and 14
the Riccati-like differential equation
R = SR 2 RR R SRAQ.
(18)
where from [14 we conchide that SRRSRR SR. = Ae. with R(O) Rg. However, in general.
one cannot back from (18) to [13 unless the conservation law (14) is explicitly imposed
for all times 4. The natural question then, when are they equivalent? Our next. result
provides the answer and additionally reveals interesting relation between
18 and
matrix factorization problem with explicit regularization (proof in Appendix F
is equivalent to Proposition 4 (Explicit regularization). The differential equation
18
t E UW --000- - VTV -Ag)):
(19)
À = -= + -VIUTO VTV - -gg))
This system corresponds to the dynamics of gradient. flow applied to the regularized problem VTV
min 0,V
[ UV - )
(20)
Morcover, if 2()) = (T()(()))))))))C
(L) obeys @(ta) AQo at some t tu- then
@(t) = AQ, for all t 2 to- In particular, if we initialize 19 -or equivalently (18) such
-1ch
that 2(0) 2 SRg Ag.. then the conservation law holds true for all t and both the
dynamics of 19 and 13 are the same.
Some remarks appropriate:
A solution to 13 implics a solution to 19h because when holds the 2nd terms on the RHS of (19) vanish, while the 1st terms are exactly 13 However, the converse is not necessarily ly true, unless (19 is initialized
the same way as 13 The above result relates implicit acceleration (as will be shown in Proposition 5 to an eaplicit regularization;
namely one can either select a particular initialization and solve an unregularized problem,
or start at an arbirary initialization and explicitly regularize.
Note the specific weight of 1/8 in (20) is special: If one replaces 1/8 by some constant a > 0, the gradient flow dynamics, 1.e. the analog (19), will not be equivalent to (18).
Note also problem
20 also appeared in (Du al. 2018a). but without such connections.
Equation [8] reminiscent of a Riccati differential equation due to the cubic term in R
(similar to the gradient flow in the symmetric case) but we believe that, in general, it cannot be solved exactly due to the last term. However, can be solved exactly in a particular case (proof in Appendix
Proposition 5 (Exact solution and convergence rate). If Af, = Apla for some constant Ao, then the differential equation (88 reduces to the following equation with closed form solution for R.T
i = 5R RR R,
((I)R(()) = 1 Ro (1+ 1 -))ook es
(21)
6

---

Under review as a conference paper at ICLR 2021
where $ S Xn DΣΦT, Ru R(O). Moreover, if S and I + 35- R, are invertible then R(t) converges exponentially to R,R,, defined as the projection of the matriz 25 on the PSD cone, S = + "/E|0". More precisely, if a square matria the convergence rate is O0--e((100 oo +X6), where Jmin is the smallest eigenvalue of Y, and otherwise the rote is O(e AAolt.
Note that (21 is the gradient flow for the symmetric factorization problem mina { 11 IR The particular case AQ, = AnLE is mathematically interesting because it is amenable to an analytical treatment. However, may not be realizable in practice because the conserved quantity 2o (or Ag,) must. have low rank, i.e.,
rank( Qo) rank(U Uo- Vo) S rank(Ug Un) + rank Vo Vo) S m + n.
(22)
Since rank(Ag = k, choosing Uo and Vo such that 0g Un - V Vo- = Ao is not gen- crally possible in an ovcr-paramcterized set- ting with & m On the other hand, the choice Ag. = An does not present a prob-
layer model
One lager model
(+) [1Qb]l 0.02
(b) [[Qi]]
1.83
[c) [|Qbll 196.55
Recomencctoon
lem if we consider the system (99 where we
have the freedom to choose any initialization.
The experiments in Section 5 illustrate that
[9), or equivalently (2I), is actually enough
1000
13.31
-0003
13.1
1330
to capture the general behaviour of 13
Singitr volue
NUMERICAL EXPERIMENTS
Imbalanced initialization. We now pro- vide numerical evidence supporting our the- oretical results. First, we generate ran- dom matrix Y with Y,, ~ NO, 1), and set m 5, n = 10 and k = 50. We approximate the dynamics of gradient flow for one-layer and two-layer linear models by using gra- dient descent. with a fixed small step size n 10-3 3 (smaller step sizes did not lead to discernible change). We evaluate the recon- struction error ||Y (f)lll/ [|,, where x(t) = Utt (t), and compare the evolu- tion of the singular values of (t).. We con-
2001
4331
2301 4331
23:01
4300
Iterations
Ieeations
Figure 1: Top row: reconstruction error for one- vs. two-layer linear models. Bottom row: evolution of singular values the solution. From left to right we use = 10-2, 2 a = 10-1, and a = 1 respectively.
Ao 111.22
(b) 10000,
[c] 200,
Recortruccion ena
sider Gaussian initializations i.e., Un and Vo
have entries ~ NN(0, o²), where a is varied to obtain different degrees of inbalance. In or- der for the two models to start from the same state, we choose X (0) VO for the one- layer case. The results are shown in Fig. 1 From our thcoretical analysis, a different behaviour for the convergence rate depending on whether the initializa- lion is balanced or imbalanced, i.e. whether [|]]] || P ||U] U-5VooV' Vy Volle is small or large, respectively. When it is very small (Fig. [p), the strength of the singular value dominates and we expect the components to be learned sequentially from the largest to the smallest, which agrees with the find- ings in (axx et al. 2019 Gidel et al.) 2019). As we make the weights more imbalanced (Fig. [b)) the singular values are learned faster, even the smaller ones. Finally,
(000)
.331
-311
2203
expect
2003
Figure 2: Top row: reconstruction error for the asymmetric factorization dynamics without reg- ularization in 13) and (14) and general Qo (red dashed line), versus the regularized dynamics in [9) with a diagonal Ago = doli (black solid line). Bottom row; evolution of the corresponding singu- lar values of the solutions. From left to right we set k k550, k= 100, and ==000, respectively.
7

---

Under review as a conference paper at ICLR 2021
Q becomes very large, the implicit acceleration becomes more prominent and the solution of the factorized problem converges significantly faster (Fig. [t). In other words, these numerical results are consistent with Proposition 2 and Proposition 5
Aa. = Apl is general enough. Since Proposition contains the case where an exact solution is available, want to investigate whether this is general enough to capture the qualitative behaviour of the general problem, i.c. the general dynamics of 13). To avoid confusion, let us refer to 0² and V as the variables of (19 , as well as A4. AoJL; here
stands for "identity." The variables U and V refer to the original dynamical system (13) with its conserved quantity Ag completely fixed by the initial conditions; see (1)). We want to show that it is possible to find an "optimal" A E I such that both cases have very close dynamics. Hence, we initialize Uo and Vo (and thus equivalently Uo and Vo) with entries sampled from [100, o²), and we fix a = 1. The same initial condition is also used for (19,, i.e. Uyf = Uo and Vol = Vo. We set 7j 10 5- Y ~ (0 m 5, n 10 and vary k. Thus, look for a A that minimizes the error [tx - xt) F. In Fig. 2 we illustrate that, indeed, this can be done. Note that in these three cases the evolution of both dynamical systems are nearly indistinguishable.
Extension to nonlinear networks. Our analysis so far has shown that the acceleration phenomenon is a result of imbal- ance and is induced through a conservation law for which the definition is expected to change when introducing nonlinearities. In fact, both the architecture of the network and the objective function should affect the invariances. As such, conducting a full analysis of the symmetries and conservation laws for more complex nonlinear networks is necessary to characterize the implicit bias in such cases, which we leave for future work.
[a] Balanced
One layer
Two layers
Trinkin
400
2x-2
(b) libal
Nonetheless, provide numerical evidence in more realistic setting than the one of linear networks, where we only add a nonlinearity (sigmoid) to the final layer in order to preserve some of the structure of the linear model. We train the two networks (one layer vs two layers) on synthetic data, i.e. we compare the dynamics of gradient descent for the objectives: (W)) =1l (XIV-( and (2(U, V) +|lY~(XUVT)) where e is the sigmoid function, X Raxd and Y E Raxd represent the training samples and labels respectively (n= 10 = 10), W E RRdd U, V E B are the weight matrices and the width is k- 100. We generated the matrices W* and X with entries drawn from // 1) and Y = @XW*)++ whhere 10 "N I). Our results shown in figure 3 interestingly
100
200
300
400
500
Iterations
Figure 3: Evolution of the training loss for nonlinear one- layer and two-layer models. Top row: [Qol12 0. Bot- tom row: ||Qollz = 4.6. Initial weights are drawn from a nor- mal distribution (00 10-1).
suggest that our conclusions about the role of imbalance hold in this case as well.
6 CONCLUSION
We considered the gradient flow dynamics for two-layer lincar neural networks, providing an analytical treatment to great level of detail. Our results establish a detailed characterization of the "implicit acceleration" phenomenon in this case, without assuming balanced or vanishingly small initializations, which so far have been present in all prior work in this vein. More specifically, our analysis shows that the implicit acceleration of gradient flow is a consequence of an emerging rotational symmetry induced by vee--araaeterizattinn and giving rise to several conservation laws that constrain the dynamics to follow specific trajectories. Moreover, conserved quantities are completely fixed by the initialization, which has profound effects on the convergence rate. Although our analysis focuses on the simple case of linear networks, it reveals a potential key to understanding implicit bias which lies in the conservation laws that arise from the symmetries of the problem. These symmetries depend on the network architecture, objective, optimization algorithm, and they constrain the dynamics to an invariant manifold that encapsulates the implicit regularization and acceleration effects. Understanding this in more complex models may thus be reduced to finding dynamical invariants, for which our results provide foundational starting point.
8

---

Under review as a conference paper at ICLR 2021
REFERENCES
Sanjeev Arora, Nadav Cohen, and Elad Hazan. On the optimization of deep networks: Implicit acceleration by overparameterization. In International Conference on Machine Learning (ICML), 2018.
Sanjeev Arora, Nadav Cohen, Wei Hu, and Yuping Luo. Implicit regularization in deep matrix factorization. In Advances in Neural Information Processing Systems (NeurIPS). 2019a.
Sanjeev Arora, Noah Golowich, Nadav Cohen, and Wei Hu, A convergence analysis of gradient descent for deep linear neural networks. In 7th International Conference on Learning Representations, ICLR 2019, 2019b.
F. M. Callier, J. Winkin, and J. L. Willems. On the exponential convergence of the time- invariant matrix Riccati differential equation. In Proceedings of the 31st IEEE Conference Decision and Control, pp. 1536 1537 vol.2, 1992.
F. M. Callier, J.
Winkin, and J. L.. Willems. Convergence the time-invariant Riccati differential equation and le-problem: mechanisms of attraction. International Journal of Control, (9): 1000, 1994.
Jacopo Cavazza, Benjamin D Haeffele, Connor Lane, Pietro Morerio, Vittorio Murino, and Rene Vidal. Dropout as low-rank regularizer for matrix factorization. In Artificial Intelligence and Statistics (AISTATS), volume 84, pp. 435-444, 2018.
Lenaic Chizat and Francis Bach. Implicit bias of gradient descent for wide two-layer neural networks trained with the logistic loss. In arXiv preprint WxXiu:200.004480 2020.
Simon Du and Wei Hu. Width provably matters in optimization for deep linear neural networks. In International Conference on Machine Learning, pp. 1655-1664, 2019.
Simon S Du, Wei Hu, and Jason D Lee. Algorithmic regularization in learning deep homo- geneous models: Layers are automatically balanced. In Advances in Neural Information Processing Systems (NewrIPS), 2018a.
Simon Du, Xiyu Zhai, Barnabas Poczos, and Aarti Singh. Gradient descent provably optimizes weerpparaeterrized neural networks. In International Conference on Learning Representations, 2018b.
G. Franca, M. Jordan, and R. Vidal. On dissipative symplectic integration with applications to gradient-based optimization. orXiv:2004.065 06844 /math. OC/, 2020.
Kenji Fukumizu. Effect of batch learning in multilayer neural networks, 1998.
Gauthier Gidel, Francis Bach, and Simon Lacoste-Julien. Implicit regularization of discrete gradient dynamics in linear neural networks. In Advances in Neural Information Processing Systems (NeurIPS), 2019.
Suriya Gunasekar, Blake Woodworth, Srinadh Bhojanapalli, Behnam Neyshabur, and Nathan Srebro. Implicit regularization in matrix factorization. In Proceedings of the 31sl. Interna- lional Conference on Newral Information Processing Systems (NeurIPS), 2017.
Suriya Gunasekar, Jason Lee, Daniel Soudry, and Nathan Srebro. Characterizing implicit bias in terms of optimization geometry. In Proceedings of the 35th International Conference on Machine Learning (TCML), pp. 1832-1841, 201 20.8a.
Suriya Gunasekar, Jason D Lee, Daniel Soudry, and Nati Srebro. Implieit bias of gradient descent on linear convolutional networks. In Advances in Neural Information Processing Systems 31 (NeurIPS), pp. 9461 14719 2018b.
Ziwei Ji and Matus Telgarsky. Gradient descent aligns the layers deep linear networks. In International Conference Learning Representations (ICLR), 2019a.

---

Under review as a conference paper at ICLR 2021
Ziwei Ji and Matus Telgarsky. The implicit bias of gradient descent nonseparable data. In Proceedings of the Thirtt-Second Conference on Learning Theory (COLT), pp. 1772-1798, 2019b.
Yuanzhi Li, Tengyu Ma,, and Hongyang Zhang. Algorithmic regularization in over- parameterized matrix sensing and neural networks with quadratic activations. In Proceed- ings of the 31st Conference On Learning Theory (COLT), 2018.
B. P. Molinari. The time-invariant linaar-uadrratc optimal control problem. Automatica, 13 (4):347 357, July 1977.
Ambar Pal, Connor Lane, René Vidal, and Benjamin D Hacffele. On the regularization propertics of structured dropout. In IEEE Conference on Computer Vision and Pattern Recognition, pp. 7671-7679, 2020.
Takao Sasagawa. On the finite escape phenomena for matrix Riccati equations. IEEE Transactions on Automatic Control, 1982.
Andrew M. Saxe, James L. McClelland, and Surya Ganguli. Exact solutions to the nonlinear dynamics of learning deep linear neural networks. In International Conference on Learning Representations (ICLR), 2014.
Andrew M. Saxe, James L. McClelland, and Surya Ganguli. A mathematical theory of semantic development in deep neural networks. Proceedings of the National Academy of Sciences, 116(23):1 153377111577-11466 2019. ISSN 0027-8424. doi: 10.1073/ /.a00/...s.... ppappppsss... mmas.1820226116
Daniel Soudry, Elad Hoffer, and Nathan Srebro. The implicit bias of gradient descent on separable data, In International Conference on Learning Representations (ICLR), 2018.
10

---

Under review as a conference paper at ICLR 2021
A SOLUTION TO THE MATRIX RICCATI DIFFERENTIAL EQUATION
Here we prove Proposition
First, we note that general matrix Riccati differential equation takes the form
P(t) AP(t) P(t)4" P())R(()) QQ,
P(0) Po.
(23)
where P(t) E R"*", and A, R. Q, Po E R are constant matrices. Associated to 23 one has the linear system
(t)
-
A
AT
R
(1
X1 [2
I Po
X2(()
Q
A
[2(1)
-
(24)
The closed form solution of (3 follows as consequence of lemma G below.
Lemma 6 (Sasagawa 1982)
Consider the two initial value problems 23 and 24
We
have:
The invitial value problem 23 has a solution in the interval [0, 4 if and only the matric (t) in the solution of the linear differential equation (24 is inverlåble for
all t [0,tt] ti). Moreover, the solution to (23 unique and given by
P(() X
(25)
Let P be a solution to the algebraic Riccati equation (ARE)
AP PA - PRP++=0.
(26)
Then the solution of 24 is given by 27 below, where À A- - PR. and À
A -RP:
X1(1) X2(O)
e-A (S Us (Po
(Po-P) - P)
(27)
Pe Pe- LA + P (fo (tee (Po PPo- ( P)
Now, we apply Lemma 6 to the Riccati equation induced by symmetric factorization, namely,
x((t) = 2x + 2Y 4
(28)
The associated linear system is
(t)
-Y
41
Xi(t) X2(1)
X1 X2
In
(29)
(t)
0
2Y
X
and the algebraic Riccati equation is given by
XY + YX 2x =
(30)
This equation admits the trivial solution X = 0. Thus, one can verify that as long the given matrix Y is invertible we have
X1 X2
21Y
+ e 20Y (- 1) 221 X
(31)
Therefore, the unique solution to
5 is explicitly given by
X(t) = (t) 2 Xo (e
21Y
(32)
as long as the matrix between parenthesis is invertible. Note that ( - -I) is always positive semidefinite for 2 O. Thus, for all Xn ž= 0, the matrix + (c - I) is
invertible and the solution of the Riccati equation is well defined for all 1, 0.
11

---

Under review as a conference paper at ICLR 2021
B PROOF OF COROLLARY J
The first part follows trivially by substituting (6) and verifying the invertibility of the matrix in parenthesis. For the second part, note from 7 that if a, > 0, then a,(t) - as t + oo at a rate O(e D(e-4tr)) and if oi < O, then (i(t)
at rate O(e4t.)). Therefore, E() + max(2,0) 0) and (1) + max(E, 0)6" at a rate O(e O(e-Homm(()) ())m(())) as claimed.
C RATE OF CONVERGENCE IN THE SYMMETRIC CASE
In this section, we show that the solution of the Riccati equation (t) 2Y (t) (1)² converges exponentially to at rate O(e (Y))m))))) where 0 min (Y) min; la, is the smallest eigenvalue of Y in magnitude. We already know that if ( converges to X,, then
X, positive semidefinite if Xo ł and rank(X.)) S raan((oo)
X, is solution to the algebraic Riccati equation 2Y 2XY 4X 0.
Our proof of proposition 2 is inspired by the proofs in (Molinari 1977 Callier et al. 1992
1994) which studied the solution of the Riccati equation under more general setting. The stralegy will be to show that. the algebraic Riceati equation has a unique PSD solution X such that the cigenvalues
Á = Å = 2Y 4X+ have negative real parts. Such solution is usually referred to the strong solution stabilizing solution in optimal control literature because it is the only solution of the algebraic Riccati equation such that the matrix A is exponentially stable, iee exp(tA) converges to O. The stability of the matrix Á is important because appears in the solution of the Riccati equation as can be observed in (27).
We start by proving that y solution to the algebraic Riccati equation, iee it is a critical point of the problem.
Due to the symmetric and positive semidefinite nature of the matrices and Y, in our case the algebraic equation (30p can be reduced to
x(X - Y) 0.
(33)
For X y, we thus have
Y(Ý -Y- = Σ maxfo, 07
Σ (maxfor, 0} o:)0107
(34)
m
max
min min.or,0)0:0 030:0!
The first sum contains only the positive eigenvalues while the second sum contains only the negative ones. Therefore, no eigenvalue will appear in both sums and using the orthogonality of the vectors o; we conclude that YY) =
Next, we prove that P = V is the unique symmetric positive semidefinite solution of the algebraic Riccati equation such that the eigenvalues of Á = A - RP have negative real parts. Note that Á = 2Y - 4 = 2Ý where Y = -mmlaaoo
> 0.
We proceed by contradiction. Let X2 be a PSD solution of the (ARE) such that the cigenvalucs of Å 2Y 4X, have negative real parts and X,
Now consider A y X2. By straightforward calculation, we can show that A is
solution of
AA + AA + 4A =.
(35)
Since A is not necessarily invertible, we consider a basis such that
0 0
A Z
n" ZÂZF.
(36)
D
12

---

Under review as a conference paper at ICLR 2021
where D is invertible. Note that our proof holds and is more trivial if A is invertible. We write
W1
W2
n" = ZÄZ
(37)
W3
W4
After the change of basis, equation 35 becomes
Aâ + AA + 4 = 0.
(38)
From
35
can deduce the following for the block matrices:
W₂ 0,
W 0, 0,
DW --4 0.
Note that the last equation is similar to equation (35) for the invertible block D. Morcover, in the new basis, Á is a block diagonal matrix. Using the change of variable T D-1, we obtain the Lyapunov equation
TW₂ + WGT =/
(39)
Since A= -2 < 0 is invertible, the block Wa is also invertible and its eigenvalues are a subset of the eigenvalues of -2Y. As a result, the Lyapunov equation has the unique trivial solution T = 2 0, therefore D - > 0.
Using the solution of 35 we can obtain a new derivation for Åς as follows
À 2Y
4X2 = 2Y 48
+4A
W 0
Ä 44 z
ZT
W4
with W1 < 0 and - W4 > 0. This contradicts with the initial assumption on the eigenvalues
of Ag having negative real parts.
Now we can use Lemma 6 with P= V to obtain a new expression of the solution. We have
Á A 28 and
(1) [2(()
20 - [T c - 2Y] + 218
P]]
(40)
(Xo -
Therefore, the solution of the Riccati equation is also given by
(()) = (1))) +e 218 (Xo
No 8)
-2
(41)
Note that the inverse exists for l when X, because we have proven the existence of the solution using the previous expression detailed proof).
32 (see lemma 2 in (Sasagawa 1982 for
We introduce the function H(t) = (Xo - + ( 8)]
Thus,
- - = 2 l(tee 218
(42)
The function H has the following properties for t 2 0:
It is decreasing:
dH(t) - <0
(43)
If I - Y) is invertible then
lim (tt) (Xo 1 - - H.
(44)
H is bounded on R+;
f S H() S H(0) = Xo V
(45)
Therefore, can conclude that there exists constant C> 0 such that [|X(t) r < Ce mim
(46)
where Fmin min(lax]} is the smallest eigenvalue of Y in absolute value.
13

---

Under review as a conference paper at ICLR 2021
D CLOSED FORM SOLUTION FOR THE ASYMMETRIC CASE UNDER SPECTRAL INITIALIZATION
Under the spectral initialization, Uo and Vo diagonal, hence are (() and ( As result, (t) and (t) remain diagonal for all t 2 O, because the components of [3) can be decoupled and the evolution of (13) will induce no change in the off-diagonal elements. To
see this, observe that
(o,
== Un
- U.)
(1sism) m).
(47)
whereas the off- off-diagnaal terms obey U, and V. = 0 Thus 44 describes the evolution of the singular values of the solution. This decouples the problem into a set of independent one-dimensional equations, therefore it suffices to consider the scalar system
û = (a ūo)e,
ē= (a ut
(48)
where we drop the index i = 1,..., m for simplicity. In (Saxe et al. 2019) there is a strong assumption, i.e. (0) (0) for all which balanced intaalztttin (Definition 2 Here we solve 477 without such an assumption. From (48), it is immediate that the conservation law [1 becomes d (a² - 52) 0. Trajectories (ū(t). B(t)) are thus constrained to lie on hyperbolas defined by
dt
(1) o²(t)
16 -t)
A const.
(49)
Since we are mostly interested in the behavior of the product x(() a())((() by making explicit use of the conservation law (49). specifically XB ãt 22 = ãt +5 2xx,, and
(ū 72 AB + 4x², obtain the tollowing first-order differential equation for z()):
b= (a - ne 2² + (a - (a-
= (a ūt) u2)2 2(0 x) x² + X3/4.
(50)
Even though this is nonlinear differential equation, it is separable, thus integrating both sides yields precisely (15) (we restore i and - a; represents the corresponding component associated with singular value o, and conserved quantity An,i), where C'> > is constant:
V 20 M 2602z 2 + XØ 4 40.o - A
==
(51)
427(0 20)
Above, only m out of k 2 conserved quantitics are used. Hence, there is degeneracy in the solution and only m effective degrees freedom regardless how large is. Note that l < m sunder--arameterized case) then (47) becomes under-determinnd.
E RATE OF CONVERGENCE IN THE ASYMMETRIC CASE
Here we prove Proposition 5 First, note that equation -= ŠR 1 RRF R follows directly
from R SR SRR R SRAO, when Ag, =XI.. This leads to Riccati differential equation for P(t) (t))
dP(t) RRT + RĀT = (ŠR- HRR R) R² + (SR-ARR"RR)
dt
(52)
SP(() + -(1)5-(PP))
with initial conditions P(O) Ro By the arguments used in Appendix A the exact solution of the above equation is given by;
(RR' (() e Ro Ro 1 - 18
(53)
One can show that RR converges exponentially to the matrix R, RT, defined as the projection 25 on the positive semidefinte cone; this is derived by the same arguments leading to Proposition 2 (see Appendix C The convergence rate depends on the eigenvalues
14

---

Under review as a conference paper at ICLR 2021
of 25 which can be determined using Schur's complement formula. For block matrix M
one has
det = det(D) det()) -0-0-AAA--CC BD)
(54)
C
Thus, using the above formula for S - AI = (A0/2-A)1m
Σ
we have
8"
-(An/2+2+)
det(S - X1) (-1) (A0 (-1mttt) \)mー"
(o] + X - A2).
(55)
i=1
Therefore, the cigenvalues of Š are:
2m eigenvalues 13, where 8, 2 + X (i - m).
2. The eigenvalue -2//2 of multiplicity least m.
In general the smallest magnitude eigenvalue will be |ool/2. However, if Y is square matrix then only the eigenvalues =5, above will be present.
F FROM IMPLICIT ACCELERATION TO EXPLICIT REGULARIZATION
We now prove Proposition Replacing 17 into immediately yields 19 It is also straightforward to verify that applying the gradient flow, _0- - and V -ve, with f being the objective function in (20), yields 19 Next, consider the formal Taylor series
((t) 2(to) + (t to) Q(to) + (t to) ((to)) +....
(56)
Define P(() 0 (1) (1)) From 19] one obtains
0() = - Ae.) PP()) - (()(((()) Aqa)
(57)
Since this is "linear" in (인 ). higher order derivatives take the form
2(1) Σ (((1 A.))
(58)
where the functions Z 's and Wt's contain sum of powers and time derivatives of P(t). For instance, the second order derivative yields
-= = - +2P) - (P² + 2P) (Q - -
- 1 AQ.))(() 1 -())() -/000---0 Ag.) (00(0)=-10000)
Therefore, if @(4)) Aoo at 1 = to then all derivatives 58 vanish identically. As consequence, the expansion 56 implies Q(t) Q(to) Ac for any other 2 to as well.
G ON MATCHING RATES BETWEEN GRADIENT FLOW AND GRADIENT DESCENT
We provide an explicit example to illustrate why studying the continuous- time dynamics of the gradient flow is expected to reproduce the behaviour of its discretization, i.c. gradient descent. For the sake of simplicity let us limit the discussion to the case considered in Section 2 and in the scalar case, Y different. algorithms, namely gradient descent applied to problem (2) versus GD applied to the factorized problem (3). We thus have
R. What we would really like to do is compare two
X = Nx + n100 Xx)
(59)
versus
Ust1 Un 2nfl U.)
(60)
15

---

Under review as a conference paper at ICLR 2021
The respective contimuous-timg time limits of these algorithms
= (=
(61)
versus
2 - U2)U,
(62)
for X Xt) and U (t). The solution of (61) is given by X(t) 5 (Xp Y)e
yielding rate O Let us introduce the perturbed variable X2 X, - a. Hence (59
readily gives Xk+1 = (1 - n7) Xx, thus a matching rate of O(e Oa.-*)). (This example trivial because both systems linear.) Now let us move to the more interesting nonlinear case.
Consider 61] and let ((t) (t). Thus X X 4X2, whose solution is
((t) 1 -e 4at
cc co 4at
(63)
which implies a continuous-tinn rate of ( (et)). (Compare this with the last phrase in
Corrolary L.) Now let see what happens for 60). Note that this is complicated nonlinear recurrence relation. Fortunately, we can solve it approximately. By introducing X = U this recurrence becomes
Xk+1 X5 4n/(0 XR) 412( Xx)
(64)
Now we do the following trick. Let X2 = X -a a be a perturbation from equilibrium. For a small enough 1) we can neglect. terms of O hence Net1 s Xx - 4n + X&) s
(1 - 4njo Xx. This implies that & + 0, or equivalently X + a, at a discrete-time rate of e- e4aank This matches the contimouss-time rate.
It is not hard to see how such nonlinear recurrence relations quickly become intractable for more complicated problems. On the other hand, even though the contimuous-ting time limit provided by the gradient flow consists a nonlinear ODE, the analysis much more feasible, besides introducing connections with interesting areas of mathematics such as Riccati differential equations.
16


Reviews:
[{"review_id": "b7S7YqsqdUO", "reviewer": "AnonReviewer2", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "Summary of review:\n\nThis paper provides a detailed analysis of gradient flow in (over-parametrized) two-layer linear neural networks. The main results state the precise dynamics of gradient flow for both symmetric and asymmetric matrix factorization, starting from certain spectral initialization. One novel insight that stems from the analysis is that for asymmetric matrix factorization, \"imbalanced initializations\", where the left and right singular values of the iterates differ, converges faster than \"balanced initializations\". Simulations further validate this insight.\n\nSetting:\n\n(i) Symmetric matrix factorization: Given a symmetric matrix Y, the problem is to solve a mean squared loss between Y and UU^T to factorize Y, where U is a (possibly over-parametrized) variable matrix.\n\n(ii) Asymmetric matrix factorization: The asymmetric setting considers an asymmetric matrix Y and the asymmetric factorization of UV^T for factorizing Y.\n\nResults\n\n(i) This paper focuses on the convergence of the gradient flow of U for minimizing the mean squared loss, starting from spectral initializations. Informally, a spectral initialization has the same eigenspace as Y.\n\nFor these spectral initializations, the gradient flow paths essentially become coordinate-wise updates over every singular value. Then, the authors went on to state the precise dynamics of gradient flow for every singular value. The results imply that \"large\" singular values (of Y) converge faster than \"small\" singular values (of Y) in gradient flow.\n\n(ii) This paper begins by studying spectral initializations, where $||U^TU - V^TV||_F$ is small, then discusses how to generalize their result to non-spectral initializations.\n\nFor spectral initializations, this paper observes crucially that $U^TU - V^TV$ is preserved throughout gradient flow. Hence if this quantity starts small, it will remain small throughout gradient flow. Furthermore, for \"large\" singular values of $U^TU - V^TV$, the results imply that these singular values converge faster than \"small\" singular values of $U^TU - V^TV$.\n\nFor non-spectral initializations, this paper observes that $U^TU - V^TV$ is still preserved during gradient flow, but this quantity now depends on how \"balanced\" the initializations of U, V are.\n\nCriticism:\n\nIt would help improve my understanding of this paper if the authors clarify the following questions.\n\n(a) The \"acceleration\" claim of this paper comes from comparing the results to a linear model baseline. However, I am not completely sold on this comparison. For example, what would the results imply if compared to properly parametrized U (and V)? Is there any provable advantage of over-parametrization in this setting?\n\n(b) The result of Corollary 1 for symmetric matrix factorization, where larger eigenvalues converge faster than smaller eigenvalues, also appears in linear regression. In particular, the gradient flow of linear regression also shows similar patterns, where larger eigenvalues (of the sample covariance matrix) converge faster than smaller eigenvalues. Therefore, in my opinion, the results in Section 3 and 4 seem more novel compared to the results in Section 2.\n\n(c) While the results state fairly precise dynamics of gradient flow, how well would they extend to settings with sampling errors? For example, what about matrix completion? Would your claim regarding \"imbalanced initializations\" still hold?\n\nWriting\n\nOverall, this paper is well-written and easy to follow. Although the paper would be easier to read if it is less dense. I do not quite understand the claims in Figure 2. You explained in Section 4 that you need $K \\le m + n$ for identity initializations and Proposition 5, but here you say that k ranges from 50 to 200? Please clarify.", "extra_scores": {}, "extra_text": {}}, {"review_id": "EzZ0Ro8mEN7", "reviewer": "AnonReviewer3", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "This paper considers the gradient flow dynamic for two-player linear neural networks. In detail, it studies the implicit acceleration of gradient flow brought by overparameterization and shows the reason for implicit acceleration is the existence of conservation law. It studies the convergence for gradient flow under both balanced or imbalanced linear networks, and with spectral or non-spectral initialization. Compared with previous work, this work is the first to provide an explicit characterization of the gradient flow with respect to their eigenvalues. Experiment results suggest that such an implicit acceleration indeed exists. \n\nHere are my detailed comments.\n\n- Page 2, (2): why call (2) ‘a symmetric one-layer linear model’? Does that suggest $m = n$? The same issue holds for (3). Is $Y$ in (3) the same as that in (2)?\n\n- An interesting observation is that under the spectral initialization, for the symmetric case, the convergence rate of the eigenvalues does not depend on the initial value of $X_0$ (it is $e^{-4t|\\sigma_i|}$). However, for the asymmetric case, the convergence rate is $e^{-t\\sqrt{4\\sigma_i^2 + \\lambda_0^2}}$, which explicitly depends on the initial matrix $X_0$. Can the authors make more comments about that?\n\n- The ‘acceleration’ compared with original gradient flow over $X$ suggests that by a matrix factorization, the convergence of eigenvalues varies may be accelerated according to the eigenvalues over the data matrix $Y$. Is it true that in the worst case, say for some $i$, $\\sigma_i \\approx 0$, then the convergence speed of $\\sigma_i(t)$ will decrease to 0? In that case, can we still call it ‘acceleration’?\n\n- In Proposition 5, the authors build the convergence results for the case $\\Lambda_{Q_0} = \\lambda_0 I_k$ for the general non-spectral initialization case. However, I hardly can find a non-spectral initialization case satisfying such a condition. Can the authors provide some examples?\n\n\n", "extra_scores": {}, "extra_text": {}}, {"review_id": "9bhtVOMNk4Q", "reviewer": "AnonReviewer1", "summary": "", "questions": "", "limitations": "", "rating": 7, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "1.  Paper Summary \n\nThis work analyzes the implicit acceleration of gradient flow for over-parameterized 2 layer networks (i.e. 1 hidden layer networks) used for matrix factorization.  By presenting a novel analysis connecting the gradient flow to Riccati type differential equations, this work demonstrates that imbalanced initializations can lead to acceleration.  The authors present convergence rates for symmetric and asymmetric matrix factorization under both spectral and non-spectral initializations.  The convergence rates for these settings are indeed faster than those in the linear model.  The authors lastly provide empirical results to support their theory.  \n \n\n######################################################################\n\n2. Strengths\n\n2.1.  The connection between gradient flow and Riccati type differential equations is novel to the best of my knowledge and provides a simpler and clearer means of understanding implicit acceleration in over-parameterized models than prior works.  \n\n2.2. The paper is written clearly and is easy to follow.  The authors present examples of acceleration under the simpler setting of spectral initialization before extending to the more nontrivial case of non-spectral initializations.  \n \n\n######################################################################\n\n\n3. Limitations\n\n3.1.  I believe the authors are missing some references to related work.  In particular, below are some related works:\n\n(1) https://arxiv.org/pdf/1810.02281.pdf - This work extends the notion of zero balancedness considered in the work by Arora et al. 2018 to the case of approximate balancedness.  The analysis presented in this work is for deep linear networks (more than 2 layers) and for gradient descent.  Technically, I believe this work does fall under the imbalanced case and also yields a fast convergence rate.  \n\n(2) https://arxiv.org/abs/2003.06340 - This work analyzes spectral initialization under gradient descent in deep linear networks (of arbitrary layer structure).  In particular, this work also demonstrates linear (fast) convergence for gradient descent under spectral initialization (Proposition 2).  \n\nI believe it is important for the authors to position their work with respect to the above works, but I do not feel that missing these references diminishes the novelty of the submitted work. \n \n3.2.  While the connection to Riccati type differential equations is interesting, it would be nice if the authors could discuss how similar analysis could be used for deep linear networks (at the moment there does not appear to be an obvious extension).  \n\n3.3.  (Minor) The authors briefly mention that discrete time convergence rates can be derived from continuous time counterparts using symplectic integrators in the introduction, but it would be nice if there were a more rigorous connection between the related work and the rates presented in the submitted work.  \n\n\n######################################################################\n\n4. Score and Rationale\n\nMy recommendation is to accept the paper.  I believe the main strength of the paper was in providing a well-presented, rigorous, and novel analysis for understanding acceleration in over-parameterized matrix factorization.  The connection to Riccati type differential equations and identifying dynamical invariants presents an interesting alternative means of gaining intuition around implicit acceleration in over-parameterized networks.  \n\n\n######################################################################\n\n5. Comments\n\n5.1.  I believe the Y in equation 3 should be m x m instead of m x n for the symmetric case. \n\n5.2.  I feel that the notation could be made a bit easier to follow in Section 3 & 4. In particular, I believe bar(U), bar(V) are overloaded to represent diagonal matrices in definition 3, but these quantities are assumed to be non-diagonal for the remainder of the work.  As this is an important point, I feel that it could be emphasized a bit more.  \n\n5.3.  I feel that there could be more intermittent references to prior results throughout the work.  For example, the symmetric matrix factorization problem is discussed extensively in https://papers.nips.cc/paper/7195-implicit-regularization-in-matrix-factorization.pdf.  Similarly, invariance and convergence rate for gradient descent and under spectral initialization is discussed in https://arxiv.org/abs/2003.06340.  \n\n5.4.  There is an important distinction between the matrix factorization problems considered in this work and prior work.  Namely, the matrix Y has fully observed entries whereas in some prior works, Y is not completely observed.  I think the analysis may get a bit more tricky for the case of unobserved entries, but the authors could maybe point this out.  ", "extra_scores": {}, "extra_text": {}}, {"review_id": "3XfTUXgqATh", "reviewer": "AnonReviewer4", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "This paper studies the implicit acceleration of gradient flow for training a two-layer linear model. Compared with the one-layer linear model, the authors show that gradient flow over an overparameterized two-layer linear model may achieve a faster convergence rate, given a nice data spectrum and proper initialization. Moreover, the authors investigate the convergence of gradient flow with an arbitrary initialization and show its connection to Riccati differential equations as well as the explicit regularization. Overall the idea is clearly presented, the experimental results also well back up the theory.\n\nMy main concern is that it may not be fair to compare the convergence rate in terms of gradient flow. For example, you can simply reparameterize the parameter by X -> 2X, and the acceleration can also be achieved from the perspective of gradient flow. I think in order to fairly compare the convergence between different parameterizations/initializations, gradient descent is a better choice. Back to the example of X->2X, in this case, the smoothness parameter will become larger, finally one can observe that the convergence rate of GD under this parameterization will remain unchanged. \n\nRegarding the convergence results, the authors still require strong conditions (parameter and data can be diagonalized simultaneously) on the initialization to prove the convergence of gradient flow. What happens if considering more general assumptions on the initialization, such as the random/orthogonal initialization used in the following papers? The authors may also need to compare the convergence rates of the derived results and those in the following papers.\n\n[1] Du, Simon S., and Wei Hu. \"Width provably matters in optimization for deep linear neural networks.\" arXiv preprint arXiv:1901.08572 (2019). \n\n[2] Wu, Lei, Qingcan Wang, and Chao Ma. \"Global convergence of gradient descent for deep linear residual networks.\" Advances in Neural Information Processing Systems. 2019. \n\n[3] Zou, Difan, Philip M. Long, and Quanquan Gu. \"On the Global Convergence of Training Deep Linear ResNets\". International Conference on Learning Representations.\n\n[4] Hu, Wei, Lechao Xiao, and Jeffrey Pennington. \"Provable Benefit of Orthogonal Initialization in Optimizing Deep Linear Networks.\" International Conference on Learning Representations.\n\nThe regularization term in Eq. (20) seems pretty similar to the commonly-used regularization in solving low-rank matrix problems [5,6] for balancing (but they use 1/16 weight), could you please briefly discuss their connections?\n\n[5] Tu, Stephen, et al. \"Low-rank solutions of linear matrix equations via Procrustes flow.\" International Conference on Machine Learning. PMLR, 2016.\n\n[6] Wang, Lingxiao, Xiao Zhang, and Quanquan Gu. \"A unified computational and statistical framework for nonconvex low-rank matrix estimation.\" Artificial Intelligence and Statistics. 2017.\n\nIs the overparameterization necessary? From the presented theorems, I do not see whether the derived results have a dependency on the dimension of the matrix U or V. The authors may clearly specify why one needs to consider overparameterization linear models. For example, assume the data matrix has rank r, can the theory in this paper still hold if set the dimension of U or V as d*r?\n\n=========== after reading rebuttal ===============\n\nI do not agree with the author's response to my first comment. Considering parameterization $X = 2U$ (here I use the notation U to avoid the confusion), then the loss function becomes $L(U) = \\frac{1}{2}\\\\|Y-2U\\\\| _2^2$ and the gradient flow with respect to $Z$ writes $\\dot {U}= -2(2Z-Y)$. Then we have $\\dot {X} = 2\\dot{U} = -4(2U-Y) = -4(X-Y)$, which gives a rate $O(e^{-4t})$, which is faster than the $O(e^{-t})$ rate achieved without using this parameterization. Therefore, comparing the convergence rate in terms of gradient flow may still not be fair and valid. \n\nBased on this issue I would like to keep my score.", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: This paper studies the implicit acceleration of gradient flow in over-parameterized two-layer linear models. The authors show that the amount of acceleration depends on the spectrum of the data without assuming small, balanced, or spectral initialization for the weights, and establish interesting connections between matrix factorization and Riccati differential equations.  While this paper provides some interesting results regarding implicit acceleration in training linear neural networks, the reviewers raised quite a few questions and concerns about some claims made in the paper, as well as an inadequate comparison with previous work. Even after the author's response and reviewer discussion,  the reviewers' doubts are still not completely cleared away. I feel the current form of the paper is slightly below the bar of acceptance, and encourage the authors to carefully address reviewers' comments in the revision.
Decision: Reject