Provable Efficient Online Matrix Completion via Non-convex Stochastic Gradient Descent

Introduction

Low rank matrix completion refers to the problem of recovering a low rank matrix by observing the values of only a tiny fraction of its entries. This problem arises in several applications such as video denoising , phase retrieval and most famously in movie recommendation engines . In the context of recommendation engines for instance, the matrix we wish to recover would be user-item rating matrix where each row corresponds to a user and each column corresponds to an item. Each entry of the matrix is the rating given by a user to an item. Low rank assumption on the matrix is inspired by the intuition that rating of an item by a user depends on only a few hidden factors, which are much fewer than the number of users or items. The goal is to estimate the ratings of all items by users given only partial ratings of items by users, which would then be helpful in recommending new items to users.

The seminal works of Candès and Recht first identified regularity conditions under which low rank matrix completion can be solved in polynomial time using convex relaxation – low rank matrix completion could be ill-posed and NP-hard in general without such regularity assumptions . Since then, a number of works have studied various algorithms under different settings for matrix completion: weighted and noisy matrix completion, fast convex solvers, fast iterative non-convex solvers, parallel and distributed algorithms and so on.

Most of this work however deals only with the offline setting where all the observed entries are revealed at once and the recovery procedure does computation using all these observations simultaneously. However in several applications , we encounter the online setting where observations are only revealed sequentially and at each step the recovery algorithm is required to maintain an estimate of the low rank matrix based on the observations so far. Consider for instance recommendation engines, where the low rank matrix we are interested in is the user-item rating matrix. While we make an observation only when a user rates an item, at any point of time, we should have an estimate of the user-item rating matrix based on all prior observations so as to be able to continuously recommend items to users. Moreover, this estimate should get better as we observe more ratings.

Algorithms for offline matrix completion can be used to solve the online version by rerunning the algorithm after every additional observation. However, performing so much computation for every observation seems wasteful and is also impractical. For instance, using alternating minimization, which is among the fastest known algorithms for the offline problem, would mean that we take several passes of the entire data for every additional observation. This is simply not feasible in most settings. Another natural approach is to group observations into batches and do an update only once for each batch. This however induces a lag between observations and estimates which is undesirable. To the best of our knowledge, there is no known provable, efficient, online algorithm for matrix completion.

On the other hand, in order to deal with the online matrix completion scenario in practical applications, several heuristics (with no convergence guarantees) have been proposed in literature . Most of these approaches are based on starting with an estimate of the matrix and doing fast updates of this estimate whenever a new observation is presented. One of the update procedures used in this context is that of stochastic gradient descent (SGD) applied to the following non-convex optimization problem

where $\mathbf{M}$ is the unknown matrix of size $d_{1}\times d_{2}$, $k$ is the rank of $\mathbf{M}$ and $\mathbf{U}\mathbf{V}^{\top}$ is a low rank factorization of $\mathbf{M}$ we wish to obtain. The algorithm starts with some $\mathbf{U}_{0}$ and $\mathbf{V}_{0}$, and given a new observation $(\mathbf{M})_{ij}$, SGD updates the $i^{\textrm{th}}$-row and the $j^{\textrm{th}}$-row of the current iterates $\mathbf{U}_{t}$ and $\mathbf{V}_{t}$ respectively by

where $\eta$ is an appropriately chosen stepsize, and $\mathbf{U}^{(i)}$ denote the $i^{\textrm{th}}$ row of matrix $\mathbf{U}$. Note that each update modifies only one row of the factor matrices $\mathbf{U}$ and $\mathbf{V}$, and the computation only involves one row of $\mathbf{U},\mathbf{V}$ and the new observed entry $(\mathbf{M})_{ij}$ and hence are extremely fast. These fast updates make SGD extremely appealing in practice. Moreover, SGD, in the context of matrix completion, is also useful for parallelization and distributed implementation .

In this work we present the first provable efficient algorithm for online matrix completion by showing that SGD (2) with a good initialization converges to a true factorization of $\mathbf{M}$ at a geometric rate. Our main contributions are as follows.

We provide the first provable, efficient, online algorithm for matrix completion. Starting with a good initialization, after each observation, the algorithm makes quick updates each taking time $O(k^{3})$ and requires $O(\mu dk\kappa^{4}(k+\log\frac{\|{\mathbf{M}}\|_{\text{F}}}{\epsilon})\log d)$ observations to reach $\epsilon$ accuracy, where $\mu$ is the incoherence parameter, $d=\max(d_{1},d_{2})$, $k$ is the rank and $\kappa$ is the condition number of $\mathbf{M}$.

Moreover, our result features both sample complexity and total runtime linear in $d$, and is competitive to even the best existing offline results for matrix completion. (either improve over or is incomparable, i.e., better in some parameters and worse in others, to these results). See Table 1 for the comparison.

To obtain our results, we introduce a general framework to show SGD updates tend to stay away from saddle surfaces. In order to do so, we consider distances from saddle surfaces, show that they behave like sub-martingales under SGD updates and use martingale convergence techniques to conclude that the iterates stay away from saddle surfaces. While shows that SGD updates stay away from saddle surfaces, the stepsizes they can handle are quite small (scaling as $1/\text{poly}(d_{1},d_{2})$), leading to suboptimal computational complexity. Our framework makes it possible to establish the same statement for much larger step sizes, giving us near-optimal runtime. We believe these techniques may be applicable in other non-convex settings as well.

2 Related Work

In this section we will mention some more related work.

Offline matrix completion: There has been a lot of work on designing offline algorithms for matrix completion, we provide the detailed comparison with our algorithm in Table 1. The nuclear norm relaxation algorithm has near-optimal sample complexity for this problem but is computationally expensive. Motivated by the empirical success of non-convex heuristics, a long line of works, and so on, has obtained convergence guarantees for alternating minimization, gradient descent, projected gradient descent etc. Even the best of these are suboptimal in sample complexity by $\text{poly}(k,\kappa)$ factors. Our sample complexity is better than that of and is incomparable to those of . To the best of our knowledge, the only provable online algorithm for this problem is that of Sun and Luo . However the stepsizes they suggest are quite small, leading to suboptimal computational complexity by factors of $\text{poly}(d_{1},d_{2})$. The runtime of our algorithm is linear in $d$, which makes $\text{poly}(d)$ improvements over it.

