Asymmetry Helps: Eigenvalue and Eigenvector Analyses of Asymmetrically Perturbed Low-Rank Matrices

Introduction

with $\bm{H}$ denoting a noise matrix. A classical problem is concerned with estimating the leading eigenvalues and eigenspace of $\bm{M}^{\star}$ given observation $\bm{M}$.

The current paper concentrates on a scenario where the noise matrix $\bm{H}$ (and hence $\bm{M}$) consists of independently generated random entries and is hence asymmetric in general. This might arise, for example, when we have available multiple (e.g. two) samples for each entry of $\bm{M}^{\star}$ and arrange the samples in an asymmetric fashion. A natural approach that immediately comes to mind is based on singular value decomposition (SVD), which employs the leading singular values (resp. subspace) of $\bm{M}$ to approximately estimate the eigenvalues (resp. eigenspace) of $\bm{M}^{\star}$. By contrast, a much less popular alternative is based on eigen-decomposition of the asymmetric data matrix $\bm{M}$, which attempts approximation using the leading eigenvalues and eigenspace of $\bm{M}$. Given that eigen-decomposition of an asymmetric matrix is in general not as numerically stable as SVD, conventional wisdom often favors the SVD-based approach, unless certain symmetrization step is implemented prior to eigen-decomposition.

When comparing these two approaches numerically, however, a curious phenomenon arises, which largely motivates the research in this paper. Let us generate $\bm{M}^{\star}$ as a random rank-1 matrix with leading eigenvalue $\lambda^{\star}=1$, and let $\bm{H}$ be a Gaussian random matrix whose entries are i.i.d. $\mathcal{N}(0,\sigma^{2})$ with $\sigma=1/\sqrt{n\log n}$. Fig. 1(a) compares the empirical accuracy of estimating the 1st eigenvalue of $\bm{M}^{\star}$ via the leading eigenvalue (the blue line) and via the leading singular value of $\bm{M}$ (the red line). As it turns out, eigen-decomposition significantly outperforms vanilla SVD in estimating $\lambda^{\star}$, and the advantage seems increasingly more remarkable as the dimensionality $n$ grows. To facilitate comparison, we include an additional green line in Fig. 1(a), obtained by rescaling the red line by $2.5/\sqrt{n}$. Interestingly, this green line coincides almost perfectly with the blue line, thus suggesting orderwise gain of eigen-decomposition compared to SVD. What is more, this phenomenon does not merely happen under i.i.d. noise. Similar numerical behaviors are observed in the problem of matrix completion — as displayed in Fig. 1(b) — even though the components of the equivalent perturbation matrix are apparently far from identically distributed or homoscedastic.

The goal of the current paper is thus to develop a systematic understanding of this phenomenon, that is, why statistical asymmetry empowers eigen-decomposition and how to exploit this feature in statistical estimation. Informally, our findings suggest that: when $\bm{M}^{\star}$ is rank-1 and $\bm{H}$ is composed of zero-mean and independent (but not necessarily identically distributed or homoscedastic) entries,

the leading eigenvalue of $\bm{M}$ could be $O(\sqrt{n})$ times (up to some logarithmic factor) more accurate than the (unadjusted) leading singular value of $\bm{M}$ when estimating the 1st eigenvalue of $\bm{M}^{\star}$;More precisely, this gain is possible when $\|\bm{H}\|$ is nearly as large as $\|\bm{M}^{\star}\|$ (up to some logarithmic factor).

the perturbation of the leading eigenvector is well-controlled along an arbitrary deterministic direction; for example, the eigenvector perturbation is well-controlled in any coordinate, indicating that the eigenvector estimation error is spread out across all coordinates.

We will further provide partial theory to accommodate the rank-$r$ case. As an important application, such a theory allows us to estimate the leading singular value and singular vectors of an asymmetric rank-1 matrix via eigen-decomposition of a certain dilation matrix, which also outperforms the vanilla SVD approach.

We would like to immediately remark that: for some scenarios (e.g. the case with i.i.d. Gaussian noise), it is possible to adjust the leading singular value of $\bm{M}$ to obtain the same accuracy as the leading eigenvalue of $\bm{M}$. As it turns out, the advantages of the eigen-decomposition approach may become more evident in the presence of heteroscedasticity — the case where the noise has location-varying and unknown variance. We shall elaborate on this point in Section 4.1.2.

All in all, when it comes to low-rank matrix estimation, arranging the observed matrix samples in an asymmetric manner and invoking eigen-decomposition properly could sometimes be statistically beneficial.

Problem formulation

where $\bm{H}=[H_{ij}]_{1\leq i,j\leq n}$ is a random noise matrix.

The current paper concentrates on independent — but not necessarily identically distributed or homoscedastic — noise. Specifically, we impose the following assumptions on $\bm{H}$ throughout this paper.

(Independent entries) The entries $\{H_{ij}\}_{1\leq i,j\leq n}$ are independently generated;

(Magnitude) Each $H_{ij}$ ($1\leq i,j\leq n$) satisfies either of the following conditions:

In what follows, the dependency of $\sigma_{n}$ and $B_{n}$ on $n$ shall often be suppressed whenever it is clear from the context, so as to simplify notation.

Note that we do not enforce the constraint $H_{ij}=H_{ji}$, and hence $\bm{H}$ and $\bm{M}$ are in general asymmetric matrices. Also, Condition 3 does not require the $H_{ij}$’s to have equal variance across different locations; in fact, they can be heteroscedastic. In addition, while Condition 4(a) covers the class of bounded random variables, Condition 4(b) allows us to accommodate a large family of heavy-tailed distributions (e.g. sub-exponential distributions). An immediate consequence of Assumption 1 is the following bound on the spectral norm $\|\bm{H}\|$ of $\bm{H}$.

Under Assumption 1, there exist some universal constants $c_{0},C_{0}>0$ such that with probability exceeding $1-C_{0}n^{-10}$,

This is a standard non-asymptotic result that follows immediately from the matrix Bernstein inequality Tropp (2015) and the union bound (for Assumption 4(b)). We omit the details for conciseness. ∎

2 Our goal

The aim is to develop non-asymptotic eigenvalue and eigenvector perturbation bounds under this family of random and asymmetric noise matrices. Our theoretical development is divided into two parts. Below we introduce our goal as well as some notation used throughout.

For the rank-1 case, we assume the eigen-decomposition of $\bm{M}^{\star}$ to be

with $\lambda^{\star}$ and $\bm{u}^{\star}$ being its leading eigenvalue and eigenvector, respectively. We also denote by $\lambda$ and $\bm{u}$ the leading eigenvalue and eigenvector of $\bm{M}$, respectively. The following quantities are the focal points of this paper (see Section 4):

Eigenvalue perturbation: $|\lambda-\lambda^{\star}|$;

Entrywise eigenvector perturbation: $\min\{\|\bm{u}-\bm{u}^{\star}\|_{\infty},\|\bm{u}+\bm{u}^{\star}\|_{\infty}\}$.

For the general rank-$r$ case, we let the eigen-decomposition of $\bm{M}^{\star}$ be

As is well-known, eigen-decomposition can be applied to estimate the singular values and singular vectors of an asymmetric matrix $\bm{M}^{\star}$ via the standard dilation trick Tropp (2015). As a consequence, our results are also applicable for singular value and singular vector estimation. See Section 5.2 for details.

3 Incoherence conditions

Finally, we single out an incoherence parameter that plays an important role in our theory, which captures how well the energy of the eigenvectors is spread out across all entries.

The incoherence parameter of a rank-$r$ symmetric matrix $\bm{M}^{\star}$ with eigen-decomposition $\bm{M}^{\star}=\bm{U}^{\star}\bm{\Sigma}^{\star}\bm{U}^{\star\top}$ is defined to be the smallest quantity $\mu$ obeying

An alternative definition of the incoherence parameter (Candès and Recht, 2009; Keshavan et al., 2010; Chi et al., 2019; Chen et al., 2019a) is the smallest quantity $\mu_{0}$ satisfying $\left\|\bm{U}^{\star}\right\|_{2,\infty}\leq\sqrt{\mu_{0}r/n}$. This is a weaker assumption than Definition 1, as it only requires the energy of $\bm{U}^{\star}$ to be spread out across all of its rows rather than all of its entries. Note that these two incoherent parameters are consistent in the rank-$1$ case; in the rank-$r$ case one has $\mu_{0}\leq\mu\leq\mu_{0}r$.

4 Notation

Preliminaries

Before continuing, we gather several preliminary facts that will be useful throughout. The readers familiar with matrix perturbation theory may proceed directly to the main theoretical development in Section 4.

We begin with a standard result concerning eigenvalue perturbation of a diagonalizable matrix Bauer and Fike (1960). Note that the matrices under study might be asymmetric.

In addition, if $\bm{A}$ is symmetric, then there exists an eigenvalue $\lambda$ of $\bm{A}$ such that

