Understanding Incremental Learning with Closed-form Solution to Gradient Flow on Overparamerterized Matrix Factorization

TL;DR

This paper provides a theoretical analysis of incremental learning in gradient flow on overparameterized matrix factorization, revealing how small initialization scales induce a sequential learning process of singular values.

Contribution

It introduces a closed-form solution for gradient flow dynamics, explaining the emergence of incremental learning through time-scale separation in symmetric matrix factorization.

Findings

Incremental learning arises from time-scale separation among different components.

Decreasing initialization scale enhances low-rank approximation capabilities.

Analysis suggests potential extensions to asymmetric matrix factorization.

Abstract

Many theoretical studies on neural networks attribute their excellent empirical performance to the implicit bias or regularization induced by first-order optimization algorithms when training networks under certain initialization assumptions. One example is the incremental learning phenomenon in gradient flow (GF) on an overparamerterized matrix factorization problem with small initialization: GF learns a target matrix by sequentially learning its singular values in decreasing order of magnitude over time. In this paper, we develop a quantitative understanding of this incremental learning behavior for GF on the symmetric matrix factorization problem, using its closed-form solution obtained by solving a Riccati-like matrix differential equation. We show that incremental learning emerges from some time-scale separation among dynamics corresponding to learning different components in the target matrix. By decreasing the initialization scale, these time-scale separations become more prominent, allowing one to find low-rank approximations of the target matrix. Lastly, we discuss the possible avenues for extending this analysis to asymmetric matrix factorization problems.

Equations

$$\underset{U\in {\mathbb{R}}^{n\times r}}{\min }\mathcal{L}(U)=\frac{1}{4}\Vert Y-U{U}^{\top }{\Vert }_F^2,$$

$$\dot{U}=-{\nabla }_U\mathcal{L}(U)=(Y-U{U}^{\top })U,U(0)=U_0,$$

$$\dot{W}=\dot{U}{U}^{\top }+U{\dot{U}}^{\top }=WY+YW-2W^2,W(0)=U_0{U}_0^{\top }.$$

$$\dot{W}=WY+YW-2W^2,W(0)=W_0⪰0$$

$$W(t)=\Phi S(t){\tilde{W}}_0{\left(I_n+G(t){\tilde{W}}_0\right)}^{-1}{S}^{\top }(t){\Phi }^{\top },$$

$$G(t)=\left[\begin{array}{cc}{\Sigma }_Y^{-1}(e^{2{\Sigma }_{Y}t}-I_K)& 0\\ 0& 2I_{n-K}t\end{array}\right],S(t)=\left[\begin{array}{cc}{e}^{{\Sigma }_{Y}t}& 0\\ 0& I_{n-K}\end{array}\right].$$

$$\dot{P}=AP+P{A}^{\top }-PRP+Q,P(0)=P_0,$$

$$\left[\begin{array}{c}{\dot{X}}_1\\ {\dot{X}}_2\end{array}\right]=\left[\begin{array}{cc}-A^{\top }& R\\ Q& A\end{array}\right]\left[\begin{array}{c}{X}_1\\ X_2\end{array}\right],X_1(0)=I_n,X_2(0)=P_0,$$

$$\underbrace{\Phi^{\top}\!W(t)\Phi}_{P(t)}\!=\!\underbrace{S(t)\Phi^{\top}\!W_{0}\Phi}_{X_{2}(t)}\!\cdot\!\bigg{(}\underbrace{\!S^{-1\!}(t)\!\!\left(I_{n}\!+\!G(t)\Phi^{\top}W_{0}\Phi\right)}_{X_{1}(t)}\!\!\bigg{)}^{\!\!\!\!-1}\!\!\!\!,\!\!$$

$$\dot{W}=WY+YW-2W^2,W(0)=\alpha \Phi {\Sigma }_{U_0}^2{\Phi }^{\top },$$

$$\dot{\tilde{W}}=\tilde{W}{\Sigma }_Y+{\Sigma }_Y\tilde{W}-2{\tilde{W}}^2,\tilde{W}(0)=\alpha {\Sigma }_{U_0}^2,$$

$${\sigma }_{i,W}(t)$$

$${\sigma }_{i,W}\left(\frac{1}{2{\sigma }_{i,Y}}\log \frac{c}{\alpha }\right)$$

$${\sigma }_{i,W}\left(\frac{1}{2{\sigma }_{i,Y}}\log \frac{C}{\alpha }\right)$$

$${\mathcal{I}}_k:=\left[\frac{1}{2{\sigma }_{k,Y}}\log \frac{C_{\epsilon }}{\alpha },\frac{1}{2{\sigma }_{k+1,Y}}\log \frac{c_{\epsilon }}{\alpha }\right],1\le k\le K$$