Other models for online matrix completion: Another variant of online matrix completion studied in the literature is where observations are made on a column by column basis e.g., . These models can give improved offline performance in terms of space and could potentially work under relaxed regularity conditions. However, they do not tackle the version where only entries (as opposed to columns) are observed.

Non-convex optimization: Over the last few years, there has also been a significant amount of work in designing other efficient algorithms for solving non-convex problems. Examples include eigenvector computation , sparse coding etc. For general non-convex optimization, an interesting line of recent work is that of , which proves gradient descent with noise can also escape saddle point, but they only provide polynomial rate without explicit dependence. Later show that without noise, the space of points from where gradient descent converges to a saddle point is a measure zero set. However, they do not provide a rate of convergence. Another related piece of work to ours is , proves global convergence along with rates of convergence, for the special case of computing matrix squareroot. During the preparation of this draft, the recent work was announced which proves the global convergence of SGD for matrix completion and can also be applied to the online setting. However, their result only deals with the case where $\mathbf{M}$ is positive semidefinite (PSD) and their rate is still suboptimal by factors of $\text{poly}(d_{1},d_{2})$.

3 Outline

The rest of the paper is organized as follows. In Section 2 we formally describe the problem and all relevant parameters. In Section 3, we present our algorithms, results and some of the key intuition behind our results. In Section 4 we give proof outline for our main results. We conclude in Section 5. All formal proofs are deferred to the Appendix.

Preliminaries

In this section, we introduce our notation, formally define the matrix completion problem and regularity assumptions that make the problem tractable.

2 Problem statement and assumptions

Low rank matrix completion is the task of recovering $\mathbf{M}$ by only observing $\mathcal{P}_{\Omega}(\mathbf{M})$. This task is ill-posed and NP-hard in general . In order to make this tractable, we make by now standard assumptions about the structure of $\mathbf{M}$.

Main Results

In this section, we present our main result. We will first state result for a special case where $\mathbf{M}$ is a symmetric positive semi-definite (PSD) matrix, where the algorithm and analysis are much simpler. We will then discuss the general case.

The algorithm uses an initial set of observations $\Omega_{\textrm{init}}$ to produce a warm start iterate $\mathbf{U}_{0}$, then enters the online stage, where it performs SGD.

The sample complexity of the warm start phase is $O(\mu dk^{2}\kappa^{2}(\mathbf{M})\log d)$. The initialization consists of a top-$k$ SVD on a sparse matrix, whose runtime is $O(\mu dk^{3}\kappa^{2}(\mathbf{M})\log d)$.

For the online phase (SGD), if we choose $\eta=\frac{c}{\mu dk\kappa^{3}(\mathbf{M})\left\|{\mathbf{M}}\right\|\log d}$, the number of observations $T$ required for the error $\|{\mathbf{U}_{T}\mathbf{U}_{T}^{\top}-\mathbf{M}}\|_{\text{F}}$ to be smaller than $\epsilon$ is $O(\mu dk\kappa(\mathbf{M})^{4}\log d\log\frac{\sigma_{\min}(\mathbf{M})}{\epsilon})$.

Since each SGD step modifies two rows of $\mathbf{U}_{t}$, its runtime is $O(k)$ with a total runtime for online phase of $O(kT)$.

Our proof approach is to essentially show that the objective function is well-behaved (i.e., is smooth and strongly convex) in a local neighborhood of the warm start region, and then use standard techniques to show that SGD obtains geometric convergence in this setting. The most challenging and novel part of our analysis comprises of showing that the iterate does not leave this local neighborhood while performing SGD updates. Refer Section 4 for more details on the proof outline.

2 General Case

Algorithm 2 and Algorithm 3 are equivalent in the sense that: given same observations from $\mathbf{M}$ and other inputs, the outputs of Algorithm 2, $\mathbf{U},\mathbf{V}$ and those of Algorithm 3, $\mathbf{U}^{\prime},\mathbf{V}^{\prime}$ satisfy $\mathbf{U}\mathbf{V}^{\top}=\mathbf{U}^{\prime}\mathbf{V}^{\prime}{}^{\top}$.

Since the output of both algorithms is the same, we can analyze Algorithm 2 (which is easier than that of Algorithm 3), while implementing Algorithm 3 in practice. The following theorem is the main result of our paper which presents guarantees on the performance of Algorithm 2.

Just as in the case of PSD matrix completion (Theorem 3.1), Algorithm 2 needs a an initial set of observations $\Omega_{\textrm{init}}$ to provide a warm start $\mathbf{U}_{0}$ and $\mathbf{V}_{0}$ after which it performs SGD.

The sample complexity and runtime of the warm start phase are the same as in symmetric PSD case. The stepsize $\eta$ and the number of observations $T$ to achieve $\epsilon$ error in online phase (SGD) are also the same as in symmetric PSD case.

However, runtime of each update step in online phase is $O(k^{3})$ with total runtime for online phase $O(k^{3}T)$.

The proof of this theorem again follows a similar line of reasoning as that of Theorem 3.1 by first showing that the local neighborhood of warm start iterate has good smoothness and strong convexity properties and then use them to show geometric convergence of SGD. Proof of the fact that iterates do not move away from this local neighborhood however is significantly more challenging due to renormalization steps in the algorithm. Please see Appendix C for the full proof.

Proof Sketch

In this section we will provide the intuition and proof sketch for our main results. For simplicity and highlighting the most essential ideas, we will mostly focus on the symmetric PSD case (Theorem 3.1). For the asymmetric case, though the high-level ideas are still valid, a lot of additional effort is required to address the renormalization step in Algorithm 2. This makes the proof more involved.

First, note that our algorithm for the PSD case consists of an initialization and then stochastic descent steps. The following lemma provides guarantees on the error achieved by the initial iterate $\mathbf{U}_{0}$.

By Lemma 4.1, we know the initialization algorithm already gives $\mathbf{U}_{0}$ in the local region given by Eq.(3). Intuitively, stochastic descent steps should keep doing local search within this local region.