However, caution needs to be exercised as the Bauer-Fike Theorem does not specify which eigenvalue of $\bm{A}$ is close to an eigenvalue of $\bm{A}+\bm{H}$. Encouragingly, in the low-rank case of interest, the Bauer-Fike Theorem together with certain continuity of the spectrum allows one to localize the leading eigenvalues of the perturbed matrix.

Suppose $\bm{M}^{\star}$ is a rank-$r$ symmetric matrix whose top-$r$ eigenvalues obey \big{|}\lambda_{1}^{\star}\big{|}\geq\cdots\geq\big{|}\lambda_{r}^{\star}\big{|}>0. If \|\bm{H}\|<\big{|}\lambda^{\star}_{r}\big{|}/2, then the top-$r$ eigenvalues $\lambda_{1},\cdots,\lambda_{r}$ of $\bm{M}=\bm{M}^{\star}+\bm{H}$, sorted by modulus, obey that: for any $1\leq l\leq r$,

In addition, if $r=1$, then both the leading eigenvalue and the leading eigenvector of $\bm{M}$ are real-valued.

This result, which we establish in Appendix A.1, parallels Weyl’s inequality for symmetric matrices. Note, however, that the above bound (9) might be quite loose for specific settings. We will establish much sharper perturbation bounds when $\bm{H}$ contains independent random entries (see, e.g. Corollary 1).

2 The Neumann trick and eigenvector perturbation

Next, we introduce a classical result dubbed as the “Neumann trick” Eldridge et al. (2018). This theorem, which is derived based on the Neumann series for a matrix inverse, has been applied to analyze eigenvectors in various settings Erdős et al. (2013); Jain and Netrapalli (2015); Eldridge et al. (2018).

Consider the matrices $\bm{M}^{\star}$ and $\bm{M}$ (see (5) and (2)). Suppose $\left\|\bm{H}\right\|<|\lambda_{l}|$ for some $1\leq l\leq n$. Then

We supply the proof in Appendix A.2 for self-containedness. ∎

In particular, if $\bm{M}^{\star}$ is a rank-1 matrix and $\|\bm{H}\|<|\lambda_{1}|$, then

An immediate consequence of the Neumann trick is the following lemma, which asserts that each of the top-$r$ eigenvectors of $\bm{M}$ resides almost within the top-$r$ eigen-subspace of $\bm{M}^{\star}$, provided that $\|\bm{H}\|$ is sufficiently small. The proof is deferred to Appendix A.3.

In addition, if $r=1$, then one further has

Perturbation analysis for the rank-1 case

This section presents perturbation analysis results when the truth $\bm{M}^{\star}$ is a symmetric rank-1 matrix. We shall start by presenting a master bound which, as we will see, immediately leads to our main findings.

Our master bound is concerned with the perturbation of linear forms of eigenvectors, as stated below.

We would like to remark on the range of the noise size covered by our theory. If the incoherence parameter of the truth $\bm{M}^{\star}$ obeys $\mu\asymp 1$, then even the magnitude of the largest entry of $\bm{M}^{\star}$ cannot exceed the order of $|\lambda^{\star}|/n$. One can thus interpret the condition (14) in this case as

In other words, the standard deviation $\sigma$ of each noise component is allowed to be substantially larger (i.e. $\sqrt{n/\log n}$ times larger) than the magnitude of any of the true entries. In fact, this condition (14) matches, up to some log factor, the one required for spectral methods to perform noticeably better than random guessing.

In words, Theorem 3 tells us that: the quantity $\frac{\bm{u}^{\star\top}\bm{u}}{\lambda/\lambda^{\star}}\bm{a}^{\top}\bm{u}^{\star}$ serves as a remarkably accurate approximation of the linear form $\bm{a}^{\top}\bm{u}$. In particular, the approximation error is at most $O(1/\sqrt{n})$ under the condition (14) for incoherent matrices. Encouragingly, this approximation accuracy holds true for an arbitrary deterministic direction (reflected by $\bm{a}$). As a consequence, one can roughly interpret Theorem 3 as

where such an approximation is fairly accurate along any fixed direction. Compared with the identity $\bm{u}=\frac{1}{\lambda}\bm{M}\bm{u}=\frac{1}{\lambda}(\bm{M}^{\star}+\bm{H})\bm{u}$, our results imply that $\bm{H}\bm{u}$ is exceedingly small along any fixed direction, even though $\bm{H}$ and $\bm{u}$ are highly dependent. As we shall explain in Section 4.3, this observation usually cannot happen when $\bm{H}$ is a symmetric random matrix or when one uses the leading singular vector instead, due to the significant bias resulting from symmetry.

This master theorem has several interesting implications, as we shall elucidate momentarily.

1.2 Eigenvalue perturbation

To begin with, Theorem 3 immediately yields a much sharper non-asymptotic perturbation bound regarding the leading eigenvalue $\lambda$ of $\bm{M}$.

Under the assumptions of Theorem 3, with probability at least $1-O(n^{-10})$ we have

Without loss of generality, assume that $\lambda^{\star}=1$. Taking $\bm{a}=\bm{u}^{\star}$ in Theorem 3, we get

From Lemma 1 and the condition (14), we know $\|\bm{H}\|<1/4$, which combines with Lemma 2 and Lemma 3 yields $\lambda\asymp|\bm{u}^{\star\top}\bm{u}|\asymp 1$. Substitution into (18) yields

For the vast majority of applications we encounter, the maximum possible noise magnitude $B$ (cf. Assumption 1) obeys $B\lesssim\sigma\sqrt{n/\log n}$, in which case the bound in Corollary 1 simplifies to

This means that the eigenvalue estimation error is not much larger than the variability of each noise component. In addition, we remind the reader that for a fairly broad class of noise (see Remark 4), the leading eigenvalue $\lambda$ of $\bm{M}$ is guaranteed to be real-valued, an observation that has been made in Lemma 2. In practice, however, one might still encounter some scenarios where $\lambda$ is complex-valued. As a result, we recommend the practitioner to use the real part of $\lambda$ as the eigenvalue estimate, which clearly enjoys the same statistical guarantee as in Corollary 1.

In order to facilitate comparison, we denote by $\lambda_{\mathsf{svd}}$ the largest singular value of $\bm{M}$, and look at $|\lambda_{\mathsf{svd}}-\lambda^{\star}|$. Combining Weyl’s inequality, Lemma 1 and the condition (14), we arrive at

When $\mu\asymp 1$, this error bound w.r.t. this (unadjusted) singular value could be $\sqrt{n}$ times larger than the perturbation bound (17) derived for the leading eigenvalue. This corroborates our motivating experiments in Fig. 1.

The reader might naturally wonder what would happen if we symmetrize the data matrix before performing eigen-decomposition. Consider, for example, the i.i.d. Gaussian noise case where $H_{ij}\overset{\text{i.i.d.}}{\sim}\mathcal{N}(0,\sigma^{2})$, and assume $\lambda^{\star}>0$ for simplicity. The leading eigenvalue $\lambda_{\mathsf{sym}}$ of the symmetrized matrix $(\bm{M}+\bm{M}^{\top})/2$ has been extensively studied in the literature Füredi and Komlós (1981); Yin et al. (1988); Péché (2006); Féral and Péché (2007); Benaych-Georges and Nadakuditi (2011); Renfrew and Soshnikov (2013); Knowles and Yin (2013). In particular, it has been shown (e.g. Capitaine et al. (2009)) that, with probability approaching one,

If $\sigma=|\lambda^{\star}|\sqrt{1/(n\log n)}$ (which is the setting in our numerical experiment), then this can be translated into

This implies that $\lambda_{\mathsf{sym}}$ suffers from a substantially larger bias than the leading eigenvalue $\lambda$ obtained without symmetrization, since in this case we have (cf. Corollary 1)

Armed with the approximation (21) in the i.i.d. Gaussian noise case, the careful reader might naturally suggest a properly corrected estimate $\lambda_{\mathsf{sym},\mathsf{c}}$ as follows (again assuming $\lambda>0$)

which is a shrinkage-type estimate chosen to satisfy $\lambda_{\mathsf{sym}}=\lambda_{\mathsf{sym},\mathsf{c}}+\frac{n\sigma^{2}}{2\lambda_{\mathsf{sym},\mathsf{c}}}$. A little algebra reveals that: if $\sigma=1/\sqrt{n\log n}$, then

thus matching the estimation accuracy of $\lambda$ (cf. (22)). In addition, some sort of universality results has been established as well in the literature (Capitaine et al., 2009), implying that the same approximation and correction are applicable to a broad family of zero-mean noise with identical variance. As we shall illustrate numerically in Section 4.2, this approach (i.e. $\lambda_{\mathsf{sym},\mathsf{c}}$) performs almost identically to the one using vanilla eigen-decomposition without symmetrization. In addition, very similar observations have been made for the SVD-based approach (Silverstein, 1994; Yin et al., 1988; Péché, 2006; Féral and Péché, 2007; Benaych-Georges and Nadakuditi, 2012; Bryc and Silverstein, 2018); for the sake of brevity, we do not repeat the arguments here.

