Fast matrix completion without the condition number

Introduction

Matrix Completion is the problem of recovering an unknown real-valued low-rank matrix from a possibly noisy subsample of its entries. The problem has received a tremendous amount of attention in signal processing and machine learning partly due to its wide applicability to recommender systems. A beautiful line of work showed that a particular convex program—known as nuclear norm minimization—achieves strong recovery guarantees under certain reasonable feasibility assumptions [CR09, CT10, RFP10, Rec11]. Nuclear norm minimization boils down to solving a semidefinite program and therefore can be solved in polynomial time in the dimension of the matrix. Unfortunately, the approach is not immediately practical due to the large polynomial dependence on the dimension of the matrix. An ongoing research effort aims to design large-scale algorithms for nuclear norm minimization [JY09, MHT10, JS10, AKKS12, HO14]. Such fast solvers, generally speaking, involve heuristics that improve empirical performance but may no longer preserve the strong theoretical guarantees of the nuclear norm approach.

A successful scalable algorithmic alternative to Nuclear Norm Minimization is based on Alternating Minimization [BK07, HH09, KBV09]. Alternating Minimization aims to recover the unknown low-rank matrix by alternatingly optimizing over one of two factors in a purported low-rank decomposition. Each update is a simple least squares regression problem that can be solved very efficiently. As pointed out in [HO14], even state of the art nuclear norm solvers often cannot compete with Alternating Minimization with regards to scalability. A shortcoming of Alternating Minimization is that formal guarantees are less developed than for Nuclear Norm Minimization. Only recently has there been progress in this direction [Kes12, JNS13, GAGG13, Har13a].

Unfortunately, despite this recent progress all known convergence bounds for Alternating Minimization have at least a quadratic dependence on the condition number of the matrix. Here, the condition number refers to the ratio of the first to the $k$-th singular value of the matrix, where $k$ is the target rank of the decomposition. This dependence on the condition number can be a serious shortcoming. After all, Matrix Completion rests on the assumption that the unknown matrix is approximately low-rank and hence we should expect its singular values to decay rapidly. Indeed, strongly decaying singular values are a typical feature of large real-world matrices.

The dependence on the condition number in Alternating Minimization is not a mere artifact of the analysis. It arises naturally with the use of the Singular Value Decomposition (SVD). Alternating Minimization is typically intialized with a decomposition based on a truncated SVD of the partial input matrix. Such an approach must incur a polynomial dependence on the condition number. Many other approaches also crucially rely on the SVD as a sub-routine, e.g., [JMD10, KMO10a, KMO10b], as well as most fast solvers for the nuclear norm. In fact, there appears to be a kind of dichotomy in the current literature on Matrix Completion: either the algorithm is not fast and has at least a quadratic dependence on the dimension of the matrix in its running time, or it is not well-conditioned and has at least a quadratic dependence on the condition number in the sample complexity. We emphasize that here we focus on formal guarantees rather than observed empirical performance which may be better on certain instances. This situation leads us to the following problem.

Main Problem: Is there a sub-quadratic time algorithm for Matrix Completion

with a sub-linear dependence on the condition number?

In fact, eliminating the polynomial dependence on the condition number was posed explicitly as an open problem in the context of Alternating Minimization by Jain, Netrapalli and Sanghavi [JNS13].

In this work, we resolve the question in the affirmative. Specifically, we design a new variant of Alternating Minimization that achieves a logarithmic dependence on the condition number while retaining the fast running time of the standard Alternating Minimization framework. This is an exponential improvement in the condition number compared with all subquadratic time algorithms for Matrix Completion that we are aware of. Our algorithm works even in the noisy Matrix Completion setting and under standard assumptions—specifically, the same assumptions that support theoretical results for the nuclear norm. That is, we assume that the first $k$ singular vector of the matrix span an incoherent subspace and that each entry of the matrix is revealed independently with a certain probability. While strong, these assumptions led to an interesting theory of Matrix Completion and have become a de facto standard when comparing theoretical guarantees.

For the sake of exposition we begin by explaining our results in the exact Matrix Completion setting, even though our results here are a direct consequence of our theorem for the noisy case. In the exact problem the goal is to recover an unknown rank $k$ matrix $M$ from a subsample $\Omega\subset[n]\times[n]$ of its entries where each entry is included independently with probability $p.$ We assume that the unknown matrix $M=U\Lambda U^{T}$ is a symmetric $n\times n$ matrix with nonzero singular values $\sigma_{1}\geq\dots\geq\sigma_{k}>0.$ Following [Har13a], our result generalizes straightforwardly to rectangular matrices. To state our result we need to define the coherence of the subspace spanned by $U.$ Intuitively, the coherence controls how large the projection is of any standard basis vector onto the space spanned by $U.$ Formally, for a $n\times k$ matrix $U$ with orthonormal columns, we define the coherence of $U$ to be

We show that our algorithm outputs a low-rank factorization $XY^{T}$ such that with high probability $\|M-XY^{T}\|^{2}\leq\varepsilon\left\|M\right\|$ provided that the expected size of $\Omega$ satisfies

Here, the exponent $c>0$ is bounded by an absolute constant. While we did not focus on minimizing the exponent, our results imply that the value of $c$ can be chosen smaller if the singular values of $M$ are well-separated. The formal statement follows from \hyperref[thm:main]Theorem 1. A notable advantage of our algorithm compared to several fast algorithms for Matrix Completion is that the dependence on the error $\varepsilon$ is only poly-logarithmic. This linear convergence rate makes near exact recovery feasible with a small number of steps.

We now discuss our more general result that applies to the noisy or robust Matrix Completion problem. Here, the unknown matrix is only close to low-rank, typically in Frobenius norm. Our results apply to any matrix of the form

where $M=U\Lambda U^{T}$ is a matrix of rank $k$ as before and $N=(I-UU^{T})A$ is the part of $A$ not captured by the dominant singular vectors. We note that $N$ can be an arbitrary deterministic matrix. The assumption that we will make is that $N$ satisfies the following incoherence conditions:

Recall that $e_{i}$ denotes the $i$-th standard basis vector so that $\left\|e_{i}^{T}N\right\|_{2}$ is the Euclidean norm of the $i$-th row of $N.$ The conditions state no entry of $N$ should be too large compared to the norm of the corresponding row in $N,$ and no row of $N$ should be too large compared to $\sigma_{k}.$ Our bounds will be in terms of a combined coherence parameter $\mu^{*}$ satisfying

We show that our algorithm outputs a rank $k$ factorization $XY^{T}$ such that with high probability