$$\Vert W(t)-{\hat{Y}}_k\Vert \le \epsilon ,\forall t\in {\mathcal{I}}_k,$$

$$|{\sigma }_{l,W}(t)-{\sigma }_{l,Y}|\le \epsilon ,\forall l\le k,\text{and }|{\sigma }_{l,W}(t)|\le \epsilon ,\forall l>k.$$

$${\sigma }_{l,Y}\stackrel{(\ast )}{>}{\sigma }_{l,W}(t)$$

$$\stackrel{\text{eqref}e{q}_s{v}_{t}ra{j}_s{i}_{e}va{l}_b}{=}{\sigma }_{l,Y}-\left({\sigma }_{l,Y}-\frac{C_{\epsilon }{\sigma }_{l,Y}{\sigma }_{l,0}}{{\sigma }_{l,Y}+(C_{\epsilon }-\alpha ){\sigma }_{l,0}}\right)$$

$$={\sigma }_{l,Y}-\frac{{\sigma }_{l,Y}^2}{({\sigma }_{l,Y}-\alpha {\sigma }_{l,0})+C_{\epsilon }{\sigma }_{l,0}}$$

$$\stackrel{(\star )}{\ge }{\sigma }_{l,Y}-\frac{{\sigma }_{1,Y}^2}{C_{\epsilon }{\sigma }_{l,0}}\ge {\sigma }_{l,Y}-\epsilon ,$$

$$0\stackrel{(\ast )}{<}{\sigma }_{l,W}(t)$$

$$\stackrel{\text{eqref}e{q}_s{v}_{t}ra{j}_s{i}_{e}va{l}_s}{=}\frac{c_{\epsilon }{\sigma }_{l,Y}{\sigma }_{l,0}}{({\sigma }_{l,Y}-\alpha {\sigma }_{l,0})+c_{\epsilon }{\sigma }_{l,0}}\stackrel{(\star )}{\le }{c}_{\epsilon }{\sigma }_{l,0}\le \epsilon ,$$

$$0<{\sigma }_{l,W}(t)=\frac{\alpha {\sigma }_{l,0}}{1+2\alpha {\sigma }_{l,0}t}\le \alpha {\sigma }_{l,0}\le \epsilon .$$

$$\dot{\tilde{W}}=\tilde{W}{\Sigma }_Y+{\Sigma }_Y\tilde{W}-2{\tilde{W}}^2,\tilde{W}(0)=\alpha {\tilde{W}}_0,$$

$$\Vert W(t)-{\hat{Y}}_k\Vert \le \epsilon ,\forall t\in {\mathcal{I}}_k,$$

$$\tilde{W}(t_k)$$

$$\stackrel{({\tilde{W}}_0^{-1}=V)}{=}\alpha S(t_k){\left(V+\alpha G(t_k)\right)}^{-1}{S}^{\top }(t_k),$$

$$V=\left[\begin{array}{cc}{V}_{11}& V_{12}\\ V_{12}^{\top }& V_{22}\end{array}\right],S(t_k)=\left[\begin{array}{cc}{S}_1& 0\\ 0& S_2\end{array}\right],G(t_k)=\left[\begin{array}{cc}{G}_1& 0\\ 0& G_2\end{array}\right],$$

$$S_2=\left[\begin{array}{cc}{e}^{{\Sigma }_k^c{t}_k}& 0\\ 0& I_{n-K}\end{array}\right],G_2=\left[\begin{array}{cc}({\Sigma }_k^c{)}^{-1}(e^{2{\Sigma }_k^c{t}_k}-I_{K-k})& 0\\ 0& 2I_{n-K}{t}_k\end{array}\right].$$

$$\tilde{W}(t_k)$$

$$=\left[\begin{array}{cc}{H}_{11}& H_{12}\\ H_{12}^{\top }& H_{22}\end{array}\right],$$

$$H_{11}$$

Full text

I Introduction

Why do deep neural networks [1, 2, 3] work so well in practice? At the center of this mystery lies the fact that practical neural networks are typically overparameterized, i.e., their number of parameters is several orders of magnitude larger than the number of available training examples. As a consequence, there are infinitely many parameter choices that can perfectly fit the data. However, different parameters have different generalization properties, and only a few choices lead to networks that make correct predictions for new samples. Behind the success of deep learning lies the ability of training algorithms, such as gradient descent, to find those networks that generalize well, which over past years has motivated many theoretical studies on such an implicit bias [4] of learning algorithms.