We would nevertheless like to single out a few statistical advantages of the eigen-decomposition approach without symmetrization. To begin with, $\lambda$ is obtained via vanilla eigen-decomposition, and computing it does not rely on any kind of noise statistics. This is in stark contrast to the bias correction (23) in the presence of symmetric data, which requires prior knowledge about (or a very precise estimate of) the noise variance $\sigma^{2}$. Leaving out this prior knowledge matter, a more important issue is that the approximation formula (21) assumes identical variance of noise components across all entries (i.e. homoscedasticity). While an approximation of this kind has been found for more general cases beyond homoscedastic noise (e.g. Bryc and Silverstein (2018)), the approximation formula (e.g. (Bryc and Silverstein, 2018, Theorem 1.1)) becomes fairly complicated, requires prior knowledge about all variance parameters, and is thus difficult to implement in practice. In comparison, the vanilla eigen-decomposition approach analyzed in Corollary 1 imposes no restriction on the noise statistics and is fully adaptive to heteroscedastic noise.

To complete the picture, we provide a simple information-theoretic lower bound for the i.i.d. Gaussian noise case, which will be established in Appendix B.

Fix any small constant $\varepsilon>0$. Suppose that $H_{ij}\overset{\text{i.i.d.}}{\sim}\mathcal{N}(0,\sigma^{2})$. Consider three matrices

In short, Lemma 4 asserts that one cannot possibly locate an eigenvalue to within a precision of $\Delta$ much better than $\sigma$, which reveals a fundamental limit that cannot be broken by any algorithm. In comparison, the vanilla eigen-decomposition method based on asymmetric data achieves an accuracy of $|\lambda-\lambda^{\star}|\lesssim\sigma\sqrt{\log n}$ (cf. Corollary 1 and (19)) for the incoherent case, thus matching the information-theoretic lower bound up to some log factor. In fact, the extra $\sqrt{\log n}$ factor arises simply because we are aiming for a high-probability guarantee.

1.3 Perturbation of linear forms of eigenvectors

The master bound in Theorem 3 admits a more convenient form when controlling linear functions of the eigenvectors. The result is this:

Under the same setting of Theorem 3, with probability at least $1-O(n^{-10})$ we have

Without loss of generality, assume that $\bm{u}^{\star\top}\bm{u}\geq 0$ and that $\lambda^{\star}=1$. Then one has

where the last inequality arises from Theorem 3 as well as the definition of $\mu$. In addition, apply Lemma 2 and Lemma 3 to obtain

Putting the above bounds together concludes the proof. ∎

The perturbation of linear forms of eigenvectors (or singular vectors) has not yet been well explored even for the symmetric case. One scenario that has been studied is linear forms of singular vectors under i.i.d. Gaussian noise (Koltchinskii and Xia, 2016; Xia, 2016). Our analysis — which is certainly different from Koltchinskii and Xia (2016) as our emphasis is eigen-decomposition — does not rely on the Gaussianality assumption, and accommodates a much broader class of random noise. Another work that has looked at linear forms of the leading singular vector is Ma et al. (2019) for phase retrieval and blind deconvolution, although the vector $\bm{a}$ therein is specific to the problems (i.e. the design vectors) and cannot be made general.

The perturbation theory for linear forms of eigenvectors has been substantially extended in our follow-up work; the interested reader is referred to Cheng et al. (2020) for details.

1.4 Entrywise eigenvector perturbation

A straightforward consequence of Corollary 2 that is worth emphasizing is sharp entrywise control of the leading eigenvector as follows.

Under the same setting of Theorem 3, with probability at least $1-O(n^{-9})$ we have

Recognizing that $\|\bm{u}-\bm{u}^{\star}\|_{\infty}=\max_{i}|\bm{e}_{i}^{\top}\bm{u}-\bm{e}_{i}^{\top}\bm{u}^{\star}|$ and recalling our assumption $|\bm{e}_{i}^{\top}\bm{u}|\leq\sqrt{{\mu}/{n}}$, we can invoke Corollary 2 and the union bound to establish this entrywise bound. ∎

2 Applications

We apply our main results to two concrete matrix estimation problems and examine the effectiveness of these bounds. As before, $\bm{M}^{\star}$ is a rank-1 matrix with incoherence parameter $\mu$ and leading eigenvalue $\lambda^{\star}$.

Low-rank matrix estimation from Gaussian noise. Suppose that $\bm{H}$ is composed of i.i.d. Gaussian random variables $\mathcal{N}(0,\sigma^{2})$.In this case, one can take $B\asymp\sigma\sqrt{\log n}$, which clearly satisfies $B\log n\ll\sqrt{n\sigma^{2}\log n}$. If $\sigma\lesssim\frac{1}{\sqrt{n\log n}}$, applying Corollaries 1-3 reveals that with high probability,

Low-rank matrix completion. Suppose that $\bm{M}$ is generated using random partial entries of $\bm{M}^{\star}$ as follows

where $p$ denotes the fraction of the entries of $\bm{M}^{\star}$ being revealed. It is straightforward to verify that $\bm{H}=\bm{M}-\bm{M}^{\star}$ is zero-mean and obeys $|H_{ij}|\leq\frac{\mu}{np}:=B$ and $\mathsf{Var}(H_{ij})\leq\frac{\mu^{2}}{pn^{2}}$. Consequently, if $p\gtrsim\frac{\mu^{2}\log n}{n}$, then invoking Corollaries 1-3 yields

Finally, we remark that all the above applications assume the availability of an asymmetric data matrix $\bm{M}$. One might naturally wonder whether there is anything useful we can say if only a symmetric matrix $\bm{M}$ is available. While this is in general difficult, our theory does have direct implications for both matrix completion and the case with i.i.d. Gaussian noise in the presence of symmetric data matrices; that is, it is possible to first asymmetrize the data matrix followed by eigen-decomposition. The interested reader is referred to Appendix J for details.

3 Why asymmetry helps?

We take a moment to develop some intuition underlying Theorem 3, focusing on the case with $\lambda^{\star}=1$ for simplicity. The key ingredient is the Neumann trick stated in Theorem 2. Specifically, in the rank-1 case we can expand

where the last inequality holds since (i) $|\bm{u}^{\star\top}\bm{u}|\leq 1$, and (ii) $\lambda$ is real-valued and obeys $\lambda\approx 1$ if $\|\bm{H}\|\ll 1$ (in view of Lemma 2). As a result, the perturbation can be well-controlled as long as $|\bm{a}^{\top}\bm{H}^{s}\bm{u}^{\star}|$ is small for every $s\geq 1$.

As it turns out, $\bm{a}^{\top}\bm{H}^{s}\bm{u}^{\star}$ might be much better controlled when $\bm{H}$ is random and asymmetric, in comparison to the case where $\bm{H}$ is random and symmetric. To illustrate this point, it is perhaps the easiest to inspect the second-order term.

Asymmetric case: when $\bm{H}$ is composed of independent zero-mean entries each with variance $\sigma_{n}^{2}$, one has

Symmetric case: when $\bm{H}$ is symmetric and its upper trangular part consists of independent zero-mean entries with variance $\sigma_{n}^{2}$, it holds that

In words, the term $\bm{a}^{\top}\bm{H}^{2}\bm{u}^{\star}$ in the symmetric case might have a significantly larger bias compared to the asymmetric case. This bias effect is substantial when $\bm{a}^{\top}\bm{u}^{\star}$ is large (e.g. when $\bm{a}=\bm{u}^{\star}$), which plays a crucial role in determining the size of eigenvalue perturbation.

The vanilla SVD-based approach can be interpreted in a similar manner. Specifically, we recognize that the leading singular value (resp. left singular vector) can be computed via the leading eigenvalue (resp. eigenvector) of the symmetric matrix $\bm{M}\bm{M}^{\top}$. Given that $\bm{M}\bm{M}^{\top}-\bm{M}^{\star}\bm{M}^{\star\top}$ is also symmetric, the aforementioned bias issue arises as well. This explains why vanilla eigen-decomposition might have an advantage over vanilla SVD when dealing with asymmetric matrices.

4 Proof outline of Theorem 3

This subsection outlines the main steps for establishing Theorem 3. To simplify presentation, we shall assume without loss of generality that

Throughout this paper, all the proofs are provided for the case when Conditions 1-3, 4(a) in Assumption 1 are valid. Otherwise, if Condition 4(b) is valid, then we can invoke the union bound to show that

As already mentioned in Section 4.3, everything boils down to controlling $|\bm{a}^{\top}\bm{H}^{s}\bm{u}^{\star}|$ for $s\geq 1$. This is accomplished via the following lemma.

The proof of Lemma 5 is combinatorial in nature, which we defer to Appendix C. ∎