where $\left\|\cdot\right\|$ denotes the spectral norm. It follows from our argument that we can have the same guaranteee in Frobenius norm as well. To achieve the above bound we show that it is sufficient to have an expected sample size

For an extended discussion of related work see \hyperref[sec:related]Section 2.2. We proceed in the next section with a detailed proof overview and a description of our notation.

Preliminaries

In this section, we will give an overview of our proof, give a more in-depth survey of previous work, and set notation.

As the proof of our main theorem is somewhat complex we will begin with an extensive informal overview of the argument. In order to understand our main algorithm, it is necessary to understand the basic Alternating Minimization algorithm first.

A natural idea ot fix these problems is the so-called deflation approach. If it so happens that $\sigma_{1}\gg\sigma_{k},$ then there must be an $r<k$ such that $\sigma_{1}\approx\sigma_{r}\gg\sigma_{k}.$ In this case, we can try to first run Alternating Minimization with $r$ vectors instead of $k$ vectors. This results in a rank $r$ factorization $XY^{T}.$ We then subtract this matrix off of the original matrix and continue with $A^{\prime}=A-XY^{T}.$ This approach was in particular suggested by Jain et al. [JNS13] to eliminate the condition number dependence. Unfortunately, as we will see next, this approach runs into serious issues.

Given any algorithm NoisyMC for noisy matrix completion, whose performance depends on the condition number of $A$, we may hope to use NoisyMC in a black-box way to obtain a deflation-based algorithm which does not depend on the condition number, as follows. Suppose that we know that the spectrum of $A$ comes in blocks,

and so on. We could imagine running NoisyMC on $P_{\Omega}(A)$ with target rank $r_{1}$, to obtain an estimate $M^{(1)}$. Then we may run NoisyMC again on $P_{\Omega}(A-M^{(1)})=P_{\Omega}(A)-P_{\Omega}(M^{(1)})$ with target rank $r_{2}-r_{1}$, to obtain $M^{(2)}$, and so on. At the end of the day, we would hope to approximate $A\approx M^{(1)}+M^{(2)}+\cdots$. Because we are focusing only on a given “flat" part of the spectrum at a time, the dependence of NoisyMC on the condition number should not matter. A major problem with this approach is that the error builds up rather quickly. More precisely, any matrix completion algorithm run on $A$ with target rank $r_{1}$ must have error on the order of $\sigma_{r_{1}+1}$ since this is the spectral norm of the “noise part” that prevents the algorithm from converging further. Therefore, the matrix $A-M^{(1)}$ might now have $2r_{1}$ problematic singular vectors corresponding to relatively large singular values, namely those vectors arising from the residuals of the first $r_{1}$ singular vectors, as well as those arising from the approximation error. This multiplicative blow-up makes it difficult to ensure convergence.

The above intuition may make a “deflation”-based argument seem hopeless. We instead use an approach that looks similar to deflation but makes an important departure from it. Intuitively, our algorithm is a single execution of Alternating Minimization. However, we dynamically grow the number of vectors that Alternating Minimization maintains until we’ve reached $k$ vectors. At that point we let the algorithm run to convergence. More precisely, the algorithm proceeds in at most $k$ epochs. Each epoch roughly proceeds as follows:

At the beginning of epoch $t,$ the algorithm has a rank $r_{t-1}$ factorization $X_{t-1}Y_{t-1}^{T}$ that has converged to within error $\sigma_{r_{t-1}+1}/100.$ At this point, the $(r_{t-1}+1)$-th singular vector prevents further convergence.

What can we say about the matrix $A_{t}=A-X_{t-1}Y_{t-1}^{T}$ at this point? We know that the first $r_{t-1}$ singular vectors of $A$ are removed from the top of the spectrum of $A_{t}.$ Moreover, each of the remaining singular vectors in $A$ is preserverd so long as the corresponding singular value is greater than $\sigma_{r_{t-1}+1}/10.$ This follows from perturbation bounds and we ignore a polynomial loss in $k$ at this point. Importantly, the top of the spectrum of $A_{t}$ corresponds is correlated with the next block of singular vectors in $A.$ This motivates the next step in epoch $t,$ which is to compute the top $k-r_{t-1}$ singular vectors of $A_{t}$ up to an approximation error of $\sigma_{r_{t-1}+1}/10.$ Among these singular vectors we now identify a gap in singular values, that is we look for a number $d_{t}$ such that $\sigma_{r_{t-1}+d_{t}}\leq\sigma_{r_{t-1}+1}/2.$

We call this algorithm SoftDeflate. The crucial difference to the deflation approach is that we always run Alternating Minimization on a subsampling $P_{\Omega}(A)$ of the original matrix $A$. We only ever compute a deflated matrix $P_{\Omega}(A-XY^{T})$ for the purpose of initializing the next epoch of the algorithm. This prevents the error accumulation present in the basic deflation approach.

This simple description glosses over many details and there are a few challenges to be overcome in order to make the idea work. For example, we have not said how to determine the appropriate “gaps" $d_{t}$. This requires a little bit of care. Indeed, these gaps might be quite small: if the (additive) gap between $\sigma_{r}$ and $\sigma_{r+1}$ is on the order of, say, $\frac{\log^{2}(k)}{k}\sigma_{r}$, for all $r\leq k$, then the condition number of the matrix may be super-polynomial in $k$, a price we are not willing to pay. Thus, we need to be able to identify gaps between $\sigma_{r}$ and $\sigma_{r+1}$ which are on the order of $\sigma_{r}/k$. To do this, we must make sure that our estimates of the singular values of $A-X_{t-1}Y_{t-1}^{T}$ are sufficiently precise.

Another major issue that such an algorithm faces is that of coherence. As mentioned above, incoherence is a standard (and necessary) requirement of matrix completion algorithms, and so in order to pursue the strategy outlined above, we need to be sure that the estimates $X_{t-1}$ stay incoherent. For our first “rough estimation" step, our algorithm carefully truncates (entrywise) its estimates, in order to preserve the incoherence conditions, without introducing too much error. In particular, we cannot reuse the truncation analysis of Jain et al. [JNS13] which incurred a dependence on the condition number. Coherence in the Alternating Minimization step is handled by the algorithm and analysis of [Har13a], upon which we build. Specifically, Hardt used a form of regularization by noise addition called SmoothQR, as well as an extra step which involves taking medians, which ensures that various iterates of Alternating Minimization remain incoherent.

2 Further Discussion of Related Work

