Connectivity Shapes Implicit Regularization in Matrix Factorization Models for Matrix Completion
Authors: Zhiwei Bai
Jiajie Zhao
Yaoyu Zhang
Venue: NeurIPS 2024 poster (2024)
Keywords: matrix completion; implicit regularization; training dynamics; low nuclear norm; low rank;
Primary area: learning_theory
TL;DR: We provide a unified understanding of when, how, and why matrix factorization models achieve different implicit regularization effects.

Abstract: Matrix factorization models have been extensively studied as a valuable test-bed for understanding the implicit biases of overparameterized models. Although both low nuclear norm and low rank regularization have been studied for these models, a unified understanding of when, how, and why they achieve different implicit regularization effects remains elusive. In this work, we systematically investigate the implicit regularization of matrix factorization for solving matrix completion problems. We empirically discover that the connectivity of observed data plays a key role in the implicit bias, with a transition from low nuclear norm to low rank as data shifts from disconnected to connected with increased observations. We identify a hierarchy of intrinsic invariant manifolds in the loss landscape that guide the training trajectory to evolve from low-rank to higher-rank solutions. Based on this finding, we theoretically characterize the training trajectory as following the hierarchical invariant manifold traversal process, generalizing the characterization of Li et al.(2020)  to include the disconnected case. Furthermore, we establish conditions that guarantee minimum nuclear norm, closely aligning with our experimental findings, and we provide a dynamics characterization condition for ensuring minimum rank. Our work reveals the intricate interplay between data connectivity, training dynamics, and implicit regularization in matrix factorization models.

Connectivity Shapes Implicit Regularization in Matrix Factorization Models for Matrix Completion
Zhiwei Bai', Jiajie Zhao', Yaoyu Zhang'*
School of Mathematical Sciences, Institute of Natural Sciences, MOE-LSC, Shanghai Jiao Tong University, Shanghai 200240, P.R. China. (bai299, zjj0216, zhyy.sjtu) @sjwueduuccn
Abstract
Matrix factorization models have been extensively studied as a valuable test-bed for understanding the implicit biases of overparameterized models. Although both low nuclear norm and low rank regularization have been studied for these models, a unified understanding of when, how, and why they achieve different implicit regularization effects remains elusive. In this work, we systematically investigate the implicit regularization of matrix factorization for solving matrix completion problems. We empirically discover that the connectivity of observed data plays a crucial role in the implicit bias, with a transition from low nuclear norm to low rank as data shifts from disconnected to connected with increased observations. We identify a hierarchy of intrinsic invariant manifolds in the loss landscape that guide the training trajectory to evolve from low-rank to higher-rank solutions. Based on this finding, we theoretically characterize the training trajectory as following the hierarchical invariant manifold traversal process, generalizing the characterization Li et al. 2020) to include the disconnected case. Furthermore, we establish conditions that guarantee minimum nuclear norm, closely aligning with experi- mental findings, and we provide a dynamics characterization condition for ensuring minimum rank. Our work reveals the intricate interplay between data connectivity, training dynamics, and implicit regularization in matrix factorization models.
1 Introduction
Overparameterized models have the capacity to easily fit data with random labels (Zhang et al.) 2017
2021). However, in real-world applications, models with more parameters than training samples still generalize well. This has led researchers to hypothesize that overparameterizee models undergo implicit regularization, favoring certain functions as outputs. Overparameterizeed matrix factorization models, fo AB with 0 = (A, B), A, B E ddxd have served as a simplified test-bed for studying this implicit regularization. In the context of matrix completion problems like the Netflix challenge, these models aim to find low-rat low-rank completion of a partially observed matrix M E [ Prior works have offered seemingly conflicting perspectives with some claiming it promotes low nuclear norm (Gunasekar et al. 2017) and others arguing for low rank (Arora al. 2019 Li et al. 2020 Razin and Cohen 2020). However, unified understanding of when, how, and why they achieve different implicit regularization effects remains elusive.
the implicit regularization play,
Unlike previous works that focus on either low rank or low nuclear norm regularization, we sys- tematically investigate the training dynamics and implicit regularization of matrix factorization for matrix completion. Through extensive experiments, of observed data plays a key role in the implicit regularization effects. Data connectivity, in the context of this paper, refers to the way observed entries in the matrix are linked through shared rows
found that certain connectivity property
*Corresponding author: mnyyyjju@sstt..dd..o
38th Conference on Neural Information Processing Systems (NcurIPS 2024).

---

Rank vs. Connectivity
Nuclear Norm vs, Connectivity
Expertmental Mucloaa Norm
Experimental rank
25-
Connected
2
Disconeetted with complete bipartite components
Disconeccedtotheen
2
7.5 10.0 12.5
15.0 17.5
Ground Truthrrakk)
Ground TuthhNuclear Norm)
Figure 1: The connectivity of observed data affects the implicit regularization. The ground truth matrix M E R4x4 has rank ranging from to 3. The sample size n covers settings where n is equal to, smaller than, and larger than the 2rd - 22 threshold required for exact reconstruction. Darker scatter points indicate a greater number of samples, while lighter points indicate fewer samples. The positions of observed entries are randomly chosen, and the experiment is repeated 10 times for each sample size. (Please refer to Appendix B for additional experiments and detailed methodology.)
or columns. set of observations is considered connected if there's path between any two observed entries via other observed entries in the same rows or columns. This concept plays a crucial role in determining the behavior of matrix factorization models, as we will demonstrate throughout this paper. As shown in Fig. [1
we sample observations randomly from a ground truth matrix M E IRdxd with rank(MM) < and train models fe AB, A, E Radxd from small random initialization without any rank constraints. For each observation set, we calculate the solutions with the minimum nuclear norm and minimum rank, which serve as the ground truth benchmarks. These are then compared with the completion matrix obtained by the model. From Fig. D
we observe that:
(i) Low rank bias in connected case: When the observed entries are connected, the model consis- tently learns the lowest-rank solution.
(ii) Low nuclear norm bias in certain disconnected case:
When the observed entries are dis- connected, the model generally does not find the minimum nuclear norm or lowest-rank solution. However, in the special case where each connected component is complete bipartite subgraph, the model consistently finds the minimum nuclear norm solution.
To understand how data connectivity modulates the implicit bias, analyze the loss landscape and optimization dynamics. We find a hierarchy intrinsic invariant manifolds S2 of different ranks in the loss landscape. These manifolds constrain the optimization trajectory, causing the model to learn by incrementally ascending through higher ranks. the disconnected case, additional sub-S2; invariant manifolds emerge within the 22 invariant manifold, preventing the model from reaching the global lowest-rank solution. However, we prove that the minimum nuclear norm solution is guaranteed in the disconnected with complete bipartite subgraph case.
The contributions of our work are summarized as follows:
(i) We systematically investigate the influence of data connectivity on the implicit regularization. Our empirical findings indicate that the connectivity of observed data plays key role in the implicit bias, leading to a transition from favoring solutions with a low nuclear norm to those with a low rank as the data becomes more connected with an increase in observations (refer to Sec. 4
(ii) We characterize the training dynamics of matrix factorization theoretically, showing that the optimization trajectory follows Hierarchical Invariant Manifold Traversal (HIMT) process. This gen- eralizes the characterization of et al. 2200)). whose proposed Greedy Low-Rank Learning GLRL) algorithm equivalence only corresponding to the connected case (refer to Sec. 3 and Sec. 6.1)
(iii) Regarding the minimum nuclear norm regularization, we establish conditions that provide guarantees closely aligned with our empirical findings, which complement the results of Gunasekar jet al. (2077). For the minimum rank regularization, we present a dynamic characterization condition that assures the attainment of the minimum rank solution (refer to Sec. 6.2

---

2
Related works
Norm minimization and rank minimization. Extensive research has been conducted on the implicit regularization of matrix factorization models, focusing norm minimization and rank mini- mization. For norm minimization, Gunasekar et al. (2017) proved that gradient flow with infinitesimal initialization converges to the minimum nuclear norm solution in the special case of commutative observations. Ji and Telgarsky p: /unasekar et al. (2018) studied norm minimization regular- ization in deep linear networks. For rank minimization, numerous works have shown that matrix factorization models favor low-rank solutions. Arora et al. 20 Gidel et al. (2019); Gissin et al.
(2019); Razin and Cohen (2020); Jiang et al. (2023); Belabbbas (2020) investigated how infinitesimal initialization of gradient flow encourages low rank in specific settings. Li et al. (2020) showed that under certain assumptions, matrix factorization dynamics are equivalent
greedy low-rank
learning heuristic. Li el al. (2018); Stöger and Soltanolkotabi (2021); Jin et al. (2023) established low-rank recovery guarantees for matrix sensing problems under the Restricted Isometry Property
(RIP) condition. Zhang et al. [2022] 2023) studied a broader class of model rank minimization for nonlinear models, of which the matrix factorization model a special case.
Nonlinear dynamics. The initialization scale can significantly influence the implicit regularization of neural networks. Large initialization typically leads to linear dynamics (Jacot et al. 2018 and poor generalization (Chizat al. 2019), while small initialization induces nonlinear dynamics (Luo) et al. 202T). In this work, we focus on the case of infinitesimal initialization, which corresponds to highly nonlinear dynamics. An important characteristic of nonlinear neural network dynamics is the
phenomenon of condensation (Luo et al. [2221] Zhou et al., (022). where the network's 's effective complexity small. The low-rank 12 invariant manifolds we propose are essentially a manifestation of condensation. Zhang et al. (2021) 2022); Bai et al. (2022) Fukumizu et al. (2019); Simssk et al. (202T) established the embedding principle of the loss landscape of neural networks and empirically demonstrated that the training process traverses critical points embedded from smaller subnetworks. Jacot al.) (202T) conjectured saddle to saddle dynamics for deep linear networks, which is conceptually analogous to the dynamics characterization in this work.
3 Preliminaries
Matrix completion problem. This study focuses on the matrix completion problem, which involves estimating missing entries within a partially observed matrix. Given an incomplete matrix M E R the goal is to predict the entirety of M based its observed elements. The set of observed entries is represented as S = {(k,, jk), M Ik=1- where (fx) indicates the row and column indices, and Mic is the corresponding value assumed non-zero in the matrix. The set of observed indices j&)} Entries that are not observed, denoted by *, are considered
is defined as S₂ = missing or unknown. The positions of observed elements in the matrix M defined by binary observation matrix P, where P. = 1 indicates that M,y is observed, and Ps 0 indicates that My is unobserved.
Matrix factorization model. Matrix factorization is prevalent approach for addressing the matrix completion problem. It reconstructs the matrix W E Rd through the product W = ABB where A E Rdd and B Brxd This work studies the overparameterized scenario with r d, aiming to understand the implicit regularization effect in the absence of explicit rank restrictions, paralleling prior research (Gunasekar et al. 2017 Arora et al. 2019 Li et al. 2020 Jín al., (023). In this work, we focus on the asymmetric factorization, which can be represented as a parametric model:
fo = AB AB, A, E Rdd
(1)
The matrix factorization model parameters are denoted by 0 (A, B), identified with its vectorized
form vec(O) E 2d2 The augmented matrix is W = T B Rddx2d and row(A) and col(B) denote the row and column spaces of A and B. respectively. The augmented matrix Wa plays a crucial role in our subsequent analysis, particularly in characterizing the intrinsic invariant manifolds Ω2 of the optimization process. Specifically, it allows us to establish the relationship rank(A) A) = rank( B = rank W, which is important to understanding the invariance property under gradient flow.
aug
3

---

Loss function. The learning process for the parameters 0 = (A, B) involves minimizing loss function that measures the difference between observed and estimated entries. In this work, we focus on the mean squared error, and the empirical risk is thus formulated as
Rs(0)= = II M),,
ai b. M1z.se)
(2)
where a, and b. represent the i-th row and j-th column of matrix and B, respectively. The residual matrix 8M (AB M)s,, has elements 6M, (AB);; M, for (i.j) S² and 5M = 0 for (f.j) $ Sx. The training dynamics follow the gradient flow of (0)
d0 dt
= (o ()),
00) 8g.
(3)
In all experiments, 0, o2) initialized from Gaussian distribution with mean 0 and small variance o². We use gradient descent with a small learning rate to approximate the gradient flow dynamics (Please refer to Appendix B. for the detailed experiment setup).
4 Connectivity affects implicit regularization
In this section, we define connectivity and present experimental results on implicit regularization for connected and disconnected observational data.
Definition (Associated Observation Graph). Given incomplete matrix M to be completed and its observation matrix P, the associated observation graph GM is the bipartite graph with adjacency [o
PT] with isolated vertices removed.
matrix P PO.
Definition 2 (Connectivity). Given a incomplete matrix M be completed, is considered connected if its associated observation graph Gm is connected; otherwise. it is disconnected. The connected components of M are defined as the connected components of G
The connectivity of the graph, as defined above, reflects the connectivity of the observed data. Appendix A Sec. A.2 provides detailed discussion on the equivalent definition of connectivity.
In the case disconnectivity, there a special case where each connected component has full observations, characterized by disconnectivity with complete bipartite components.
Definition (Disconnectivity with Complete Bipartite Components). incomplete matrix M is considered disconnected with complete bipartite components if its associated observation graph Gm is disconnected and each connected component forms complete bipartite subgraph.
We present examples to demonstrate how connectivity influences the characteristics of the learned solutions. Consider three matrices be completed, each obtained by adding one more observation the previous matrix: M, (disconnected). M, (disconnected with complete bipartite components), and M3 (connected). Fig. Al of Appendix B illustrates the associated graphs GM.
1
2
*
[1
2 3
1
2
*
M =
3
* M2 =
*
* M₁ =
3
4
*
(4)
6
5
Figs. 2 a-bb compare the learned matrices with the ground truth (GT) solutions having the smallest nuclear norm and rank. For disconnected M1 (blue bars), the learned solution achieves neither the smallest nuclear norm rank. For disconnected M2 with complete bipartite components (green bars), the learned matrix has the smallest nuclear norm but not rank. For connected M3 (red bars), the lowest rank-2 solution is not unique; the model identifies particular lowest rank-2 solution, but does not correspond to the one with the minimum nuclear norm.
To thoroughly study all possible cases, we examine all sampling patterns of the 3x3 matrix completion. Fig. c)) shows that the model consistently learns the lowest-rank solution for connected sampling patterns but fails to do so for disconnected patterns. Fig. 2 d) further verifies the impact of connectivity on low-rank matrix recovery by comparing the reconstruction error for 100 randomly sampled rank-1 matrices using two connected sampling patterns (red and blue dots) and one disconnected sampling

---

12.5
G7:1
612
Nodbar horm
20.0
165
Lowest cams
Connected
10'
Cornectiina Decernuevily
20m1
2.5
0.9
10m
t0 Tatmandonly gerented
Samp dbt
(a) Nuclear norm
(b) Rank
(c) 3 x 3 matrix
(d) 10 X 10 matrix
Figure 2: (a) Nuclear norms of the learned solutions for M.. M2. and M3. Dashed lines represent
theoretically computed smallest nuclear norms. (b) Singular values of the learned matrices for
matrix. The thick vertical lines partition significantly nonzero singular values, which serves the empirical rank. The text (GT) shows the ground truth minimum rank. Mcan and standard deviation are recorded over 100 repetitions. (c) All equivalent sampling patterns of the 3 3 matrix completion problem (see Appendix B for details). Cyan stars marked the case learning the lowest-rank solution. (d) Reconstruction error the solutions for 10 10 matrix reconstruction problem with M' randomly sampled at rank 1 and sample size set to the minimum reconstruction setting = 2rd r².
M1, M2. M3. Each set of three bars represents the singular values
pattern (green dots). The model consistently achieves small reconstruction errors under connected sampling patterns, while the error is significantly larger for the disconnected pattern.
These empirical results demonstrate an implicit preference for low rank induced by connectivity and a preference for low nuclear norm in particular kind of disconnection. In the following section, we will investigate the training dynamics under both connected and disconnected scenarios.
Training dynamics in connected and disconnected cases
5.1 Connected case
This section empirically demonstrates the detailed dynamics of connected observed data. Fig. 3 shows the connected target matrix M with a single unknown element denoted by *. The rank of M is at least three and equals three if and only if 1.2.
Learning lowest-rank solution. We initialize A and B with different scales and record the singular values of the learned matrix. As depicted in Fig. ()), when starting with larger initialization, the learned solutions are almost always rank-4. Conversely, as the initialization scale decreases, the first three singular values of the learned solution are consistently maintained in magnitude, but the fourth singular value keeps decreasing, resulting in the model learning the lowest rank-3 solution.
Traversing progressive optima at each rank. For a small random initialization (Gaussian distribu- tion with mean and variance 10-16) (0-), the loss curves exhibit steady, stepwise decline (Fig. 3 c)).
The flat periods correspond small gradient norms, indicating potential saddle points (Fig. 3 d)).
We compare the matrices learned at these saddle points with the optimal approximation of each rank and plot their difference Fig. ((d), which is very small. These findings suggest that the model starts near 0 (rank-0) and progressively finds optimal approximations within rank-1, rank-2, and higher-rank manifolds until reaching a global minimum.
Alignment of the row space of and the column space B. Starting with small initialization. we track the rank (number of significantly non-zero singular values) of W = AB, A, B, and the augmented matrix aug during the training process. We observe that the rank gradually increases. with singular values growing rapidly one after another (Fig. (e)-h)) Throughout the entire process, we consistently find that rank( A) rank = rank W aug)) which implies that the row space of A and the column space of B remain aligned at all times. This alignment corresponds to special structure that we refer as the "Hierarchical Intrinsic Invariant Manifold" in Sec. 6. which plays crucial role in the overall dynamics of the system.
5

---

4
6
0.6
0.9
2.7
1.8
0.8
10
g- value
16
11-
Tork
X
Suddl point
Offercnod Neam
2.2
2.6
1.6
8
2.7
5.1
3.6
10
1244
10-4
10-1
(a) Matrix completion
(h) Initialization scale
(c) Training loss
(d) Gradient norm
30
30
30
(c) Singular values of W (f) Singular values of A Singular values of B (h) Singular values of W Y
Figure 3: (a) The matrix M to be completed, with the position unknown. (b) The four singular values of the learned solution different initialization scale (Gaussian distribution, mean 0, variance from 10 to 10-16). (c) Training loss for 16 connected sampling patterns in a 4 x 4 matrix, each covering element and observing the remaining 15 in fixed rank-3 matrix. (d) Evolution of the 12-norm of the gradients throughout the training process. The cyan crosses represent the difference between the matrix corresponding the saddle point and the optimal approximation each rank. (e-h) Evolution of singular values for matrices W, A, B, and Waug during training.
The dynamics increasing ranks step by step aligns with the description of Greedy Low Rank Learning (GLRL) (Li et al. 2020). However, we will show next that when the observed data disconnected, the learning process is not equivalent to GLRL.
5.2 Disconnected case
In this section, we present typical experiment in the disconnected situation. As depicted in Fig, 4aa), the target matrix M contains four unknown elements denoted by * and is disconnected. The rank of M is at least one, and there are infinitely many rank-1 solutions.
Alignment of the row space of A and the column space of B. As shown in Fig. 4 b-e). the learning process the disconnected case similar the previous experiment: the model naturally evolves from low-rank to high-rank, with each step increasing singular value and satisfying rank(AA rank B rank W "mug). Fig. 4( illustrates that the initialization scale decreases. the model tends to learn symmetric solutions. However, unlike the connected case, the output does approach a particular solution as the initialization decreases. For this specific disconnected M. we will show that every symmetric solution learned is a minimal nuclear norm solution(see Sec. 6.2 Thm.
For fewer observations, the experimental phenomena are similar (see Appendix Fig. B5)
Lowest-rank solution is not learned. Despite the adaptive learning behavior, the final learned solution has rank 2, as evidenced by the two significantly non-zero singular values Fig. b-d). Examining the dynamics B3. we find that they decouple into two independent systems: one for the 1st and 3rd rows of A and columns of B, and another for the 2nd row of A and column of B. Fig. shows that the model first learns the surrounding elements 1, 3, 3, 9 (rank-1 saddle point), then learns the middle element 5 in the next stage. The decoupling of dynamics is equivalent to the definition of disconnection (see Appendix explored all nine disconnected sampling patterns with 5 observations. For each pattern, we conducted experiments with small initializations. The loss curves consistently indicate that in disconnected cases, the model learns a sub-optimal solution in the rank-1 manifold, ultimately resulting in rank-2 solution. This demonstrates that regardless of the specific disconnected sampling pattern, the model fails to achieve the optimal low-rank solution.
proof). Fig. {e), we fixed rank-1 matrix and

---

*
3
10c
12-
127
*
5
*
3
*
9
10-
10
10"
Epsch
Egocch
Epech
(a) Matrix completion (b) Singular values of
(c) Singular values B (d) Singular values of W Vong
lea
1,0000
0. 0.0000
1,0000
1,0000
0.0000
3,5000
Curpon sontion
rank-1 subroptinal
0.0000
0.
0.0000
0,0000
5,0000
0. 0.0000
rank-2
1092
5,0000
0.0000
9,0000
Euoch
position
(c) Training loss
(f) Output at
(g) Saddle point
(h) GLRL solution
Figure 4: (a) The matrix to be completed, with unknown entries marked by *. (b-d) Evolution of singular values for A, B, and Wang during training. (e) Training loss for 9 disconnected sampling patterns in 3 3 matrix, each covering 4 elements and observing the remaining 5 in fixed rank- matrix. (f) Learned values at symmetric positions (1, 2) and (2, 1) under varying initialization scales (zero mean, varying variance). Each point represents one of ten random experiments per variance; labels show initialization variance. Other symmetric positions exhibit similar behavior. (g) Learned output the saddle point corresponding to the red dot (e). (h) Final learned solution the GLRL algorithm (Li et al. 2020).
Not equivalent GLRL in disconnected case. We compare the GLRL algorithm (Li al. 2020 with the matrix factorization model for solving the same matrix completion problem (Fig. al. 2020) claim that the matrix factorization dynamics mathematically equivalent to the RL.I. algorithm under reasonable assumptions. While GLRL learns the same rank-1 saddle point shown in Fig. 4 the first stage, then fills unobserved elements with 0, resulting unique rank-2 solution 4h)). In contrast, the matrix factorization model learns symmetric solutions with some degree freedom depending on the random seed (Fig. (t). The key difference is that the first critical point (Fig. g)) reached by the trajectory is a sub-optimal and not second-order stationary point of the rank- manifold as assumed by Li et al. (2020). Therefore, the equivalence assumption between GLRL and matrix factorization does not hold in the disconnected case.
6 Theoretical analysis of training dynamics and implicit regularization
6.1 Characterizatiion of training dynamics
Matrix factorization models exhibit distinctive adaptive learning behavior, progressively evolving from low rank to high rank. Understanding this phenomenon is rooted in grasping the global dynamics of matrix factorization models, where the role of intrinsic invariant manifolds becomes critical.
Proposition (Hierarchical Intrinsic Invariant Manifold (HIIM)). (see Appendix A Prop. A.7 Proof) Lei fo = AB be a matrix factorization model and {a),.. ,0x) be k linearly independent vectors. Define the manifold Zy as Ωr := Pe(o...... .ra) spa. ,,x,). The manifold Sr possesses the following properties:
(0 (A, B) row( A) col =
(i) Invariance under Gradient Flow: Given data S and the gradient flow dynamics θ - (R)
if the initial point 0g E Sik then ((t) E Ωy for all t 2 0.
(ii) Intrinsic Property: Sy is a data-independent invariant manifold, meaning that for any data S.
Sr remains invariant under the gradient flow dynamics.
(iii) Hierarchical Structure: The manifolds L2 form hierarchy: Ro S U2 _ ... 5. 5Pa-1 & Sk.
7

---

Figs. f and Figs. (b-d) show that the training process with small initialization consistently satishies rank(A) r rank( ang)). aligning with the Z2 invariant manifold. Since non-zero initialization in practice, the training trajectory is close to the Q2 invariant manifold, approaches a critical point, and transitions the next level invariant manifold without getting trapped.
In both connected and disconnected scenarios, we observe step-by-ster hierarchical 22 invariant manifold traversal. the connected case, each level observe that the model reaches an optimal solution (Fig. 3). However, in the disconnected case, we can prove that each connected component induces a sub-52; S2g invariant manifold, leading the experimentally observed sub-optimal solution (see Fig. 4..
Proposition 2 (Intrinsic Sub-S2; Invariant Manifold). (see Appendix [ Prop. A.2 for Proof) Let
fe AB be matrix factorization model. M be an incomplete matrix and S2y, be an invariant manifold defined in Prop. [ If M is disconnected with m connected components, then there exist m sub-St, manifolds w.e such that wi & S4 each possessing the following properties:
(i) Invariance under Gradient Flow: Given data S and the gradient flow dynamics θ = - if the initial point 0g E Wk, then (t) E we for all 2 O.
(ii) Intrinsic Property: Wh a datta-vauu-iin calue-independent invariant manifold, meaning that for fixed sampling pattern in M and any observed values S, wP. remains invariant under the gradient flow.
(iii) Strict Subset Relation:
The output set {fo O w/} is proper subset of {fe 0 12x).
namely, {fe 0 E (w) 2 (fo |0 3.
Optimal within each rank
Initialization
Sub-optimal within each level
135
Initi Final
Sub-O, Inveriant Manifold
Evalad
220 :20 200-200
Align to one deoctoon
Hierarchical Intrinsic Invariant Manifold ft2
225
315'
Hierarchical Intrinsic Invariant Manifold 0z
(a) Illustration of training trajectories in disconnected case
(b) Alignment of row (A) and col(
Figure 5: (a) Illustrated trajectories for the experiment in Fig. The blue line represents the trajectory converging to the lowest-rank solution, and the red line represents the actual trajectory experienced by the model. (b) The parameter trajectory escaping from second-order stationary point to reach the next critical point for the experiment in Fig- The scatter points represent the 4 row vectors of matrix A and the 4 column vectors of matrix B. For ease visualization, we randomly project them onto dimensions and plot them in polar coordinates.
Fig. 5(a) illustrates the trajectory of the experiment in Fig. In the disconnected case, sub-12g invariant manifolds exist and attract the dynamics, leading the model to learn sub-optimal solutions on the entire S2 invariant manifold. In fact, can prove that these sub-optimal solutions necessarily strict saddle points. This loss landscape result extends Theorem 5.10 from Li et al. (2020), which established the findings for the specific case of symmetric matrix factorization models (see Appendix A Sec. [A.3 a detailed discussion).
Theorem (Loss Landscape). (see Appendix A Thm. A.3 for Proof) Given any data S, the critical points of (0) are either strict saddle points global minima.
Gradient descent easily escapes saddle points (Lee et al. 2016 2019). Fig. 5]b) shows that when the model escapes a saddle point, the parameters initially appear chaotic but align in one direction after some time, consistent with the "condensation" phenomenon in neural networks (Luo et al. 2021 Zhou et al. 2022). For matrix factorization models, by meticulously analyzing the Hessian matrix structure (see Appendix A.5) we find that this alignment corresponds to an Ω2 invariant manifold,

---

resulting in a rank increase of one at time. Under reasonable assumptions, we prove that the training trajectory follows the Stg invariant manifold step by step.
Assumption (Unique Top Singular Value). Let öM (A.B. M) be the residual matrix at the critical point 0g = (Ac, B.). Assume that the largest singular value of öM is unique.
Assumption 2 (Second-order Stationary Point). Let St be an Qu invariant manifold or sub-S1k invariant manifold defined in Prop. 0r Assume 0. is second-order stationary point within Z2 Le., V Rf2) = O and 0 2 Rs( (0.) 0 2 O for all E Ω.
Theorem 2 (Transition to the Next Rank-lever Invariant Manifold), (see Appendix A Thm. A.4 for proof) Consider the dynamics i - (s)) Let (p 1) denote the value of (t) when (()) = Og.
Let 22 be an S2g sub-Sty invariant manifold. Let O. S2 be critical point satisfying Assump.
and 2 Then, for randomly selected Og. with probability with respect to Oo, the limit
log 2) A1
1
!
@(0c. 1) == liim vp a-0
(o. (Pa+aloy a00, +
(5)
a
exists and falls into an invariant manifold S2+11. Here X1 is the top eigenvalue of - Rss
Proof sketch. The main idea to analyze the local dynamics near the critical point 0.. The nonlinear dynamics can be approximated linearly in the vicinity of 0g: de dt ( 0z), where H
-V² R.s (f.) is the negative Hessian matrix. For exact linear approximation, the solution is: O(t)
HT (0,-0.) 0.)+ -e. Let X1 > X2 > > X be the eigenvalues of H, with corresponding eigenvectors
Zí-- We can express (t) as: (()) = Ei= Ei e
(0o Oc. 9:g) 0.. For sufficiently large
to, the dynamics follows dominant eigenvalue dynamics: 0 (to) Σラ= e- (Oo Oc. 4
O(e2zo)) Through detailed analysis of the eigenvalues and eigenvectors of the Hessian matrix (please refer to Lems A.2 A.4 of Appendix [A). we show that if the largest singular value of residual matrix 5M at 0g is unique and O. is a second-orde stationary point within 52, the first principal component 4
J= e eAAioo (0o- Oc, ) will happen to be an 22 invariant manifold. Consequently, escaping Or increases the rank by 1, entering 22+11
Remark. Assump. Mensures that upon departing from a critical point p.. the trajectory is constrained to escape along single dominant eigendirection corresponding to the largest singular value. This assumplion holds for randomly generated matrix with probability 1, making it a reasonable condition in most practical scenarios. In Sec Bec|A.7off Appendix A
we provide an special example illustrate the
situation where Assump. 7 does not hold.
Remark. To ensure the escape direction falls within the 24x+1 invariant manifold, the Hessian's 's top eigenvectors must satisfy rank( A) rank( B ) rank (ale(Wang Wagg).. The condition that 0g is a second- order stationary point within 52 in Assump. [ guarantees this Hessian structure. Our Assump. P is more general than conditions proposed by Li et al. (2020), as it remains valid across both connected and disconnected configurations. Empirical findings finding (Figs. 3 and 4) indicate that this assumption consistently holds in practical scenarios.
Thm. 2 provides a characterization of the escape trajectory. It shows that as the point approaches a second-order stationary point O. S2g, the trajectory generically converges to well-defined limit within 2 Since the origin O is always a second-order stationary point of Z. the theorem implies that the trajectory escaping from a small initialization will be close to G2g This iterative process gives rise to the phenomenon of Hierarchical Invariant Manifold Traversal (HIMT), which involves a sequential progression through these Ω2 manifolds.
6.2 Implicit regularization analysis
Rank minimization is a challenging non-convex optimization problem. Li et al. (0118) : Jin et al. (2023) proved that the Restricted Isometry Property (RIP) condition ensures a minimal rank solution. However, the RIP condition is often too stringent for practical matrix completion. For instance, the matrix M3 in Eq. (4) does not satisfy the RIP criteria, yet the model still finds the minimum rank solution. Our empirical findings (Figs. 22 3 suggest that more lenient condition, specifically the connectivity of the observed data, trequently leads to convergence towards the minimal rank solution. Proving this result directly, however, would necessitate comprehensive examination of

---

the convergence characteristies within each Q2 invariant manifold, which is an endeavor we leave for future work. Despite this, our insights into the system's dynamics, i.e., hierarchical invariant manifold traversal, allow us to assert that if trajectory successfully navigates through the optimal on each rank-level invariant manifold S2g, a solution of minimal rank can be achieved naturally.
Theorem 3 (Minimum Rank). (see Appendix A Thm. A.5 for proof) Consider the dynamics 0 RR)(O), where O(t) (A(t), B())), and denote W, A(t))((), Assume W. achieves an optimal within each invariant manifold Sy For a full rank initialization Wo. if the limit
W = -0 W (o exists and is a global optimum with W
Aj M2 for all (i,j) E Sw. then
W E argminw (ank(W) S.L. Wij M,j. V(,jj) E S
(6)
For disconnected matrix M, theoretical results (Prop. 2 and experiments (Fig. 4 confirm the existence of sub-Ω2, invariant manifolds. These manifolds attract the training trajectory, leading to sub-optimal solutions and preventing convergence to the lowest-rank solution.
However, in a specific disconnected case, such as disconnection with complete bipartite components, as illustrated in Figs. and 2] the minimum nuclear norm may still serve as characterization. Gunasekar et al. (2017) proved à special case: the observations are commutative, then the symmetric model will learn the minimum nuclear norm solution. Intriguingly, for the example M2 in Eq. 4. even though the observations are not commutative, the model still learns a minimum nuclear norm prove the following result, which aligns well with practical experiments.
solution. In fact,
Theorem 4 (Minimum Nuclear Norm Guarantee). (see Appendix Thm. A.olfor proof) Consider the dynamics Ô - (O),,, where O()) (A(I), B()),, and let W, A()) If the observation graph ussociated with the incomplete matrix M disconnected with complete bipartite components, and if for a full rank initialization Wo, the limit W = lima-10 (a: exists and is global optimum with Wij Miyoor all (i,j) Sw, then
W E argminw [|Il-
s.1. Wi Ma, v(i.j) E S..
7 Conclusion and future work
This study presents a comprehensive experimental and theoretical investigation of matrix factorization models. The primary objective was to develop a cohesive framework for understanding the conditions, mechanisms, and reasons behind the diverse implicit regularization effects exhibited by matrix factorization models. A key finding of this research is the pivotal role of the connectivity of observed data in shaping the implicit regularization behavior. To elucidate this phenomenon, we identified the significance of hierarchical invariant manifold traversal within the training dynamics.
Our experiments (Figs. I 3 provide strong evidence that connected observed data leads to minimum-rank solutions, as the model learns the optimal of the 22 invariant manifold. However, further investigation needed to uncover the underlying mechanisms by which connectivity facilitates optimal attainment across different f2 invariant manifolds. Additionally, the trade-offs between initialization scale and training efficiency warrant further research, as certain cases may require extremely small initialization, potentially impacting training speed (see Appendix B Sec. B.4)
Generalizing the insights gained from matrix factorization models to other architectures is also important avenue for future work. Our preliminary experiments indicate that the learning phenomenon from low rank high rank persists in deep multi-layer matrix factorization and the query-key factorization model in Transformer attention mechanisms (see Appendix B Figs. B9 B 111 These
findings suggest that the hierarchical invariant manifold traversal process uncovered in our study may have broader implications and merit further exploration.
Acknowledgments and Disclosure of Funding
This work sponsored by the National Natural Science Foundation China Grant No. 12101402, the National Key R&D Program of China Grant No. 2022YFA 1008200, the Lingang Laboratory Grant No. L.--QS-202202-08, Shanghai Municipal of Science and Technology Major Project No. 2021SHZDZX0102.
10

---

References
Z. Li, Y. Luo, K. Lyu, Towards resolving the implicit bias of gradient descent for matrix factorization: Greedy low-rank learning, in: International Conference Learning Representations, 2020.
C. Zhang, S. Bengio, M. Hardt, B. Recht, O. Vinyals, Understanding deep learning requires rethinking generalization. iclr 2017, arXiv preprint arXiv: 01116611133000 (2017).
C. Zhang, S. Bengio, M. Hardt, B. Recht, O. Vinyals, Understanding deep learning (still) requires rethinking generalization, Communications of the ACM 64 (2021) 107-115.
S. Gunasekar, B. E. Woodworth, S. Bhojanapalli, B. Neyshabur, N. Srebro, Implicit regularization in matrix factorization, Advances in neural information processing systems 30 (2017).
S. Arora, N. Cohen, W. Hu, Luo, Implicit regularization in deep matrix factorization, Advances in Neural Information Processing Systems 32 (2019).
N. Razin, N. Cohen, Implicit regularization in deep learning may not be explainable by norms, Advances in neural information processing systems 33 (2020) 21174-21187.
Z. Ji, M. Telgarsky, Gradient descent aligns the layers of deep linear networks, in: 7th International Conference Learning Representations, ICLR 2019, 2019.
Gunasekar, J. D. Lee, D. Soudry, N. Srebro, Implicit bias of gradient descent on linear convolutional networks, Advances in neural information processing systems 31 (2018).
G. Gidel, F. Bach, S. Lacoste-Julien Implicit regularization of discrete gradient dynamics in linear neural networks, Advances in Neural Information Processing Systems 32 (2019).
D. Gissin, S. Shalev-Shwartz, A. Daniely, The implicit bias of depth: How incremental learning drives generalization, in: International Conference Learning Representations, 2019.
L. Jiang, Y. Chen, L. Ding, Algorithmic regularization in model-free overparametrized asymmetric matrix factorization, SIAM Journal on Mathematics of Data Science 5 (2023) 723-744.
M. A. Belabbas, On implicit regularization: Morse functions and applications to matrix factorization, arXiv preprint arXiv:2001.04266 (2020).
Y. Li, T. Ma, H. Zhang, Algorithmic regularization in over wer-parameterized matrix sensing and neural networks with quadratic activations, in: Conference On Learning Theory, PMLR, 2018, pp. 2-47.
D. Stöger, M. Soltanolkotabi, Small random initialization akin to spectral learning: Optimization and generalization guarantees for overparameterized low-rank matrix reconstruction, Advances in Neural Information Processing Systems 34 (2021) 23831-23843.
J. Jin, Z. Li, K. Lyu, S. S. Du, J. D. Lee, Understanding incremental learning of gradient descent: A linengrained analysis of matrix sensing, in: International Conference on Machine Learning, PMLR, 2023, pp. 15200-15234
Y. Zhang, Z. Zhang. L. Zhang, Z. Bai, T. Luo, Z.-Q. J. Xu, Linear stability hypothesis and rank stratification for nonlinear models, arXiv preprint arXiv:2211 1.12211111612 1623 (2022).
Y. Zhang, Z. Zhang, L. Zhang, Z. Bai, T. Luo, Z.-Q. J. Xu, Optimistic estimate uncovers the potential of nonlinear models, arXiv preprint arXiv:2307.0800 (2023).
A. Jacot, F. Gabriel, C. Hongler, Neural tangent kernel: Convergence and generalization in neural networks, Advances in neural information processing systems 31 (2018).
I. Chizat, F. Oyallon, Bach, On lazy training in differentiable programming. Advances in neural information processing systems 32 (2019).
Luo, Z.-Q. Xu, Z. Ma, Y. Zhang, Phase diagram for two-layer relu neural networks infinite- width limit, Journal of Machine Learning Research 22 (2021) 1-47.

---

H. Zhou, Z. Qixuan, T. Luo, Y. Zhang, Z.-Q. Xu. Towards understanding the condensation of neural networks at initial training, Advances Neural Information Processing Systems 35 (2022) 21 2184-2196.
Y. Zhang, Z. Zhang, Luo. Z. J. Xu, Embedding principle of loss landscape of deep neural networks, Advances in Neural Information Processing Systems 34 (2021)
Y. Zhang, Y. Li, Z. Zhang, T. Luo, Z. J. Xu, Embedding principle: A hierarchical structure of loss landscape of deep neural networks, Journal of Machine Learning (2022) 113.
Z. Bai, T. Luo, Z.-Q. J. Xu, Y. Zhang, Embedding principle in depth for the loss landscape analysis of deep neural networks, arXiv preprint arXiv:2205, 13283 (2022).
K. Fukumizu, S. Yamaguchi, Y.-i. Mototake, M. Tanaka, Semi-flat minima and saddle points by embedding neural networks to overparameterization,, Advances in Neural Information Processing Systems 32 (2019) 13868-13876.
B. Simsek, F. Ged, A. Jacot, F. Spadaro, C. Hongler, W. Gerstner, J. Brea, Geometry of the loss landscape in overparameterized neural networks: Symmetries and invariances, in: Proceedings of the 38th International Conference Machine Learning. PMLR, 2021. pp. 9722-9732.
A. Jacot, Ged, B. Şimşek, C. Hongler, F. Gabriel, Saddle-to-saddle dynamics in deep linear networks: Small initialization training, symmetry, and sparsity, arXiv preprint arXiv:2106. 159333333 (2021).
Lee, M. Simehowitz, M. Jordan, Recht, Gradient descent only converges minimizers, in: Conference on learning theory, PMLR, 2016, pp.
J. D. Lee, I. Panageas, G. Piliouras, M. Simehowitz, M. I. Jordan, B. Recht, First-order methods almost always avoid strict saddle points, Mathematical programming 176 (2019) 311-337.
Du, J. Lee, On the power of over over--aramaan parametrization in neural networks with quadratic activation, in: International conference machine learning, PMLR, 2018, pp. 1329-1338.
12

---

A Proofs of Theoretical Results
In this section, we give all proofs for theoretical results mentioned in the main text.
A.1 Hierarchical Intrinsic Invariant Manifold and Sub Invariant Manifold
Proposition A.1 (Hierarchical Intrinsic Invariant Manifold (HIM)). Let fa= AB be matrix fuctorization model and {ar, 00k be k linearly independent veclors. Define the munifold Sh as S2 = Ag(or,,"" (Ok. ) = (0 (A, B) row(A) col(B) = span (t1," ,. The
manifold S2 possesses the following properties:
(i) Invariance under Gradient Flow: Given data S and the gradient flow dynamics 0 - RRs(O, the initial point Og- E Sia, then O(t) St for all t 0.
(ii) Intrinsic Property: Si is a data data-indepen independent invariant manifold, meaning that for any data S, Sy remains invariant under the gradient flow dynamics.
(iii) Hierarchical Structure: The manifolds S2y form a hierarchy: So S S2 - ... 52-1 1 f Szk.
Proof. (i) Invariance under Gradient Flow.
By definition, S2k Du(a1, ax)
{o
(A. B)
row(A)
col(B)
span(or, , ) Consider the gradient flow dynamics in
à
Σ
(a; b.j - M
jEl:
(8)
b.
Σ
(a: b. - My
ieh
where 1, = f/33: (1,5) E S²), 1 {+]]] (i.5) E Sa), a; and b. represent the i-th row and j-th column of A and B, respectively.
For any (i,j) E Sa, the evolution of a; is coupled with b.5 for j E I, The condition row(A) col (B) = span (an, .ok) ensures the existence of k linearly independent vectors "I,""", Org Rd such that a, b.j span(ar, ) for all 1 S i.j < d.
Consequently, if a; and b.j are initially in span .2...... ; within this subspace under the gradient flow dynamics. Additionally, for (t,j) & Sz. the gradients for the corresponding a, and b will be zero, provided their initial values are zero, maintaining this state throughout the evolution.
ax they will continue evolve
(ii) Intrinsic Property.
As demonstrated in part (i), S2 is invariant under gradient flow dynamics for any dataset S, confirming its status data-independent invariant manifold.
(iii) Hierarchical Structure.
Th invariant manifold 2g encompasses matrices of rank up to k, including those of lower ranks. Consequently, the manifolds exhibit the following hierarchical nesting:
Zg C S C ... ... S C Ωr.
Proposition A.2 (Intrinsic Sub-St Invariant Manifold), Let fe = AB be matrix factorization model, M be an incomplete matrix and Sr be an invariant manifold defined in Prop. [ If M is disconnected with connected components. then there exist m sub-St manifolds wy. such that We Dok. each possessing the following properties:
(i) Invariance under Gradient Flow: Given data S and the gradient flow dynamics θ = (Rs(O), the initial point O, E Wy, then O(I) wh. for all 2 0.
(ii) Intrinsic Property: wi a data-valee value- inlluee-dddeendentt invariant manifold, meaning that for a fixed sampling pattern in M and any observed values S. wy remains invariant under the gradient flow.
13

---

(iii) Strict Subset Relation: The output set {fe 0 wE is proper subset of (fe 0 E lsh. namely, (fo 0 wa) & (fo 9g
Proof. Existence.
Let consider an incomplete matrix M whose associated observational graph is divided into m connected components, denoted by L2, L2,.... Lm. For each component Ly, we define Sh" as the subset of observed indices within Ly. where 1 S S m and Sn is the set of all observed indices.
For each Lp. we can identify row indices Rp and column indices C₂ corresponding to the observed entries in Lp as follows:
Rp {i3 E E Cp = {j]i: (i,j) E "²").
Here, Rp includes the row indices and Cp includes the column indices of the entries observed in Ty
Define A Ue and _ as the submatrices A and B corresponding to Rp and Cp. respectively, and let A 1 and B be the remaining not in Rp and Cp-
Let 9eg := (01 aa.) (0 (A, row (A) col B) span (a1, 0ck be the given f2 invariant manifold.
The sub-S2; invariant manifold associated with the connected component 4, can be defined
w = ( (A, (ow(A'>) col col((B)>)) = span .a4.: A B = 0}. (9)
It is easy to check " is proper subset of S2g. (i) Invariance under Gradient Flow.
The condition (x) col((B span/ar, aa3 along with A B guarantees that a;, b.j E span (oi...,, 0 for all (i.j) E shr,, and a, b.g 0 for all (i,j) $ st
In other words, the sub-S2 invariant manifold " is the set of all pairs (A, B) where, for each observed position (i.j) in the connected component Lp. the vectors a, and b.j lie within the span of fai, .ax). ) and for any position not in Sh." the vectors zero.
Considering the dynamics expressed in equation (8). it is evident that the evolution of a; influenced by b. for (i,j) E shr Hence, if a, and b.j are initially in the span of (a11 022, ... oxs ) they will continue evolve within this span under the gradient flow dynamics. Moreover, for positions (i,j) & sh", we consider the following scenarios:
For (i & Sz, since the matrix entry Mij does not contribute to the loss Rs(0), the gradients for corresponding a, and b., will perpetually be zero. Thus, if their initial values are zero, they will remain zero throughout the evolution.
For (i.j) E Sg but sh,, the dynamics corresponding to different connected compo- nents are decoupled. Therefore, if the initial values for a; and b. are zero, they will stay zero during the evolution.
(ii) Intrinsic Property.
As established in (i), the manifold whr is invariant under gradient flow for any data S with a fixed sampling pattern, qualifying it data--valu-inddeeeneen invariant manifold.
(iii) Striet Subset Relation.
The output set (So 0 E 22 encompasses all matrices rank k whereas {So E whr is limited to rank-k matrices with specific row and column indices confined to Rd and Cp. Consequently, (fo E wx) forms strict subset of {fo }2 as stated by (fo E w (fo 0 )
0
14

---

A.2 Connectivity
Definition A.1 (Associated Observation Graph), Given a incomplete matrix M to be completed and its observation matrix P. the associated observation graph Gm is the bipartite graph with
0
PT]
adjacency matrix
0 0 with isolated vertices removed.
P
Definition A.2 (onneccivvtty). A matrix M to be completed is considered connected its associated observation graph GM is connected, otherwise, we call it disconnected. The connected components of M are defined as the connected components of this graph.
Definition A.3 (Disconnectivity with Complete Bipartite Components). A matrix M to be completed is considered disconnected with complete bipartite components if associated observation graph GM is disconnected and each connected component forms a complete bipartite subg subgraph.
Remark. In the bipartite graph representation of the observed data, isolated vertices correspond to entire rows or columns of the matrix M that are not observed. These rows columns do not contribute to the loss calculation and have no influence on the dynamics of the matrix factorization under infinitesimal initialization. Consequently, when analyzing the connectivity of the observed data and its impact on the learning dynamics, these isolated vertices can be safely disregarded.
Remark. The disconnectivity of the bipartite graph representing the observed data is equivalent
to the reducibility of the adjacency matrix O P
where P is the binury observation matrix indicating the positions of the observed entries in M.
P
0
In the context of matrix completion problems, such as the Netflix problem, connectivity has practical interpretation. Connected components in the bipartite graph indicate groups of users and movies that are linked by the users' viewing history. Users within the same connected component are related through the movies they have watched common. Due to this practical significance, we prefer to use the term "connectivity" instead of "reducibility" when discussing the structure of the observed data in matrix completion problems.
Definition A.4 (Connectivity of Observed Data). Given matrix M to be completed, undirected simple graph G can be induced from it: the nodes of the graph are the observed elements in the matrix, and two nodes are adjacent if and only if they are in the same or column of the matrix M. A matrix M to be completed is considered connected if its induced graph G is connected, otherwise, we call it disconnected.
Lemma A.1. For any simple graph G. if we remove all isolated vertices from G to obtain new graph G', then G' is connected if and only if the line graph of G', denoted as L(G"), is connected.
Proof. = Assume G' is connected. Consider any two nodes in L(G"), which correspond to two edges in G², say 61 and 62 Since G' is connected, there exists a path connecting the endpoints of er and 62. This path corresponds to a sequence of edges in G', which in turn corresponds to path connecting the nodes representing 61 and e2 in L(G"). Therefore, L(G') is connected.
= Conversely, assume L(G') is connected. Consider any two vertices ei and 02 in G'. Since G' has no isolated vertices, each of vn and 1/ is incident e2, respectively. Since I(G")) is connected, there exists a path connecting the nodes representing e and e2 in L(G'.. This path corresponds to a sequence of edges in G², which turn corresponds to path connecting U1 and U2 in G'. Therefore, G' is connected.
least edge. Let these edges be e1 and
In conclusion, we have proven that for any simple graph G, if we remove all isolated vertices from G to obtain a new graph G', then G' connected if and only if the line graph of G', denoted as L(G"). connected.
Proposition A.3. Given incomplete matrix M, the connectivity of M defined in Def. A.2 and Def. A.4 is equivalens.
Proof. By definition, each edge of a bipartite graph corresponds to an observed data item, and two edges in a bipartite graph are adjacent and only the two corresponding observed data items are in the same row or column. Therefore, the connectivity of the observed data is equivalent to the connectivity of the edges of the bipartite graph, which is, turn, equivalent to the connectivity of the line graph of the bipartite graph.
15

---

According to Lem. A.1 for any graph G, if we remove all isolated vertices from G to obtain new graph G'. then G' is connected if and only if the line graph G'. denoted as L(G"), is connected.
In the context of the bipartite graph representation of the observed data, removing isolated vertices corresponds to removing
and columns that contain observed entries. Thus, the connectivity of the bipartite graph after removing isolated vertices is equivalent to the connectivity of the observed data as defined in Def. A.4
Consequently, the connectivity the observed data as defined in Def. A.2 (based on the line graph of the bipartite graph) is equivalent the connectivity of the observed data as defined in Def. A.4
(based on the connectivity of observed data).
Definition A.S (Decoupling of Dynamics), Given an incomplete matrix M, consider the gradient flow dynamics of matrix factorization models, V1 i, S d,
à;
Σ
(a b.j Mij)b
jEL:
(10)
b.
2
Σ iEI,
(a: b. M.)aaa
The dynamics are said to be decoupled if there exist disjoint subsets of indices R1, Ro..... RR S {1, 2....... d} for the rows of A and C1, C2...., Cj 21. for the columns of B, such that for each E (1,2..... k) the dynamics of {a, : i Ri and {b jj E C1} form an independent system of equations. In other words, the dynamics can be divided into k(kk 1) independent subsystems, each involving subset of rows of A and a subset of columns of B. If such division is not possible, the dynamics are said be coupled.
Proposition A.4. Given incomplete matrix M, if is disconnected as defined by Def. A.2 then
the dynamics are decoupled as defined by Def. A.5 if it connected as defined by Def. A.2 then the dynamics are coupled as defined by Def. A.5
Proof. Consider matrix M be completed, with its associated observation graph comprising m connected components, denoted as L2, Lm. Let Shr C Se represent the subset of observed indices corresponding the connected component Lp. where 1 S S m and See denotes the complete set of observed indices. If the incomplete matrix M is disconnected, then for each connected component Lp. the subset 5 can be partitioned into two subsets R, and Cp. s m,
such that
R₂ Ei (i,j) str}. Cp {il3i (i,j) s
(11)
In other words, Rp contains the row indices and C₂ contains the column indices of the observed entries in the connected component Lp.
It can be easily verified that the dynamics are decoupled in this case, as the subsets Rpp, Cp Jm satisfy the conditions in Def. A.5 Each connected component Lp corresponds to an independent subsystem involving the rows of A indexed by Rd and the columns of B indexed by Cp.
If M is connected, then its associated observation graph consists of single connected component, and the entire dynamics are coupled.
Examples of connectivity and disconnectivity. Consider three matrices to be completed, each ob- tained by adding one more observation to the previous matrix: Mi (disconnected), M2 (disconnected with complete bipartite components), and M3 (connected).
2
*
1
M
3
M2
3
4
M2
3
4
(12)
*
*
5
6
*
The observation matrix Pis:
0
1
P
0
0
P3
0
P3
1
1
(13)
0
1
16

---

And the adjacency matrix is:
O A = Pi 0
PT
P2 P2
0
0
P T
A =
A3
(14)
P3
Given the adjaceney matrix A, can obtain graph GM- Fig. AT illustrates the associated graphs Gm. from which we can see that M1 disconnected, with its associated observation graph consisting of two connected components. M2 is also disconnected, but each connected component of its associated observation graph forms a complete bipartite subgraph. In contrast, M3 is connected, and its associated observation graph consists of single connected component.
(a) M
(b) M2
(c) M2
Figure Al: The associated observation graphs Gm of the incomplete matrices M1, M2, and M3 in Eq. 4. M, is disconnected, with its associated observation graph consisting of two connected components. M2 is also disconnected, but each connected component of its associated observation graph forms complete bipartite subgraph. In contrast, M3 is connected, and its associated observa- tion graph consists of a single connected component.
A.3 Loss Landscape
In this paper, we focus on the problem of asymmetric matrix factorization. Previous literature (Gu- nasekar al., 2017: Li et al. 2018 2020 Jin et al., 2023) has predominantly concentrated
symmetric matrix factorization problems. Although asymmetric matrix factorization models can be transformed into symmetric cases, studying symmetric matrix factorization does not necessarily cover all aspects of the asymmetric scenarios.
Generally, asymmetric matrix factorization model W situation by setting
AB can be transformed into symmetric
U
R 2ddd
B
We then consider the model W UU problem:
which corresponds to the following matrix completion
AAT B A
AB
B
B
We define the loss as
((A A
+ 1
L'
D
((CB))+;((()))
Li
2020
established the following results:
Theorem A.I (Theorem 5.10 in Li et al. (2020)). Let f Rdxd_ R be a convex -mmoohh function. (1). All stationary points of L: Raxd IR, L(U)= f (UUT) either strict saddles or global minimizers; (2). For any random initialization, GF converges to strict saddles of L(U) with probability
The proof this theorem relies heavily on Theorem A.2 of Du Du and Lee 2018), which requires the parameter matrix U E Bd satisfy the condition that & > d. In the case of symmetric matrix factorization, where U E gddd, this condition is naturally met. However, for asymmetric matrix
E p2dxd, this condition is not satisfied, and thus the proof of Thcorem A.I is only applicable to the symmetric case of matrix factorization.
factorization, where
B
17

---

Theorem A.2 (Theorem 3.1 in Du and Lee (018))). Let f: R - be a C² convex function. Then C Bd + R. L(U) f f(UT)) k 2 satisfies that (1). Every local minimizer of C is also global minimizer; (2). All saddles are strict. Here saddles denote those stationary points whose hessian are not positive semi-definite (thus including local maximizers).
Below we give a direct proof of the loss landscape of an asymmetric matrix factorization model.
Theorem A.3 (Loss Landscape). For any data S, the critical points of Rs(0) are either strict saddle points or global minima.
Proof. We start by recalling the definition of the loss function:
Rs(0) = C(A,B) AB M 2
Σ ((AB), - M.j)².
2
(i.JES
Let (A, B) denote critical point. We define matrix, 8M, the difference between the product of A and B and the matrix M. with this difference being computed only over the indices in the set Sn. More formally, 5M (AB M) where the elements of 6M are given by:
For (i,5) E S2., we have 5 MM = (AB); kJj Mig-
For (i,j) & Sn, we have MM =
This definition of 5M ensures that only consider the differences in the entries that belong the set Sg, while all other entries are set to zero.
Consider the function:
1
C(A + B B)
8M + EB + An + En
2
2 18M S. + (0M, + An) |I + Anli'_
(15)
+ (6M, Em) s= - =[|==| ||
where the inner product of two matrices A, B is defined as (A, B) = Tr(AB At the critical point, the first order term (6M, eB An)s. equals The Hessian operator, representing the second order term, is given by:
hA,B(E, 7) = || Anlis +
-
Our goal is to demonstrate that if 5M بf O, there always exists e, n such that h A.,((, n) To this end, we consider the ranks of matrices A and B in two cases:
(i) rank(A)
rank(B) d:
Without loss of generality, we assume 8M, := MM * 0 for some (i. E Sa, and rank < Under these conditions, there exists non-zero vector such that Av O.
We set 7 v and v7=- for s + j. where n" denotes the j-th column of the matrix n. Let
E = w E d and e, = 0 for s + i, where E denotes the i-th row of the matrix e.
We then have:
1
hA.BlE ,n")
|s B
+ w
2
1
<
w
8M T
2
1
BB
+ 5M, w To.
2 2
18

---

We define g(w, v)
BB" w + ww"BB'
+MM wiv and consider:
a( a(-aôMa;,0)
7a a a TBB
BB v a
1 2-0 v). ve)
z 2 M
BB
For a AB.max g(-aoM,, v v) < O.
where m,m represents the top eigenvalue of BBT, we find
Therefore, when
M obtain h ,n²) 0. This immediately implies the critical point = (A, B) is a strict saddle point.
(ii (aaakaA) rank (B) d:
Let adMB-1 and 7 0. In this scenario, the first order term (6M, An) s. in Eq. (15) simplifies to a MI At a critical point, this quantity equals zero, it implies that 8M = 0. This in turn implies that the critical point = (A, B) is a global minimum.
This concludes the proof, establishing that the critical points of Rg(0) either strict saddle points or global minima.
A.4 Escaping from Top Eigendirection
In this section, we focus on the dynamics of escaping from a critical point. According to Prop. A.3 the loss landscape consists solely of strict saddle points and global minimum. Consequently, gradient-based methods can readily escape from a critical point that is not a global minimum.
In the following, will demonstrate that the escaping dynamics approximated by a linearized version of these dynamics. For this, consider the following:
critical point be
Ô (O).
(16)
Assume 0, is saddle point for which (R, O. We can apply first-order Taylor expansion to the right-hand side of Eq. (16), yielding:
-Rs - VRs((.) v Rs 0.) + O 0.ll²).
(17)
where ( represents the Hessian matrix. Given that (s.( = 0, the gradient flow dynamics around O. be approximated
H - 0c),
(18)
where H = -- Eq. 18 is a classic linear ordinary differential equation, with the solution:
@(t) = "H - 02) + 0..
(19)
The dynamics near a critical point can be approximated by linearized version. Hence, in the vicinity of critical point, we can analyze the linearized dynamics understand the escape mechanism. In the following, we will show that during this escape process, the dynamics follow a pattern referred to as dominant eigenvalue dynamics.
Eq. 199 clucidates that the dynamics near the critical point 0g are predominantly dictated by the x2d Its eigendecompposition is given by:
properties of H, real symmetric matrix in x2d2
H == -== =AAQ
(20)
where A is a diagonal matrix and Q is an orthogonal matrix. Let A A2 > ... X E IR denote the eigenvalues of H. and let qu, 4:2 qu, represent the eigenvectors corresponding A,
Given that A > X2, the ratio /e for > 1 grows exponentially fast. Consequently, near O.. the evolution of the system primarily driven by eigenvectors Qı1, "12 qu1, associated with the largest eigenvalue \1-
19

---

This following proposition formalizes the intuitive idea that in the vicinity of a saddle point, the dynamies primarily follow the direction associated with the largest eigenvalue. This leading eigendi- rection becomes increasingly dominant as time evolves, allowing for an escape from the saddle point and facilitating specific structured transition.
Let's consider Og to be a saddle point. Consider:
= 2 R0(0) 0c), 0 - 0o
(21)
Here, 0g and 0. are in close enough proximity for the linearized dynamics be valid over sufficiently long period. We can then establish the following proposition:
Proposition A.5 (Escape from Saddle Points Following Dominant Eigenvalue Dynamics). Consider the linearized dynamics given by I 2 R, (0))( Qz), denote H := -=s Let d E TR be the largest eigenvalue of H, with corresponding eigenvectors q11, Q12, ... qu..- Denote Cj (6o- oz, quz), V1 S I h. Assume there exists such that cj # 0, then, given any
& > 0, there exists to > O such that for all t 2 to, the following holds:
O(t)
O.
(22)
11
edit Σ cjq1i
This means that, as time progresses, the direction of the parameter evolution increasingly aligns with the dominant eigenvectors of H.
This proposition is a consequence of the fact that the solution of the differential equation is given by (tt) (0)), and as t tends to infinity, the term corresponding to the dominant eigenvalue in the matrix exponential Ht becomes dominant. Therefore, near the saddle point, we have ((t) m
C9p + O..
(00-00.) 0.) O.. Here, Proof. The solution of the ordinary differential equation 21 is e() (00
H - (0.) a real symmetric matrix, which can be diagonalized. Let X > A > > X be the cigenvalues of H, and let qı1 9:2 -r qqt, be the eigenvectors corresponding to A₁. We can then express (t) as:
4
(tt)
e A Pc, qu + O
(23)
i= 1=1
Next, we analyze the norm of the relative difference between (()) and term dominating its growth:
4
Σ
Σ det (Oo Oo. q:s)
((1) O.
1=2 2=
Ax (0. tc, q1))
Σ (00 0c. 91))
5 4
(0. O., qiy)
(24)
jml
Σ 10 O.r 91a 4
5 4
00 0., 4:g 4
Se -0 x2)
Σ 00 0 c: 91 41
20

---

(0o- c: 9:5 9ij
log C
We define C =
Σ
By choosing to
ensure that for all
A1 22
2=
Σ (00 0., 9 q q
1=1
t > to the following condition is met:
(()) 0.
1 <.
(25)
c 0c: q1))911
Prop. A.5
describes that under the linearized dynamics, the parameters will escape from the saddle point along specific direction. However, when considering the original nonlinear dynamics @(t) -(),(0), we encounter a trade-off: we should choose la sufficiently large so that the trajectory can align well with the dominant eigendirection while escaping the saddle point, but if to is too large, the linearization approximation will fail as (to) moves away from 0.. Li et al. (2020) (Theorem 5.3) proved general dynamical result through careful analysis and error control: assuming the eigenvector corresponding to the maximum eigenvalue is unique and the initialization is sufficiently
close to the saddle point, there always exists suitable to such that the linear dynamics can align with the dominant eigendirection before the linearization breaks down. We can generalize this result to the case where the eigenvector corresponding to the largest eigenvalue is not unique:
Proposition A.6. Consider the dynamics given by @(t) - - Rs))), we use (0,,, 1) to denote the value of (t)) )) the case of (0) 0g- At a critical point Oz, we denote the negative Hessian as H = - (o.).). Let X1 E IR be the largest eigenvalue of H, with corresponding eigenvectors
711: 912 ... -hh Denote Cj (0g 0,, (1,)), VV S S 4, and V1 there exists j such that Cj + O. Let (t) = (0. + avi, Ai log 2)) for every a > 0, then xtt) == lima-+0 za(t) exists and is also a solution of the given dynamics, i.e., x(t) = (0), Furthermore, Vt E I. there exists constant C' >> O such that
=Ei= Cy91j- Assume
| 06 abo
(0.
A1-d2 1 log 2) 10-20 z(t) 5 Ca 221
for every sufficiently small where d1 22 > 0 is the eigenvalue gap.
Proof. By Theorem 5.3 in Section 5.1 of Li al. (2020), know that if the eigenspace corresponding to Aj is one-dimensional, i.e., 2 = 1, then the escaping direction will be the top eigenvector direction, and the convergence rate is (a AT- Therefore, Proposition A.6 holds in this case.
Now, if the eigenspace corresponding to Å₁ is not one-dimensional, denote cj (0o O.r (.)))), < < li and e1 = Σj Cy will the escaping direction. Following the same technique as in Li al. (2020),
A1-Aa easily verify that the convergence rate remains O(a -x Therefore, Proposition A.6 holds in this case as well.
A.5 Eigenvalues and Eigenvectors of Hessian
Suppose the dominant directions fulfill specific conditions, such any combination c1 C2 + c 9 leading to rank 1 model parameters (A, B). In such scenarios, we may observe a phenomenon where the rank of the matrix increases incrementally.
Firstly, we analyze the eigenvector structure of the Hessian matrix the critical point O_ (Ac, to understand why the parameter will enter the rank-1 invariant manifold.
21

---

Computation of the Hessian Matrix at Critical Point. To compute the Hessian matrix, we first consider the gradient:
Rs(0) = (fx,
(o (f Vefo(#).
d
Σ Es(fe Vo)
where 0,e (fe. f*) is the i-th element of (f(e, 0). (x))) and (fo); is the i-th element of the vectorization of fe.
For the Hessian matrix H0), we have
H(0) = VoVoRs(O)
EsVo (z) (o) + Σ (fo
d² Esdyy (So. 8 Vol (Vo(Jo),)
L
Σ
(fo); (Vo (fo),,
E Es(fo -
:==1
where &ey (So, f*) is the j.,j)tth element of (fl 0), (z))
We define matrices H (1) (0) and H (0) as follows:
( = Σ Vo (fo); (Vo (Fo)))
&
(0) =
Es(fe ft VoVo
We further denote that H(0) H ( (0) + (0).
For matrix factorization model, the eigenvectors of (2) has special structure, characterized by
(2)
Lem.
A.2
Lemma A.2 (Data-Independen) Interleaved Structure Eigenvectors (2))), Let O. (Ac, B.) be any critical point of the matrix factorization model. If X is eigenvalue of (0.) RRad x2d then there exist at least eigenvectors associated with These d eigenvectors take the form v t9 e1, u 3 62. u & eg E 2R2 where v E 2R2d is vector to be determined and e, the unit vector representing the i-th column of the identity matrix Is E Bgdd
Proof. Let's denote the residual matrix the critical point 5M (A.B. M) where (Ac. B.) a critical point. For the vectorized parameter 0.,, by direct calculation the matrix (2) = - V²Rg can be formulated as a block matrix, with the diagonal blocks being 0. The
specific format as follows:
-5M & la]
(2)
(26)
-MM & la
0
Next, we compute the eigenvectors of (2), Let À be an eigenvalue of (2), We need to verify that E1, e2, ... X ed 2d the eigenvectors of H corresponding to X, for particular
22

---

E 2dd yet to be determined. That is, we need to ensure that for all 1 S I d, the equation H [2) 222 (v @ es) has a non-zero solution for v. Notice that
PIK ② Is
-MM @ Id
H 2) 1222) & =
(0 2 ei)
8M X La
-XK-
0 La
-Xds
-8M
& la
(v @ @e)
(27)
5M
-Aa
Ald
- 5M
5M
-Aa
v
X C2.
Since is eigenvalue of H2),, the determinant of the matrix H (2)
222 equals zero. Hence
-XIa
-8M) det ([-MM
-VIa -MM
det 6M -Xtx
La
0.
(28)
-Ma
Consequently, from Eq. (27 we conclude that there always exists a non-zero vector RR2d such that (H) - X222) 2 e,) = O. Since v + 0, it is evident that v 8 C1. v & cg. ... , & ea E R linearly independent, and thus they represent eigenvectors corresponding X.
Proposition A.7 (Figenvectors Structure of H at the Origin). Consider the dynamics given by Eq. (2I). where we denote H - (s))..). If X is an eigenvalue of H E RRd x2d? then there exist at least d eigenvectors associated with in H. These d eigenvectors take the form @ e1, v @ "2, v & ed E g2d,, where E R2dd is a vector to be determined and e, is the unit vector representing the i-th column of the identity matrix Ia Rdxd
Proof. In the matrix factorization model, at the origin the gradient Va (fo) = 0 and thus ()) O and the Hessian matrix reduces to H (0), making H = 2 Rss (0) (2)) (0).
Let's denote the residual matrix at the origin as 5M = (A.B. M)... S.. where at the origin
(Ao: B.) = (0, 0). For the vectorized parameter Oc. the matrix H == - Rs can be formulated
block matrix, with the diagonal blocks being The specific format is as follows:
5M 09 1a
H 88M @ la
(29)
0
By Lem. A.2 the proof is completed.
Lemma A.3 (Eigenvectors Structure of H at Second-order Stationary Point). Let Ω denote an Sr invariant manifold or sub-Sty invariant manifold defined in Prop. A. and A.2 and consider a second-order stationary point 0. within Ω, i.e., VRs = 0 and Ø v2 Rss( (0.) 2 0 for all Z. Then, the eigenvectors corresponding to the negative eigenvalues of the Hessian matrix H (0,) are contained within the span of the eigenvectors corresponding to the negative eigenvalues of (2) (0z).
Proof. Recall the definitions of H and H (2)(()) (0) given by:
H (0) = Σ
Vo (fo); (Vo (Fo),)
J=
d
H
Es [fo - f"); V (fo).
t=1
The Hessian matrix (0)) at 0 is H()) = H ()) + H (2) (0).
The manifold Q is affine subspace with orthogonal complement denoted by - Let H have the following block representation in the bases of Z2 and Ω
H
H 12
(30)
H
Hz
23

---

Since 52 is an invariant subspace under the gradient flow, we have Hθ E S2 for all 0 E Z, which implies that 1 O. Since H symmetry, have H
Let < be a negative eigenvalue of = (02) with v as the corresponding eigenvector. Since H11 is positive semi-definite, v must lie 52+
At a critical point 0. (Ac, Be), by direct calculation, the gradient Vefe, can be structured as:
Be
B.
I T B A &
Vofo. =
011 I
I
a I a
(31)
012l
ag2
I
28 xd2
02a
a.d
where & denotes the Kronecker product.
Note that Vo (fo-); is the j-th column of matrix VoJe, and it falls precisely within the defined Oz invariant manifold or sub-S2; invariant manifold Ω. Therefore, we have:
(Vo (foo), v 0 v1 si
(32)
which implies that is orthogonal to the image of Vo (foz). placing it in the null space of (O.)
As result, we have:
v H (0.)v Av. H(O.) v (0.) + H H
Thus, the eigenvector the Hessian H (0.), corresponding to the negative eigenvalue also an eigenvector of (2) (0.), confirming that is within the span of the eigenvectors of (2) (0.). O
A.6 Transition to the Next Rank-level Invariant Manifold
Proposition A.8. The linear combination of the eigenvectors of (?); c1 @ e1) + c22(0 @ e2)
+ ca(o @ ea) falls within the invariant manifold (c),),, where (ci, 02, .a)
Proof. Notice that
C1 (v 2 G1) 02 ( 02) + ca(o ea)
= X |cje1 + 02e2 + .... .............. cdea]
(33)
v & c.
By Definition A.I the data data-independent invariant manifold generated by c 0 =
((A, B)lai, by E span(e), v1 í.j d). 0 (A, B) E Ri(c), then A, B must take the form
3
Ia+ic
3
Bu42C
=
B
(34)
B&c
B2dc
for some - :2..... 32d] 2g and the vectorized parameter 0 vec((A, ()) E 2d takes the form B @ c. Let B = v, and the proof is complete.
Lemma A.4. Suppose 0n1 a....... ,k-1 E Rd are linearly independent, the data- independent invariant manifold exhibits the property au(an, "2.." ai) Ω (a(+)))) 1x-1 .2..... -1)-
Proof. Assume that 0 (), E E (141((1110 02,, 0k41). Then A, B should adopt the form:
b ma + 7
3ia + B2410 B
(35)
i=1
24

---

[. Denote C; Yi
E RR" then the vectorized parameters 0 vec(O) can be expressed as:
k
G 2 az + Ck+1 @ ag+1-
(36)
Since know that Σ C. a E 12 02, ..)) and ck+1 @ ak+1 S2((a), it is straight- (=1 forward to validate that (211 022 ,x,11 = 22 02, .okk + l2 (o0++1).
Assumption A.1 (Unique Top Eigenvalue). Let 5M = (A.B. -)ss be the residual matrix at the critical point 0g = (Ac, Be). Assume that the top eigenvalue of the matrix unique.
-8M
-5M
is
Let Si be an S2 invariant manifold or sub-Sbg invariant manifold defined in Prop. A.7 or A.2 Assume 0. is second-order stationary point within
Assumption A.2 (Second-order Stationary Point).
S2 i.e. V Rs(O., (02) = O and 0 42 (s( 2 for all Ω
Theorem A.4 (Transition the Next Rank-level Invariant Manifold). Consider the dynamics (O)) Let cp(0,, t) denote the value of (t) when ((0) Og. Let be S2g invariant manifold or sub-St invariant manifold. Let 0. St be critical point satisfying Assump. A.7 and Then, for randomly selected Og with probability with respect Qu, the limit
(0 t) lim sp (0. (0. ato, t+ a-
- log a.) À
1
1
(37)
exists and falls into an invariant manifold S2k+1. Here A is the top eigenvalue of negative Hessian 2 Rg((
Proof. At the critical point o.., we denote the negative Hessian as - Rs Let A E be the largest eigenvalue of H, with corresponding eigenvectors 911, Q12, -911-
Denote Cj (0.- c.. qıy) v1 S S 41, and vı E/= cjq1j. For randomly selected 0b. with probability there exists at least one j such that Cj # O.
Consider the path (t)
( + ovi, t + 2 log _ for every > 0. By Prop. A.6 the limit =tt) == lima-s0 =a (t) exists and satisfies the dynamics xft) = (0),
Furthermore, Vt E IR, there exists constant C such that
p 0g + ato, t + A log _ x(t)
I Ca 2t -1-22
for every sufficiently small a, where \1 \2 > the eigenvalue gap.
This implies that the limit lima-s0 0. + tab, t 1 log _ exists and
@(0., t) = lim &o (0. a0,, -
1
log 2) lim -0
X log 2
0, + av
Assuming 0, E S2 and satisfies Assumps, A.l and A.2 invariant manifold 2x+1 containing 0. av.
we aim to show the existence a (ank-(k+1)
Define the following matrices:
H H =
Σ
Vo. (fo.); (Vo. (fo.)))
+,1=1
d
H (0) = Σ
Es [fo. - ff), Vg. (fo.).)
im
25

---

The Hessian matrix (0.) at 0g can be expressed as H(O)) H (0.) H (0,), and we have H (O.))
Assump. A.1 and Lem. A.2 imply that there exist exactly d eigenvectors associated with the top eigenvalue A of - H (0.). These eigenvectors are of the form @ e, for i .,.... d, where RR2d a vector to be determined and e; the i-th standard basis vector R By Assump. A.2 and Lem. A.3 the eigenvectors corresponding to λ of H are contained within the span of the eigenvectors associated with the negative cigenvalues of H (0.).
Prop. A.8 ensures that the escaping direction ei lies within rank-1 invariant manifold Qt Lem. A.4 then guarantees the existence of an invariant manifold 22++1 that includes 0. + a09 Since 2 is invariant under the gradient flow, the trajectory yp o01, + Xi log !) remains within 12
Finally, since 1 is closed subspace, the limit (0 8) lies in S2 concluding the proof.
A.7 Example of Coincident Top Eigenvalues
2 * IIn this case, the two numbers on the
Consider the x 2 matrix completion problem: M
2
diagonal are identical, which causes the maximum singular value of the residual matrix at the origin to be non-unique, violating Assump. [] Consequently, the training process will jump directly from the rank the rank invariant manifold, thereby missing the lowest rank solution of rank This behavior is demonstrated in Fig. A2 which shows experimental results for this scenario.
10
3
10-1
10
10
10
10
109
(a) Training loss
(b) Singular values of (c) Singular values of (d) Singular values of Way
Figure A2: Analysis of matrix completion for M with identical diagonal elements. (a) Training loss under small initialization. (b-d) Singular value evolution for A, B, Wang. Simultaneous growth of singular values results in direct convergence to rank-2 invariant manifold.
A.8 Minimum Rank
Theorem A.5 (Minimum Rank). Let Ω denote an invariant as defined previously. Assume W achieves a global minimum within each invariant manifold S2k- If the limit W lim (o) exists and is global optimum with (i.j) (i.,)) ()r all (i.j) Sw, then
W E argminw rank(W) s.l. W M M(..)). V(..)) E
(38)
Proof. Consider the invariant manifold 24, which is defined as follows:
Ω = f4 aa) {(A, B),a, b.j E spanfaa..... ax). V1 i,j S d}. where denotes the i-th row of A, b.y denotes the j-th column of B, and ai...., Ot are indepen- dent vectors that span the invariant subspace associated with 2t
According to Thm. A.4 the training trajectory adheres to a hierarchical traversal across invariant manifolds. For any matrix C with rank(C)) < k will show that there always exists Ωr such that AB = C. Therefore, Ωr contains all matrices of rank k.
(A, B)
In fact, Since rank (C) S k. we can express as a sum of rank-one matrices: C Σ ii u.n where u, and 2s column vectors. By the definition of Zg. for any 0 (A, B) E Sk. each
26

---

Tow of A and each column of B can be expressed as linear combination of (a (
a; Dj Cij(i b E. dija; ai where Gjj and diy scalars. We can write:
A =
B = Ex d1a Ct-idaack
Σ a
Now, we can express the product AB as: AB Σ j== EEiz ardy)) a; By choosing appropriate values for Cij and dij. we can make AB = C. This is possible because the outer products a span the same subspace as the rank-one matrices u in the expression of C.
Therefore, for any matrix with rank( S there always exists (A, B) E Sz such that AB = C.
If the output matrix W, attains optimums within each Stg it suggests that the optimization process selecting the best approximation from the set of all possible rank-k matrices. Provided that each step in the optimization is optimal, the resulting solution will naturally be the matrix with the lowest feasible rank that satisfies the matrix completion criteria, thereby completing the proof.
A.9 Minimum Nuclear Norm Guarantee
Lemma A.S (Minimal Nuclear Norm Computation). Given matrix M to be completed with observed diagonal entries, ie., tiag(M) e, the minimal nuclear norm solution among all possible completions is ||||
Proof. The nuclear norm of a matrix is the dual of the spectral norm || 1122, defined as:
|||| = max X).
XX 2251
Given that diag(sign(+)) )) < 1, for any matrix with diag( A) v, it follows that:
||| 2 (A, diag(sign(()))) = (v, sign (v)) = ||||1_
Specifically, the nuclear norm of the diagonal matrix with its diagonal is (diag(v)ll ||0| which establishes that the diagonal matrix with v indeed a minimizer for the nuclear norm.
Diagonal Observations
Proposition A.9 (Minimum Nuclear Norm Guarantee in Diagonal Case). Consider the dynamics
Rs where @(t) (A(t), B(t)) and denote W, B(t). If the observation data diagonal, and if for a full rank initialization Wo. the limit W = lim Wo (oW) exists and is global optimum with W Myj for all (i.j) E Sz. then
W E argminw [|
s.t.
Wy My, v(i.j) E Sz.
(39)
Proof. Without loss of generality, assume that M a diagonal matrix given by:
f1
12
M =
P
By Lem. A.5 the minimal nuclear norm among all possible completions is |uy] +|uaa Jua/++0 +.. | When the matrix be completed diagonal, evolution the i-th row A is influenced only by the i-th column of B. Hence, the dynamics decouple into d independent parts, each equivalent to
27

---

exists and is a global optimum with learning scalar 14. The learning process thus unfolds in d stages, with each stage passing through critical point to learn a respective pę-
By Lem.
A.2 the second term of the Hessian matrix can be expressed as:
0
-MM @ ld
H
6M ② la
0
According Lem. A.2 the eigenvectors of H (2) take the form a & C1, C2,.... ...,, ",, R2d where v E Rdd is a vector to be determined and e; is the unit vector corresponding to the i-th column of the identity matrix Ia Raxd
(i) Suppose > >
Pa 0.
With infinitesimal initialization, the training dynamics first focus on the element with the largest singular value, then proceed sequentially to blocks with smaller singular values. This pattern of learning is consistent with the concept of "sequential learning" as reported in the literature (Gidel ct al. 2019 Gissin et al. 2019 Jiang et al., (202)).
For diagonal observation matrix, the residual matrix 8M at any critical point remains diagonal matrix. While starting learn H4, from a critical point, direct calculation confirms that vector = [er,e]] E I22d. Escaping from each saddle point d.. the trajectory ( 0g approximates E C e.), which satisfies B, (A)),- Thus, learning diagonal matrix M using the asymmetric model AB is equivalent using symmetric model AAT. The final outcome ensures that diag( AB) diag( AA diag(MM
The nuclear norm AAT equals the sum of its eigenvalues, which precisely the trace of the matrix, and tr(AAT tr(M) = la + #2 +... Mda Therefore, the nuclear norm of the learned matrix W = AB = AAT remains | + [u2] ........ Ipal.
(ii) If some μ4 0, assume without loss of generality that |al |2a| >...> |un1 0.
While starting
learn H4; from a critical point, direct calculation confirms that
[e,, sign((;;);] E 2ddd Escaping from each saddle point 0., the trajectory 0(t) O approxi-
mates Σ Eis1 (v @ e,), satisfying B sign(u)((A)). Hence, AB AAT Q, where Q is an orthogonal matrix given by:
(i))(Hn)
sign(/2)
Q
(pan
The final result ensures that diag( (AB) diag( AATQ diag( [M) meaning diag(AA) (AAT)
diag (QM). The nuclear norm of AA equals the sum of its eigenvalues, which is the trace of the matrix, and
tr(QM) |يr |pa] +... + |pal-
Since an orthogonal transformation does not change the nuclear norm of a matrix, the nuclear norm of the final learned matrix W = AB = AATO Q is still |ps| + |12| + ... + [na]-
Disconnected with Complete Bipartite Components
Theorem A.6 (Minimum Nuclear Norm Guarantee). Consider the dynamics Ô -Vs(((())
where (tt) = (A(t), B()), and let W, B()) If the observation graph associated with the matrix M to be completed is disconnected with complete bipartite components, and if for full rank initialization Wo. the limit W lim Wao(aWo)
Wi = Mir all (i.j) E Sz. then
W E argminw [ s.t.
Wy Maj, V(i.j) E Sa.
(40)
Proof. Consider a matrix M E Rdxd composed of m connected components, with each component forming complete bipartite subgraph. Since M is disconnected, it can be represented in block
28

---

diagonal form without loss of generality:
Mi
M2
M
M
where each block M, E Rdxd, and Ei=1d d. = d. Ei=1 d' = d, representing the sum of the dimensions of the blocks.
Each block M, corresponds some singular values of the corresponding Hessian matrix at a critical point. With infinitesimal initialization, the training dynamics first focus on the block with the largest singular value, then proceed sequentially to blocks with smaller singular values. This pattern of learning is consistent with the concept of "sequential learning" as reported in the literature (Gidel ct al. 2019) Gissin al. 2019; Jiang et al., 2023).
Since each connected component of M forms a complete bipartite subgraph, the block M, is fully observed. We do singular value decomposition (SVD) on each sub-block M, M, U,E, V,". where U, and Vi are orthogonal matrices, and Σ, is a diagonal matrix with the singular values of M..
Construct block diagonal matrices U and V as follows:
U1
V1
U2
V2
U =
V =
Um
V m_
This leads to the diagonal matrix:
Σ1
M1
Σ²
P2
UMV T
Em
pa
Orthogonal transformations preserve the nuclear norm, so by Lem. A.5 the minimal nuclear norm among all possible completions is the sum of the absolute values of the diagonal entries, i.e., |P1] |22| +...+ |pa|-
Consider an incomplete matrix M whose associated observational graph is divided into m connected components, denoted by L1. Lz..... Lm. For each component Ly. we define s as the subset of observed indices within Ly. where p S m and S. is the set of all observed indices.
For each Lpp, we can identify row indices Rg and column indices Cp corresponding to the observed entries in Lpp as follows:
Rp Ei (i,5) S'r), Cp ( (1,5) sh").
Here, Rp includes the indices and C. includes the column indices of the entries observed in Lp-
Define A l. and 1 as the submatrices of A and B corresponding to Rp and Cp. respectively. The evolution of A 4 is influenced only by Bd,, Thus, the dynamics decouple into m independent parts, each equivalent learning fully observed matrix M..
Accordingly, we partition and B into m blocks:
A
A
B = B
Bnl,
A
where A, AL. and B; B4j
Denote L; = I; M; 12. The overall loss function = Σm L, can be decomposed into m
independent parts. Performing orthogonal transformations A, = U, A, and B; = B, we obtain
diagonal loss L, || A. E.1l for sisc.
29

---

Since gradient descent is the steepest descent in the l2 norm and orthogonal transformations preserve this norm, the dynamics of optimizing L; are equivalent to those of optimizing L,.
Without loss of generality, assume μ1 l2 orthogonal transformation Q can be applied without changing the nuclear norm.
pa > O. Otherwise, with Prop. A.9 sign
By Prop. A.9 the learning result for a diagonal matrix implies B, AT, which means B,V, U,T A,
The final learning result
A
=
= B1
Bm]
A
satisfies B = A
The result ensures that diag(AB) diag (AAT) diag( Σ). The nuclear norm of ÁA equals
the sum of its eigenvalues. which is the trace of the matrix, and Tr(AA Tr(E) |p1l + |
+.+ Ipal-
Since orthogonal transformations do not alter the nuclear norm of matrix, the nuclear norm of W = AB is also | + |pz| + +... [pal, concluding the proof.
B Experimental Setup and Supplementary Experiments
In this section, we present the supplementary experiments mentioned the main text and the details of experiments.
B.1 Experimental Setup
For all our experiments, we employ gradient descent with carefully chosen small learning rate. A learning rate is deemed suitable when yields a smooth, monotonically decreasing training trajectory for the loss function, free from any abrupt fluctuations or oscillations. We initialize all model parameters using a Gaussian distribution with mean of zero and variance that detailed for each specific experiment. Because of the small size of the experiment, the experiment can be completed on a single CPU.
The criterion for the sufficiency of training in all cases is a training loss that falls below 10 10-10 To ascertain the rank of the matrix produced by the learning process, we utilize a technique of extrapolation with an infinitesimally small initialization. As depicted in Fig. (3),,, if a singular value persistently diminishes in response decreasing initialization magnitudes, it then inferred that such a singular value will not contribute to the rank in the context of an infinitesimal initialization.
We have included code in the Supplementary Material that determines the connectivity of partially observed matrix and provides specific examples illustrating the implicit regularization effects. This code can be used to reproduce our results and explore the relationship between data connectivity and the implicit biases of matrix factorization models in various matrix completion scenarios.
B.2 Connectivity Experiments
In the connectivity experiments corresponding to Fig. [ we explore the behavior of randomly generated 4 x matrices with intrinsic ranks of 1, 2, and 3. To investigate the impact of sampling density on matrix reconstruction, we sample matrices at three different levels: 2rd - r² which meets the threshold for exact reconstruction, 2rd - 2rd ? + 1, which exceeds the threshold.
which is just below the threshold, and
For each sampling size, we randomly generate 10 sets of sampling positions. We then assess the connectivity of the sampled positions and compute both the rank and the nuclear norm the solutions obtained through gradient descent. As an illustration, in Fig. B1 panel (a) presents scenario with connected sampling positions, panel (b) shows disconnected sampling positions, and panel (c) depicts disconnected sampling with each disconnected component forming complete bipartite graph.
30

---

X
x
X
X
X
X
X
X
x
X
x X X
X
x
X
x
x
x
X
x
X
X
X
X
X
X
x
X
x
x
X
(a) Connected
(b) Disconnected
(c) Disconnected with complete bipartite components
Figure Bl: Examples of connected sampling and disconnected sampling patterns in Fig.
In the connectivity experiments depicted in Figs. 2 c)d,,, examine the behavior of randomly generated matrices of size x 4 and 10 x 10 with a rank 1. The matrices are sampled at size of 2rd r² which corresponds the threshold for exact reconstruction. We evaluate two connected and one disconnected sampling patterns.
Fig. B2 a) displays the first connected sampling pattern, where all entries in the first row and the first column are sampled. Fig. B2 illustrates the second connected sampling pattern, which forms a "Z" shape across the matrix. Fig. B2(c) shows the disconnected sampling pattern, where the samples split into two unconnected blocks, one in the top-left and the other in the bottom-right of the matrix. A similar approach is taken for the 10 X 10 matrices.
X
(a) Connected sampling 1
(b) Connected sampling 2
(c) Disconnected sampling
Figure B2: Examples of connected sampling and disconnected sampling patterns Fig. c).
For Figs. (a-b), we performed 100 random initializations for each initialization scale, recorded the mean and standard deviation, and plotted them on the figure. For sequence æ,, 22, ... Fya the mean is calculated by T Dn T4, and the standard deviation is calculated by:
n
I 2)2
72
Figs. ccdd demonstrate that when the target matrix has a rank of and the number samples meets the minimum requirement for reconstruction with connected sampling positions, the matrix factorization model is capable of accurately reconstructing the original target matrix.
In scenarios where the target matrix has higher rank, we have extended experiments accordingly. For randomly chosen x 4 matrix with rank 2, we selected sample count less than the threshold of 2rd , 12, specifically 10 samples, while ensuring that the sampling pattern is connected, as shown in Fig. B3 The resulting solution from the matrix completion has rank of 2, which the minimal rank that fits the sampled data.
Fig. B4 reveals distinct behaviors of the matrix completion depending on the scale of initialization. With larger initialization, the third and fourth singular values of the completed matrix remain
31

---

Figure B3: The sample pattern in Fig. B4
relatively significant, suggesting that the model does converge to the lowest rank solution. On the other hand, with smaller initialization, the third and fourth singular values are uniformly small, indicating that the model successfully converges to the lowest rank solution.
107
Singulle values
20
Singular vakess
107
Reconstruction entor
10"
Squib vilue
10
Sigjor vaue
Coneccivittl
10"1
Siguirr vauus
10
Connechvitt2
10-5
Sigu
Disconnectivity
10-4
10-7
10-4
10"4
20
Trailtrandmly generated matri))
Taillr hallirandomly generated matrix)
Trall Talltrandomly generated matrixt
(a) Initialization scale 10"
(b) Initialization scale a 10
(c) Reconstruction error
Figure B4: For randomly selected 4 x 4 matrix with rank 2, we chose 10 samples, fewer than the threshold 2rd - x² 12, ensuring connected sampling positions, as shown in Fig. B3] The figures show the experimental results for the four singular values of the matrix learned under Gaussian initialization with mean 0 and standard deviations of 10-1 (a) and 10 (b), respectively. (c) Reconstruction error of the solutions for 4 x matrix reconstruction problem with M randomly sampled at rank r = 1 and sample size set to the minimum reconstruction setting = 2rd r Red and blue scatter points represent two connected sampling patterns, while green points represent disconnected pattern.
B.3 Equivalent Sampling Patterns
For given sample size, there are different sampling models corresponding to connected or discon- nected sampling. As shown in Fig. BT 7 observations sampled, but different sampling positions affect connectivity or disconnection. To thoroughly study all possible cases, we examine all sampling cases 3 3 matrix completion, as illustrated Figure (
For a matrix, the sample size varies from 1 to 9. When using the matrix decomposition model fe = AB for matrix completion, the dynamics obtained by exchanging transposing the matrix to be completed are equivalent. These three operations allow us to divide all sampling patterns equally.
columns
For sample size there only one sampling pattern the equivalent sense, and the observation matrix P is:
0
0
P
0
0
where 1 indicates that the position observed and non-zero, and 0 means that the position is not observed or is O.
32

---

Lon
*
3
107
100
5
9
12-
16-1
10
Ayar
Fpech
(a) Matrix completion
(b) Singular values
(c) Singular values B (d) Singular values of Wang Ving
10
o at tnaee th, 01
Ie-16 Ie-12
1,0000
0,0000
3 3,0000
1,0000
0.0000
1,0000
sub-optimal
0.0000
0.0000
6.0000
5,0000
0.0000
12.4
rank-2
5.0000
0
(c) Training loss
(f) Output at
(g) Saddle point
(h) GLRL solution
Figure B5: (a) The matrix be completed, with unknown entries marked by (b-d) Evolution of singular values for A, B, and Waug during training. (e) Training loss for the disconnected sampling pattern. (f) Learned values symmetric positions (0, 1) and (1,0) under varying initialization scales (zero mean, varying variance). Each point represents one of ten random experiments per variance; labels show initialization variance. Other symmetric positions exhibit similar behavior. (g) Learned output at the saddle point corresponding to the red dot in (e). (h) Final learned solution of the GLRL algorithm (Li al. 2020).
For sample size 2, there are only sampling patterns in the equivalent sense:
1
01
0 0
P
0
0
0
P2
0 0
0
0
0
0
0 0
For sample size = 3, there are only sampling patterns in the equivalent sense:
1
1
1 07
1
1
0
1
0
07 0
0000
P
=
0 0
0
P2
= 1 0
P3
0 0
1
P1
=
0
1
0
0
0
0
0
0
0
0
0
0
1
For sample size 4, there are only sampling patterns in the equivalent sense:
1
1
07
1 0
I
07
0
P =
0
P
=
P3
=
P1
=
Ps
=
0
0
0
0
0
0
0
0
0
0
0
1
For sample size 5, there are only sampling patterns in the equivalent sense:
1
17
1
17
[
0
P
P.
=
P
0
0
0
1
For sample size 6, there only sampling patterns in the equivalent sense:
1
11
[
1
11
1
17
1 07
P
1
1
Pı
1
1
0
P3
1
0
P4
0
1
1
0
0
0
0
0
1
1
0
0
1
0
1
For sample size 7, there only sampling patterns in the equivalent sense:
1
1
1 1
P
1
1
1
Pa
1
0
0
0
1 1
33

---

For sample size 8, there is only sampling pattern in the equivalent sense:
1 1 0 B.4 Initialization Scale Analysis
[ P =
For sample size =,, there is only 1 sampling pattern in the equivalent sense:
1 P
1
Our experimental findings indicate that when the observational data connected, matrix factorization models often learn the lowest-rank solution starting from small initialization. However, the required scale of initialization not constant across different instances. We empirically observed that if the magnitude of the numerical values in the matrix to be completed varies significantly, an extremely small initialization is necessary, which, in some cases, can exceed machine precision.
Consider the following two simple 2 x 2 matrix completion problems, with the only difference being that the number 3 in the first row is replaced by 20. When training begins from small initialization, for M4 the fourth element only needs to learn the value 6 to be a rank-1 solution. However, for Ms, the fourth element needs to learn the value 40 to achieve rank-1.
1
2
2
M
M;
20
X
Fig. Billsstrrs the difficulty in learning these two examples. For M,, an initialization variance of approximately 10-7 is sufficient to learn the lowest-rank solution. In contrast, for Mj. extremely small initialization variance required, making it challenging to learn rank-1 solution. Yet, with exceedingly small initialization variance of 10 10-83
still observe the second singular value plummeting to zero. The origin is a saddle point. The smaller the initialization, the longer it will stay the origin. If initialization continues to decrease, training will stagnate. Therefore, if the magnitude difference is even greater, such as replacing 20 with 100 in Ms, then with the small initialization allowed by machine precision, is nearly impossible to learn the value 200 completely.
For the matrix factorization model fe = AB, the Hessian matrix at has strictly negative eigenvalues, making the origin a strict saddle point. Under small random initialization, gradient descent escapes this saddle at an exponential rate. Our Theorem ensures that only strict saddle points and global minima exist as critical points in the loss landscape. Subsequent saddle points invariant manifolds are also strict, with exponential escape speeds, maintaining acceptable optimization process.
Reaching the lowest rank solution requires parameters to escape the saddle point along the unique top eigen-direction. Fig. [6 illustrates the relationship between required initialization scale and obser- vation magnitude differences. For lowest possible rank with large numerical magnitude differences observations, extremely small initialization necessary. However, for approximately low rank solutions, relatively small initialization suffices.
Trail frandom seed!
10
Singlar
Dursue
10-1
frandom
1
10-17
10:30
(a) Output: M.
(b) Singular values: M4
(c) Output: Ms
(d) Singular values: Ms
Figure B6: (a, c) The value of the (2, 2) element learned by the model as the initialization decreases. (b. d) The singular values the matrix learned by the model with decreasing initialization.
34

---

B.5 High-dimensional Experiments
To validate the scalability of our findings, we extended our experiments higher-dimensional matri- ces. We conducted tests on 20 x 20 matrices, employing both connected (Fig. [B7] and disconnected (Fig. [B8) sampling patterns, while monitoring rank evolution during training. Our results consistently corroborated the main findings:
(i) Connected observations converged to optimal low-rank solutions.
i Disconnected observations yielded higher-rank solutions.
(ii) The Hierarchical Invariant Manifold Traversal (HIMT) process was observed in both connected and disconnected scenarios.
101
13-
121
101
101
10
10
20
10
10
300
300
30
Spach
(a) Training loss
(b) Singular values of A
(c) Singular values of B (d) Singular values of W
Figure B7: Connected sampling pattern analysis for 20 x 20 random rank-2 matrix completion
problem. (a) Training loss under small initialization. (b-d) Singular value evolution for A, B, W
mug.
10-1
121
101
101
10-4
191
191
191
10
12"
10
10
10-1
Epoee
(a) Training loss
(b) Singular values of A
(c) Singular values of B (d) Singular values of Wang
Figure B8: Disconnected sampling pattern analysis for a 20 20 random rank-2 matrix com- pletion problem. (a) Training loss under small initialization. (b-d) Singular value evolution for A, B, W augg-
B.6 Dynamics of Deep Matrix Factorization
In the context of depth-3 matrix factorization models, we consider the functional form:
fo ABC, where B, C gdxd
Figs. B9 and B B10 suggest that even for a depth-3 model, the learning process exhibits a progression from low rank to high rank structures.
B.7
Incorporating Attention Mechanisms
Within the Transformer architecture, the matrix factorization component retains its significance. The attention mechanism is formalized as follows:
X XWQ, W K.
fx - Σ softmaxrow
XWv. Wo.
Vdle
where the row-wise softmax operation applied the attention scores, and the sum is over the different attention heads, with WQ.. WK., Wv, Wo, representing the learnable weight matrices
35

---

4
0.6
1.8
0.8
10
6
0.9
2.7
(Px1.2)
8
2.2
2.6
1.6
rank-2
8
2.7
5.1
3.6
10-1
10-1
Wnlialaatance
(a) Matrix completion
Output
Training loss
(d) Adaptive rank
Figure B9: Deep Matrix factorization models learn adaptively from low rank to high rank with small initialization. (a) The matrix to completed, the position unknown. (b) The value at learned under different initialization scales (mean is zero, variance changes), 10 random experiments were done under each variance, and each blue point represents an experiment. (c) The training loss curve with an initial variance of 10-14, (d) Evolution of singular values for fe = ABC during training. The count significantly non-zero singular values indicative rank.
50
30
12
Gradinn 16111
Excth
(a) Gradient norm
(b) Singular values of A (c) Singular values of B (d) Singular values C'
Figure B10: Deep Matrix factorization models learn adaptively from low rank to high rank with small initialization. (a) Evolution of the 12-norm of the gradients for all parameters throughout the training process. (b-d) Evolution of singular values for matrices A, B, and C during training. The count of non-zero singular values is indicative of the rank.
for queries, keys, values, and output transformations, respectively, and d2 is the dimensionality the key vectors.
The attention module's 's ability to capture low-rank representations reflected the depth-2 -2 matrix factorization model. As illustrated in Fig, BII the attention models consistently learns representations that evolve from lower higher ranks.
10
32
10
12-4
33
10a
00007
Singular values of WQ(b) Singular values of Wx(C) Singular values of Wy(d) Singular values of Wo
Figure B11: The attention modules in Transformer learn adaptively from low rank to high rank with small initialization. (a-d) The evolution of singular values for the matrices Wo. WK. Wv. and Wo throughout the training process. The number of significantly non-zero singular values suggests the effective rank of each matrix.
36

---

NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification: Please the Abstract.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No or
NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings.
is fine to include aspirational goals as motivation long is clear that these goals are not attained by the paper.
2. Limitations
Question:
Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification: Please see Section 7 for the limitation discussion and see Section and Section 6.2 for the discussion of assumptions.
Guidelines:
The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors are encouraged to create a separate "Limitations" section their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution is low images are taken in low lighting. Or speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness.
While the authors might fear that complete honesty about limitations might used by reviewers as grounds for rejection, worse outcome might be that reviewers discover limitations that aren't 't acknowledged judgment and recognize that individual actions in favor of transparency play impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
the paper. The authors should use their best
3. Theory Assumptions and Proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
Answer: [Yes]
37

---

Justification:
Please see Section 6.1 ann and Section 6.2 for assumptions, Please see Appendix A for proofs of all theoretical results.
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide short proof sketch provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental Result Reproducibility
Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that it affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?
Question:
Answer: [Yes]
Justification: Please refer to Appendix for comprehensive details the experiments. The design ensures that all experiments can be readily replicated.
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, a No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is novel architecture, describing the architecture fully might suffice, or if the contribution a specific model and empirical evaluation, it may be necessary either make it possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g., in the case of large language model), releasing of a model checkpoint, or other means that appropriate to the research performed.
While NeurIPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) the contribution is primarily a new algorithm, the paper should make clear how reproduce that algorithm. (b) If the contribution is primarily new model architecture, the paper should describe the architecture clearly and fully.
(c) the contribution is a new model (e.g., a large language model), then there should either be way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome deseribe the particular way they provide for reproducibility. In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but it should possible for other researchers to have some path to reproducing or verifying the results.
5. Open access to data and code
38

---

Does the paper provide open access to the data and code, with sufficient instruc- tions faithfully reproduce the main experimental results, as described in supplemental material?
Question:
Answer: [Yes)
Justification: We have included example code the supplementary material to facilitate replication of the matrix completion experiment. Given the simplicity of both the model and the data in the matrix factorization context, this code enables straightforward verification of our experimental results.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurlPS code and data submission guidelines (https://nips.. cc/
public/ gwwi///gd//////eeee SrrodeSuumissinme
for more details.
While we encourage the release of code and data, we understand that this might not be possible, so "No" is an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurlPS code and data submission guidelines (https c/nips.ccpppullcc.. httpppplllppgggpp des desd Cet///oodSSubb ss onPol iiyyy for more details. //nips.co cc/
The authors should provide instructions on data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible in supplemental material (appended to the paper) recommended, but including URLs to data and code permitted.
6. Experimental Setting/Details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: Please refer to Appendix for comprehensive details of the experiments. The design ensures that all experiments can be readily replicated.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core of the paper to a level of detail that is necessary appreciate the results and make sense of them.
The full details can provided either with the code, appendix, or as supplemental material.
7. Experiment Statistical Significance
Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Answer: [Yes]
Justification: For information on error bars, please consult Section and for a detailed caleulation, refer to Appendix B
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by emoor bars, confi- dence intervals, or statistical significance tests, least for the experiments that support the main claims of the paper.
39

---

The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.)
The assumptions made should be given (e.g.. Normally distributed errors).
It should be clear whether the error bar is the standard deviation or the standard error of the mean.
is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CL, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that out of range (e.g. negative error rates).
If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text.
8. Experiments Compute Resources
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time execution) needed to reproduce the experiments?
Question:
Answer: [Yes]
Justification: Please see Appendix B for the setup of experiments.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs as well estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary failed experiments that didn't 't make into the paper).
9. Code Of Ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurlPS Code of Ethics http::////.eaaresse.ccoc/o/lle/a
Answer:
[Yes]
Justification:
We have reviewed the NeurIPS Code of Ethies.
Guidelines:
The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.
If the authors answer No. they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., there a special consid- eration due to laws regulations in their jurisdiction).
10. Broader Impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer:[Yor Yes]
Justification: Please see Section 7
Guidelines:
The answer NA means that there no societal impact of the work performed.
If the authors answer NA or No, they should explain why their work has no societal impact why the paper does not address societal impact.
40

---

Examples of negative societal impacts include potential malicious unintended uses (e.g., disinformation, generating fake profiles, surveillance), faireess considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional unintentional) misuse of the technology.
If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms monitor how system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Does the paper deseribe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification: The paper poses no such risks.
Guidelines:
The answer NA means that the paper poses no such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Question:
Answer: [NA]
Justification:
The paper does not use existing assets.
Guidelines:
The answer NA means that the paper does use existing assets. The authors should cite the original paper that produced the code package dataset.
The authors should state which version of the asset used and, possible, include URL.
The name of the license (e.g., CC-BY 4.0) should be included for each asset.
For scraped data from a particular source (e.g., website). the copyright and terms of service of that source should be provided.
assets are released, the license, copyright information, and terms of use package should be provided. For popular datasets, baaassttssppprrssesee papersswss thoode. com datasstt has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
41

---

For existing datasets that are re-packaged, both the original license and the license of the derived asset (if has changed) should be provided.
If this information is not available online, the authors are encouraged to reach out to the asset's s creators.
13. New Assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
Justification:
The paper does not release new assets.
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details the dataset/code/moodl as part their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset used.
At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsoureing and Research with Human Subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well details about compensation (if any)?
Answer: [NA]
Justification:
The paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Including this information the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included the main paper.
According the NeurIPS Code of Ethics, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional Review Board (IRB) Approvals or Equivalent for Research with Human Subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or equivalent approval/review based the requirements of your country or institution) were obtained?
Answer: [NA]
Justification: The paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Depending the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review.
42


Reviews:
[{"review_id": "3J4xR86835", "reviewer": "Reviewer_CGEg", "summary": "This paper systematically investigates the implicit regularization of matrix factorization for solving matrix completion problems. The authors find by empirical results that the connectivity of observed data plays a crucial role in the implicit bias, with a transition from low nuclear norm to low rank as data becomes more connected with increased observations. They identify a hierarchy of intrinsic invariant manifolds in the loss landscape that guide the training trajectory to evolve from low-rank to higher-rank solutions. The optimization trajectory follows a Hierarchical Invariant Manifold Traversal (HIMT) process, generalizing the characterization of Li et al. (2020), whose proposed Greedy Low-Rank Learning (GLRL) algorithm equivalence only corresponding to the connected case. Regarding the minimum nuclear norm regularization, the authors establish conditions that provide guarantees closely aligned with the empirical findings, and they present a dynamic characterization condition that assures the attainment of the minimum rank solution.", "questions": "No", "limitations": "N/A", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "- The paper systematically studies the training dynamics and implicit regularization of matrix factorization for matrix completion, considering both connected and disconnected cases, which provides a more comprehensive understanding of the model. \n- The paper identifies the hierarchical invariant manifolds in the loss landscape and characterizes the training trajectory, providing a theoretical basis for understanding the model's behavior. \n- The authors conduct extensive experiments to support their findings, demonstrating the influence of data connectivity on the implicit regularization and the training dynamics of the matrix factorization model.", "weaknesses": "- In some cases, an extremely small initialization is required, which may potentially impact the training speed. \n- The proofs of the theoretical results are quite complex and may be difficult to follow for readers who are not familiar with the mathematical background. A more intuitive explanation or proof sketch could help readers better understand the key ideas.\n- It is suggested to provide justifications to demonstrate the reasonableness of Assumption 1.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "* **Point 1: In some cases, an extremely small initialization is required, which may potentially impact the training speed.**\n\n* Reply: Our empirical analysis in Appendix B.4 (Fig. B6) demonstrates the relationship between observation magnitude differences and the required initialization scale. To elucidate the trade-off between initialization scale and training speed, we have expanded Appendix B.4 with the following details:\n\n  1. Theoretical optimization insights: \"For the matrix factorization model $f_\\theta = AB$, the Hessian matrix at 0 has strictly negative eigenvalues, making the origin a strict saddle point. Under small random initialization, gradient descent escapes this saddle at an exponential rate. Our Theorem 1 ensures that only strict saddle points and global minima exist as critical points in the loss landscape. Subsequent saddle points on invariant manifolds are also strict, facilitating exponential escape speeds and maintaining an efficient optimization process.\"\n\n  2. Practical low-rank attainment: \"Achieving the lowest rank solution requires parameters to escape the saddle point along the unique top eigen-direction and keep evolving near the invariant manifold $\\Omega_k$. Fig. B6 in Appendix B.4 illustrates the relationship between required initialization scale and observation magnitude differences. For the lowest possible rank solution with large numerical magnitude differences in observations, an extremely small initialization is necessary. However, for approximately low rank solutions (some singular values ​​are relatively small), which may be sufficient in practice, a relatively small initialization suffices without significantly impacting training speed.\"\n\n* **Point 2: The proofs of the theoretical results are quite complex and may be difficult to follow for readers who are not familiar with the mathematical background. A more intuitive explanation or proof sketch could help readers better understand the key ideas.**\n\n* Reply: We acknowledge that Theorem 2 is a comprehensive result with a relatively complex proof. To enhance understanding, we have added a more detailed proof sketch for Theorem 2:\n\n  \"Proof sketch: We analyze the local dynamics in the vicinity of the critical point $\\theta_c$. The nonlinear dynamics can be approximated linearly near $\\theta_c$: $\\frac{d\\theta}{dt} \\approx H(\\theta_0 - \\theta_c)$, where $H = -\\nabla^2R_S(\\theta_c)$ is the negative Hessian matrix. For exact linear approximation, the solution is: $\\theta(t) = e^{tH}(\\theta_0 - \\theta_c) + \\theta_c$. Let $\\lambda_1 > \\lambda_2 > ... > \\lambda_s$ be the eigenvalues of $H$, with corresponding eigenvectors $q_{ij}$. We can express $\\theta(t)$ as: $\\theta(t) = \\sum_{i=1}^s \\sum_{j=1}^{l_i} e^{\\lambda_i t}\\langle\\theta_0 - \\theta_c, q_{ij}\\rangle q_{ij} + \\theta_c$. For sufficiently large $t_0$, the dynamics follows a dominant eigenvalue trajectory: $\\theta(t_0) = \\sum_{j=1}^{l_1} e^{\\lambda_1 t_0}\\langle\\theta_0 - \\theta_c, q_{1j}\\rangle q_{1j} + O(e^{\\lambda_2t_0}).$ Through detailed analysis of the eigenvalues and eigenvectors of the Hessian matrix (Lemmas A.2-A.4), we demonstrate that if the largest singular value of residual matrix $\\delta M$ at $\\theta_c$ is unique and $\\theta_c$ is a second-order stationary point within $\\Omega$, the first principal component $\\sum_{j=1}^{l_1} e^{\\lambda_1 t_0}\\langle\\theta_0 - \\theta_c, q_{1j}\\rangle q_{1j}$ will correspond to an $\\Omega_1$ invariant manifold. Consequently, escaping $\\theta_c$ increases the rank by 1, entering $\\Omega_{k+1}$. We defer the details to Appendix A.\"\n\n* **Point 3: It is suggested to provide justifications to demonstrate the reasonableness of Assumption 1.**\n\n* Reply: 1. We appreciate this suggestion and have revised Assumption 1 and its accompanying remark to provide better justification:\n\n  \"Assumption 1 (Unique Top Singular Value): Let $\\delta M=\\left(A_c B_c-M\\right)_{S_x}$ be the residual matrix at the critical point $\\theta_c=\\left(A_c, B_c\\right)$. Assume that the largest singular value of $\\delta M$ is unique.\"\n\n  \"Remark: Assumption 1 ensures that upon departing from a critical point $\\theta_c$, the trajectory is constrained to escape along a single dominant eigendirection corresponding to the largest singular value. This assumption holds for randomly generated matrix with probability 1, making it a reasonable condition in most practical scenarios.\"\n\n  2. To further illustrate the implications when this assumption does not hold, we have added a simple example in the appendix:\n\n  \"Consider the $2\\times 2$ matrix completion problem: $M = \\begin{bmatrix}2 & \\star\\\\\\\\ \\star & 2\\end{bmatrix}$. In this case, the two numbers on the diagonal are identical, which causes the maximum singular value of the residual matrix at the origin to be non-unique. Consequently, the training process will jump directly from the rank 0 to the rank 2 invariant manifold, thereby missing the lowest rank solution of rank 1. This behavior is demonstrated in the attached PDF (Figure 1), which shows experimental results for this scenario.\"\n\nWe sincerely appreciate the reviewer's thorough examination of our work, which has significantly contributed to improving the clarity and precision of our presentation.", "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": "", "extra_scores": {}}, {"review_id": "9K4pGudJHn", "reviewer": "Reviewer_FMau", "summary": "This paper attempts to present a unified understanding of when and how matrix factorization models have different implicit regularization effects. Their key finding is that connectivity of the observed data plays an important role: (i) in the connected cases, the model learns the lowest-ranked solution and (ii) in the disconnected case, it seems to sort of depend, but generally does not find the minimum nuclear norm solution. They identify invariant manifolds in the loss landscape that guides the training trajectory from these lower ranked to higher ranked solutions. They support a few of their findings with theoretical results.", "questions": "- I am trying to understand Figure 4(e): do the different dashed lines correspond to different variances of the initialization? In that case, I am assuming that the figure is trying to show that despite the initialization scales, the model always learns a sub-optimal solution in this disconnected case.", "limitations": "I have addressed a few limitations in the weaknesses section.", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "- Overall, I think this paper presents some neat findings while addressing an interesting question.\n- The paper is well-written and supported with examples to help elucidate the definitions / settings. \n- The paper fits well into the literature of implicit regularization, and helps unify some of the observations in the literature as well, such as greedy low-rank learning.", "weaknesses": "- I am unsure of how restrictive the assumptions used for the theoretical results are which are used to prove Theorem 2. The authors do make a remark about the assumptions, but it is still not immediately clear to me that these are fair assumptions to make.\n- As the authors mention, it would be interesting (and important) future work to show how the connectivity of the observed data holds attainment across the invariant manifolds.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "* **Point 1: I am unsure of how restrictive the assumptions used for the theoretical results are which are used to prove Theorem 2.**\n\n* Reply: We appreciate the reviewer's concern regarding the restrictiveness of our assumptions. We have taken steps to clarify and justify these assumptions:\n\n  1. Regarding Assumption 1, as other reviewers have noted, it can be more concisely stated as the uniqueness of the largest singular value of the residual matrix. We have revised it as follows:\n\n    \"Assumption 1 (Unique Top Singular Value): Let $\\delta M=\\left(A_c B_c-M\\right)_{S_x}$ be the residual matrix at the critical point $\\theta_c=\\left(A_c, B_c\\right)$. Assume that the largest singular value of $\\delta M$ is unique.\"\n\n  2. We have also updated the accompanying remark:\n\n    \"Remark: Assumption 1 ensures that upon departing from a critical point $\\theta_c$, the trajectory is constrained to escape along a single dominant eigendirection corresponding to the largest singular value. This assumption holds for randomly generated matrix with probability 1, making it a reasonable condition in most practical scenarios.\"\n\n  3. To illustrate the implications when Assumption 1 does not hold, we've added an example in the appendix:\n\n    \"Consider the $2\\times 2$ matrix completion problem: $M = \\begin{bmatrix}2 & \\star\\\\\\\\ \\star & 2\\end{bmatrix}$. In this case, the two numbers on the diagonal are identical, which causes the maximum singular value of the residual matrix at the origin to be non-unique. Consequently, the training process will jump directly from the rank 0 to the rank 2 invariant manifold, thereby missing the lowest rank solution of rank 1. This behavior is demonstrated in the attached PDF (Figure 1), which shows experimental results for this scenario.\"\n\n  4. To enhance understanding of Theorem 2, we've included a detailed proof sketch:\n\n   \"Proof sketch: We analyze the local dynamics in the vicinity of the critical point $\\theta_c$. The nonlinear dynamics can be approximated linearly near $\\theta_c$: $\\frac{d\\theta}{dt} \\approx H(\\theta_0 - \\theta_c)$, where $H = -\\nabla^2R_S(\\theta_c)$ is the negative Hessian matrix. For exact linear approximation, the solution is: $\\theta(t) = e^{tH}(\\theta_0 - \\theta_c) + \\theta_c$. Let $\\lambda_1 > \\lambda_2 > ... > \\lambda_s$ be the eigenvalues of $H$, with corresponding eigenvectors $q_{ij}$. We can express $\\theta(t)$ as: $\\theta(t) = \\sum_{i=1}^s \\sum_{j=1}^{l_i} e^{\\lambda_i t}\\langle\\theta_0 - \\theta_c, q_{ij}\\rangle q_{ij} + \\theta_c$. For sufficiently large $t_0$, the dynamics follows a dominant eigenvalue trajectory: $\\theta(t_0) = \\sum_{j=1}^{l_1} e^{\\lambda_1 t_0}\\langle\\theta_0 - \\theta_c, q_{1j}\\rangle q_{1j} + O(e^{\\lambda_2t_0}).$ Through detailed analysis of the eigenvalues and eigenvectors of the Hessian matrix (Lemmas A.2-A.4), we demonstrate that if the largest singular value of residual matrix $\\delta M$ at $\\theta_c$ is unique and $\\theta_c$ is a second-order stationary point within $\\Omega$, the first principal component $\\sum_{j=1}^{l_1} e^{\\lambda_1 t_0}\\langle\\theta_0 - \\theta_c, q_{1j}\\rangle q_{1j}$ will correspond to an $\\Omega_1$ invariant manifold. Consequently, escaping $\\theta_c$ increases the rank by 1, entering $\\Omega_{k+1}$.\"\n\n  5. Regarding Assumption 2, we've revised the remark to better explain its role:\n\n  \"To ensure the escape direction falls within the $\\Omega_{k+1}$ invariant manifold, the Hessian's top eigenvectors must satisfy $\\text{rank}(A) = \\text{rank}(B^\\top) = \\text{rank}(W_{\\text{aug}})$. The condition that $\\theta_c$ is a second-order stationary point within $\\Omega$ in Assumption 2 guarantees this Hessian structure. Our Assumption 2 is more general than conditions proposed by Li et al. (2020), as it remains valid across both connected and disconnected configurations. Empirical findings (Fig. 3 and Fig. 4) indicate that this assumption consistently holds in practical scenarios.\"\n\n* **Point 2: It would be interesting (and important) future work to show how the connectivity of the observed data holds attainment across the invariant manifolds.**\n\n* Reply: We concur that investigating the relationship between data connectivity and attainment across invariant manifolds is a crucial direction for future research. Our planned approach includes a comprehensive analysis of the loss landscape structure at each invariant manifold level, with a particular focus on the characteristics of critical points. This investigation will provide valuable insights into how data connectivity influences the optimization trajectory and the ultimate attainment of low-rank solutions.\n\n* **Point 3: Do the different dashed lines correspond to different variances of the initialization?**\n\n* Reply: We appreciate the reviewer's attention to detail regarding Figure 4(e). To clarify:\n\n  1. The dashed lines in Figure 4(e) represent different disconnected sampling patterns, each covering 4 elements and observing the remaining 5 elements in a fixed rank-1 3x3 matrix. We comprehensively explored all nine possible disconnected sampling patterns with 5 observations. To enhance clarity, we have updated the caption of Figure 4(e):\n\n  \"Figure 4(e): Training loss for all nine disconnected sampling patterns in a 3x3 matrix, each covering 4 elements and observing the remaining 5 in a fixed rank-1 matrix.\"\n\n  2. Furthermore, we have added a brief discussion in Section 5.2, Line 199:\n\n  \"In Figure 4(e), we fixed a rank-1 matrix and explored all nine disconnected sampling patterns with 5 observations. For each pattern, we conducted experiments with small initializations. The loss curves consistently indicate that in disconnected cases, the model learns a sub-optimal solution in the rank-1 manifold, ultimately resulting in a rank-2 solution. This demonstrates that regardless of the specific disconnected sampling pattern, the model fails to achieve the optimal low-rank solution.\"\n\nWe sincerely appreciate the reviewer's thorough examination of our work.", "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": "", "extra_scores": {}}, {"review_id": "HqD4VIFcdP", "reviewer": "Reviewer_AjkV", "summary": "This paper studies the training dynamics of matrix factorisation for matrix completion, optimising the vanilla mean-squared-error loss function via gradient descent with small initialisation. The authors characterise the observation pattern of the underlying matrix via the connectivity of its associated bipartite graph, and show both empirically and theoretically that this plays an important role in the implicit regularisation of the learned solution.\n\nThe paper provides empirical evidence that in both the connected and disconnected cases, the row/column spaces of the factor matrices remain aligned at all times, and that the loss curves exhibit steady step-wise decline as the solutions traverse solutions of increasing rank. In the connected case, it is observed that the learned solutions of each rank are near optimal, however in the disconnected case, suboptimal solutions are learned.\n\nThese empirical observations and explained by theoretical results which examine the training dynamics of the continuous dynamical system following the gradient flow of the loss function with infinitesimal initialisation. The authors begin by defining a Hierarchical Inrinsic Invariant Manifold and Proposition 1 describes how the gradient flow behaves on these manifolds and that they form a hierarchy. Proposition 2 shows that in the special case that when the observation graph is disconnected, the gradient flow dynamics are restricted to sub-manifolds. Theorem 1 (which extends an analgous result for symmetric matrix factorization models in Li et al. (2020)) shows that the critical points of the gradient flow are either strict saddle points or global minima and Theorem 2 shows that the gradient flow traverses the hierarchy of Intrinsic Invariant Manifolds (or sub-manifolds in the disconnected case), a phenomenon coined \"Hierarchical Invariant Manifold Traversal\". The flagship results of the paper are Theorems 3 and 4 which show that in the connected case, the gradient flow converges to the minimum rank solution, and in the a special disconnected case (disconnected with complete bipartite components), the gradient flow converges to the minimum nuclear norm solution.", "questions": "- What the three bars for each matrix in Fig 2b? What is the thick vertical line?\n- In Fig 3d, is there reason to expect the norms of the gradients to be on the same scale as the difference between the matrix at the saddle points and the optimal approximation for the corresponding rank? They seem to be close to the gradient norms at the saddle points. Is this a coincidence or is there a reason to expect this?\n- In Fig 3c and 4e, loss trajectories for different sampling patterns are plotted. If I understand correctly, in Fig 3a for example, each line corresponds to a different entry of the matrix being unobserved. In this case is $\\star = 1.2$, so the underlying rank is three, or does it take a different value? Is it the case that the matrix has rank three if and only if the unobserved entry takes a specific value?\n- In Fig 4f, are you counting position indexes from zero? In Section 3 (lines 97 and 99), you count position indexes from one. \n- The red dot in Fig 4e would perhaps be best a little larger and plotted on top of the other lines. At present it is quite hard to see.\n- Assumption 1 is equivalent to the top singular value of $-\\delta \\mathbf{M}$​ being unique.  This seems like a more intuitive way to write this assumption. Is there a reason you chose to write the assumption in terms of the eigenvalues of the symmetric dilation matrix?\n- In 231, it is claimed that \"in the connected case, at each level [the] model reaches an optimal solution\". Is this implied by one of the theorems, by previous work, or is this an empirical observation?", "limitations": "The authors are upfront about the limitations of their work in the conclusions, which present exciting directions for future research.", "rating": 9, "confidence": 5, "soundness": 4, "presentation": 4, "contribution": 4, "strengths": "This paper makes a significant and original contribution to the understanding of implicit bias of gradient descent for matrix factorisation. The key strength of this paper is showing, both empirically and theoretically, how the connectivity of the observation graph impacts the training dynamics, and how this affects the final solution. Distinguishing the connected and disconnected cases is novel and enlightening, and the detail with which the training dynamics are theoretically characterised is eye-opening.\n\nThe paper is beautifully written: the problem is set out very clearly, the relevant existing literature is summarises comprehensively and concisely. The reader is guided through simple toy experiments which  precisely demonstrate the phenomenon which is described. The theoretical results then provide a satisfying explanation for the empirical observations. I really felt like I had gained a deep understanding of the training dynamics for this problem by the time I had finished reading this paper. While the theoretical results are technical, the intuition I gained in the first half of the paper really helped me to understand the details. In addition, the results and the proofs are written very clearly and accurately.", "weaknesses": "I find this work hard to fault. While the theoretical results could be criticised for being fairly complex, having spend some time understanding them, the payoff makes it absolutely worthwhile. I don't see that they could be made simpler.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "* **Point 1: I find this work hard to fault. While the theoretical results could be criticised for being fairly complex, having spend some time understanding them, the payoff makes it absolutely worthwhile.**\n\n* Reply: We are deeply grateful for the reviewer's careful examination and understanding of our theoretical results. Your recognition of the value and payoff of our complex theoretical framework is greatly appreciated. \n\n* **Point 2: What the three bars for each matrix in Fig 2b? What is the thick vertical line?**\n\n* Reply: We appreciate the reviewer's attention to detail regarding Figure 2b. To clarify:\n\n  - The three bars for each matrix represent the singular values of the learned matrix.\n  - The thick vertical line partitions significantly nonzero singular values, which serves as the empirical rank.\n\n  To enhance clarity, we have updated the caption of Figure 2b as follows:\n\n  \"Figure 2(b): Singular values of the learned matrices for $M_1, M_2, M_3$. Each set of three bars represents the singular values of a matrix. The thick vertical lines partition significantly nonzero singular values, which serves as the empirical rank.\"\n\n* **Point 3: In Fig 3d, is there reason to expect the norms of the gradients to be on the same scale as the difference between the matrix at the saddle points and the optimal approximation for the corresponding rank?**\n\n* Reply: We appreciate the reviewer's insightful observation about the similarity between gradient norms and matrix differences at saddle points. This relationship indeed stems from the structure of optimization and can be more precisely characterized.\n\n  As $R_S(\\theta)$ approaches $R_S(\\theta^*)$, we observe that $\\nabla R_S(\\theta)$ also approaches zero. This relationship can be formalized under certain conditions:\n\n  1. For $L$-smooth functions, we have the inequality: $\\\\|\\nabla R_S(\\theta)\\\\|^2 \\leq 2L(R_S(\\theta) - R_S(\\theta^*))$, where $L$ is the Lipschitz constant of the gradient.\n\n  2. For $\\mu$-strongly convex functions, we also have the inequality (the Polyak-Łojasiewicz condition): $\\\\|\\nabla R_S(\\theta)\\\\|^2 \\geq 2\\mu(R_S(\\theta) - R_S(\\theta^*))$.\n\n  In our matrix factorization setting, while the problem isn't globally strongly convex, it may exhibit local strong convexity restricted to invariant manifold $\\Omega_k$ near optimal solutions. This local behavior possibly contributes to the observed similarity. \n\n* **Point 4: Clarification on Figures 3c and 4e. For Figure 3c, is it the case that the matrix has rank three if and only if the unobserved entry takes a specific value?**\n\n* Reply: We appreciate the reviewer's attention to detail regarding Figures 3c and 4e. To clarify:\n\n  1. The dashed lines in both figures correspond to different sampling patterns. For Figure 3(c), the matrix indeed has rank three if and only if the unobserved entry takes a specific value. This is because, given 15 observations, the rank-3 matrix is uniquely determined.\n\n  2. To improve clarity, we have updated the captions of Figures 3(c) and 4(e):\n\n  \"Figure 3(c): Training loss for 16 connected sampling patterns in a $4\\times 4$ matrix, each covering 1 element and observing the remaining 15 in a fixed rank-3 matrix.\"\n\n  \"Figure 4(e): Training loss for 9 disconnected sampling patterns in a $3\\times 3$ matrix, each covering 4 elements and observing the remaining 5 in a fixed rank-1 matrix.\"\n\n* **Point 5: In Fig 4f, are you counting position indexes from zero? In Section 3 (lines 97 and 99), you count position indexes from one.**\n\n* Reply: We sincerely thank the reviewer for pointing out this inconsistency. Indeed, in Figure 4f, we counted position indexes from zero, while in Section 3 we counted from one. To maintain consistency throughout the paper, we have updated Figure 4f to count from one.\n\n  We've revised the caption of Figure 4f as follows:\n\n  \"Figure 4f: Learned values at symmetric positions $(1, 2)$ and $(2, 1)$ under varying initialization scales (zero mean, varying variance).\"\n\n* **Point 6: The red dot in Fig 4e would perhaps be best a little larger and plotted on top of the other lines.**\n\n* Reply: We appreciate this suggestion for improving the figure's clarity. We have implemented this change by increasing the size of the red dot and ensuring it's plotted on top of the other lines in Figure 4e. This modification enhances the visibility of this key data point in the figure.\n\n* **Point 7: Assumption 1 is equivalent to the top singular value of $-\\delta M$ being unique. This seems like a more intuitive way to write this assumption.**\n\n* Reply: We thank the reviewer for this insightful suggestion. We agree that reformulating Assumption 1 in terms of the top singular value of $-\\delta M$ is indeed more intuitive. We have revised Assumption 1 accordingly:\n\n  \"Assumption 1 (Unique Top Singular Value): Let $\\delta M=\\left(A_c B_c-M\\right)_{S_x}$ be the residual matrix at the critical point $\\theta_c=\\left(A_c, B_c\\right)$. Assume that the largest singular value of $\\delta M$ is unique.\"\n\n* **Point 8: In 231, it is claimed that \"in the connected case, at each level [the] model reaches an optimal solution\". Is this implied by one of the theorems, by previous work, or is this an empirical observation?**\n\n* Reply:  The statement \"in the connected case, at each level the model reaches an optimal solution\" is based on our empirical observations. Our experiments consistently show this behavior in connected cases, as illustrated in Figure 3.\n\n  To improve precision, we have revised the statement in line 231 to:\n\n  \"In the connected case, at each level we observe that the model reaches an optimal solution (Figure 3).\"\n   \n   This revision clarifies that our claim is based on empirical evidence rather than a theoretical guarantee.\n\nWe sincerely appreciate the reviewer's thorough examination of our work, which has significantly contributed to improving the clarity and precision of our presentation.", "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": "", "extra_scores": {}}, {"review_id": "Q1Uac1voH0", "reviewer": "Reviewer_PxmR", "summary": "This study empirical reveals that the connectivity of observed data significantly influences the implicit bias and identifies a hierarchy of intrinsic invariant manifolds in the loss landscape, providing a preliminary framework for understanding the mechanisms behind implicit regularisation in matrix factorisation models.", "questions": "1, Have you verified your findings on high dimensional matrix with experiments?\\\n2, What's the intuition of the augmented matrix? See $\\textit{Line 110}$?\\\n3, What exactly \"small initialisation\" you defined in the experiments? See $\\textit{Line 164}$\\\nMinor suggestion(s):\\\n$\\textbf{A.}$ $Eq.3$ after $\\textit{Line 117}$: You have mentioned $\\theta(0) = \\theta_0$, which indicates the initialisation (if I understand correctly). However, the initialisation method was mentioned far from here at $\\textit{Line 36-37}$. It could be easier to follow if you can optimise this.\\\n$\\textbf{B.}$ I expect a detailed example of connectivity and disconnectivity in the appendix (see $\\textit{Line 561 - 564 and Fig. A1}$), which will make it easier for a wider audience to understand.", "limitations": "Yes.", "rating": 6, "confidence": 2, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "1, The paper is clearly written.\\\n2, Great presented figures illustrate experimental results and findings, making them easier for readers to understand.\\\n3, Differences from prior similar works have been properly addressed.", "weaknesses": "Please refer to the \"Questions\".", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "* **Point 1: Have you verified your findings on high dimensional matrix with experiments?**\n\n* Reply: We sincerely appreciate the reviewer's important question regarding the scalability of our findings. While our main results focus on smaller matrices for clarity of presentation, we have indeed conducted additional experiments on higher-dimensional matrices to verify the broader applicability of our findings.\n\n  We've added a new subsection in the Appendix titled \"B.7 High-Dimensional Experiments\" with the following content:\n\n  \"To validate the scalability of our findings, we extended our experiments to higher-dimensional matrices. We conducted tests on  randomly generated $20\\times 20$ matrices, employing both connected (Fig. 2) and disconnected (Fig. 3) sampling patterns, while monitoring rank evolution during training. Our results consistently corroborated the main findings (please refer to high-dimensional experiments in the attached PDF):\n\n  - (i) Connected observations converged to optimal low-rank solutions.\n  - (ii) Disconnected observations yielded higher-rank solutions.\n  - (iii) The Hierarchical Invariant Manifold Traversal (HIMT) process was observed in both connected and disconnected scenarios.\"\n\n* **Point 2: What's the intuition of the augmented matrix? See Line 110?**\n\n* Reply: We appreciate the request for clarification on the augmented matrix. The augmented matrix $W_{\\text{aug}}$ is indeed a key construct in our analysis, and we recognize the need for a more intuitive explanation.\n\n  We've added the following explanation after Line 110:\n\n  \"The augmented matrix $W_\\text{aug}$ plays a crucial role in our subsequent analysis, particularly in characterizing the intrinsic invariant manifolds $\\Omega_k$ of the optimization process. Specifically, it allows us to establish the relationship $\\text{rank}(A) = \\text{rank}(B^\\top) = \\text{rank}(W_{\\text{aug}})$, which is important to understanding the invariance property under gradient flow.\"\n\n* **Point 3: What exactly \"small initialization\" you defined in the experiments? See Line 164.**\n\n* Reply: We thank the reviewer for pointing out the need for a clearer definition of \"small initialization\". To address this, we've revised the description in Line 164 as follows:\n\n  \"For a small random initialization (Gaussian distribution with mean 0 and variance $10^{-16}$), the loss curves exhibit a steady, stepwise decline\"\n\n* **Point 4: You have mentioned $\\theta(0) = \\theta_0$, which indicates the initialization (if I understand correctly). However, the initialization method was mentioned far from here at Line 36-37. It could be easier to follow if you can optimise this.**\n\n* Reply: We appreciate the reviewer's suggestion to improve the flow and clarity of our presentation regarding initialization. We have revised the relevant section as follows:\n\n  \"The training dynamics follow the gradient flow of $R_S(\\boldsymbol{\\theta})$:\n  $$ \\frac{\\mathrm{d} \\boldsymbol{\\theta}}{\\mathrm{d} t}=-\\nabla_{\\boldsymbol{\\theta}} R_S(\\boldsymbol{\\theta}), \\quad \\boldsymbol{\\theta}(0)=\\boldsymbol{\\theta}_0 . $$\n  In all experiments, $\\theta_0 \\sim N(0, \\sigma^2)$ is initialized from a Gaussian distribution with mean 0 and small variance $\\sigma^2$. We use gradient descent with a small learning rate to approximate the gradient flow dynamics (please refer to Appendix B.1 for the detailed experimental setup).\"\n\n* **Point 5: I expect a detailed example of connectivity and disconnectivity in the appendix (see Line 561 - 564 and Fig. A1), which will make it easier for a wider audience to understand.**\n\n* Reply: We appreciate the reviewer's suggestion to provide a more detailed example of connectivity and disconnectivity. We've expanded the examples in the appendix (Lines 561-564 and Fig. A1) with a step-by-step explanation:\n\n  \"Examples of connectivity and disconnectivity. \n   - Consider three matrices to be completed, each obtained by adding one more observation to the previous matrix:\n    $$ \\boldsymbol{M}_1=\\left[\\begin{array}{ccc} 1 & 2 & \\star \\\\\\\\ 3 & \\star & \\star \\\\\\\\ \\star & \\star & 5 \\end{array}\\right], \\boldsymbol{M}_2=\\left[\\begin{array}{ccc} 1 & 2 & \\star \\\\\\\\ 3 & 4 & \\star \\\\\\\\ \\star & \\star & 5 \\end{array}\\right], \\boldsymbol{M}_3=\\left[\\begin{array}{ccc} 1 & 2 & \\star \\\\\\\\ 3 & 4 & \\star \\\\\\\\ 6 & \\star & 5 \\end{array}\\right]. $$\n\n  - The corresponding observation matrices $\\boldsymbol{P}$ are:\n  $$ \\boldsymbol{P}_1=\\left[\\begin{array}{ccc} 1 & 1 & 0 \\\\\\\\ 1 & 0 & 0 \\\\\\\\ 0 & 0 & 1 \\end{array}\\right], \\boldsymbol{P}_2=\\left[\\begin{array}{ccc} 1 & 1 & 0 \\\\\\\\ 1 & 1 & 0 \\\\\\\\ 0 & 0 & 1 \\end{array}\\right], \\boldsymbol{P}_3=\\left[\\begin{array}{ccc} 1 & 1 & 0 \\\\\\\\ 1 & 1 & 0 \\\\\\\\ 1 & 0 & 1 \\end{array}\\right]. $$\n\n  - And by Definition 1, the corresponding adjacency matrices are:\n  $$ \\boldsymbol{A}_1=\\left[\\begin{array}{ccc} \\boldsymbol{0} & \\boldsymbol{P}_1^\\top \\\\\\\\ \\boldsymbol{P}_1 & \\boldsymbol{0} \\end{array}\\right], \\boldsymbol{A}_2=\\left[\\begin{array}{ccc} \\boldsymbol{0} & \\boldsymbol{P}_2^\\top \\\\\\\\ \\boldsymbol{P}_2 & \\boldsymbol{0} \\end{array}\\right], \\boldsymbol{A}_3=\\left[\\begin{array}{ccc} \\boldsymbol{0} & \\boldsymbol{P}_3^\\top \\\\\\\\ \\boldsymbol{P}_3 & \\boldsymbol{0} \\end{array}\\right]. $$\n\n  - Given the adjacency matrix $A_i$, according to Definition 1, we can obtain a bipartite graph $G_{M_i}$, which we refer to as the associated observation graph. Fig. A1 in Appendix A illustrates the associated graphs $G_{M_i}$, from which we can see that $M_1$ is disconnected, with its associated observation graph consisting of two connected components. $M_2$ is also disconnected, but each connected component of its associated observation graph forms a complete bipartite subgraph. In contrast, $M_3$ is connected, and its associated observation graph consists of a single connected component.\"\n\nWe sincerely appreciate the reviewer's thorough examination of our work, which has significantly contributed to improving the clarity and precision of our presentation.", "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": "", "extra_scores": {}}, {"review_id": "I5K9pJsUnC", "reviewer": "Reviewer_1kCy", "summary": "The paper studies implicit regularization in matrix completion, a problem that estimates missing entries from a partially observed matrix. In particular, in this paper, matrix completion is formulated as an optimization problem of the form: $\\min_{A, B} f(A,B) = \\|M - AB\\|_S^2$ where $S$ is the subset of observed entries and the function $\\|\\cdot\\|_S$ is the Frobenius norm of a matrix after masking entries outside $S$ to zero. Using this formulation, the authors study the implicit regularization of the gradient flow associated with $f$ whose initial condition is close to zero. Their main results characterize the influence of connectivity - a property defined based on $S$ - on the dynamics and solutions obtained by the gradient flow. The results are supported both theoretically and numerically.", "questions": "1. Major questions/suggestions:\n- I wonder under which regime of $n$ (the size of the matrix) and $m$ (the number of observed entries) do we have more than/exactly one connected component? I think this question is important to justify the motivation of this work.\n2. Minor questions/suggestions:\n- The experiments in Section 4 need to be provided with more details: What is the algorithm used? If it is gradient descent, what is the learning rate? Or the solution is obtained by solving an ODE.\n- Throughout the introduction, the term \"connectivity\" is repeatedly used but its description or definition is not provided (it is delayed until Section 3). The authors might consider describing/defining this term earlier to improve the readability.\n- Line 612: The second inequality should be $\\|w^T B\\|\\_F^2$ instead of $\\|w^T B\\|\\_{S_x}^2$.", "limitations": "There are no particular limitations or potential negative impacts of this paper.", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 4, "contribution": 2, "strengths": "- The results seem to be correct. Although I did not check all the details, there are certain lemmas/theorems that I checked carefully and I did not find any serious problem.\n- The presentation is excellent. Many definitions are followed by examples and/or remarks to provide intuitions/insightful discussion. The experiments are quite convincing.", "weaknesses": "- The main concern is the interest of the problem: why should we care about data connectivity (in the sense of this paper)? I think the motivation of this study needs to be further justified.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "* **Point 1: Why should we care about data connectivity (in the sense of this paper)? I think the motivation of this study needs to be further justified.**\n\n* Reply: We appreciate the reviewer's question regarding the importance of data connectivity in our study. As the title of our paper suggests, \"Connectivity Shapes Implicit Regularization in Matrix Factorization Models for Matrix Completion\", connectivity plays a crucial role in understanding the implicit regularization behavior of matrix factorization models. This concept is closely related to recent research on nuclear norm minimization and rank minimization mentioned in the Introduction and Related Works.\n\n  - From an experimental perspective, our initial motivation stemmed from extensive experiments exploring the conditions under which low-rank tendencies emerge or fail to emerge. Through a series of simple experiments and analysis of their dynamics, we discovered that the connectivity of the observed data played a pivotal role, as illustrated in the Introduction.\n\n  - Theoretically, the ongoing debate about whether the implicit regularization of matrix factorization models tends towards low nuclear norm or low rank necessitates a clear definition of the conditions that determine these outcomes. Our experiments provided the initial motivation for developing a theoretical framework. In fact, connectivity is closely related to whether the dynamic system can be decoupled into several parts (please refer to Proposition A.4 in Appendix A), which inspired our subsequent theoretical investigations.\n\n  This study thus bridges a crucial gap in our understanding of how data structure influences the behavior of matrix factorization models, providing insights that are valuable for both theoretical analysis and practical applications in matrix completion tasks.\n\n* **Point 2: I wonder under which regime of $n$ (the size of the matrix) and $m$ (the number of observed entries) do we have more than/exactly one connected component? I think this question is important to justify the motivation of this work.**\n\n* Reply: We thank the reviewer for this insightful question about the relationship between matrix size, number of observations, and connectivity. To address this, we've added the following discussion in Appendix B:\n\n  \"Threshold function: For a square $n\\times n$ matrix, the threshold for connectivity is $m = (n-1)^2 + 1$. Below this threshold, disconnected components are likely; above it, a single connected component is guaranteed. The intuition is that for a square $n\\times n$ matrix, if the upper left $(n-1)\\times (n-1)$ sub-block is fully observed, and the last $(n,n)$ position is observed, this represents the case with the maximum number of observations while still being disconnected. If one more entry is observed, the entire observation set becomes connected. However, if there are fewer than $m = (n-1)^2 + 1$ observations, both connected and disconnected situations may occur, depending on the specific sampling pattern.\"\n\n* **Point 3: The experiments in Section 4 need to be provided with more details: What is the algorithm used? If it is gradient descent, what is the learning rate? Or the solution is obtained by solving an ODE.**\n\n* Reply: We appreciate the reviewer's request for more experimental details. We have added the following information in Section 3 to briefly introduce the experimental setup, with further details provided in Appendix B.1:\n\n  \"The training dynamics follow the gradient flow of $R_S(\\boldsymbol{\\theta})$:\n  $$ \\frac{\\mathrm{d} \\boldsymbol{\\theta}}{\\mathrm{d} t}=-\\nabla_{\\boldsymbol{\\theta}} R_S(\\boldsymbol{\\theta}), \\quad \\boldsymbol{\\theta}(0)=\\boldsymbol{\\theta}_0 . $$\n  In all experiments, $\\theta_0 \\sim N(0, \\sigma^2)$ is initialized from a Gaussian distribution with mean 0 and small variance $\\sigma^2$. We use gradient descent with a small learning rate to approximate the gradient flow dynamics (please refer to Appendix B.1 for the detailed experimental setup).\"\n\n* **Point 4: Throughout the introduction, the term \"connectivity\" is repeatedly used but its description or definition is not provided (it is delayed until Section 3). The authors might consider describing/defining this term earlier to improve the readability.**\n\n* Reply: We thank the reviewer for this suggestion to improve readability. We've added a brief description of connectivity in the Introduction:\n\n  \"Data connectivity, in the context of this paper, refers to the way observed entries in the matrix are linked through shared rows or columns. A set of observations is considered connected if there's a path between any two observed entries via other observed entries in the same rows or columns. This concept plays a crucial role in determining the behavior of matrix factorization models, as we will demonstrate throughout this paper.\"\n\n- **Point 5: Line 612: The second inequality should be $\\\\|w^\\top B\\\\|^2_F$ instead of $\\\\|w^\\top B\\\\|^2_{S_x}$.**\n\n* Reply: We sincerely thank the reviewer for catching this typographical error. We have corrected the inequality to $\\\\|w^\\top B\\\\|^2_F$ as suggested. This correction ensures the mathematical accuracy of our presentation.\n\nWe greatly appreciate the reviewer's thorough examination of our work. Your insightful comments and suggestions have significantly contributed to improving the clarity, precision, and overall quality of our presentation.", "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": "", "extra_scores": {}}]

Author rebuttal:
Dear Reviewers,

We sincerely thank all reviewers for their thoughtful and insightful comments. We have carefully addressed every comment, and we believe that the reviewers' collective feedback has significantly improved the manuscript. To address the common concerns raised, we have made the following key improvements:

1. **Clarified Definitions:** We have refined the definition of connectivity and provided a more intuitive explanation of the augmented matrix, enhancing the overall readability of the paper.

2. **Enhanced Experimental Details:** We have added comprehensive information about our experimental setup, including initialization methods and training dynamics, in both the main text and Appendix B.1.

3. **Expanded Theoretical Justifications:** We have provided more detailed proof sketches and intuitive explanations for our theoretical results, particularly for Theorem 2.

4. **Improved Figures and Captions:** We have updated several figures (e.g., Fig. 2b, 3c, 4e, 4f) and their captions to provide clearer visualizations and explanations of our results.

5. **Added High-Dimensional Experiments:** To address concerns about scalability, we have included new experiments on 20x20 matrices, demonstrating the applicability of our findings to higher-dimensional problems.

6. **Strengthened Motivation:** We have elaborated on the importance of data connectivity in shaping implicit regularization, providing a stronger justification for our study.

7. **Corrected Mathematical Notations:** We have fixed minor typographical errors in mathematical expressions to ensure accuracy.

8. **Extended Appendix:** We have added more detailed examples and explanations in the appendix, particularly regarding connectivity and disconnectivity.

Additionally, we have included a one-page PDF in the attachment presenting key experimental results, including the example of coincident top eigenvalues and the high-dimensional experiments.

We believe these revisions comprehensively address the reviewers' concerns while maintaining the core contributions of our work. The revised manuscript now offers a clearer, more rigorous, and more insightful exploration of how connectivity shapes implicit regularization in matrix factorization models for matrix completion.

We sincerely hope that the revised manuscript now satisfies the reviewers' requirements, and we hereby submit it for publication consideration. We are grateful for the opportunity to improve our work based on such valuable feedback.

Best regards,

The Authors

Meta-review: The paper investigates the implicit regularization of matrix factorization for solving matrix completion problems, revealing that the connectivity of observed data has a strong influence on the impact of implicit regularization. It received a positive response from the reviewers, who unanimously found it novel and original, with sound numerical and theoretical analysis and a high-quality presentation.
Decision: Accept (poster)