For function $f(\mathbf{U})=\|{\mathbf{M}-\mathbf{U}\mathbf{U}^{\top}}\|_{\text{F}}^{2}$ and any $\mathbf{U}_{1},\mathbf{U}_{2}\in\{\mathbf{U}|\left\|{\mathbf{U}}\right\|\leq\Gamma\}$, we have:

For function $f(\mathbf{U})=\|{\mathbf{M}-\mathbf{U}\mathbf{U}^{\top}}\|_{\text{F}}^{2}$ and any $\mathbf{U}\in\{\mathbf{U}|\sigma_{\min}(\mathbf{X}^{\top}\mathbf{U})\geq\gamma\}$, we have:

Lemma 4.2 tells function $f$ is smooth if spectral norm of $\mathbf{U}$ is not very large. On the other hand, $\sigma_{\min}(\mathbf{X}^{\top}\mathbf{U})$ not too small requires both $\sigma_{\min}(\mathbf{U}^{\top}\mathbf{U})$ and $\sigma_{\min}(\mathbf{X}^{\top}\mathbf{W})$ are not too small, where $\mathbf{W}$ is top-k eigenspace of $\mathbf{U}\mathbf{U}^{\top}$. That is, Lemma 4.3 tells function $f$ has a property similar to strongly convex in standard optimization literature, if $\mathbf{U}$ is rank k in a robust sense ($\sigma_{k}(\mathbf{U})$ is not too small), and the angle between the top k eigenspace of $\mathbf{U}\mathbf{U}^{\top}$ and the top k eigenspace $\mathbf{M}$ is not large.

Within the region $\mathcal{D}=\{\mathbf{U}|\left\|{\mathbf{M}-\mathbf{U}\mathbf{U}^{\top}}\right\|_{F}\leq\frac{1}{10}\sigma_{k}(\mathbf{M})\}$, we have:

Lemma 4.4 tells inside region $\{\mathbf{U}|\left\|{\mathbf{M}-\mathbf{U}\mathbf{U}^{\top}}\right\|_{F}\leq\frac{1}{10}\sigma_{k}(\mathbf{M})\}$, matrix $\mathbf{U}$ always has a good spectral property which gives preconditions for both Lemma 4.2 and 4.3, where $f(\mathbf{U})$ is both smooth and has a property very similar to strongly convex.

With above three lemmas, we already been able to see the intuition behind linear convergence in Theorem 3.1. Denote stochastic gradient

where $SG(\mathbf{U})$ is a random matrix depends on the randomness of sample $(i,j)$ of matrix $\mathbf{M}$. Then, the stochastic update step in Algorithm 1 can be rewritten as:

Let’s suppose ideally, we always have $\mathbf{U}_{0},\ldots,\mathbf{U}_{t}$ inside region $\mathcal{D}$, this directly gives:

One interesting aspect of our main result is that we actually show linear convergence under the presence of noise in gradient. This is true because for the second-order ($\eta^{2}$) term above, we can roughly see from Eq.(4) that $\|{SG(\mathbf{U})}\|_{\text{F}}^{2}\leq h(\mathbf{U})\cdot f(\mathbf{U})$, where $h(\mathbf{U})$ is a factor depends on $\mathbf{U}$ and always bounded. That is, $SG(\mathbf{U})$ enjoys self-bounded property — $\|{SG(\mathbf{U})}\|_{\text{F}}^{2}$ will goes to zero, as objective function $f(\mathbf{U})$ goes to zero. Therefore, by choosing learning rate $\eta$ appropriately small, we can have the first-order term always dominate the second-order term, which establish the linear convergence.

Now, the only remaining issue is to prove that “$\mathbf{U}_{0},\ldots,\mathbf{U}_{t}$ always stay inside local region $\mathcal{D}$”. In reality, we can only prove this statement with high probability due to the stochastic nature of the update. This is also the most challenging part in our proof, which makes our analysis different from standard convex analysis, and uniquely required due to non-convex setting.

Let $f(\mathbf{U})=\left\|{\mathbf{U}\mathbf{U}^{\top}-\mathbf{M}}\right\|_{F}^{2}$ and $g_{i}(\mathbf{U})=\left\|{\mathbf{e}_{i}^{\top}\mathbf{U}}\right\|^{2}$. Suppose initial $\mathbf{U}_{0}$ satisfying:

Then, there exist some absolute constant $c$ such that for any learning rate $\eta<\frac{c}{\mu dk\kappa^{3}(\mathbf{M})\left\|{\mathbf{M}}\right\|\log d}$, with at least $1-\frac{T}{d^{10}}$ probability, we will have for all $t\leq T$ that:

Note function $\max_{i}g_{i}(\mathbf{U})$ indicates the incoherence of matrix $\mathbf{U}$. Theorem 4.5 guarantees if inital $\mathbf{U}_{0}$ is in the local region which is incoherent and $\mathbf{U}_{0}\mathbf{U}_{0}^{\top}$ is close to $\mathbf{M}$, then with high probability for all steps $t\leq T$, $\mathbf{U}_{t}$, $\mathbf{U}_{t}$ will always stay in a slightly relaxed local region, and $f(\mathbf{U}_{t})$ has linear convergence.

It is not hard to show that all saddle point of $f(\mathbf{U})$ satisfies $\sigma_{k}(\mathbf{U})=0$, and all local minima are global minima. Since $\mathbf{U}_{0},\ldots,\mathbf{U}_{t}$ automatically stay in region $f(\mathbf{U})\leq(\frac{\sigma_{\min}(\mathbf{M})}{10})^{2}$ with high probability, we know $\mathbf{U}_{t}$ also stay away from all saddle points. The claim that $\mathbf{U}_{0},\ldots,\mathbf{U}_{t}$ stays incoherent is essential to better control the variance and probability 1 bound of $SG(\mathbf{U}_{t})$, so that we can have large step size and tight convergence rate.

The major challenging in proving Theorem 4.5 is to both prove $\mathbf{U}_{t}$ stays in the local region, and achieve good sample complexity and running time (linear in $d$) in the same time. This also requires the learning rate $\eta$ in Algorithm 1 to be relatively large. Let the event $\mathfrak{E}_{t}$ denote the good event where $\mathbf{U}_{0},\ldots,\mathbf{U}_{t}$ satisfies Eq.(5). Theorem 4.5 is claiming that $P(\mathfrak{E}_{T})$ is large. The essential steps in the proof is contructing two supermartingles related to $f(\mathbf{U}_{t})1_{\mathfrak{E}_{t}}$ and $g_{i}(\mathbf{U}_{t})1_{\mathfrak{E}_{t}}$ (where $1_{(\cdot)}$ denote indicator function), and use Bernstein inequalty to show the concentration of supermartingales. The $1_{\mathfrak{E}_{t}}$term allow us the claim all previous $\mathbf{U}_{0},\ldots,\mathbf{U}_{t}$ have all desired properties inside local region.