Our work is most closely related to recent works on convergence bound for Alternating Minimization [Kes12, JNS13, GAGG13, Har13b]. Our bounds are in general incomparable. We achieve an exponential improvement in the condition number compared to all previous works, while losing polynomial factors in $k$. Our algorithm and analysis crucially builds on [Har13a]. In particular we use the version and analysis of Alternating Minimization derived in that work more or less as a black box. We note that the analyses of Alternating Minimization in other previous works would not be sufficiently strong to be used in our algorithm. In particular, the use of noise addition to ensure coherence already gets rid of one source of the condition number that all previous papers incur.

We are not aware of a fast nuclear norm solver that has theoretical guarantees that do not depend polynomially on the condition number. The work of Keshavan et al. [KMO10a, KMO10b] gives another alternative to nuclear norm minimization that has theoretical guarantees. However, these bounds have a quartic dependence on the condition number. We are not aware of any fast nuclear norm solver with theoretical guarantees that do not depend polynomially on the condition number. The work of Keshavan et al. [KMO10a, KMO10b] gives another alternative to nuclear norm minimization that has theoretical guarantees. However, these bounds have a quartic dependence on the condition number. There are a number of fast algorithms for matrix completion: for example, based on (Stochastic) Gradient Descent [RR13]; (Online) Frank-Wolfe [JS10, HK12]; or CoSAMP [LB10]. However, the theoretical guarantees for these algorithms are typically in terms of the error on the observed entries, rather than on the error between the recovered matrix and the unknown matrix itself. For the matrix completion problem, convergence on observations does not imply convergence on the entire matrix. For some matrix recovery problems—in particular, those where the observations obey a rank-restricted isometry property—convergence on the observations is enough to imply convergence on the entire matrix. However, for matrix completion, the relevant operator does not satisfy this condition [CR09]. Further, these algorithms typically have polynomial, rather than logarithmic, dependence on the accuracy parameter $\varepsilon$. Since setting $\varepsilon\approx\sigma_{k}/\sigma_{1}$ is required in order to accurately recover the first $k$ singular vectors of $A$, a polynomial dependence in $\varepsilon$ implies a polynomial dependence on the condition number.

3 Notation

to the be matrix $A$, restricted to the entries indexed by $\Omega$ and renormalized.

Our algorithm, and its proof, will involve choosing a sequence of integers $r_{1}<\cdots<r_{t}\leq k$, which will mark the significant “gaps” in the spectrum of $A$. Given such a sequence, we will decompose $A$ as

where $M^{(\leq t)}$ has the spectral decomposition $M^{(\leq t)}={U}^{(\leq t)}\Lambda_{(\leq t)}({U}^{(\leq t)})^{T}$ and $\Lambda_{(\leq t)}$ contains the eigenvalues corresponding to singular values $\sigma_{1}\geq\cdots\geq\sigma_{r_{t}}$. We may decompose $M^{(\leq t)}$ as the sum of $M^{(j)}$ for $j=1\dots t,$ where each $M^{(j)}$ has the spectral decomposition $M^{(j)}=U^{(j)}\Lambda_{j}\left(U^{(j)}\right)^{T}$ corresponding to the singular values $\sigma_{r_{j-1}+1},\ldots,\sigma_{r_{j}}$. Similarly, the matrix $N_{t}$ may be written as $N_{t}=(V_{t})\Lambda_{(>t)}(V_{t})^{T},$ and contains the singular values $\sigma_{r_{t}+1},\ldots,\sigma_{n}$. Eventually, our algorithm will stop at some maximum $t=T$, for which $r_{t}=k$, and we will have $A=M+N=M^{(\leq T)}+N_{T}$ as in (2). We will use the notation $U^{(\leq j)}$ to denote the concatenation

For an index $r\leq n$, we quantify the gap between $\sigma_{r}$ and $\sigma_{r+1}$ by

By definition, we always have $\gamma\geq 1/4k$; for some matrices $A$, it may be much larger, and this will lead to improved bounds. Our analysis will also depend on the “final" gap quantified by $\gamma_{k}$, whether or not it is larger than $1/4k$. To this end, we define

Algorithms and Results

In Algorithm 1 we present our main algorithm SoftDeflate. It uses several subroutines that are presented in \hyperref[sec:subroutines]Section 3.1.

In the Matrix Completion literature, the most common assumption on the distribution of the set $\Omega$ of observed entries is that each index $(i,j)$ is included independently with some probability $p$. Call this distribution $\mathcal{D}(p)$. In order for our results to be comparable with existing results, this is the model we adopt as well. However, for our analysis, it is much more convenient to imagine that $\Omega$ is the union of several subsets $\Omega_{t}$, so that the $\Omega_{t}$ themselves follow the distribution $\mathcal{D}(p_{t})$ (for some probability $p_{t}$, where $\sum_{t}p_{t}=p$), and so that all of the $\Omega_{t}$ are independent. Algorithmically, the easiest thing to do to obtain subsets $\Omega_{t}$ from $\Omega$ is to partition $\Omega$ into random subsets of equal size. However, if we do this, the subsets $\Omega_{t}$ will not follow the right distribution; in particular they will not be independent. For theoretical completeness, we show in Appendix A (Algorithm 6) how to split up the set $\Omega$ in the correct way. More precisely, given $p_{t}$ and $p$ so that $\sum_{t}p_{t}=p$, we show how to break $\Omega\sim\mathcal{D}(p)$ into (possibly overlapping) subsets $\Omega_{t}$, so that the $\Omega_{t}$ are independent and each $\Omega_{t}\sim\mathcal{D}(p_{t})$.

SoftDeflate uses a number of subroutines that we outline here before explicitly presenting them:

S-M-AltLS (Algorithm 2) is the main Alternating Least Squares procedure that was given and analyzed in [Har13a]. We use this algorithm and its analysis. S-M-AltLS by itself has a quadratic dependence on the condition number which is why we can only use it as a subroutine.

SmoothQR (Algorithm 3) is a subroutine of S-M-AltLS which is used to control the coherence of intermediate solutions arising in S-M-AltLS. Again, we reuse the analysis of SmoothQR from [Har13a]. SmoothQR orthonormalizes its input matrix after adding a Gaussian noise matrix. This step allows tight control of the coherence of the resulting matrix. We defer the description of SmoothQR to Section 6 where we need it for the first time.

SubsIt is a standard textbook version of the Subspace Iteration algorithm (Power Method). We use this algorithm as a fast way to approximate the top singular vectors of a matrix arising in SoftDeflate. We use only standard properties of SubsIt in our analysis. For this reason we defer the description and analysis of SubsIt to \hyperref[sec:subsit]Section B.3.

2 Statement of the main theorem

$\mu^{*}$ satisfies (4) and $\mu_{0}\geq\frac{C}{(\gamma^{*})^{2}}\left(\mu^{*}\left(k+\left(\frac{k^{4}\Delta}{\varepsilon\sigma_{1}}\right)^{2}\right)+\log(n)\right)$

$L_{t}\geq\frac{C}{\gamma^{*}}\log\left(\frac{k\sigma_{r_{t}}}{\sigma_{r_{t}+1}+\varepsilon\sigma_{1}}\right),$ and $L\geq Ck^{7/2}\log(n)$

There is a choice of $p_{t},p_{t}^{\prime}$ (given in the proof below) so that

Suppose that each element of $[n]\times[n]$ is included in $\Omega$ independently with probability $p$. Then the matrices $X,Y$ returned by SoftDeflate satisfy with probability at least $1-1/n,$

The guarantee of $\left\|A-XY^{T}\right\|\leq\left(1+o(1)\right)\left\|N\right\|+\varepsilon\left\|M\right\|$ is what naturally falls out of our analysis: the natural value for the $o(1)$ term is polynomially small in $k$. It is not hard to see in the proof that we may make this term as small as we like, say, $(1+\alpha)\left\|N\right\|$, by paying a logarithmic penalty $\log(1/\alpha)$ in the choice of $p$. It is also not hard to see that we may have a similar conclusion for the Frobenius norm.

As written, then algorithm requires the user to know several parameters which depend on the unknown matrix $A$. For some parameters, these requirements are innocuous. For example, to obtain $p_{t}^{\prime}$ or $L_{t}$ (whose values are given in Section 4.1), the user is required to have a bound on $\log(\sigma_{r_{t}}/\sigma_{r_{t}+1})$. Clearly, a bound on the condition number of $A$ will suffice, but more importantly, the estimates $s_{t}$ which appear in Algorithm 1 may be used as proxies for $\sigma_{r_{t}}$, and so the parameters $p_{t}^{\prime}$ can actually be determined relatively precisely on the fly. For other parameters, like $\mu^{*}$ or $k$, we assume that the user has a good estimate from other sources. While this is standard in the Matrix Completion literature, we acknowledge that these values may be difficult to come by.

3 Running Time

The running time of SoftDeflate is linear in $n$, polynomial in $k$, and logarithmic in the condition number $\sigma_{1}/\sigma_{k}$ of $A$. Indeed, the outer loop performs at most $k$ epochs, and the nontrivial operations in each epoch are S-M-AltLS, QR, and SubsIt. All of the other operations (truncation, concatenation) are done on matrices which are either $n\times k$ (requiring at most $nk$ operations) or on the subsampled matrices $P_{\Omega_{t}}(A)$, requiring on the order of $pn^{2}$ operations.

Running SubsIt requires $L=O(k^{7/2}\log(n))$ iterations; each iteration includes multiplication by a sparse matrix, followed by QR. The matrix multiplication takes time on the order of

where the $\widetilde{O}$ hides logarithmic factors in $n$.

Proof of Main Theorem

In this section, we prove Theorem 1. The proof proceeds by maintaining a few inductive hypotheses, given below, at each epoch. When the algorithm terminates, we will show that the fact that these hypotheses still hold imply the desired results. Suppose that at the beginning of step $t$ of Algorithm 1, we have identified some indices $r_{1},\ldots,r_{t-1}$, and recovered estimates $X_{t-1},Y_{t-1}$ which capture the singular values $\sigma_{1},\ldots,\sigma_{r_{t-1}}$ and the corresponding singular vectors. The goals of the current step of Algorithm 1 are then to (a) identify the next index $r_{t}$ which exhibits a large “gap" in the spectrum, and (b) estimate the singular values $\sigma_{r_{t-1}+1},\ldots,\sigma_{r_{t}}$ and the corresponding singular vectors.

Letting $r_{t}$ be the index obtained by Algorithm 1, we will decompose $A=M^{(<t)}+N_{t-1}=M^{(\leq t)}+N_{t}$ as in (6). To help keep the notation straight, we include a diagram below, which indicates which singular values of $A$ are included in which matrix.

We will maintain the following inductive hypotheses. At the beginning of epoch $t$ of SoftDeflate, we assert

for some sufficiently large constant $C_{0}$ determined by the proof. We also maintain that the current estimate $X_{t-1}$ is incoherent:

for a constant $C_{5}$. Above, equation (H3) defines $\mu_{t-1}$. Observe that when $t=1$, everything in sight is zero and the hypotheses (H1), (H2),(H3) are satisfied. Finally, we assume that the estimate $s_{t-1}$ of $\sigma_{r_{t-1}+1}$ is good.

where we used the incoherence bounds (33) and (34) in the appendix to bound $\left\|A\right\|_{\infty}$ and $\left\|e_{i}^{T}A\right\|_{2}$. Thus, as long as

Now, suppose that the inductive hypotheses (H1), (H2), (H3), and (H4) hold. We break up the inner loop of SoftDeflate into two main steps. In the first step, lines 1 to 1 in Algorithm 1, the goal is to obtain an estimate $r_{t}$ of the next “gap," as well as an estimate $W_{t}$ of the subspace $U^{(\leq t)}$. We analyze this step in Lemma 2 below.

There exists a constants $C,C_{1}$ so that the following holds. Suppose that

where $\varepsilon_{0}\leq\frac{1}{4C_{1}k^{5/2}}.$ Further assume that the inductive hypotheses (H1), (H2), (H3), and (H4) hold. Then with probability at least $1-1/n^{2}$ over the choice of $\Omega_{t}$ and the randomness in SubsIt, one of the following statements must hold:

Algorithm 1 terminates at line 1, and returns $X_{t-1},Y_{t-1}$ so that $\left\|A-X_{t-1}Y_{t-1}^{T}\right\|\leq C\varepsilon\left\|M\right\|$; or

Algorithm 1 does not terminate at line 1, and the following conditions hold:

The error level $\varepsilon$ has not yet been reached:

The matrix $W_{t}$ has orthonormal columns, and satisfies

The proof of Lemma 2 is given in Section 5. In the second part of SoftDeflate, lines 1 to 1 in Algorithm 1, we run S-M-AltLS, initialized with the subspace $W_{t}$ returned by the first part of the algorithm. Lemma 3 below shows that S-M-AltLS improves the estimate $W_{t}$ to the desired accuracy, so that we may move on to the next iteration of SoftDeflate.

Assume that the conclusion (b) of Lemma 2 holds, as well as the inductive hypotheses (H1), (H2), (H3), and (H4) . There is a constant $C$ so that the following holds. Let $\gamma^{*}$ be as in (10). Suppose that

Then after $L_{t}$ steps of S-M-AltLS with the initial matrix $W_{t}$, and parameters $\mu_{t},\varepsilon$, the following hold with probability at least $1-1/n^{2}$, over the choice of $\Omega_{t}^{\prime}$.

The inductive hypothesis (H1) holds for the next round:

The inductive hypothesis (H2) holds for the next round:

The inductive hypothesis (H3) holds for the next round: $\mu(X_{t})\leq\mu_{t}.$

The proof of Lemma 3 is addressed in Section 6.

Theorem 1 now follows using 2 and 3. First, we choose $\mu_{0}$ as in the statement of Theorem 1. Because $\mu_{t}\geq\mu_{0}$ for all $t=1,\ldots,T$, this implies that $\mu_{t}$ satisfies the requirements of Lemma 3. Then, the hypotheses of Lemma 3 are implied by the conclusions of the favorable case of Lemma 2. Now, a union bound over at most $k$ epochs of SoftDeflate ensures that with probability at least $1-\nicefrac{{2k}}{{n^{2}}}\geq 1-1/n$, the conclusions of both lemmas hold every round that their hypotheses hold.

If SoftDeflate terminates with the guarantees (a) of Lemma 2, then $\left\|A-X_{T}Y^{T}_{T}\right\|\leq C\varepsilon\left\|M\right\|.$ On the other hand, if (b) holds, then Lemma 2 implies (H4) and the hypotheses of Lemma 3, and then Lemma 3 implies that with probability $1-1/n^{2}$, the remaining inductive hypotheses (H1), (H2), and (H3) for the next round.

Thus, if the situation (a) above never occurs, then the hypotheses of Lemma 3 hold until SoftDeflate terminates because $r_{t}=k$. In this case, Lemma 3 implies that

Finally, we tally up the number of samples. The base case (11) required

For Lemma 3, we required, for a sufficiently large constant $C$,

From the definition of $\mu_{t}$ (in (H3)), we may bound $\mu_{t}$ for all $t\leq k$ by

for some constant $C$. Summing over $t$ gives the result.

Proof of Lemma 2

In this section, we prove Lemma 2, which shows that either Algorithm 1 hits the precision parameter $\varepsilon$ and returns, or else produces an estimate $W_{t}$ for $U^{(\leq t)}$ that is close enough to run S-M-AltLS on. There are several rounds of approximations between the beginning of iteration $t$ and the output $W_{t}$. For the reader’s convenience, we include an informal synopsis of the notation in Figure 1.

We will first argue that the matrix $N_{t-1}$ is close to the truncated, subsampled, noisy estimate $T_{t}$.

Let $T_{t}$ be as in Algorithm 1, and choose any constant $C_{1}>0$. Suppose that the inductive hypotheses (H2) and (H4) hold. Suppose that $p_{t}$ is as in the statement of Lemma 2. Then, for a sufficiently large choice of $C_{0}$ in the hypothesis (H2) (depending only on $C_{1}$), with probability at least $1-1/n^{2}$,

Let $\mathcal{T}$ denote the Truncate operator. As in Algorithm 1, consider

where as in Line 1, $\tau_{t}=\frac{\mu^{*}}{np_{t}}\left(2ks_{t-1}+\Delta\right).$ Above, use used that the sampling operation $P_{\Omega_{t}}$ and the truncate operator $\mathcal{T}$ commute after adjusting for the normalization factor $p_{t}^{-1}$ in the definition of $P_{\Omega_{t}}$. Because $N_{t-1}$ is incoherent, each of its entries is small. More precisely, by the incoherence implication (35) along with the guarantee (H4) on $s_{t-1}$, we have

Thus, each entry of $\widetilde{N_{t-1}}=N_{t-1}+E_{t-1}$ is the sum of something smaller than $p_{t}\tau_{t}$ from $N_{t-1}$, and an error term from $E_{t-1}$, and so truncating entrywise to $p_{t}\tau_{t}$ can only remove mass from the contribution of $E_{t-1}$. This implies that for all $i,j$,

by (H4). Thus, our choice of $p_{t}$ implies that

The choice of $\varepsilon_{0}$ and a sufficient choice of $C_{0}$ (depending only on $C_{1}$) completes the proof. ∎

Suppose for the rest of the proof that the conclusion of Lemma 4 holds. The first thing SoftDeflate does after computing $T_{t}$ is to obtain estimates $\widetilde{U}_{t}$ and $\widetilde{\sigma}_{1},\ldots,\widetilde{\sigma}_{k-r_{t}}$ for the top singular values and vectors of $T_{t}$. These estimates are recovered by SubsIt in Line 1 of Algorithm 1. We first wish to show that the estimated singular values are close to the actual singular values of $T_{t}$. For this, we will invoke Theorem 16 in the appendix, which implies that as long as the number of iterations $L$ of SubsIt satisfies

Above, we took a union bound over all $j$. Again, we condition on this event occuring. Thus, with our choice of $L$, the estimates $\widetilde{\sigma}_{j}$ are indeed close to the singular values $\sigma_{j}(T_{t})$, which by Lemma 4 are with high probability close to the singular values $\sigma_{r_{t-1}+j}$ of $N_{t-1}$ itself.

Before we consider the next step (to $\widetilde{Q}_{t}$) in Figure 1, consider the case when Algorithm 1 returns at line 1. Then $\widetilde{\sigma}_{1}\leq 10\varepsilon s_{0}\leq 20\varepsilon\sigma_{1}$, and so using (17) above we find that $\left\|T_{t}\right\|\leq 21\varepsilon\sigma_{1}.$ Then by Lemma 4,

Thus, for sufficiently large $C_{1}$, we conclude $\sigma_{r_{t-1}+1}\leq 22\varepsilon\sigma_{1}$. In this case, we are done:

and case (a) of the conclusion holds, as long as Lemma 4 does.

On the other hand, suppose that Algorithm 1 does not return at line 1 (and continue to assume that Lemma 4 holds). As above, (17) implies that since $\widetilde{\sigma}_{1}\geq 10\varepsilon$, we must have

This establishes the conclusion (12). With (18), Lemma 4 and (17) together imply that

Above, we use Lemma 13 in the appendix in the first inequality.

We now show that the choice of $d_{t}$ in Line 1 of Algorithm 1 accurately identifies a “gap" in the spectrum.

Suppose that the hypotheses and conclusions of Lemma 4 hold, and in particular that (19) holds. Then the value $r_{t}=r_{t-1}+d_{t}$ obtained in Line 1 of Algorithm 1 satisfies:

Let $d_{t}^{*}$ be the “correct" choice of $d_{t}$; that is, $d_{t}^{*}$ be the smallest positive integer $d\leq k-r_{t-1}$ so that

or let $d_{t}^{*}=d-r_{t-1}$ if such an index does not exist. Write $r_{t}^{*}=r_{t-1}+d_{t}^{*}$. By definition, because $d_{t}^{*}$ is the smallest such $d$ (or smaller than any such $d$ in the case that $r_{t}^{*}=k$), we have

Suppose that, for some $j\leq d_{t}^{*}$, we have

assuming $C_{1}$ is sufficiently large. In Algorithm 1, we choose $d_{t}$ in Line 1 so that there is no $j<d_{t}$ with

Thus, if there were a big gap, the algorithm would have found it: more precisely, using the definition of $\gamma$, we have

This establishes the first conclusion of the lemma. Now, a similar analysis as above shows that if for any $j\leq d_{t}^{*}$ we have

assuming $C_{1}$ is sufficiently large. That is, our algorithm will always find a small gap, if it exists. In particular, if $r_{t}^{*}<k$, we have

and hence $d_{t}\leq d_{t}^{*}$. On the other hand, if $r_{t}^{*}=k$, then we must have $d_{t}=d_{t}^{*}$. In either case, $d_{t}\leq d_{t}^{*},$ and so

Now, we are in a position to verify the inductive hypothesis (H4) for the next round, in the favorable case that Lemma 4 holds. By definition, we have $s_{t}=\widetilde{\sigma}_{d_{t}}$, and (19), followed by Lemma 5 implies that

Suppose that the conclusions of Lemma 4 and Lemma 5 hold, and that (18) holds. Then

We will use a $\sin\theta$ theorem (Theorem 14, due to Wedin, in the appendix) to control the perturbation of the subspaces. Theorem 14 implies

Now, we show that $Q_{t}$ is close to $\widetilde{Q}_{t}$.

Suppose that the conclusions of Lemma 4 and Lemma 5 hold, and that (18) holds. Then with probability $1-1/n^{2}$,

By (17), Lemma 4, and Lemma 5, a similar computation as in the proof of Lemma 6 shows that

Together, Lemmas 6 and 7 imply that, when Lemma 4 and the favorable case for SubsIt hold,

and using the fact that $U^{(t)}$ and $\widetilde{Q}_{t}$ have rank at most $k$, we have that

As in Algorithm 1, let $B$ be a random orthogonal matrix, and let $\overline{Q}_{t}$ be the truncation

The reason for the random rotation is that while $U^{(t)}O$ is reasonably incoherent (because $U^{(t)}$ is), $U^{(t)}OB$ is, with high probability, even more incoherent. More precisely, as in [Har13a], we have

where we define $Z_{t}:=(I-X_{t-1}X_{t-1}^{T})\overline{Q}_{t}$ to be the projection of $\overline{Q}_{t}$ onto $\mathcal{R}(X_{t-1})_{\perp}$. Because $\overline{Q}_{t}$ is close to $U^{(t)}OB$, and $X_{t-1}$ is close to $U^{(<t)}$, $Z_{t}$ is close to $U^{(t)}OB$. More precisely,

Further, the Gram-Schmidt process gives a decomposition

where the triangular matrix $R$ has the same spectrum as $Z_{t}$. In particular,

where above we used that $(U^{(\leq t)}_{\perp})^{T}U^{(t)}=0$. Next,

where we have used the definition of $\overline{Q}_{t}$, the incoherence of $X_{t-1}$, and the computations above in the final line. Thus,

for some constant $C_{5}$. Thus, when the conclusions of Lemma 4 hold, $P_{t}$ is both close to $U^{(t)}$ and incoherent. By induction, the same is true for $W_{t}$. Indeed, if $t=1$, then $P_{t}=W_{t}$, and we are done. If $t\geq 2$, then we have

Then, the inductive hypothesis (H1) and our conclusion (25) imply that

for suitably large $C_{0},C_{1}$. Finally, (26), along with the inductive hypothesis (H3) implies that

We remark that this last computation is the only reason we need $\sin\theta(P_{t},U^{(t)})\lesssim 1/k$, rather than bounded by $1/4$; eventually, we will iterate and have

and we need that $\left(1+\frac{C_{5}}{k}\right)^{T}\leq e^{C_{5}}$ is bounded by a constant (rather than exponential in $T$).

Finally, we have shown that with probability $1-1/n^{2}$ (that is, in the case that Lemma 4 holds and SubsIt works), all of the conclusions of Lemma 2 hold as well. This completes the proof of Lemma 2.

Proof of Lemma 3

In the proof of Lemma 3 we will need an explicit description of the subroutine SmoothQR that we include in Algorithm 3.

We will maintain the following inductive hypothesis:

Above, the tangent of the principal angle obeys

To establish the base case of (J1) for $j=t$, we have

by conclusion (14) of Lemma 2, and hence by (27),

If $t=1$, then $W_{t}=R_{0}$, and we are done with the base case for (J1); if $t\geq 2$, then for $j\leq t-1$, we have

Thus, for $j\leq t-1$, (J1) is implied by (27) again, along with the fact that

which is the (outer) inductive hypothesis (H1), followed by the conclusions (12) and (13) from Lemma 2. This establishes the base case for (J1). The base case for (J2) follows from the conclusion (14) of Lemma 2 directly.

The proof of Lemma 9 is similar to the analysis in [Har13a]. For completeness, we include the proof in Appendix C. Using the inductive hypothesis (J1), and the fact that $\left\|M^{(j)}\right\|_{F}\leq\sqrt{k}\sigma_{r_{j}}$ ,

for a constant $C_{3}$ to be chosen sufficiently large. Observe that with this choice of $\delta$, the requirement on $p_{t}^{\prime}$ in Lemma 9 is implied by the requirement on $p_{t}^{\prime}$ in the statement in Lemma 3. Then the choice of $\delta$ implies

Then, for every $\zeta\leq\tau\nu$ satisfying $\log(n/\zeta)\leq n$, we have with probability at least $1-1/n^{4}$ that the algorithm SmoothQR $(AR+G,\zeta,\mu_{t})$ terminates in $\log(n/\zeta)$ iterations, and the output $R^{\prime}$ satisfies $\mu(R^{\prime})\leq\mu_{t}$. Further, the final noise matrix $H$ added by SmoothQR satisfies $\left\|H\right\|\leq\tau\nu$.