An interesting one is the simplicity bias of gradient descent or gradient flow under certain initial conditions, typically random initialization with small variance, through which the learned network model possesses certain notions of low complexity. Depending on the network architecture, such notions range from classic ones, such as low-rankness [5, 6, 7, 8] or sparsity [9, 10, 11] of the network weights, to complexity measures in some function space [12]. In many scenarios, the analyses for such simplicity biases often call for theoretical understandings of the incremental learning [13, 6, 7] phenomenon, which describes the evolution of the network model during training as sequential learning of components of the underlying target ground truth function (suppose there is one) in some canonical decompositions, without learning any spurious components that overfit to the training data.

One important testbed for understanding incremental learning is the matrix factorization [14] problem, studied in the context of linear regression [13, 6, 15, 16, 17], matrix sensing [5, 18, 19, 7], and recently the low-rank adaptation [20, 21]. Under random initialization with small variance, training with gradient descent/flow on the factorized model finds a stationary point that corresponds exactly (noiseless) or approximately (noisy) to the underlying ground truth matrix, following a particular dynamic trajectory. During training, the factorized model grows its singular components to fit/learn the ground truth matrix in a sequential manner such that the growth happens for one singular component at a time to fit one component in the ground truth, and the ground truth component with the largest singular values get learned first, then followed by the second largest and so on so forth. As such, the training algorithm finds the ground truth by every time incrementally learning one of its components.

A body of theoretical work on incremental learning in matrix factorization is motivated by some of the early works [13, 5] that investigate the implicit bias of deep learning algorithms. One initial line of attempts studies closed-form solutions: If one carefully initializes the singular vectors of the factorized model to align with those of the ground truth, then the training dynamics are decoupled into multiple scalar dynamics, each corresponding to one principal component in the ground truth [13, 6, 15]. In this case, understanding incremental learning reduces to studying the time-scale separation among these scalar dynamics. However, these results are achieved under very restrictive conditions on initialization.

To analyze this phenomenon under general initial conditions, recent work focuses on a more fine-grained trajectory analysis[19, 7] that introduces some orthogonal decomposition of the parameter space into signal and error subspaces based on the ground truth. This analysis shows that, at certain time epoch, the components within the signal subspace learn some principal components of the ground truth, together with upper bounds on the growth within the error subspace. While the trajectory analysis produces the most general results [7] for matrix factorization, it does not fully answer many quantitative questions. For example, how small does the variance of the random initialization need to be in order for incremental learning to happen? If one seeks one of the low-rank approximations of the ground truth, within what time interval does the algorithm need to stop? These questions have become increasingly important alongside the growing interest (for example, [22]) in utilizing the incremental mechanism to obtain low-complexity models, as opposed to the classic regularized approach.

In this paper, we extend the closed-form solution-based line of work to gradient flow on matrix factorization under general initial conditions, facilitating a rich quantitative understanding of the incremental phenomenon. Our theorems aim to exactly address several questions previously mentioned by showing the precise dependence of incremental learning on the variance of the random initialization and providing time intervals within which the factorized model serves as a good low-rank approximation of the ground truth. Although our results are primarily for symmetric matrix factorization, we outline the path toward extending them to the asymmetric matrix factorization problems, opening up new opportunities to utilize closed-form solutions to study a broader class of problems.

The organization of this paper is as follows. In Section II, we introduce the gradient flow dynamics on the symmetric matrix factorization problem, the primary focus of our theoretical results, and we discuss the closed-form solution of the factorized model to the gradient flow. This forms the basis of our main results in Section III on the incremental learning in the symmetric matrix factorization problem. Then, in Section IV, we conclude by briefly discussing the path toward extending current results to the asymmetric matrix factorization problems.

Notations: For an $n\times m$ matrix $A$, we let $\|A\|$ and $\|A\|_{F}$ denote the spectral and Frobenius norm of $A$, respectively. We write $A\succeq 0$ ($A\preceq 0$) when $A$ is positive semi-definite (negative semi-definite). Then, we let $\mathrm{diag}\{a_{i}\}_{i=1}^{n}$ be a diagonal matrix with $a_{i}$ being its $i$-th diagonal entry.

II Preliminaries: Gradient Flow on Symmetric Matrix Factorization

The problem of our main interest in this paper is the following symmetric matrix factorization:

[TABLE]

where $Y\in\mathbb{R}^{n\times n}$ with $Y\succeq 0$ is the target matrix to be factorized. Although one may find solving (1) easy by simply taking the SVD of $Y$, the goal of this paper is rather to understand training trajectories when one seeks to solve (1) by gradient flow (GF) on $\mathcal{L}(U)$, i.e.,

[TABLE]