Finally, we see Theorem 3.1 as a immediate corollary of Theorem 4.5.

Conclusion

In this paper, we presented the first provable, efficient online algorithm for matrix completion, based on nonconvex SGD. In addition to the online setting, our results are also competitive with state of the art results in the offline setting. We obtain our results by introducing a general framework that helps us show how SGD updates self-regulate to stay away from saddle points. We hope our paper and results help generate interest in online matrix completion, and our techniques and framework prompt tighter analysis for other nonconvex problems.

References

Appendix A Proof of Initialization

In this section, we will prove Lemma 4.1 and a corresponding lemma for asymmetric case as follows (which will be used to prove Theorem 3.3):

We will focus mostly on Lemma A.1, and prove Lemma 4.1 as a special case. Most of the argument of this section follows from . We include here for completeness. The remaining of this section can be viewed as proving both the Frobenius norm claim and incoherence claim of Lemma A.1 seperately.

In this section, We always denote $d=\max\{d_{1},d_{2}\}$. For simplicity, WLOG, we also assume $\left\|{\mathbf{M}}\right\|=1$ in all proof. Also, when it’s clear from the context, we use $\kappa$ to specifically to represent $\kappa(\mathbf{M})$. Then $\sigma_{\min}(\mathbf{M})=\frac{1}{\kappa}$. Also in the proof, we always denote $\text{SVD}(\mathbf{M})=\mathbf{X}\mathbf{S}\mathbf{Y}^{\top}$, and $\text{SVD}(\mathbf{U}\mathbf{V}^{\top})=\mathbf{W}_{\mathbf{U}}\mathbf{D}\mathbf{W}_{\mathbf{V}}^{\top}$, where $\mathbf{S}$ and $\mathbf{D}$ are $k\times k$ diagonal matrix.

A finite sequence $\{\mathbf{X}_{t}\}$ of independent, random matrices with dimension $d_{!}\times d_{2}$. Assume that each matrix satisfies:

Let $|\Omega|=m$, then there exists universal constant $C,c_{0}$, for any $m\geq c_{0}\mu dk\log d$, with probability at least $1-\frac{1}{d^{10}}$, we have:

where $Z_{ij}$ are independence Bernoulli$(m/d_{1}d_{2})$ random variables. Let matrix

Then, by matrix Bernstein (Theorem A.2), we have:

That is, with probability at least $1-\frac{1}{d^{10}}$, for some universal constant $C$, we have:

For $m\geq\mu dk\log d$, we finishes the proof. ∎

Let $\mathbf{U}_{0}\mathbf{V}_{0}^{\top}$ be the top-$k$ SVD of $\frac{d_{1}d_{2}}{m}\mathcal{P}_{\Omega}(\mathbf{M})$, where $|\Omega|=m$ then there exists universal constant $c_{0}$, for any $m\geq c_{0}\mu dk^{2}\kappa^{2}\log d$, with probability at least $1-\frac{1}{d^{10}}$, we have:

Since $\mathbf{M}$ is a rank $k$ matrix, we know $\sigma_{k+1}(\mathbf{M})=0$, thus

Meanwhile, since $\text{rank}(\mathbf{M})=k$, $\text{rank}(\mathbf{U}_{0}\mathbf{V}_{0}^{\top})=k$, we know: $\text{rank}(\mathbf{M}-\mathbf{U}_{0}\mathbf{V}_{0}^{\top})\leq 2k$, and therefore:

by choosing $m\geq c_{0}\mu dk^{2}\log d\cdot\kappa^{2}$ for large enough constant $c_{0}$ and apply Lemma A.3, we finishes the proof. ∎

A.2 Incoherence of Initialization

Let $\mathbf{U}\mathbf{V}^{\top}$ be the top-$k$ SVD of $\frac{d_{1}d_{2}}{m}\mathcal{P}_{\Omega}(\mathbf{M})$, where $|\Omega|=m$. then there exists universal constant $c_{0}$, for any $m\geq c_{0}\mu dk\kappa^{2}\log d$, with probability at least $1-\frac{1}{d^{10}}$, we have:

For the first term, since $\mathbf{W}_{\mathbf{U}}^{\top}\mathbf{W}_{\mathbf{U},\perp}=0$, we have:

The last step is due to sample $m\geq\mu dk\kappa^{2}\log d$, and theorem A.4.

Therefore, with $m\geq\mu dk\kappa^{2}\log d$, by matrix Bernstein, we have with probability at least $1-\frac{1}{d^{10}}$, we know that for all $j\in[d_{2}]$, there exists some absolute constant $C^{\prime}$ so that:

Combine above two equations and choose $m\geq c_{0}\mu dk\kappa^{2}\log d$ for some large enough $c_{0}$. We have:

Let $\mathbf{U}_{0}\mathbf{V}_{0}^{\top}$ be the top-$k$ SVD of $\frac{d_{1}d_{2}}{m}\mathcal{P}_{\Omega}(\mathbf{M})$, where $|\Omega|=m$. then there exists universal constant $c_{0}$, for any $m\geq c_{0}\mu dk\kappa^{2}\log d$, with probability at least $1-\frac{1}{d^{10}}$, we have:

By Theorem A.5, we know for any $j\in[d_{2}]$:

By symmetry, we also know for any $i\in[d_{1}]$

Let $\mathbf{U}_{0}\mathbf{U}_{0}^{\top}$ be the top-$k$ SVD of $\frac{d^{2}}{m}\mathcal{P}_{\Omega}(\mathbf{M})$, where $|\Omega|=m$. then there exists universal constant $c_{0}$, for any $m\geq c_{0}\mu dk\kappa^{2}\log d$, with probability at least $1-\frac{1}{d^{10}}$, we have:

On the other hand, by Theorem A.4, we have:

Finally, Lemma A.1 can be easily concluded from Theorem A.4 and Theorem A.6, while Lemma 4.1 is also directly proved by Theorem A.4 and Corollary A.7.

Appendix B Proof of Symmetric PSD Case

In this section, we prove Theorem 3.1. WLOG, we continue to assume $\left\|{\mathbf{M}}\right\|=1$ in all proof. Also, when it’s clear from the context, we use $\kappa$ to specifically to represent $\kappa(\mathbf{M})$. Then $\sigma_{\min}(\mathbf{M})=\frac{1}{\kappa}$. Also in this section, we always denote $\text{SVD}(\mathbf{M})=\mathbf{X}\mathbf{S}\mathbf{X}^{\top}$, and $\text{SVD}(\mathbf{U}\mathbf{U}^{\top})=\mathbf{W}\mathbf{D}\mathbf{W}^{\top}$.

The most essential part to prove Theorem 3.1 is proving following Theorem:

Let $f(\mathbf{U})=\left\|{\mathbf{U}\mathbf{U}^{\top}-\mathbf{M}}\right\|_{F}^{2}$ and $g_{i}(\mathbf{U})=\left\|{\mathbf{e}_{i}^{\top}\mathbf{U}}\right\|^{2}$. Suppose after initialization, we have:

Then, there exist some absolute constant $c$ such that for any learning rate $\eta<\frac{c}{\mu dk\kappa^{3}\log d}$, with at least $1-\frac{T}{d^{10}}$ probability, we will have for all $t\leq T$ that:

Theorem B.1 says once initialization algorithm provides $\mathbf{U}_{0}$ in good local region, with high probability $\mathbf{U}_{t}$ will always stay in this good region and $f(\mathbf{U}_{t})$ is linear converging to 0. With this theorem, we can then immediately conclude Theorem 3.1 from Theorem B.1 and Lemma 4.1.

The rest of this section all focus on proving Theorem B.1. First, we prepare with a few lemmas about the property of objective function, and the spectral property of $\mathbf{U}$ in a local Frobenius ball around optimal. Then, we prove Theorem B.1 by constructing two supermartingales related to $f(\mathbf{U}_{t}),g_{i}(\mathbf{U}_{t})$ each, and applying concentration argument.

For symmetric PSD case, we denote the stochastic gradient as:

The update in Algorithm 1 can be now written as:

First, we prove two lemmas w.r.t the smoothness and property similar to strongly convex for objective function:

(restatement of Lemma 4.2) Within the region $\mathcal{D}=\{\mathbf{U}|\left\|{\mathbf{U}}\right\|\leq\Gamma\}$, we have function $f(\mathbf{U})=\|{\mathbf{M}-\mathbf{U}\mathbf{U}^{\top}}\|_{\text{F}}^{2}$ satisfying for any $\mathbf{U}_{1},\mathbf{U}_{2}\in\mathcal{D}$:

where smoothness parameter $\beta=16\max\{\Gamma^{2},\left\|{\mathbf{M}}\right\|\}$.

(restatement of Lemma 4.3) Within the region $\mathcal{D}=\{\mathbf{U}|\sigma_{\min}(\mathbf{X}^{\top}\mathbf{U})\geq\gamma\}$, then we have function $f(\mathbf{U})=\|{\mathbf{M}-\mathbf{U}\mathbf{U}^{\top}}\|_{\text{F}}^{2}$ satisfying:

Inside region $\mathcal{D}$, recall we denote $\mathbf{W}\mathbf{D}\mathbf{W}^{\top}=\text{SVD}(\mathbf{U}\mathbf{U}^{\top})$, thus we have:

Next, we show as long as we are in some Frobenious ball around optimum, then we have good spectral property over $\mathbf{U}$ which guarantees the preconditions for Lemma B.2 and Lemma B.3.

(restatement of Lemma 4.4) Within the region $\mathcal{D}=\{\mathbf{U}|\left\|{\mathbf{M}-\mathbf{U}\mathbf{U}^{\top}}\right\|_{F}\leq\frac{1}{10}\sigma_{k}(\mathbf{M})\}$, we have:

For spectral norm of $\mathbf{U}$, we have:

For the minimum singular value of $\mathbf{U}^{\top}\mathbf{U}$, we have:

Let the principal angle between $\mathbf{X}$ and $\mathbf{W}$ to be $\theta$. This gives $\sin^{2}\theta=\left\|{\mathbf{X}_{\perp}^{\top}\mathbf{W}}\right\|^{2}\leq\frac{1}{9}$. Thus $\cos^{2}\theta=\sigma_{\min}^{2}(\mathbf{X}^{\top}\mathbf{W})\geq\frac{8}{9}$. Therefore:

B.2 Proof of Theorem B.1

Now, we are ready for our key theorem. By Lemma B.2, Lemma B.3, and Lemma B.4, we already know the function has good property locally in the region $\mathcal{D}=\{\mathbf{U}|\left\|{\mathbf{M}-\mathbf{U}\mathbf{U}^{\top}}\right\|_{F}\leq\frac{1}{10}\sigma_{k}(\mathbf{M})\}$ which alludes linear convergence. Then, the work remains and also the most challenging part is to prove that once we initialize inside this region, our algorithm will guarantee $\mathbf{U}$ never leave this region with high probability even with relatively large stepsize. The requirement for tight sample complexity and near optimal runtime makes it more challenging, and require us to further control the incoherence of $\mathbf{U}_{t}$ over all iterates in addition to the distance $\left\|{\mathbf{M}-\mathbf{U}\mathbf{U}^{\top}}\right\|_{F}$.

Define event $\mathfrak{E}_{t}=\{\forall\tau\leq t,f(\mathbf{U}_{\tau})\leq(1-\frac{\eta}{2\kappa})^{t}(\frac{1}{10\kappa})^{2},\max_{i}g_{i}(\mathbf{U}_{\tau})\leq\frac{20\mu k\kappa^{2}}{d}\}$. Theorem B.1 is equivalent to prove event $\mathfrak{E}_{T}$ happens with high probability. The proof achieves this by contructing two supermartingales for $f(\mathbf{U}_{t})1_{\mathfrak{E}_{t}}$ and $g_{i}(\mathbf{U}_{t})1_{\mathfrak{E}_{t}}$ (where $1_{(\cdot)}$ denote indicator function), applies concentration argument.