A similar result in (Tao, 2013, Lemma 2.3) has studied the bilinear forms of the high order terms of an i.i.d. random matrix, with a few distinctions. First of all, Tao (2013) assumes that each entry of the noise matrix is i.i.d. and has finite fourth moment (if the noise variance is rescaled to be 1); these assumptions break in examples like matrix completion. Moreover, Tao (2013) focuses on the case with $k=2$, and does not lead to high-probability bounds (which are crucial for, e.g. entrywise error control).

Using Markov’s inequality and the union bound, we can translate Lemma 5 into a high probability bound as follows.

Under the assumptions of Lemma 5, there exists some universal constant $c_{2}>0$ such that

In addition, in view of Lemma 1 and the condition (14), one has

with probability $1-O(n^{-10})$, which together with Lemma 2 implies $\lambda\geq 3\|\bm{H}\|$. This further leads to

Putting the above bounds together and using the fact that $\lambda$ is real-valued and $\lambda\geq 1/2$ (cf. Lemma 2), we have

as long as \max\big{\{}B\log n,\sqrt{n\sigma^{2}\log n}\big{\}} is sufficiently small. Here, the last line also uses the fact that $\mu\geq 1$ (and hence $\sqrt{\mu/n}\gg n^{-10}$). This concludes the proof.

Extension: perturbation analysis for the rank-r𝑟r case

The eigenvalue perturbation analysis in Section 4 can be extended to accommodate the case where $\bm{M}^{\star}$ is symmetric and rank-$r$, as detailed in this section. As before, assume that the $r$ non-zero eigenvalues of $\bm{M}^{\star}$ obey \lambda_{\max}^{\star}=\big{|}\lambda_{1}^{\star}\big{|}\geq\cdots\geq\big{|}\lambda_{r}^{\star}\big{|}=\lambda_{\min}^{\star}. Once again, we start with a master bound.

This result allows us to control the perturbation of the linear form of eigenvectors. The perturbation upper bound grows as either the rank $r$ or the condition number $\kappa$ increases.

One of the most important consequences of Theorem 4 is a refinement of the Bauer-Fike theorem concerning eigenvalue perturbations as follows.

Consider the $l$th ($1\leq l\leq r$) eigenvalue $\lambda_{l}$ of $\bm{M}$. Under the assumptions of Theorem 4, with probability at least $1-O(n^{-10})$, there exists $1\leq j\leq r$ such that

for some sufficiently small constant $c_{1}>0$.

In comparison, the Bauer-Fike theorem (Lemma 2) together with Lemma 1 gives a perturbation bound

For the low-rank case where $r\ll\sqrt{n}$, the eigenvalue perturbation bound derived in Corollary 5 can be much sharper than the Bauer-Fike theorem.

Another result that comes from Theorem 4 is the following bound that concerns linear forms of the eigen-subspace.

Under the same setting of Theorem 4, with probability $1-O(n^{-9})$ we have

Consequently, by taking $\bm{a}=\bm{e}_{i}$ ($1\leq i\leq n$) in Corollary 6, we arrive at the following statement regarding the alternative definition of the incoherence of the eigenvector matrix $\bm{U}$ (see Remark 2).

Under the same setting of Theorem 4, with probability $1-O(n^{-8})$ we have

Given that $\left\|\bm{U}\right\|_{2,\infty}=\max_{1\leq i\leq n}\left\|\bm{e}_{i}^{\top}\bm{U}\right\|_{2}$ and recalling our assumption implies $\left\|\bm{U}^{\star}\right\|_{2,\infty}\leq\sqrt{\mu r/n}$, we can invoke Corollary 6 and the union bound to derive the advertised entrywise bounds. ∎

The eigenvector matrix is often employed to form a reasonably good initial guess for several nonconvex statistical estimation problems Keshavan et al. (2010), and the above kind of incoherence property is crucial in guaranteeing fast convergence of the subsequent nonconvex iterative refinement procedures Ma et al. (2019).

Unfortunately, these results fall short of providing simple perturbation bounds for the eigenvectors; in other words, the above-mentioned bounds do not imply the size of the difference between $\bm{U}$ and $\bm{U}^{\star}$. The challenge arises in part due to the lack of orthonormality of the eigenvectors when dealing with asymmetric matrices. Analyzing the eigenspace perturbation for the general rank-$r$ case will likely require new analysis techniques, which we leave for future work. There is, however, some special case in which we can develop eigenvector perturbation theory, as detailed in the next subsection.

The theory for the rank-$r$ case has recently been significantly improved; see our follow-up work Cheng et al. (2020) for details.

where $\bm{H}_{1}$ and $\bm{H}_{2}$ are independent noise matrices. The goal is to estimate the singular value and singular vectors of $\bm{M}^{\star}$ from $\bm{M}_{1}$ and $\bm{M}_{2}$.

We attempt estimation via the standard dilation trick (e.g. Tao (2012)). This consists of embedding the matrices of interest within a larger block matrix

Here, we place $\bm{M}_{1}$ and $\bm{M}_{2}$ in two different subblocks, in order to “asymmetrize” the dilation matrix. The rationale is that $\bm{M}_{\mathsf{d}}^{\star}$ is a rank-2 symmetric matrix with exactly two nonzero eigenvalues

whose corresponding eigenvectors are given by

respectively. This motivates us to perform eigen-decomposition of $\bm{M}_{\mathsf{d}}$, and use the top-2 eigenvalues and eigenvectors to estimate $\lambda^{\star}$, $\bm{u}^{\star}$ and $\bm{v}^{\star}$, respectively.

Eigenvalue perturbation analysis. As an immediate consequence of Corollary 5, the two leading eigenvalues of $\bm{M}_{\mathsf{d}}$ provide fairly accurate estimates of the leading singular value $\lambda^{\star}$ of $\bm{M}^{\star}$, as stated below.

for some sufficiently small constant $c_{1}>0$.

To begin with, it follows from Corollary 5 that both $\lambda_{1}^{\mathsf{d}}$ and $\lambda_{2}^{\mathsf{d}}$ are close to either $\lambda^{\star}$ or $-\lambda^{\star}$. Repeating similar arguments as in the proof of Lemma 2 (which we omit here), we can immediately show the separation between these two eigenvalues, namely, $\lambda_{1}^{\mathsf{d}}$ (resp. $\lambda_{2}^{\mathsf{d}}$) is close to $\lambda^{\star}$ (resp. $-\lambda^{\star}$). ∎

Eigenvector perturbation analysis. We then move on to studying the eigenvector perturbation bounds. Specifically, denote by $\bm{u}_{1}^{\mathsf{d}}$ and $\bm{u}_{2}^{\mathsf{d}}$ the eigenvectors of $\bm{M}_{\mathsf{d}}$ associated with its two leading eigenvalues $\lambda_{1}^{\mathsf{d}}$ and $\lambda_{2}^{\mathsf{d}}$, respectively. Without loss of generality, we assume that $\lambda_{1}^{\mathsf{d}}\geq\lambda_{2}^{\mathsf{d}}$. If we write

then we can employ $\bm{u}_{1,1}^{\mathsf{d}}$ and $\bm{u}_{1,2}^{\mathsf{d}}$ to estimate $\bm{u}^{\star}$ and $\bm{v}^{\star}$ after proper normalization, namely,

The following theorem develops error bounds for both $\bm{u}$ and $\bm{v}$, which we establish in Appendix I. Here, we denote $\min\|\bm{x}\pm\bm{y}\|_{2}=\min\{\|\bm{x}-\bm{y}\|_{2},\|\bm{x}+\bm{y}\|_{2}\}$, and $\min\|\bm{x}\pm\bm{y}\|_{\infty}=\min\{\|\bm{x}-\bm{y}\|_{\infty},\|\bm{x}+\bm{y}\|_{\infty}\}$.

provided that there exists some some sufficiently small constant $c_{1}>0$ such that

Similar to the symmetric rank-1 case, the estimation errors of the estimates ${\bm{u}}$ and ${\bm{v}}$ are well-controlled in any deterministic direction (e.g. the entrywise errors are well-controlled). This allows us to complete the theory for the case when $\bm{M}^{\star}$ is a real-valued and rank-1 matrix.

Further, we conduct numerical experiments for matrix completion when $\bm{M}^{\star}$ is a rank-1 and asymmetric matrix in Fig. 4. Here, we suppose that at most 1 sample is observed for each entry, and we estimate the singular value and singular vectors of $\bm{M}^{\star}$ via the above-mentioned dilation trick, coupled with the asymmetrization procedure discussed in Section J. The numerical performance confirms that the proposed technique outperforms vanilla SVD in spectral estimation.

Finally, we remark that the asymptotic behavior of the eigenvalues of asymmetric random matrices has been extensively explored in the physics literature (e.g. (Sommers et al., 1988; Khoruzhenko, 1996; Brezin and Zee, 1998; Chalker and Mehlig, 1998; Feinberg and Zee, 1997; Lytova and Tikhomirov, 2018)). Their focus, however, has largely been to pin down the asymptotic density of the eigenvalues, similar to the semi-circle law in the symmetric case. Nevertheless, a sharp perturbation bound for the leading eigenvalue — particularly for the low-rank case — is beyond their reach. A few recent papers began to explore the locations of eigenvalue outliers that fall outside the bulk predicted by the circular law Tao (2013); Rajagopalan (2015); Benaych-Georges and Rochet (2016); Bordenave and Capitaine (2016). The results reported therein either do not focus on obtaining the right convergence rate (e.g. providing only a bound like $|\lambda-\lambda^{\star}|=o(|\lambda^{\star}|)$) or are restricted to a special family of ground truth (e.g. the one with a diagonal block equal to identity) or i.i.d. noise. As a result, these prior results are insufficient to demonstrate the power and benefits of the eigen-decomposition method in the presence of data asymmetry.

3 Proof of Theorem 4

Without loss of generality, we shall assume $\lambda_{\max}^{\star}=\lambda_{1}^{\star}=1$ throughout the proof. To begin with, Lemma 2 implies that for all $1\leq l\leq r$,

as long as $\left\|\bm{H}\right\|<1/(2\kappa)$. In view of the Neumann trick (Theorem 2), we can derive

where the third line follows since \sum_{j=1}^{r}\big{|}\bm{u}_{j}^{\star\top}\bm{u}_{l}\big{|}^{2}\leq\|\bm{u}_{l}\|_{2}^{2}=1, and the last inequality makes use of (49). Apply Corollary 4 to reach

with the proviso that $|\lambda_{l}|>1/(2\kappa)$ and $\max\left\{B\log n,\sqrt{n\sigma^{2}\log n}\right\}\leq c_{1}/\kappa$ for some sufficiently small constant $c_{1}>0$. The condition $|\lambda_{l}|>1/(2\kappa)$ follows immediately by combining Lemma 2, Lemma 1 and the condition (34).

Discussions

In this paper, we demonstrate the remarkable advantage of eigen-decomposition over SVD in the presence of asymmetric noise matrices. This is in stark contrast to conventional wisdom, which is generally not in favor of eigen-decomposition for asymmetric matrices. Our results only reflect the tip of an iceberg, and there are many outstanding issues left answered. We conclude the paper with a few future directions.

Sharper eigenvalue perturbation bounds for the rank-$r$ case. Our current results in Section 5 provide an eigenvalue perturbation bound on the order of $r/\sqrt{n}$, assuming the truth is rank-$r$. However, numerical experiments suggest that the dependency on $r$ might be improvable. It would be interesting to see whether further theoretical refinement is possible, e.g. whether it is possible to improve it to $O(\sqrt{r/n})$.

Eigenvector perturbation bounds for the rank-$r$ case. As mentioned before, the current theory falls short of providing eigenvector perturbation bounds for the general rank-$r$ case. The main difficulty lies in the lack of orthogonality of the eigenvectors of the observed matrix $\bm{M}$. Nevertheless, when the size of the noise is not too large, it is possible to establish certain near-orthogonality of the eigenvectors, which might in turn lead to sharp control of eigenvector perturbation.

A challenging signal-to-noise ratio regime. Take the rank-1 case for example: the present work focuses on the regime where $\|\bm{H}\|\lesssim\|\bm{M}^{\star}\|/\sqrt{\log n}$, and it is known that spectral methods fail to yield reliable estimation if $\|\bm{H}\|\gg\|\bm{M}^{\star}\|$. There is, however, a “gray” region (which includes, for example, the case with $\|\bm{H}\|\approx\|\bm{M}^{\star}\|$) that has not been addressed. Developing non-asymptotic yet informative perturbation bounds for this regime is likely very challenging and requires new analysis techniques, which we leave for future investigation.

Correlated noise. The current theoretical development relies heavily on the assumption that the noise matrix $\bm{H}$ contains independent random entries. There is no shortage of examples where the noise matrix is asymmetric but is not composed of independent entries. For instance, in blind deconvolution Li et al. (2018), the noise matrix is a sum of independent asymmetric matrices. Can we develop eigenvalue perturbation theory for this class of noise?

Statistical inference of eigenvalues and eigenvectors. In various applications like network analysis and inference, one might be interested in determining the (asymptotic) eigenvalue and eigenvector distributions of a random data matrix, in order to produce valid confidence intervals Johnstone (2001); Bai and Yao (2008); Cai et al. (2017); Cape et al. (2018); Xia (2018); Bao et al. (2018); Chen et al. (2019c). Can we use the current framework to characterize the distributions of the leading eigenvalues as well as certain linear forms of the eigenvectors of $\bm{M}$ when the noise matrix is non-symmetric?

with $\bm{v}$ being a unit vector. This falls under the category of the spiked covariance model (Johnstone and Lu, 2009). One strategy to estimate the spectral norm $\lambda^{\star}=2$ of $\bm{\Sigma}^{\star}$ is to look at the spectrum of the sample covariance matrix $\hat{\bm{\Sigma}}=\frac{1}{n}\sum_{i=1}^{n}\bm{X}_{i}\bm{X}_{i}^{\top}$. Motivated by the results of this paper, we propose an alternative strategy by looking at the following asymmetrized sample covariance matrix

where $\mathsf{Upper}(\cdot)$ (resp. $\mathsf{Lower}(\cdot)$) extracts out the upper (resp. lower) triangular part of the matrix, including (resp. excluding) the diagonal entries. As can be seen from Fig. 5, the largest eigenvalue of the asymmetrized $\widehat{\bm{\Sigma}}_{\mathsf{asym}}$ is much closer to the true spectral norm of $\bm{\Sigma}^{\star}$, compared to the largest singular value of the sample covariance matrix $\hat{\bm{\Sigma}}$. We leave the theoretical understanding of such findings to future investigation.

Acknowledgment

Y. Chen is supported in part by the AFOSR YIP award FA9550-19-1-0030, by the ARO grant W911NF-18-1-0303, by the ONR grant N00014-19-1-2120, by the NSF grants CCF-1907661 and IIS-1900140, and by the Princeton SEAS innovation award. J. Fan is supported in part by the NSF grants DMS-1662139 and DMS-1712591, by the ONR grant N00014-19-1-2120, and by the NIH grant 2R01-GM072611-13. C. Cheng is supported in part by the Elite Undergraduate Training Program of School of Mathematical Sciences in Peking University, and by the William R. Hewlett Stanford graduate fellowship. We thank Cong Ma for helpful discussions, and Zhou Fan for telling us an example of asymmetrizing the Gaussian matrix.

References

Appendix A Proofs for preliminary facts

Without loss of generality, suppose the leading eigenvalue of $\bm{M}^{\star}$ is 1.

(1) We start by proving the rank-1 case, in which the eigenvalues of $\bm{M}^{\star}$ are either 0 or 1. Theorem 1 immediately implies that

To establish the lemma, our aim is to show that there is exactly one eigenvalue of $\bm{M}$ — denoted by $\lambda$ — lying in the disk $\mathcal{B}\left(1,\left\|\bm{H}\right\|\right)$, and it has multiplicity 1. If this is true, then the assumption $\left\|\bm{H}\right\|<1/2$ indicates that $|\lambda|\geq 1-\|\bm{H}\|>\|\bm{H}\|\geq z$ for any $z\in\mathcal{B}(0,\|\bm{H}\|)$, and hence $\lambda$ must be the leading eigenvalue. Furthermore, since $\bm{M}$ is a real-valued matrix, both $\lambda$ and its complex conjugate $\overline{\lambda}$ are eigenvalues of $\bm{M}$. However, if $\lambda\neq\overline{\lambda}$, then both of these two eigenvalues fall within $\mathcal{B}(1,\|\bm{H}\|)$, resulting in contradiction. As a result, one necessarily has $\lambda=\overline{\lambda}$ and both of them are real-valued. A similar argument demonstrates that the eigenvector $\bm{u}$ associated with $\lambda$ is also real-valued.

We then justify the existence and uniqueness of an eigenvalue in $\mathcal{B}\left(1,\left\|\bm{H}\right\|\right)$. Denote $\Lambda(\bm{M})=\left\{\lambda_{1},\cdots,\lambda_{n}\right\}$, and define a set of auxiliary matrices

Recognizing that the set of eigenvalues of $\bm{M}_{t}$ depends continuously on $t$ (e.g. (Embree and Trefethen, 2001, Theorem 6)), we can write

with each $\lambda_{j}(t),1\leq j\leq n$ being a continuous function in $t$. Meanwhile, as long as $\left\|\bm{H}\right\|<1/2$ and $0\leq t\leq 1$, the two disks $\mathcal{B}\left(1,t\left\|\bm{H}\right\|\right)$ and $\mathcal{B}\left(0,t\left\|\bm{H}\right\|\right)$ are always disjoint. Thus, the continuity of the spectrum (w.r.t. $t$) requires $\lambda_{j}(t)$ to always stay within the same disk where $\lambda_{j}(0)\in\{0,1\}$ lies, namely,

Given that $\bm{M}^{\star}$ (or $\bm{M}_{0}$) has $n-1$ eigenvalues equal to and one eigenvalue equal to $1$, we establish the lemma for the rank-1 case.

(2) We now turn to the rank-$r$ case. Repeating the above argument for the rank-1 case, we can immediately show that: if $\|\bm{H}\|<\lambda_{r}^{\star}/2$, then (i) there are exactly $n-r$ eigenvalues lying within $\mathcal{B}(0,\|\bm{H}\|)$; (ii) all other eigenvalues lie within $\cup_{1\leq j\leq r}\mathcal{B}(\lambda_{j}^{\star},\|\bm{H}\|)$, which are exactly the top-$r$ leading eigenvalues of $\bm{M}$. This concludes the proof.

A.2 Proof of Theorem 2

When $\left\|\bm{H}\right\|_{2}<|\lambda_{l}|$, one can invert $\bm{I}-\frac{1}{\lambda_{l}}\bm{H}$ and use the assumption (5) to reach

where the last line follows by rearranging terms. Finally, replacing \big{(}\bm{I}-\frac{1}{\lambda_{l}}\bm{H}\big{)}^{-1} with the Neumann series $\sum_{s=0}^{\infty}\frac{1}{\lambda_{l}^{s}}\bm{H}^{s}$, we establish the theorem.

A.3 Proof of Lemma 3

(1) We start with the rank-1 case. Towards this, we resort to the Neumann trick in Theorem 2, which in the rank-1 case reads

From Lemma 2, we know that $\lambda$ is real-valued and that $\lambda>1-\|\bm{H}\|\geq 3/4>\|\bm{H}\|$ under our assumption. This together with (57) yields

where the last inequality holds since $\lambda\geq 3/4$ and $\lambda-\|\bm{H}\|\geq 1-2\|\bm{H}\|\geq 1/2$.

Next, by decomposing $\bm{u}$ into two orthogonal components $\bm{u}=(\bm{u}^{\star\top}\bm{u})\bm{u}^{\star}+(\bm{u}-(\bm{u}^{\star\top}\bm{u})\bm{u}^{\star})$, we obtain

The inequality (59) holds since $\left(\bm{u}^{\star\top}\bm{u}\right)\bm{u}^{\star}$ is orthogonal projection of $\bm{u}$ onto the subspace spanned by $\bm{u}^{\star}$, and hence $\|\bm{u}-\left(\bm{u}^{\star\top}\bm{u}\right)\bm{u}^{\star}\|_{2}\leq\|\bm{u}-\frac{1}{\lambda}\left(\bm{u}^{\star\top}\bm{u}\right)\bm{u}^{\star}\|_{2}$.

Finally, (60) together with the fact that $\bm{u}$ is real-valued (cf. Lemma 2) leads to the advertised bound:

(2) For the rank-$r$ case, it is seen that for any $1\leq l\leq r$,

where the inequality arises since \sum\nolimits_{j=1}^{r}\big{(}\bm{u}_{j}^{\star\top}\bm{u}_{l}\big{)}\bm{u}_{j}^{\star} is the Euclidean projection of $\bm{u}_{l}$ onto the span of $\{\bm{u}_{1}^{\star},\cdots,\bm{u}_{r}^{\star}\}$. In addition, observe that

This taken collectively with Theorem 2 leads to

Appendix B Proof for the lower bound in Lemma 4

where the first identity comes from the property of KL divergence for product measures, the second identity follows since $\mathsf{KL}(\mathcal{N}(\mu_{1},\sigma^{2})\,\|\ \mathcal{N}(\mu_{2},\sigma^{2}))=\frac{(\mu_{1}-\mu_{2})^{2}}{2\sigma^{2}}$, and the third one holds since $\|\bm{u}^{\star}\|_{2}=1$. The same argument yields \mathsf{KL}\big{(}\hat{P}\,\|\,P\big{)}=\Delta^{2}/\sigma^{2}. In view of (Tsybakov, 2009, Corollary 2.6), if

Appendix C Proof of Lemma 5

To establish this lemma, we exploit entrywise independence of $\bm{H}$ and develop a combinatorial trick.

To begin with, we expand the quantity of interest as

to denote such a collection of $(s+1)k$ indices. Thus, one can write

We shall often think of this sum-product graphically by viewing each index pair \big{(}i_{t-1}^{(b)},i_{t}^{(b)}\big{)} as a directed edge over the vertex set $[n]^{(s+1)k}$. As a result, this gives us a set of $sk$ edges in total $\{\bm{e}_{t}^{(b)}\mid 1\leq t\leq s,1\leq b\leq k\}$, where $\bm{e}_{t}^{(b)}$ represents \big{(}i_{t-1}^{(b)},i_{t}^{(b)}\big{)}. See Fig. 6 for an illustration. In what follows, two edges $\bm{e}_{t_{1}}^{(b_{1})}$ and $\bm{e}_{t_{2}}^{(b_{2})}$ are said to be equivalent in $\mathcal{I}$, denoted by

if the values of the corresponding vertices are identical (namely, $i_{t_{1}-1}^{(b_{1})}=i_{t_{2}-1}^{(b_{2})}$ and $i_{t_{1}}^{(b_{1})}=i_{t_{2}}^{(b_{2})}$).

When $(s+1)k\ll n$, most of the summands in the above expansion vanish. In fact, for any summand associated with a given $\mathcal{I}$: as long as there exists a distinct edge \big{(}i_{t^{*}-1}^{(b^{*})},i_{t^{*}}^{(b^{*})}\big{)} (i.e. not equal to any other edge associated with $\mathcal{I}$), then the contribution of this summand is zero, namely,

As a result, for any term with non-zero contribution, every edge \big{(}i_{t-1}^{(b)},i_{t}^{(b)}\big{)} must appear at least twice.

To enable simple yet effective upper bounds on the non-vanishing terms, we group those terms of the same “type” and look at each “type” separately. Specifically:

We represent a collection of edges \big{\{}(i_{t-1}^{(b)},i_{t}^{(b)})\mid 1\leq t\leq s,1\leq b\leq k\big{\}} as \big{\{}\bm{e}_{t}^{(b)}\mid 1\leq t\leq s,1\leq b\leq k\big{\}}. See Fig. 6 for an illustration.

Each type encodes a set of constraints across edges, namely, which edges correspond to the same pair of vertices. To be precise, we define each type to be a partition of $\{\bm{e}_{t}^{(b)}\mid 1\leq t\leq s,1\leq b\leq k\}$ into a disjoint union of subsets, so that all edges falling within the same subset are equivalent. For instance, when $s=k=2$, one possible type is \big{\{}\{\bm{e}_{1}^{(1)},\bm{e}_{2}^{(2)}\},\{\bm{e}_{1}^{(2)},\bm{e}_{2}^{(1)}\}\big{\}}, which encodes the constraints $\bm{e}_{1}^{(1)}=\bm{e}_{2}^{(2)}$ and $\bm{e}_{1}^{(2)}=\bm{e}_{2}^{(1)}$.

For each index collection $\mathcal{I}$ (defined in (62)), we write $\mathsf{type}(\mathcal{I})=\mathcal{T}$ if the associated edge set of $\mathcal{I}$ satisfies the constraints encoded by a type $\mathcal{T}$.

With this grouping strategy and (63) in mind, we can derive

where the last line follows since all types outside $\Gamma^{+}$ have vanishing contributions (cf. (63)). To bound the right-hand side of (65), we need the following lemma. Here and throughout, $|\mathcal{T}|$ represents the number of non-empty subsets in the partition associated with $\mathcal{T}$.

Armed with this lemma, we can further obtain

where we have grouped the types based on their cardinality. The last inequality results from $|\mathcal{T}|\leq sk/2$ in Lemma 9. The following lemma bounds the number of distinct types having the same cardinality.

$\sum_{\mathcal{T}\in\Gamma^{+}}\operatorname{\mathds{1}}_{\{|\mathcal{T}|=l\}}\leq\binom{sk}{l}\,l^{sk-l}\,\leq 2^{sk}l^{sk-l}$.

The quantity of interest is the number of ways to partition $sk$ edges into $l$ disjoint subsets, where each subset contains at least 2 edges. To bound this quantity, we first pick $l$ edges and assign each of them to a distinct subset; there are $\binom{sk}{l}$ different ways to achieve it. We still have $sk-l$ edges left, and the number of ways to assign them to $l$ subsets is clearly upper bounded by $l^{sk-l}$. This concludes the proof. ∎

where (i) uses the condition $l\leq sk/2$, and (ii) relies on the fact that $a^{sk-2l}b^{l}=a^{sk-2l}(\sqrt{b})^{2l}\leq\max\{a^{sk},(\sqrt{b})^{sk}\}$ for any $a,b>0$.

Taken collectively, the preceding bounds conclude the proof of Lemma 5, provided that Lemma 6 can be established.

Appendix D Proof of Lemma 6

where the last inequality holds since each entry of $\bm{u}^{\star}$ is bounded in magnitude by $\sqrt{\mu/n}$ (see Definition 1).

According to the definition of $\Gamma^{+}$ (see (64)), for any type $\mathcal{T}\in\Gamma^{+}$, each edge \big{(}i_{t-1}^{(b)},i_{t}^{(b)}\big{)} is repeated at least twice. The total number of distinct edges is exactly $|\mathcal{T}|$. For notational simplicity, suppose the distinct edges are $\bm{e}_{1},\cdots,\bm{e}_{|\mathcal{T}|}$, where $\bm{e}_{i}$ is repeated by $l_{i}\geq 2$ times. Then we can write

It then boils down to controlling \sum_{\mathcal{I}\in[n]^{(s+1)k}:\,\mathsf{type}(\mathcal{I})=\mathcal{T}}\prod_{b=1}^{k}\big{|}a_{i_{0}^{(b)}}\big{|}. This is achieved by the following lemma, which we establish in Appendix E.

This lemma together with (69) concludes the proof.

Appendix E Proof of Lemma 8

The proof of Lemma 8 is combinatorial in nature. Before proceeding, we introduce a few graphical notions that will be useful.

To begin with, recall that the index collection we have used so far is

which involves $(s+1)k$ vertices. As it turns out, it is helpful to further augment it into $2(s+1)k$ vertices via simple duplication. Specifically, introduce the following set of $2(s+1)k$ vertices

view $\{q_{t}^{(b)}\mid 0\leq t\leq s,1\leq b\leq k\}$ as a copy of $\{i_{t}^{(b)}\mid 0\leq t\leq s,1\leq b\leq k\}$;

view $\{p_{t}^{(b)}\mid 1\leq t\leq s+1,1\leq b\leq k\}$ as another copy of $\{i_{t-1}^{(b)}\mid 1\leq t\leq s+1,1\leq b\leq k\}$.

The main incentive is that it allows us to reparametrize

Recall that we have categorized $\mathcal{I}$ into different “types”. Given that $\mathcal{V}$ is simply a redundant representation of $\mathcal{I}$, we can also associate a type with each $\mathcal{V}$ by enforcing proper constraints. These constraints can be encoded via the following graphical notion.

For any given type $\mathcal{T}$, we define the induced graph $\mathcal{G}(\mathcal{T})$ with the vertex set $\mathcal{V}$ and the (undirected) edge set $\mathcal{E}(\mathcal{T})$ as $\mathcal{E}(\mathcal{T})=\mathcal{E}_{1}\cup{\mathcal{E}_{2}(\mathcal{T})}$, where

we connect the vertices $p_{t}^{(b)}$ (which is a copy of $i_{t}^{(b)}$) and $q_{t-1}^{(b)}$ (which represents $i_{t-1}^{(r)}$) by an edge;

whenever two edges \big{(}i_{t_{1}-1}^{(b_{1})},i_{t_{1}}^{(b_{1})}\big{)} and \big{(}i_{t_{2}-1}^{(b_{2})},i_{t_{2}}^{(b_{2})}\big{)} are equivalent in $\mathcal{I}$ (or $\mathcal{T}$), we draw two edges connecting the corresponding vertices in the induced graph.

As an illustration, Fig. 7 displays the induced graph for a simple example with $s=4$ and $k=2$, where

Here, the 45-degree lines correspond to the edges in $\mathcal{E}_{1}$, and the remaining edges come from $\mathcal{E}_{2}(\mathcal{T})$. Throughout the rest of the paper, we will abuse the notation and write $\mathsf{type}(\mathcal{V})=\mathcal{T}$ if $\mathcal{V}$ is induced by an index set of type $\mathcal{T}$.

One useful feature of the above induced graph is that: each edge subset in the partition associated with $\mathcal{T}$ corresponds to a connected component in $\mathcal{G}(\mathcal{T})$. For our purpose, it is convenient to divide all connected components of $\mathcal{G}(\mathcal{T})$ into two classes. To this end, we first define \mathcal{Q}_{0}\triangleq\big{\{}q_{0}^{(b)}\mid 1\leq b\leq k\big{\}} (as illustrated in Fig. 7).

Class 1: a connected component $\mathcal{C}$ belongs to Class 1 if $|\mathcal{C}\cap\mathcal{Q}_{0}|\leq 1$;

Class 2: a connected component $\mathcal{C}$ belongs to Class 2 if $|\mathcal{C}\cap\mathcal{Q}_{0}|\geq 2$.

We denote by $m_{1}(\mathcal{T})$ (resp. $m_{2}(\mathcal{T})$) the total number of Class 1 (resp. 2) connected components in $\mathcal{G}(\mathcal{T})$.

E.2 Proof of Lemma 8

Making use of the augmented set $\mathcal{V}$ (defined in (72)) and the way of reparametrization (73), we can write

As we shall see, the benefit of this bound is to allow us sum over the connected components of $\mathcal{G}(\mathcal{T})$ in a separate manner, since all indices in the same connected component must be identical. Specifically,

Let $\mathcal{X}_{1},\mathcal{X}_{2},\cdots$ (resp. $\mathcal{Y}_{1},\mathcal{Y}_{2},\cdots$) denote the collection of Class 1 (resp. Class 2) connected components in $\mathcal{G}(\mathcal{T})$;

Denote by $x_{i}$ (resp. $y_{i}$) the value assigned to all indices in $\mathcal{X}_{i}$ (resp. $\mathcal{Y}_{i}$);

With these notations in mind, we can decompose

where $m_{1}(\mathcal{T})$ is the total number of Class 1 connected components in the induced graph $\mathcal{G}(\mathcal{T})$. Here, (i) comes from the definitions of Class 1 and Class 2 connected components, and (ii) follows since $\|\bm{a}\|_{2}^{2}=1$.

We can thus finish the proof by observing that $m_{1}(\mathcal{T})\leq|\mathcal{T}|$, as claimed in the following lemma.

For any type $\mathcal{T}\in\Gamma^{+}$, one has $m_{1}(\mathcal{T})\leq|\mathcal{T}|\leq sk/2$.

E.3 Proof of Lemma 9

Clearly, one can define a partition of $\mathcal{Q}_{\textbackslash{}0}$ — denoted by $\mathcal{T}_{\mathcal{Q}_{\textbackslash{}0}}$ — induced by $\mathcal{T}$. Specifically, we say that $q_{t_{1}}^{(b_{1})}$ and $q_{t_{2}}^{(b_{2})}$ belong to the same connected subgraph of $\mathcal{Q}_{\textbackslash{}0}$ if \big{(}i_{t_{1}-1}^{(b_{1})},i_{t_{1}}^{(b_{1})}\big{)}=\big{(}i_{t_{2}-1}^{(b_{2})},i_{t_{2}}^{(b_{2})}\big{)} in $\mathcal{T}$. Clearly, $|\mathcal{T}|=\left|\mathcal{T}_{\mathcal{Q}_{\textbackslash{}0}}\right|$.

In addition, for any connected subgraph $\mathcal{C}_{q}$ in $\mathcal{Q}_{\textbackslash{}0}$, we denote by $\mathcal{C}$ the corresponding connected component in $\mathcal{G}(\mathcal{T})$. We can thus define a mapping $\psi$ that maps $\mathcal{C}_{q}$ to $\mathcal{C}$, whose domain is the set of all connected subgraphs in $\mathcal{Q}_{\textbackslash{}0}$. We claim that the collection of Class 1 connected components — denoted by $\{\mathcal{X}_{j}\}$ — obeys

To justify the claim (78), it suffices to show that $\mathcal{X}_{j}\cap\mathcal{Q}_{\textbackslash{}0}\neq\emptyset$. Given that $\mathcal{T}\in\Gamma^{+}$ (so that each subset of the partition associated with $\mathcal{T}$ has cardinality at least 2) and that $\mathcal{X}_{j}$ is defined over the induced graph $\mathcal{G}(\mathcal{T})$ (which is a redundant representation of the original index collection), one must have

Further, from the construction of the induced graph, it is easily seen that

However, by definition of Class 1 connected component, we have $\left|\mathcal{X}_{j}\cap\mathcal{Q}_{0}\right|\leq 1$, implying that $\mathcal{X}_{j}\cap\mathcal{Q}_{\textbackslash{}0}\neq\emptyset$. This establishes the claim (78).

Finally, when $\mathcal{T}\in\Gamma^{+}$, each connected subgraph of $\mathcal{T}_{\mathcal{Q}_{\textbackslash{}0}}$ also contains at least two nodes. As a consequence,

Appendix F Proof of Corollary 4

In the sequel, we assume that $20\log n$ is an integer to avoid the clumsy notation $\lfloor 20\log n\rfloor$. But it is straightforward to extend it to the case where $20\log n$ is not an integer.

It follows from Markov’s inequality that for any even integer $k$,

For any $s\leq 20\log n$, choose $k$ such that $sk\in\left[20\log n,40\log n\right]$. It is straightforward to show — using the union bound — that with probability at least $1-O(n^{-10})$,

as long as $c_{2}>0$ is some sufficiently large constant.

Appendix G Proof of Corollary 5

Taking $\bm{a}=\bm{u}_{i}^{\star}$ for some $1\leq i\leq r$ in (35), we see that with high probability,

In addition, since $\bm{u}_{1}^{\star},\cdots,\bm{u}_{r}^{\star}$ are orthogonal to each other, we have that: for any $1\leq i\leq r$,

Therefore, combining the above two bounds yields

Finally, it comes from Lemma 3 that if $\|\bm{H}\|\ll 1/\kappa^{2}$, then \sum_{1\leq i\leq r}\big{|}\bm{u}_{i}^{\star\top}\bm{u}_{l}\big{|}^{2}\gtrsim 1, and hence

This combined with (80) establishes the claim, as long as the spectral norm condition on $\|\bm{H}\|$ can be guaranteed. In view of Lemma 1, we have $\|\bm{H}\|\ll 1/\kappa^{2}$ under the condition (38), thus concluding the proof.

Appendix H Proof of Corollary 6

Since $\|\bm{a}\|_{2}=1$, it follows immediately from (36) that

Appendix I Proof of Theorem 5

To simplify presentation, we introduce the following notation throughout this section:

In addition, we denote $\widetilde{\bm{u}}_{1,1}=\bm{u}^{\mathsf{dilation}}_{1,1}$ and $\widetilde{\bm{u}}_{1,2}=\bm{u}^{\mathsf{dilation}}_{1,2}$. Recall that we assume $\widetilde{\lambda}_{1}\geq\widetilde{\lambda}_{2}$. We also denote $\min\{\|\bm{x}\pm\bm{y}\|_{2}\}=\min\{\|\bm{x}-\bm{y}\|_{2},\|\bm{x}+\bm{y}\|_{2}\}$.

To begin with, applying Theorem 4 on $\bm{M}_{\mathsf{dilation}}$ and taking $\bm{a}=\widetilde{\bm{u}}_{2}^{\star}$, we derive

where the identity arises since $\widetilde{\bm{u}}_{2}^{\star\top}\widetilde{\bm{u}}_{1}^{\star}=0$. Given that $\widetilde{\lambda}_{2}^{\star}=-1$ and $\widetilde{\lambda}_{1}>0$ (see Corollary 8), we have $1-\widetilde{\lambda}_{2}^{\star}/\widetilde{\lambda}_{1}>1$. The near orthogonality property can then be described as follows

It then boils down to showing that $\min\{\|\bm{u}\pm\bm{u}^{\star}\|_{2}\}\lesssim\min\{\|\widetilde{\bm{u}}_{1}\pm\widetilde{\bm{u}}_{1}^{\star}\|_{2}\}$. To this end, we see that the estimate (46) satisfies

Here, (86) makes use of the triangle inequality, and (87) follows since

where the last inequality relies on the fact that $\left\|\widetilde{\bm{u}}_{1,1}-\bm{u}^{\star}/\sqrt{2}\right\|_{2}\leq\left\|\widetilde{\bm{u}}_{1}-\widetilde{\bm{u}}_{1}^{\star}\right\|_{2}$ (the first $n_{1}$ coordinates). Similarly, one can derive the above bounds for $\|\bm{u}+\bm{u}^{\star}\|_{2}$ as well. Therefore, we are left with

I.2 Perturbation bounds for linear forms of eigenvectors

Given that the first $n_{1}$ coordinates of $\widetilde{\bm{u}}_{1}^{\star}$ and $\widetilde{\bm{u}}_{2}^{\star}$ are both $\bm{u}^{\star}/\sqrt{2}$, we can invoke Corollary 8 (i.e. $\widetilde{\lambda}_{1}\asymp 1$) to obtain

In view of (84) and the fact $\left|\bm{a}^{\top}\bm{u}^{\star}\right|\leq 1$, one has

Recall that $\widetilde{\lambda}_{1}>0$ (cf. Corollary 8). If we further have

then we can use the triangle inequality to reach

as claimed. It then remains to justify (94). Towards this, it suffices to combine a series of consequences from (85), Corollary 8, and (90), namely,

The proof for the bounds on $\bm{v}$ is similar and is thus omitted.

I.3 Entrywise eigenvector perturbation bounds

Recognizing that $\left\|\bm{u}-\bm{u}^{\star}\right\|_{\infty}=\max_{i}\left|\bm{e}_{i}^{\top}\bm{u}-\bm{e}_{i}^{\top}\bm{u}^{\star}\right|$ and using the incoherence $\left|\bm{e}_{i}^{\top}\bm{u}\right|\leq\sqrt{\mu/n}$, we can prove this claim directly by invoking the results established in Appendix I.2 and taking the union bound.

Appendix J Asymmetrization of data samples: two examples

As mentioned earlier, an independent and asymmetric noise matrix arises when we collect two samples for each entry of the matrix of interest and arrange the samples in an asymmetric manner (i.e. placing 1 sample on the upper triangular part and the other on the lower triangular part). Interestingly, our results might be applicable for some cases where we only have 1 sample for each entry. In what follows, we describe two examples of this kind similar to the ones discussed in Section 4.2, but with a symmetric noise matrix. Once again, it is assumed that $\bm{M}^{\star}$ is a rank-1 matrix with leading eigenvalue 1 and incoherence parameter $\mu$.

When $\sigma$ is known, one can decouple the upper and low triangular parts of $\bm{H}$ by adding a skew-symmetric Gaussian matrix $\bm{\Delta}$. Specifically, our strategy is:

For each $1\leq i\leq j\leq n$, generate $\Delta_{ij}\sim\mathcal{N}(0,\sigma^{2})$ independently, and set ${\Delta}_{ji}=-{\Delta}_{ij}$;

Compute the leading eigenvalue and eigenvector of $\bm{M}+\bm{\Delta}$.

This is motivated by a simple observation from Gaussianality: $\bm{H}+\bm{\Delta}$ is now an asymmetric matrix whose off-diagonal entries are i.i.d. $\mathcal{N}(0,2\sigma^{2})$; in fact, it is easy to verify that $H_{ij}+\Delta_{ij}$ and $H_{ji}-\Delta_{ji}$ are independent Gaussian random variables. As a result, the orderwise bounds (26) continue to hold if $\lambda$ and $\bm{u}$ are taken to be the leading eigenvalue and eigenvector of $\bm{M}+\bm{\Delta}$, respectively.

While this asymmetrization procedure comes with the price of doubling the noise variance, the eigenvalue perturbation bound may still be significantly smaller — up to a factor of $O(\sqrt{n})$ — than the bound for the SVD approach. This is also confirmed in the numerical simulations reported in Fig. 8.

While this asymmetrization procedure achieves enhanced eigenvalue estimation accuracy compared to the SVD approach (see Fig. 8(a)(b)), it results in higher eigenvector estimation errors (see Fig. 8(c)). This is perhaps not surprising as we have added extra noise to the observed matrix. One way to mitigate this issues is to (1) generate $K$ independent copies of $\bm{\Delta}$; (2) compute the leading eigenvector of each copy of $\bm{M}+\bm{\Delta}$, denoted by $\{\bm{u}^{(l)}\mid 1\leq l\leq K\}$; and (3) aggregate these eigenvectors, namely, compute the leading eigenvector of $\frac{1}{K}\sum_{l=1}^{K}\bm{u}^{(l)}\bm{u}^{(l)\top}$. As can be seen in the green line of Fig. 8(c), this allows us to mitigate the effect of the extra noise component.

Low-rank matrix completion. Suppose that each entry $M_{ij}^{\star}$ ($i\geq j$) is observed independently with probability $p$. This is different from the settings in Section 4.2, as we do not have additional samples for $M_{ji}^{\star}$ ($i\geq j$). In order to arrange the data samples in an asymmetric and independent fashion, we employ a simple resampling technique to decouple the statistical dependency between $M_{ij}$ and $M_{ji}$:

Define $p^{\mathsf{asym}}$ so that $p=1-(1-p^{\mathsf{asym}})^{2}$ (i.e. $p^{\mathsf{asym}}=\frac{p}{1+\sqrt{1-p}}$). For any pair $i>j$, set

Compute the leading eigenvalue and eigenvector of $\bm{M}^{\mathsf{asym}}$.

As can be easily verified, this scheme is equivalent to saying that (i) with probability $p$, either $M_{ij}^{\mathsf{asym}}$ or $M_{ji}^{\mathsf{asym}}$ is taken to be a rescaled version of $M_{ij}^{\star}$; (ii) for any $i\neq j$, the entries $\{M^{\mathsf{asym}}_{ij}\}$ are independently drawn. Since $p^{\mathsf{asym}}\asymp p$, our results (28) in Section 4.2 remain valid, as long as $\lambda$ and $\bm{u}$ are set to be the leading eigenvalue and eigenvector of $\bm{M}^{\mathsf{asym}}$, respectively. Numerical simulations have been carried out in Fig. 9 to verify the effectiveness of this scheme.