from some initial condition $U_{0}$. More specifically, as stated in the introduction, we are interested in the incremental learning phenomenon when $U(t)U^{\top}(t)$ learns the singular components of $Y$ sequentially along the GF trajectory. Therefore, it suffices to study the induced dynamic evolution of the factorized model $U(t)U^{\top}(t):=W(t)$ from (2), which is

[TABLE]

(3) is a special form of matrix Riccati differential equations, thus can be analyzed from its closed-form solution.

II-A Closed-form solution

The closed-form solutions to matrix Riccati differential equations are well-studied, for example, in [23, 24]. We thus have the following result regarding the solution to (3)

Proposition 1**.**

If $Y$ has $K$ non-zero singular values and the full SVD of $Y\succeq 0$ is $\Phi\begin{bmatrix}\Sigma_{Y}&0\\ 0&0\end{bmatrix}\Phi^{\top}$, then the matrix Riccati differential equation

[TABLE]

has unique solution of the following:

[TABLE]

where $\tilde{W}_{0}=\Phi^{\top}W_{0}\Phi$ and

[TABLE]

Proof.

From [23, 24], one knows that the solution to general matrix Riccati equation

[TABLE]

is given by $P(t)=X_{2}(t)X_{1}^{-1}(t)$, where $\{X_{1}(t),X_{2}(t)\}$ are solution to an LTI system

[TABLE]

To apply this, let $P\!=\!\Phi^{\top}\!W\!\Phi$, $A\!=\!\Sigma_{Y}$, $Q\!=\!0$ and $R\!=\!I_{n}$, by solving (8) we have (an equivalent expression for (5))

[TABLE]

where $G(t),S(t)$ are defined in (6). Notice that $\forall t$, all eigenvalues of $G(t)W_{0}$ are real and non-negative since both $G(t)$ and $W_{0}$ are positive semi-definite, then $X_{1}(t)$ is invertible $\forall t$, ensuring (9) is the unique solution [24, Corollary 1]. ∎

II-B Additional settings

With the closed-form solution available, one can study this GF under any initial condition $W_{0}=U_{0}U_{0}^{\top}\succeq 0$. However, different types of initialization (for example, random v.s. deterministic initialization) may require different analyses. Due to space constraints, we opt to discuss the case of overparametrized matrix factorization under deterministic initialization with varying scales.

As we described in the introduction, the incremental learning phenomenon generally happens under initialization with a small scale. To isolate this critical effect of the scale/norm of the initialization, it is convenient to write our initial condition as $U(0)=U_{0}=\alpha^{1/2}\bar{U}_{0}$ for some initialization scale $\alpha>0$ and initialization shape $\bar{U}_{0}\succeq 0$. Our analysis

• •

assumes a fixed shape $\bar{U}_{0}$ and describes how the scale $\alpha$ affects the incremental learning phenomenon;

• •

In addition, we assume $W_{0}=\alpha\bar{U}_{0}\bar{U}_{0}^{\top}$ is invertible. This necessarily requires that the problem (1) is overparametrized [15], i.e. $r\geq n$, which imply that the optimal value of (1) is always zero (target is exactly factorized) regardless of the target $Y$.

III Main results: Incremental Learning in Symmetric Matrix Factorization

In this section, we discuss incremental learning in symmetric factorization problems. First, we revisit the GF under spectral initialization with a small scale and refine the existing asymptotic results [6] to a non-asymptotic one. Then, we show that a similar quantitative characterization of incremental learning holds for general initialization. These main results follow the problem settings stated in the previous section, and afterward, we discuss how to adapt them to other settings.

III-A Incremental learning under spectral initialization

Give the GF (2) with the SVD of the target matrix $Y=\Phi\Sigma_{Y}\Phi^{\top}$, the spectral initialization of $U$ is any $U_{0}$ that has an SVD of the form $U_{0}=\alpha^{1/2}\Phi\Sigma_{U_{0}}V_{U_{0}}^{\top}$. Recall that $\alpha>0$ is the initialization scale. Under a spectral initialization, the induced ode on $W$ (3) becomes:

[TABLE]

Through a change of variables $\tilde{W}=\Phi^{\top}W\Phi$, the dynamics on $\tilde{W}$ is given by

[TABLE]

whose solution $\tilde{W}(t)$ is diagonal $\forall t$, and the dynamics of diagonal entries are decoupled from each other. Indeed, examining the closed-form solution in Proposition 1 under a spectral initialization, we have the following:

Corollary 2**.**

Suppose $\Sigma_{Y}=\mathrm{diag}\{\sigma_{i,Y}\}_{i=1}^{n}$. Under some spectral initialization $U(0)=\alpha^{1/2}\Phi\Sigma_{U_{0}}V_{U_{0}}^{\top}$, the solution $U(t)$ to the GF dynamics (2) induces $W(t)=U(t)U^{\top}(t)$ of the form $W(t)=\Phi\mathrm{diag}\{\sigma_{i,W}(t)\}_{i=1}^{n}\Phi^{\top}$ with

[TABLE]

where $\sigma_{i,0}=[\Sigma_{U_{0}}^{2}]_{ii}\geq 0,\forall i$.

Spectral initialization offers a simplified analysis of the GF trajectory. Under the induced dynamics (3), we expect the factorized model $W(t)$ to converge to the target matrix $Y$. Indeed, for each $i$, (12) is the learning trajectory for a non-zero singular mode $\sigma_{i,Y}$ of $Y$, and it is easy to see that as long as $\sigma_{i,0}\neq 0$, we have $\lim_{t\rightarrow\infty}\sigma_{i,W}(t)=\sigma_{i,Y}$; (13) is the learning trajectory for a zero singular mode $\sigma_{i,Y}$ of $Y$, and we always have $\lim_{t\rightarrow\infty}\sigma_{i,W}(t)=0$. Therefore, if $\sigma_{i,0}\neq 0,\forall i$ with $\sigma_{i,Y}\neq 0$, one has $\lim_{t\rightarrow\infty}W(t)=Y$, thereby learning the target matrix. Note that a similar condition is also necessary in the general case. That is, learning the exact $Y$ requires $\Phi_{Y}^{\top}W(0)\Phi_{Y}$ to have full rank, where columns of $\Phi_{Y}$ form an orthonormal basis of the range space of $Y$.

As shown in the previous discussion, the asymptotic behavior of these scalar dynamics (12)(13) are easy to study, concerning whether $W(t)$ can successfully learn the target matrix $Y$. Now, we turn to the questions regarding the transient behavior under small initialization with scale $\alpha\ll 1$. Notice that when $0<\alpha\sigma_{i,0}<\sigma_{i,Y}$, $\sigma_{i,W}(t)$ is strictly monotonically increasing w.r.t. time $t$. Our following discussions assume $\alpha$ is sufficiently small so that $\sigma_{i,W}(t)$ is monotonically increasing.

For every $\sigma_{i,W}(t)$ following (12) that learns a non-zero mode $\sigma_{i,Y}$ of $Y$, we have $\sigma_{i,W}(0)=\alpha\sigma_{i,0}$, which is close to zero for sufficiently small $\alpha$, and we know asymptotically $\lim_{t\rightarrow\infty}\sigma_{i,W}(t)=\sigma_{i,Y}$. Therefore, each $\sigma_{i,W}(t)$ must learn its target $\sigma_{i,Y}$ through some transitional phase during which the value $\sigma_{i,W}(t)$ grow from a small value close to $\alpha\sigma_{i,0}$ to a value close to $\sigma_{i,Y}$. If there is a time-scale separation between these transitional phases for different $i$, then we expect that at some intermediate time along GF, a subset of $\sigma_{i,W}(t)$ has learned their corresponding target $\sigma_{i,Y}$ while the rest stays close to zero, from which the resulting $W(t)$ serves as some low-rank approximation of $Y$.

The above discussion suggests a closer look at the transient behavior in the solution (12) when the initialization scale is small: For some $c,C>0$ such that $c\sigma_{i,0}\ll\sigma_{i,Y}$ and $C\sigma_{i,0}\gg\sigma_{i,Y}$, $\forall i$. Then we see that for every $i$ (for any $\alpha\leq c$),

[TABLE]

suggesting a sharp transition in $\sigma_{i,\bm{W}}(t)$ from a small value $c\sigma_{i,0}$ to one close to its target value $\sigma_{i,Y}$ at a time that scales as $\Theta\left(\frac{1}{\sigma_{i,Y}}\log\frac{1}{\alpha}\right)$. Notably, the transition time depends inverse-proportionally on $\sigma_{i,Y}$ (See Figure 1). Thus, this time-scale separation can result in singular values of $Y$ being learned one by one, with larger ones learned first and smaller ones later, which is exactly an incremental learning phenomenon. Based on these discussions, the conditions to achieve this are:

• •

Singular values $\sigma_{i,Y}$ of $Y$ are distinct;

• •

The scale $\alpha$ is sufficiently small such that the transitional phases $\left[\frac{1}{2\sigma_{i,Y}}\log\frac{c}{\alpha},\frac{1}{2\sigma_{i,Y}}\log\frac{C}{\alpha}\right]$ are non-overlapping.

Formally, we have the following theorem:

Theorem 3** (Incremental learning under small spectral initialization).**

Suppose the $K$ non-zero singular values $\sigma_{1,Y},\sigma_{2,Y},\cdots,\sigma_{K,Y}$ of $Y$ are distinct and ordered in decreasing order. Let the spectral initialization $U(0)=\alpha^{1/2}\Phi\Sigma_{U_{0}}V_{U_{0}}^{\top}$ have a initialization shape $\Sigma_{U_{0}}^{2}=\mathrm{diag}\{\sigma_{i,0}\}_{i=1}^{d}$ with $\sigma_{i,0}>0,\forall 1\leq i\leq K$. Given any error tolerance $0<\varepsilon\leq\sigma_{K,Y}$, let $c_{\varepsilon}=\frac{\varepsilon}{\max_{i}\sigma_{i,0}}$, and $C_{\varepsilon}=\frac{\sigma_{1,Y}^{2}}{\varepsilon\min_{1\leq i\leq K}\sigma_{i,0}}$, Suppose $\alpha$ is sufficiently small such that $\alpha\leq c_{\varepsilon}$ and that $\frac{-\log\alpha+\log c_{\varepsilon}}{-\log\alpha+\log C_{\varepsilon}}>\max_{1\leq k\leq K-1}\frac{\sigma_{k+1}}{\sigma_{k}}$, then the time intervals

[TABLE]

are all non-empty (we have let $\frac{1}{0}:=\infty$), and the $W(t)$ in Corollary 2 satisfies that $\forall 1\leq k\leq K$,

[TABLE]

where $\hat{Y}_{k}:=\arg\min_{\mathrm{rank}(Z)=k}\|Y-Z\|_{F}$ is the best rank-$k$ approximation of $Y$.

Proof.

We have assumed small $\alpha$ such that $\frac{-\log\alpha+\log c_{\varepsilon}}{-\log\alpha+\log C_{\varepsilon}}>\max_{1\leq k\leq K-1}\frac{\sigma_{k+1}}{\sigma_{k}}$, which implies that for every $1\leq k\leq K$, $\mathcal{I}_{k}\neq\emptyset$. For a fixed $k$, choose any $t$ within this interval. Since $\Phi^{\top}W(t)\Phi$ and $\Phi^{\top}\hat{Y}_{k}\Phi$ are diagonal, write (17) as

[TABLE]

Given the solution in Corollary 2, we have, $\forall 1\leq l\leq k$,

[TABLE]

and $\forall k<l\leq K$,

[TABLE]

where $(*)$ uses the monotonicity of $\sigma_{i,W}(t)$, and $(\star)$ uses the fact that $\sigma_{l,Y}-\alpha\sigma_{l,0}\geq 0,\forall l\leq K$, derived from $\alpha\leq c_{\varepsilon}$. Lastly, since we have assumed $\alpha\leq c_{\varepsilon}$, we have $\forall l>K$,

[TABLE]

With (19)(20)(21), the spectral bound in (18) is verified. ∎

As shown in Theorem 3 and the preceding discussions, from a dynamical system perspective, the incremental learning phenomenon comes from some time-scale separation among dynamics corresponding to learning different components in the target matrix. By decreasing the initialization scale $\alpha$, these time-scale separations become more prominent, allowing one to find low-rank (best rank-$k$) approximations of the target matrix (See Figure 1) after one transitional phase concludes ($t=\frac{1}{2\sigma_{k,Y}}\log\frac{C_{\varepsilon}}{\alpha}$) and before another kicks in ($t=\frac{1}{2\sigma_{k+1,Y}}\log\frac{c_{\varepsilon}}{\alpha}$). Our Theorem exactly quantifies the scale needed to ensure the existence of such low-rank approximations (up to error $\varepsilon$), and the time intervals within which one can find these approximations.

III-B Incremental learning under general initialization

When one moves to a general initialization $U(0)=\alpha^{1/2}\bar{U}_{0}$ for some $U_{0}$, the dynamics

[TABLE]

can no longer be decomposed into scalar dynamics, where $\tilde{W}_{0}=\Phi^{\top}\bar{U}_{0}\bar{U}_{0}^{\top}\Phi$. Nonetheless, similar discussions in the previous section still work: Due to the small initialization scale, the dynamic evolution of the part in $\tilde{W}$ that learns top-$k$ singular components of $Y$ has a faster time scale than those that learn the rest of the singular values, and the interval $\mathcal{I}_{k}$ again separates the two dynamics. By carefully partitioning the closed-form solution obtained in (1) and analyzing the growth of each block component, we have the following theorem suggesting the same incremental learning occur under general initialization:

Theorem 4** (Incremental learning under general small initialization).**