Their probability 1 bound and variance bound in order to apply Azuma-Bernstein inequality

Final combination of concentration results to conclude the proof

First, let filtration $\mathfrak{F}_{t}=\sigma\{SG(\mathbf{U}_{0}),\cdots,SG(\mathbf{U}_{t-1})\}$ where $\sigma\{\cdot\}$ denotes the sigma field. Note by definiton of $\mathfrak{E}_{t}$, we have $\mathfrak{E}_{t}\subset\mathfrak{F}_{t}$. Also $\mathfrak{E}_{t+1}\subset\mathfrak{E}_{t}$, and thus $1_{\mathfrak{E}_{t+1}}\leq 1_{\mathfrak{E}_{t}}$. Note $\mathfrak{E}_{t}$ denotes the event which up to time $t$, $\mathbf{U}_{\tau}$ always stay in a local region which both close to $\mathbf{M}$ and incoherent.

By Lemma B.4, we immediately know that conditioned on $\mathfrak{E}_{t}$, we have $\left\|{\mathbf{U}_{t}}\right\|\leq\sqrt{2}$, $\sigma_{\min}(\mathbf{X}^{\top}\mathbf{U}_{t})\geq 1/\sqrt{2\kappa}$ and $\sigma_{\min}(\mathbf{U}_{t}^{\top}\mathbf{U}_{t})\geq 1/2\kappa$. We will use this fact throughout the proof.

Since $g_{i}(\mathbf{U})=\mathbf{e}_{i}^{\top}\mathbf{U}\mathbf{U}^{\top}\mathbf{e}_{i}$ is a quadratic function, we know for any change $\Delta\mathbf{U}$, we have:

The last step is true by choosing constant $c$ in learning rate $\eta$ to be small enough.

Let $G_{it}=(1-\frac{4\eta}{\kappa})^{-t}(g_{i}(\mathbf{U}_{t})1_{\mathfrak{E}_{t-1}}-15\frac{\mu k\kappa^{2}}{d})$. This gives:

Probability 1 bound for G𝐺G:

Since when sample $(i,j)$ entry of matrix $\mathbf{M}$, for any $l\in[d]$, we have:

Therefore, by Eq.(9), we have with probability 1:

Variance bound for G𝐺G:

Bernstein’s inequality for G𝐺G:

Since $G_{i0}=g_{i}(\mathbf{U}_{0})-15\frac{\mu k\kappa^{2}}{d}$, let $s^{\prime}=O(1)(1-\frac{4\eta}{\kappa})^{t}[\sqrt{\sigma^{2}\log d}+R\log d]$, we know

by $\eta<\frac{c}{\mu dk\kappa^{3}\log d}$ and choosing $c$ to be small enough, we have:

Since initialization gives $\max_{i}g_{i}(\mathbf{U}_{0})\leq\frac{10\mu k\kappa^{2}}{d}$, therefore:

Construction of supermartingale F𝐹F:

Let $F_{t}=(1-\frac{\eta}{\kappa})^{-t}f(\mathbf{U}_{t})1_{\mathfrak{E}_{t-1}}$, we know $F_{t}$ is also a supermartingale.

Probability 1 bound for F𝐹F:

where $\zeta_{t}$ depends on $SG(\mathbf{U}_{t})$.

First, recall we denote $\text{SVD}(\mathbf{M})=\mathbf{X}\mathbf{S}\mathbf{X}^{\top}$, and $\text{SVD}(\mathbf{U}\mathbf{U}^{\top})=\mathbf{W}\mathbf{D}\mathbf{W}^{\top}$ , and observe that:

Then, when sample $(i,j)$ entry of matrix $\mathbf{M}$, we have:

Therefore, by decomposition Eq.(13), we have with probability 1:

Variance bound for F𝐹F:

Therefore, by decomposition Eq.(13), we have:

Bernstein’s inequality for F𝐹F:

Let $s^{\prime}=O(1)(1-\frac{\eta}{\kappa})^{t}[\sqrt{\sigma^{2}\log d}+R\log d]$, this gives:

by $\eta<\frac{c}{\mu dk\kappa^{3}\log d}$ and choosing $c$ to be small enough, we have:

Since $F_{0}=f(\mathbf{U}_{0})\leq(\frac{1}{20\kappa})^{2}$, therefore:

Finally, combining the concentration result for martingale $G$ (Eq.(12)) and martingale $F$ (Eq.(16)), we conclude:

Appendix C Proof of General Asymmetric Case

In this section, we first prove Lemma 3.2, set up the equivalence between Algorithm 2 and Algorithm 3. Then we prove the main theorem for general asymmetric matrix (Theorem 3.3). WLOG, we continue to assume $\left\|{\mathbf{M}}\right\|=1$ in all proof. Also, when it’s clear from the context, we use $\kappa$ to specifically to represent $\kappa(\mathbf{M})$. Then $\sigma_{\min}(\mathbf{M})=\frac{1}{\kappa}$. Also in this section, we always use $d=\max\{d_{1},d_{2}\}$ and denote $\text{SVD}(\mathbf{M})=\mathbf{X}\mathbf{S}\mathbf{Y}^{\top}$, and $\text{SVD}(\mathbf{U}\mathbf{V}^{\top})=\mathbf{W}_{U}\mathbf{D}\mathbf{W}_{V}^{\top}$.

Clearly with same initialization algorithm, we have $\mathbf{U}_{0}\mathbf{V}_{0}^{\top}=\mathbf{U}^{\prime}_{0}\mathbf{V}^{\prime}_{0}{}^{\top}$, by induction, we finish the proof. ∎

Now we proceed to prove Theorem 3.3. Since Algorithm 2 and Algorithm 3 are equivalent, we will focus our analysis on Algorithm 2 which is more theoretical appealing. As for the symmetric PSD case, we first present the essential ingradient:

Let $f(\mathbf{U},\mathbf{V})=\left\|{\mathbf{U}\mathbf{V}^{\top}-\mathbf{M}}\right\|_{F}^{2}$, $g_{i}(\mathbf{U},\mathbf{V})=\left\|{\mathbf{e}_{i}^{\top}\mathbf{U}\mathbf{V}^{\top}}\right\|^{2}$, and $h_{j}(\mathbf{U},\mathbf{V})=\left\|{\mathbf{e}_{j}^{\top}\mathbf{V}\mathbf{U}^{\top}}\right\|^{2}$, for $i\in[d_{1}]$ and $j\in[d_{2}]$. Suppose after initialization, we have:

Then, there exist some absolute constant $c$ such that for any learning rate $\eta<\frac{c}{\mu dk\kappa^{3}\log d}$, with at least $1-\frac{T}{d^{10}}$ probability, we will have for all $t\leq T$ that:

Theorem 3.3 can easily be concluded from Theorem C.1 and Lemma A.1. Theorem C.1 also provides similar guarantees as Theorem B.1 in symmetric case. However, due to the additional invariance between $\mathbf{U}$ and $\mathbf{V}$, Theorem C.1 need to keep track of more complicated potential function $g_{i}(\mathbf{U},\mathbf{V})$ and $h_{j}(\mathbf{U},\mathbf{V})$ to control the incoherence, which makes the proof more involved.

The rest of this section all focus on proving Theorem C.1. Similar to symmetric PSD case, we also first prepare with a few lemmas about the property of objective function, and the spectral property of $\mathbf{U},\mathbf{V}$ in a local Frobenius ball around optimal. Then, we prove Theorem C.1 by constructing three supermartingales related to $f(\mathbf{U}_{t},\mathbf{V}_{t}),g_{i}(\mathbf{U}_{t},\mathbf{V}_{t}),h_{j}(\mathbf{U}_{t},\mathbf{V}_{t})$ each, and applying concentration argument.

Also denote the stochastic gradient $SG(\mathbf{U},\mathbf{V})$ by (if we sampled entry $(i,j)$ of matrix $\mathbf{M}$)

Similar to symmetric PSD case, we first prove two lemmas w.r.t the smoothness and property similar to strongly convex for objective function:

Within the region $\mathcal{D}=\{(\mathbf{U},\mathbf{V})|\left\|{\mathbf{U}}\right\|\leq\Gamma,\left\|{\mathbf{V}}\right\|\leq\Gamma\}$, we have function $f(\mathbf{U},\mathbf{V})=\|{\mathbf{M}-\mathbf{U}\mathbf{V}^{\top}}\|_{\text{F}}^{2}$ satisfying:

where smoothness parameter $\beta=8\max\{\Gamma^{2},\left\|{\mathbf{M}}\right\|\}$.

The last step is by similar technics as in the proof of Lemma B.2, by expanding

Within the region $\mathcal{D}=\{(\mathbf{U},\mathbf{V})|\sigma_{\min}(\mathbf{X}^{\top}\mathbf{U})\geq\gamma,\sigma_{\min}(\mathbf{Y}^{\top}\mathbf{V})\geq\gamma\}$, then we have function $f(\mathbf{U},\mathbf{V})=\|{\mathbf{M}-\mathbf{U}\mathbf{V}^{\top}}\|_{\text{F}}^{2}$ satisfying:

Let $\hat{\mathbf{U}},\hat{\mathbf{V}}$ be the left singular vectors of $\mathbf{U},\mathbf{V}$. Inside region $\mathcal{D}$, we have:

Next, we show as long as we are in some Frobenious ball around optimum, then we have good spectral property over $\mathbf{U},\mathbf{V}$ which guarantees the preconditions for Lemma C.2 and Lemma C.3.

Within the region $\mathcal{D}=\{(\mathbf{U},\mathbf{V})|\left\|{\mathbf{M}-\mathbf{U}\mathbf{V}^{\top}}\right\|_{F}\leq\frac{1}{10}\sigma_{k}(\mathbf{M})\}$, and for $\mathbf{U}=\mathbf{W}_{U}\mathbf{D}^{\frac{1}{2}},\mathbf{V}=\mathbf{W}_{V}\mathbf{D}^{\frac{1}{2}}$ where $\mathbf{W}_{U}\mathbf{D}\mathbf{W}_{V}=\text{SVD}(\mathbf{U}\mathbf{V}^{\top})$, we have:

For spectral norm of $\mathbf{U}$, we have:

For the minimum singular value of $\mathbf{U}^{\top}\mathbf{U}$, we have:

By symmetry, the same holds for $\mathbf{V}$. On the other hand, we have:

Let the principal angle between $\mathbf{X}$ and $\mathbf{W}_{U}$ to be $\theta$. This gives $\sin^{2}\theta=\left\|{\mathbf{X}_{\perp}^{\top}\mathbf{W}_{U}}\right\|^{2}\leq\frac{1}{9}$. Thus $\cos^{2}\theta=\sigma_{\min}^{2}(\mathbf{X}^{\top}\mathbf{W}_{U})\geq\frac{8}{9}$. Therefore:

C.2 Proof of Theorem C.1

Now, we are ready for our key theorem. By Lemma C.2, Lemma C.3, and Lemma C.4, we already know the function has good property locally in the region $\mathcal{D}=\{(\mathbf{U},\mathbf{V})|\left\|{\mathbf{M}-\mathbf{U}\mathbf{V}^{\top}}\right\|_{F}\leq\frac{1}{10}\sigma_{k}(\mathbf{M})\}$ which alludes linear convergence. Similar to the symmetric PSD case, the work remains is to prove that once we initialize inside this region, our algorithm will guarantee $\mathbf{U},\mathbf{V}$ never leave this region with high probability even with relatively large stepsize. Again, we also need to control the incoherence of $\mathbf{U}_{t},\mathbf{V}_{t}$ over all iterates additionally to achieve tight sample complexity and near optimal runtime.

For simplicity of notation, we assume $d=d_{1}=d_{2}$, and do not distinguish $d_{1}$ and $d_{2}$. However, it is easy to check our proof never use the property $\mathbf{M}$ is square matrix. The proof easily extends to $d_{1}\neq d_{2}$ case by replacing $d$ in the proof with suitable $d_{1},d_{2}$.