Suppose that $k=o(n/\log(n))$. There is a constant $C_{2}$ so that the following holds. Suppose that

by the inductive hypothesis (J1) for $j=t$.

Next, we compute the parameters that show up in Lemma 10. From Lemma 9, we have

where we have used the inductive hypothesis (J1) in the final line. Then, the requirement of Lemma 10 on $\mu_{t}$ reads

We may simplify and bound the requirement on $\mu_{t}$ as

for some constant $C_{2}$, which was the requirement in the statement of the lemma. Thus, as long as the hypotheses of the current lemma hold, Lemma 10 implies that with probability at least $1-1/n^{4}$,

using (30) in the final inequality. Now, we wish to apply Theorem 8. The hypothesis (J1), along with the conclusion (12) from Lemma 2, immediately implies that

for all $j\leq t,$ and so in particular the first requirement of Theorem 8 is satisfied. To satisfy the second requirement of Theorem 8, we must show that

provided $C_{3}$ is suitably large. A union bound over all $j$ establishes (J1) for the next iteration of S-M-AltLS. After another union bound over

steps of S-M-AltLS, for some constant $C$ depending on $C_{0}$, we conclude that with probability at least $1-1/n^{2}$, for all $j$,

To establish the second conclusion, we note that we have already conditioned on the event that (30) holds, and so we have

using (13) in the final inequality. Finally, the third conclusion, that (H3) holds, follows from the definition of SmoothQR.

Simulations

In this section, we compare the performance of SoftDeflate to that of other fast algorithms for matrix completion. In particular, we investigate the performance of SoftDeflate compared to the Frank-Wolfe (FW) algorithm analyzed in [JS10], and also compared to the naive algorithm which simply takes the SVD of the subsampled matrix $A_{\Omega}$. All of the code that generated the results in this section can be found online at \urlhttp://sites.google.com/site/marywootters.

We implemented SoftDeflate, as described in Algorithm 1, fixing $30,000$ observations per iteration; to increase the number of measurements, we increased the parameters $L_{t}$ (which were the same for all $t$). For simplicity, we used a version of S-M-AltLS which did not implement the smoothing in SmoothQR or the median. We implemented the Frank-Wolfe algorithm as per the pseudocode in Algorithm 5, with accuracy parameter $\varepsilon=0.05$. We remark that decreasing the accuracy parameter did improve the performance of the algorithm (at the cost of increasing the running time), but did not change its qualitative dependence on $m$, the number of observations. We implemented SVD via subspace iteration, as in Algorithm 8, with $L=100$.

The error was measured in two ways: the Frobenius error $\left\|A-XY^{T}\right\|_{F}$, and the error between the recovered subspaces, $\sin\Theta(U,X)$. The results are shown in Figure 2.

The experiments show that SoftDeflate significantly outperforms the other “fast" algorithms in both metrics. In particular, of the three algorithms, SoftDeflate is the only one which converges enough to reliably capture the singular vector associated with the $0.1$ eigenvalue; none of the algorithms converge enough to find the $0.01$ eigenvalue with the number of measurements allowed. The other two algorithms show basically no progress for these small values of $m$. To illustrate what happens when FW and SVD do converge, we repeated the same experiment for $n=1000$ and $k=2$; for this smaller value of $n$, we can let the number of measurements to get quite large compared to $n^{2}$. We find that even though FW and SVD do begin to converge eventually, they are still outperformed by SoftDeflate. The results of these smaller tests are shown in Figure 3.

2 Further comments on the Frank-Wolfe algorithm

As algorithms like Frank-Wolfe are often cited as viable fast algorithms for the Matrix Completion problem, the reader may be surprised by the performance of FW depicted in Figures 2 and 3. There are two reasons for this. The first reason, noted in Section 2.2, is that while FW is guaranteed to converge on the sampled entries, it may not converge so well on the actual matrix; the errors plotted above are with respect to the entire matrix. To illustrate this point, we include in Figure 4 the results of an experiment showing the convergence of Frank-Wolfe (Algorithm 5), both on the samples and off the samples. As above, we considered random $10,000\times 10,000$ matrices with a pre-specified spectrum. We fixed the number of observations at $5\times 10^{6}$, and ran the Frank-Wolfe algorithm for 40 iterations, plotting its progress both on the observed entries and on the entire matrix. While the error on the observed entries does converge as predicted, the matrix itself does not converge so quickly.

Acknowledgements

We thank the Simons Institute for Theoretical Computer Science at Berkeley, where part of this work was done.

References

Appendix A Dividing up ΩΩ\Omega

In this section, we show how to take a set $\Omega\subset[n]\times[n]$, so that each index $(i,j)$ is included in $\Omega$ with probability $p$, and return subsets $\Omega_{1},\ldots,\Omega_{L}$ which follow a distribution more convenient for our analysis. Algorithm 6 has the details. Observe that the first thing that Algorithm 6 does is throw away samples from $\Omega$. Thus, while this step is convenient for the analysis, and we include it for theoretical completeness, in practice it may be uneccessary—especially if the assumption on the distribution of $\Omega$ is an approximation to begin with.

The correctness of Algorithm 6 follows from the following lemma, about the properties of Algorithm 7.

Next, we observe that for any fixed $S$, the events $\left\{E(u,S)\right\}_{u\in\mathcal{U}}$ are independent under the distribution induced by Algorithm 7. This follows from the fact that in all of the random steps (including the generation of $\Omega$ and within Algorithm 7), the $u\in\mathcal{U}$ are treated independently. Notice that these events are also independent under $\mathcal{D}$ by definition.

Now, for any instantiation $\vec{\Omega^{\prime}}=(\Omega_{1}^{\prime},\ldots,\Omega_{L}^{\prime})$ of the random variables $(\Omega_{1},\ldots,\Omega_{L})$, consider the event

Thus the probability of any outcome $\vec{\Omega}^{\prime}$ is the same under $\mathcal{D}$ and under Algorithm 7, and this completes the proof of the lemma. ∎

Appendix B Useful statements

In this appendix, we collect a few useful statements upon which we rely.

First, we record some consequences of the bound (4) on the coherence of $A$. We always have

It will also be useful to notice that since $\left\|e_{i}^{T}U^{(>t)}\right\|_{2}\leq\left\|e_{i}^{T}U\right\|_{2},$ (4) implies that for all $t$,

B.2 Perturbation statements

Next, we will use the following lemma about perturbations of singular values, due to Weyl.

In order to compare the singular vectors of a matrix $A$ with those of a perturbed version $\widetilde{A}$, we will find the following theorem helpful. We recall that for subspaces $U,V$, $\sin\theta(U,V)$ refers to the sine of the principal angle between $U$ and $V$. (See [SS90] for more on principal angles).

Suppose that $A$ has the singular value decomposition

and let $\widetilde{A}=A+E$ be a perturbed matrix with SVD

Suppose there are numbers $\alpha,\delta>0$ so that $\sigma_{\min}(\widetilde{\Sigma_{1}})\geq\alpha+\delta$ and $\sigma_{\max}(\Sigma_{2})\leq\alpha.$ Then,

We will also use the fact that if the angle between (the subspaces spanned by) two matrices is small, then there is some unitary transformation so that the two matrices are close.

We have $V=\Pi_{U}V+\Pi_{U_{\perp}}V=U(U^{T}V)+\Pi_{U_{\perp}}V.$ Since $\sin\theta(U,V)\leq\varepsilon$, we have $\left\|\Pi_{U_{\perp}}V\right\|\leq\varepsilon,$ and $\sigma_{k}(U^{T}V)=\cos\theta(U,V)\geq\sqrt{1-\varepsilon^{2}}.$ Thus, we can write $U^{T}V=Q+E,$ where $\left\|E\right\|\leq 1-\sqrt{1-\varepsilon^{2}}.$ The claim follows from the triangle inequality. ∎

B.3 Subspace Iteration

Our algorithm uses the following standard version of the well-known Subspace Iteration algorithm—also known as Power Method.

We have the following theorem about the convergence of SubsIt.

Here, $c^{\prime}$ can be made any constant by increasing $c$ and $C$ is an absolute constant. Fix $i$ and let $x_{i}=(R_{L})_{i}$ denote the $i-th$ column of $R_{L}$. Suppose that $i\in\left\{r_{j}+1,\ldots,r_{j+1}\right\}$. Then, the estimates $\widetilde{\sigma_{i}}$ satisfy

By definition, as there are no significant gaps between $\sigma_{r_{j}+1}$ and $\sigma_{r_{j}}$, we have

and so this completes the proof after collecting terms. ∎

B.4 Matrix concentration inequalities

We will repeatedly use the Matrix Bernstein and Matrix Chernoff inequalities; we use the versions from [Tro12]:

[Matrix Bernstein [Tro12]] Consider a finite sequence $\left\{Z_{k}\right\}$ of independent, random, $d\times d$ matrices. Assume that each matrix satisfies

One corollary of Lemma 17 is the following lemma about the concentration of the matrix $P_{\Omega}(A)$.

Let $\xi_{ij}$ be independent Bernoulli-$p$ random variables, which are $1$ if $(i,j)\in\Omega$ and otherwise.

which is a sum of independent random matrices. Using the Matrix Bernstein inequality, Lemma 17, we conclude that

almost surely. This concludes the proof. ∎

Finally, we will use the Matrix Chernoff inequality.

[Matrix Chernoff [Tro12]] Consider a finite sequence $\left\{X_{k}\right\}$ of independent, self-adjoint, $d\times d$ matrices. Assume that each $X_{k}$ satisfies

B.5 Medians of vectors

Suppose that $v^{(s)}$, for $s=1,\ldots,T$ are i.i.d. random vectors, so that for all $s$,

Let $S\subset[T]$ be the set of $s$ so that $\left\|v^{(s)}\right\|_{2}^{2}\leq\alpha.$ By a Chernoff bound,

Suppose that the likely event occurs, so $|S|>3T/4$. For $j\in[k]$, let

Because $|S|>3T/4$, we have $|S_{j}|\geq T/4$. Then

Appendix C Proof of Lemma 9

First, we observe that with very high probability, $B_{i}^{(s)}$ is close to the identity.

There is a constant $C$ so that the following holds. Suppose that $p^{\prime}\geq Ck\mu_{t}\log(n)/(n\delta^{2})$. Then

where $\xi_{r}$ is $1$ with probability $p^{\prime}$ and otherwise. We apply the Matrix Chernoff bound (Lemma 19); we have

The claim follows from the choice of $p^{\prime}$. ∎

There is a constant $C$ so that the following holds. Suppose that $p^{\prime}\geq\frac{C\mu_{t}k}{n\delta^{2}}$. Then for each $s$,

Along with Markov’s inequality, this completes the proof. ∎

There is a constant $C$ so that the following holds. Suppose that $p^{\prime}\geq Ck\log(n)\mu_{t}/(\delta^{2}n)$ for a constant $C$. Then for each $s\leq T$,

We have already bounded $\left\|(B_{i}^{(s)})^{-1}\right\|$ with high probability in Claim 22, when the bound on $p^{\prime}$ holds, and so we now bound $\left\|y_{1}\right\|_{2}$ and $\left\|y_{2}\right\|_{2}$ with decent probability. As we did in Claim 23, we compute the expectation of $\left\|y_{1}\right\|_{2}^{2}$ and use Markov’s inequality.

Next, we turn our attention to the second term $\left\|y_{2}\right\|_{2}$. We have

By Claim 22, we established that with probability $1-1/n^{5}$, $\left\|I-(B_{i}^{(s)})\right\|\leq\frac{\delta}{2}$, with our choice of $p^{\prime}$. Thus, with probability at least $1-1/n^{5}$,

Altogether, we conclude that with probability at least $1-1/20-2/n^{5}$, we have

as long as $\delta\leq 1/2$. This proves the claim. ∎

Putting Claims 22, 23 and 24 together, along with the choice of $p_{t}^{\prime}=L_{t}s_{\max}p^{\prime}$, we conclude that, for each $s\in[T]$ and for any $\delta<1/2$,

is small with exponentially large probability. Indeed, by Lemma 20,

for some constant $c$. By the choice of $s_{\max}$, the failure probability is at most $1/n^{6}$, and a union bound over all $i$ shows that, with probability at least $1-1/n^{5}$,

This was the second claim in Lemma 9. Now, we show that in the favorable case that (38) holds, so does the first claim of Lemma 9, and this will complete the proof of the lemma. Suppose that (38) holds. Then

Notice that, for any real numbers $(a_{i,j})$, $i\in[n],j\in[t]$, and for any real number $b_{j}$, $j\in[t]$, we have

Thus, we may bound the second term above by

Altogether, we conclude that, in the favorable case the (38) holds,

as desired. This completes the proof of Lemma 9.