Suppose the $K$ non-zero singular values $\sigma_{1,Y},\sigma_{2,Y},\cdots,\sigma_{K,Y}$ of $Y$ are distinct and ordered in decreasing order. Suppose the initialization of $U(0)=\alpha^{1/2}\bar{U}_{0}$ satisfies that $\tilde{W}_{0}=\Phi^{\top}\bar{U}_{0}\bar{U}_{0}^{\top}\Phi$ has an inverse denoted by $V$. Given some error tolerance $0<\varepsilon\leq\min\{\sigma_{K,Y},1\}$, let $c_{\varepsilon}=\frac{\varepsilon}{16M^{2}}$, $C_{\varepsilon}=\frac{16\sigma_{1,Y}^{2}M^{2}}{\varepsilon}$, where $M:=\max\{\|V\|,\|V^{-1}\|\}$. If the initialization scale $\alpha$ is sufficiently small so that $\frac{-\log\alpha+\log c_{\varepsilon}}{-\log\alpha+\log C_{\varepsilon}}>\max_{1\leq k\leq K-1}\frac{\sigma_{k+1}}{\sigma_{k}}$ and $\alpha\leq\frac{c_{\varepsilon}}{M}$, then the time intervals $\mathcal{I}_{k}$ defined in (16) are all non-empty, and the $W(t)$ in Proposition 1 satisfies that $\forall 1\leq k\leq K$,

[TABLE]

where $\hat{Y}_{k}:=\arg\min_{\mathrm{rank}(Z)=k}\|Y-Z\|_{F}$ is the best rank-$k$ approximation of $Y$.

Proof.

Similar to the proof of Theorem 3, our assumption on $\alpha$ ensures every $\mathcal{I}_{k}$ is non-empty. To prove (23), we show that for every choice of $k$, and for any $t_{k}\in\mathcal{I}_{k}$, we have $\|W(t_{k})-\hat{Y}_{k}\|\leq\varepsilon$, the following proof proceeds given such a choice of $k$ and $t_{k}$.

Rewrite solution in blocks: From (5) in Proposition 1,

[TABLE]

Let $\Sigma_{k}=\mathrm{diag}\{\sigma_{l,Y}\}_{l=1}^{k}$, and $\Sigma_{k}^{c}=\mathrm{diag}\{\sigma_{l,Y}\}_{l=k+1}^{K}$ (note that $\Sigma_{k}^{c}$ can be an empty matrix), write $V,S(t_{k}),G(t_{k})$ with the following block forms:

[TABLE]

where $S_{1}=e^{\Sigma_{k}t_{k}}$, $G_{1}=\Sigma_{k}^{-1}(e^{2\Sigma_{k}t_{k}}-I_{k})$,

[TABLE]

Then we have

[TABLE]

where (by the inverse formula for a block matrix)

[TABLE]

Our goal is to show that

[TABLE]

Thus it suffices to show the following

[TABLE]

The rest of the proof is to show (27) holds, for which we start by deriving norm bounds on $S_{1}$ and $S_{2}$, and then move to $H_{11}-\Sigma_{k}$, $H_{12}$ and $H_{22}$.

Norm bounds on $S_{1}$ and $S_{2}$: From that $t_{k}\in\mathcal{I}_{k}$, we have

[TABLE]

and similarly since $e^{\Sigma_{k}^{c}t_{k}}\preceq e^{\Sigma_{k}^{c}\frac{1}{2\sigma_{k+1}}\log\frac{c_{\varepsilon}}{\alpha}}\preceq e^{\frac{1}{2}\log\frac{c_{\varepsilon}}{\alpha}}I_{k}$, we have $S_{2}\preceq\begin{bmatrix}\sqrt{\frac{c_{\varepsilon}}{\alpha}}I_{K-k}&0\\ 0&I_{n-K}\end{bmatrix}$. Together we have

[TABLE]

Norm bounds on $H_{11}-\Sigma_{k}$, $H_{12}$ and $H_{22}$: The primary focus is the bound on $\|H_{11}-\Sigma_{k}\|$, for the sake of simplicity, we let $(V_{11}-V_{12}(V_{22}+\alpha G_{2})^{-1}V_{12}^{\top}):=\tilde{V}_{1}$, then

[TABLE]

Notice that

[TABLE]

from which we know $\|\tilde{V}_{1}\|\leq\|V_{11}\|\overset{(\textsf{a})}{\leq}M$. Then by (28),

[TABLE]

As long as $\|\Sigma_{k}\|\|\Delta\|\overset{(\textsf{c})}{\leq}\frac{1}{2}$, we have, from [25],

[TABLE]

the desired bound on $\|H_{11}-\Sigma_{k}\|$. Since we have assumed $\varepsilon\leq\sigma_{1,Y}$, this also suggests $\|H_{11}\|\leq\|\Sigma_{k}\|+\varepsilon\leq 2\sigma_{1,Y}$.

For $H_{12},H_{22}$, we use the fact that $\|(V_{22}+\alpha G_{2})^{-1}\|\leq\|V_{22}^{-1}\|\overset{(\textsf{a})}{\leq}M$, and $\|V_{12}\|\overset{(\textsf{a})}{\leq}M$, together with (28), to obtain

[TABLE]

and

[TABLE]

This finishes the proof. For the sake of readability, we omitted several arguments $(\textsf{a})(\textsf{b})(\textsf{c})$ used in the last part of the proof, and we explain them here. ${\footnotesize(\textsf{a})}$: $\|V_{11}\|\leq M,\|V_{12}\|\leq M$ because the norm of submatrix is bounded by that of the full matrix, and $\|V_{22}^{-1}\|\leq M$ is due to the interlacing theorem [25] for principal submatrix. ${\footnotesize(\textsf{b})}$: we have $\sigma_{K,Y}^{-1}\alpha\leq\sigma_{K,Y}^{-1}\varepsilon\leq 1\leq M$. ${\footnotesize(\textsf{c})}$: with $\|\Delta\|\leq\frac{2M}{C_{\varepsilon}}$, it is easy to verify that $\|\Sigma_{k}\|\|\Delta\|\leq\sigma_{1,Y}\|\Delta\|\leq\frac{1}{2}$. ∎

The incremental learning results in Theorem 4 resemble the one in Theorem 3 for spectral initialization, except for some differences in the definitions of $c_{\varepsilon},C_{\varepsilon}$. It characterizes the quantitative effect of the initialization scale in determining the time-scale separation among the learning dynamics for each target singular component, which provides guidance on how to utilize these results to design training algorithms that can learn low-rank approximations of the ground truth by early stopping. In addition, the proof of the symmetric factorization case might serve as a basis for broader settings.

III-C Remarks on other settings

III-C1 Random initialization

A more natural setting is arguably the random initialization with small variance, which corresponds to having $U(0)=\alpha^{1/2}\bar{U}_{0}$ for some $\bar{U}_{0}$ whose entries are i.i.d. samples from a standard Gaussian (or other sub-Gaussian distributions). Applying Theorem 4 to this random initialization setting requires an extra concentration results on the matrix $\bar{W}_{0}=\bar{U}_{0}\bar{U}_{0}^{\top}$. If one can show that for $\delta\in(0,1)$, $\exists M_{\delta}>0$ such that $\mathbb{P}\left(\max\{\|\bar{W}_{0}\|,\|\bar{W}_{0}^{-1}\|\}\leq M_{\delta}\right)\geq 1-\delta$, then one can define $c_{\varepsilon},C_{\varepsilon}$ by $M_{\delta}$ and find initialization scale $\alpha$ accordingly, and Theorem 4 suggests that the incremental learning phenomenon as described by (23) happens with probability at least $1-\delta$.

III-C2 Rank-deficient initialization

Our proof assumes an initialization shape $\bar{U}_{0}$ such that $\tilde{W}_{0}=\Phi^{\top}\bar{U}_{0}\bar{U}_{0}^{\top}\Phi$ is invertible. This is satisfied if $r<n$ when $\tilde{W}_{0}$ is always rank-deficient. This requires significant changes in the proof since now we should analyze \tilde{W}(t_{k})=S(t_{k})\alpha\tilde{W}_{0}\big{(}I_{n}+\alpha G(t_{k})\tilde{W}_{0}\big{)}^{-1}S^{\top}(t_{k}) directly in (24). With Woodbury’s matrix identity, one rewrite $\tilde{W}(t_{k})$ as

[TABLE]

where $\tilde{U}=\Phi^{\top}\bar{U}_{0}$. Detailed analysis on $\tilde{U}^{\top}G(t_{k})\tilde{U}$ leads to similar results as Theorem 4, which is left to future work.

IV Conclusion

In this paper, we develop a quantitative understanding of the incremental learning phenomenon in GF on symmetric matrix factorization problems through its closed-form solution. From a dynamical system perspective, the incremental learning phenomenon comes from some time-scale separation among dynamics corresponding to learning different components in the target matrix. By decreasing the initialization scale, these time-scale separations become more prominent, allowing one to find low-rank approximations of the target matrix.

Future work includes extending current results to asymmetric matrix factorization. Notably, for asymmetric matrix factorization, [15] shows that if the factors are initialized to satisfy some balancedness conditions, then the learning dynamics can be related to GF on a symmetric factorization problem. Moreover, under a small scale, any initialization is close to one that satisfies such balancedness conditions; therefore, the analysis for general initialization can be viewed as analyzing perturbed dynamics from a nomial dynamics with closed-form solution.