Define event $\mathfrak{E}_{t}=\{\forall\tau\leq t,f(\mathbf{U}_{\tau},\mathbf{V}_{\tau})\leq(1-\frac{\eta}{2\kappa})^{t}(\frac{1}{10\kappa})^{2},\max_{i}g_{i}(\mathbf{U}_{\tau},\mathbf{V}_{\tau})\leq\frac{100\mu k\kappa^{2}}{d},\max_{j}h_{j}(\mathbf{U}_{\tau},\mathbf{V}_{\tau})\leq\frac{100\mu k\kappa^{2}}{d}\}$. Theorem C.1 is equivalent to prove event $\mathfrak{E}_{T}$ happens with high probability. The proof achieves this by contructing two supermartingales for $f(\mathbf{U}_{t},\mathbf{V}_{t})1_{\mathfrak{E}_{t}}$, $g_{i}(\mathbf{U}_{t},\mathbf{V}_{t})1_{\mathfrak{E}_{t}}$ and $h_{i}(\mathbf{U}_{t},\mathbf{V}_{t})1_{\mathfrak{E}_{t}}$ (where $1_{(\cdot)}$ denote indicator function), applies concentration argument.

The proofs also follow similar structure as symmetric PSD case:

Their probability 1 bound and variance bound in order to apply Azuma-Bernstein inequality

Final combination of concentration results to conclude the proof

Then let filtration $\mathfrak{F}_{t}=\sigma\{SG(\mathbf{U}_{0},\mathbf{V}_{0}),\cdots,SG(\mathbf{U}_{t-1},\mathbf{V}_{t-1})\}$ where $\sigma\{\cdot\}$ denotes the sigma field. Also let event , note $\mathfrak{E}_{t}\subset\mathfrak{F}_{t}$. Also $\mathfrak{E}_{t+1}\subset\mathfrak{E}_{t}$, and thus $1_{\mathfrak{E}_{t+1}}\leq 1_{\mathfrak{E}_{t}}$.

By Lemma C.4, we immediately know that conditioned on $\mathfrak{E}_{t}$, we have $\left\|{\mathbf{U}_{t}}\right\|\leq\sqrt{2}$, $\left\|{\mathbf{V}_{t}}\right\|\leq\sqrt{2}$, $\sigma_{\min}(\mathbf{X}^{\top}\mathbf{U}_{t})\geq 1/\sqrt{2\kappa}$, $\sigma_{\min}(\mathbf{Y}^{\top}\mathbf{V}_{t})\geq 1/\sqrt{2\kappa}$. We will use this fact throughout the proof.

For simplicity, when it’s clear from the context, we denote:

First, since potential function $g_{l}(\mathbf{U},\mathbf{V})$ is forth-order polynomial, we can expand:

Where we denote $R_{1}$ as the sum of first order terms and higher order terms (all second/third/forth order terms), and $R_{2}$ as the sum of second order terms and higher order terms.

We now give a proposition about properties of $R_{1}$ and $R_{2}$ which involves a lot calculation, and postpone its proof in the end of this section.

With above notations, we have following inequalities hold true.

Then by taking conditional expectation, we have:

The first order term can be calculated as:

In second last inequality, we use key observation:

The last inequality is achieved by choosing $c$ small enough.

The right inequality is true since $1_{\mathfrak{E}_{t}}\leq 1_{\mathfrak{E}_{t-1}}$. This implies $G_{it}$ is supermartingale.

Probability 1 bound for G𝐺G:

By Proposition C.5, we know with probability 1 that $|R_{1}|1_{\mathfrak{E}_{t}}\leq\eta O(\mu^{2}k^{2}\kappa^{5})1_{\mathfrak{E}_{t}}$. This gives with probability 1:

Variance bound for G𝐺G:

Bernstein’s inequality for G𝐺G:

by $\eta<\frac{c}{\mu dk\kappa^{3}\log d}$ and choosing $c$ to be small enough, we have:

Since initialization gives $\max_{i}g_{i}(\mathbf{U}_{0},\mathbf{V}_{0})\leq\frac{10\mu k\kappa^{2}}{d}$, therefore:

Construction of supermartingale F:

Where we denote $Q_{1}$ as the sum of first order terms and higher order terms (all second/third/forth order terms), and $Q_{2}$ as the sum of second order terms and higher order terms.

We also now give a proposition about properties of $Q_{1}$ and $Q_{2}$ which involves a lot calculation, and postpone its proof in the end of this section.

With above notations, we have following inequalities hold true.

Probability 1 bound:

By Proposition C.6, we know with probability 1 that $|Q_{1}|1_{\mathfrak{E}_{t}}\leq\eta O(\mu dk\kappa^{3})f(\mathbf{U}_{t},\mathbf{V}_{t})1_{\mathfrak{E}_{t}}$. This gives with probability 1:

Variance bound:

Bernstein’s inequality:

Let $s^{\prime}=O(1)(1-\frac{\eta}{\kappa})^{t}[\sqrt{\sigma^{2}\log d}+R\log d]$, this gives:

by $\eta<\frac{c}{\mu dk\kappa^{3}\log d}$ and choosing $c$ to be small enough, we have:

Since $F_{0}=f(\mathbf{U}_{0})\leq(\frac{1}{20\kappa})^{2}$, therefore:

Finally, combining the concentration result for martingale $G$ (Eq.(19)) and martingale $F$ (Eq.(22)), we conclude:

Finally we give proof for Proposition C.5 and Proposition C.6. The proof mostly consistsof expanding every term and careful calculations.

For simplicity of notation, we hide the term $1_{\mathfrak{E}_{t}}$ in all following equations. Reader should always think every term in this proof multiplied by $1_{\mathfrak{E}_{t}}$. Recall that:

We first prove first three inequality. Recall that:

By expanding the polynomial, we can write out the first order term:

and recall we choose $\eta<\frac{c}{\mu dk\kappa^{3}\log d}$, where $c$ is some universal constant, then we have:

With equation (23), (24), (25), (26), now we are ready to prove Lemma.

This gives in sum that, with probability 1:

Similarly to the proof of Proposition C.5, we hide the term $1_{\mathfrak{E}_{t}}$ in all following equations. Reader should always think every term in this proof multiplied by $1_{\mathfrak{E}_{t}}$. Recall that:

By expanding the polynomial, we can write out the first order term:

Again, in addition to equation (23), (23), (24), (25), we also need following inequality:

This gives in sum that, with probability 1: