# Sketching sparse low-rank matrices with near-optimal
sample- and time-complexity

Xiaoqi Liu and Ramji Venkataramanan Xiaoqi Liu was supported in part by a
Schlumberger Cambridge International Scholarship funded by the Cambridge
Trust. X. Liu and R. Venkataramanan are with the Department of Engineering,
University of Cambridge. Email: {xl394<EMAIL_ADDRESS>This paper will be
presented in part at the 2022 IEEE International Symposium on Information
Theory.

###### Abstract

We consider the problem of recovering an $n_{1}\times n_{2}$ low-rank matrix
with $k$-sparse singular vectors from a small number of linear measurements
(sketch). We propose a sketching scheme and an algorithm that can recover the
singular vectors with high probability, with a sample complexity and running
time that both depend only on $k$ and not on the ambient dimensions $n_{1}$
and $n_{2}$. Our sketching operator, based on a scheme for compressed sensing
by Li et al. [1] and Bakshi et al. [2], uses a combination of a sparse parity
check matrix and a partial DFT matrix. Our main contribution is the design and
analysis of a two-stage iterative algorithm which recovers the singular
vectors by exploiting the simultaneously sparse and low-rank structure of the
matrix. We derive a nonasymptotic bound on the probability of exact recovery.
We also show how the scheme can be adapted to tackle matrices that are
approximately sparse and low-rank. The theoretical results are validated by
numerical simulations.

##  1 Introduction

We consider the problem of _sketching_ a large data matrix which is low-rank
with sparse singular vectors. A sketch is a compressed representation of the
matrix, obtained via linear measurements of the matrix entries [3, 4]. The
goal is to design a sketching scheme so that the singular vectors can be
efficiently recovered from a sketch of a much smaller dimension.

Given a data matrix $\boldsymbol{X}\in{\mathbb{R}}^{n_{1}\times n_{2}}$, its
sketch $\boldsymbol{y}\in{\mathbb{C}}^{n_{1}\times n_{2}}$ is obtained as
$\boldsymbol{y}=\mathcal{A}(\boldsymbol{X})$, where $\mathcal{A}$ is a linear
operator defined via $m$ prespecified matrices
$\boldsymbol{A}_{1},\ldots,\boldsymbol{A}_{m}\in{\mathbb{C}}^{{n_{1}}\times{n_{2}}}$.
The entries of the sketch are:

$y_{i}={\sf
trace}(\boldsymbol{A}_{i}^{T}\boldsymbol{X})=\sum_{j=1}^{n_{1}}\sum_{\ell=1}^{n_{2}}(A_{i})_{j,\ell}X_{j,\ell},\quad
i=1,\ldots,m.$ (1)

We consider data matrices $\boldsymbol{X}$ of the form

$\boldsymbol{X}=\boldsymbol{X}_{0}+\boldsymbol{W},$ (2)

with

$\boldsymbol{X}_{0}=\sum_{i=1}^{r}\lambda_{i}\boldsymbol{v}_{i}\boldsymbol{v}_{i}^{T}\
(\text{symmetric case}),\ \text{ or }\
\boldsymbol{X}_{0}=\sum_{i=1}^{r}\sigma_{i}\boldsymbol{u}_{i}\boldsymbol{v}_{i}^{T}\
(\text{non-symmetric case}).$ (3)

In the symmetric case, $\lambda_{1}\geq\lambda_{2}\geq\ldots\geq\lambda_{r}$
are the eigenvalues, $\\{\boldsymbol{v}_{i}\\}_{i\in[r]}$ are the
corresponding unit-norm eigenvectors. In the non-symmetric case,
$\sigma_{1}\geq\sigma_{2}\geq\ldots\geq\sigma_{r}>0$ are the singular values,
and $\\{\boldsymbol{u}_{i},\boldsymbol{v}_{i}\\}_{i\in[r]}$ the corresponding
unit-norm singular vectors. In both cases, $\boldsymbol{W}$ is a noise matrix
(assumed to be symmetric when $\boldsymbol{X}_{0}$ is symmetric). Letting
$n:=\min(n_{1},n_{2})$, we consider the regime where $n$ grows and the rank
$r$ of $\boldsymbol{X}_{0}$ is a constant that does not scale with $n$. In the
rest of the paper, we refer to the $\\{\boldsymbol{v}_{i}\\}$ in the symmetric
case and the $\\{\boldsymbol{u}_{i},\boldsymbol{v}_{i}\\}$ in the non-
symmetric case as _signal vectors_.

##### Sparsity constraint

We assume that the signal vectors are $k$-sparse, i.e., they each have at most
$k$ nonzero entries, where $k\ll n$. The theoretical results are derived under
the assumption that the locations of the nonzeros (i.e., supports) are
uniformly random. The ‘with high probability’ statements throughout the paper
are with respect to the randomness in the support locations.

The goal is to design a sketching operator $\mathcal{A}$ and an efficient
algorithm to recover the signal vectors from the sketch $\boldsymbol{y}$. In
particular, we want the sample complexity and the running time of the
algorithm to both depend only on the sparsity level $k$, and not on the
ambient dimension $n$. (The sample complexity is the sketch size $m$.) We note
that the signal vectors can be recovered only up to a sign ambiguity as the
pair $(\boldsymbol{u}_{i},\boldsymbol{v}_{i})$ cannot be distinguished from
$(-\boldsymbol{u}_{i},-\boldsymbol{v}_{i})$. Furthermore, if there are
repeated eigenvalues or singular values, the corresponding signal vectors are
recovered up to any rotations within the subspace they span.

Simultaneously sparse and low-rank matrices arise in applications such as
sparse PCA [5, 6], sparse bilinear inverse problems [7], community detection
[8] and biclustering [9, 10]. In particular, the adjacency matrices of graphs
with community structures such as social networks and protein interaction
datasets are symmetric, sparse and block-diagonal in appropriate bases [11].
In biclustering, the sample-variable associations in high dimensional data,
like gene expression microarray datasets, often correspond to sparse non-
overlapping submatrices (or biclusters) in the data matrix.

###  1.1 Main contributions

Our sketching operator is defined via a combination of a sparse parity check
matrix and a partial DFT matrix. This operator was proposed for compressed
sensing of (vector) signals in [1, 2]. Our main contribution is the design and
analysis of a two-stage iterative algorithm to recover the signal vectors from
the sketch $\boldsymbol{y}$. The two stages use the sparsity and low-rank
properties of the signal matrix $\boldsymbol{X}_{0}$ to iteratively identify
and solve equations with a single unknown.

##### Noiseless setting

The sketching scheme and the recovery algorithm are described in Section 2 for
the noiseless case, that is, when $\boldsymbol{W}=0$ in (2). In Section 3, we
provide theoretical guarantees on the performance of the scheme. For
sufficiently large $k$, Theorem 1 shows that the scheme has the following
features when $\boldsymbol{X}_{0}$ is symmetric in (3):

  1. 1)

When the supports of the signal vectors $\\{\boldsymbol{v}_{i}\\}$ are
disjoint, fix any $\delta\in(0,1)$, with probability at least
$1\,-2r\exp(-\frac{1}{30}k^{1-\delta})$, the algorithm recovers the signal
vectors with sample complexity $m=3rk^{2}/(\delta\ln k)$ and running time
$\mathcal{O}(rk^{2}/\ln k)$. We note that with $k=\mathcal{O}(n^{\nu})$ for
$\nu\in(0,\frac{1}{2})$, and the support locations of each
$\boldsymbol{v}_{i}$ uniformly random, the supports of the
$\\{\boldsymbol{v}_{i}\\}$ will be disjoint with probability at least
$1-3r^{2}k^{2}/n$ .

  2. 2)

For the general case, where the supports of the signal vectors are not
disjoint, with probability at least $1-\mathcal{O}(k^{-2})$, the algorithm
recovers the signal vectors with sample complexity $m=2rk^{2}$ and running
time $\mathcal{O}\left((rk)^{3}\right)$.

Theorem 2 provides similar guarantees when $\boldsymbol{X}_{0}$ is non-
symmetric as defined in (3). The numerical simulations in Section 3.4 validate
the theoretical guarantees.

_Near-optimal sample complexity_ : In the symmetric case, with $r$ orthonormal
eigenvectors $\\{\boldsymbol{v}_{i}\\}$, each $k$-sparse, even if the
locations of the nonzeros in $\\{\boldsymbol{v}_{i}\\}$ are known, the degrees
of freedom in the eigendecomposition of $\boldsymbol{X}_{0}$ is close to $rk$.
In the non-symmetric case, the degrees of freedom in the singular value
decomposition of $\boldsymbol{X}_{0}$ is close to $2r(k-\frac{1}{2})$.
Therefore the sample complexity of the proposed scheme is larger than the
degrees of freedom by a factor of at most $O(k)$. Moreover, neither the sample
complexity or the running time depends on the ambient dimension $n$.

##### Noisy setting

The sketching scheme and recovery algorithm are extended to the noisy case in
Section 4 (i.e., when $\boldsymbol{W}\neq 0$ in (2)). Here, the sample
complexity is $\mathcal{O}(rk^{2}\log n/\log k)$, a factor of
$\mathcal{O}(\log n)$ greater than that in the noiseless case. This is because
additional sketches are needed to reliably identify the locations of nonzero
matrix entries in the presence of noise. This also increases the running time
of the recovery algorithm to $\mathcal{O}(\max\\{n^{2}\log(n/k),(rk)^{3}\\})$.
We do not provide theoretical performance guarantees in the noisy setting, but
show via numerical simulations that the scheme is robust to moderate levels of
noise.

###  1.2 Related work

The sketching operator we use was proposed in [2, 1] for compressed sensing of
vectors. Variations of the operator have been used for numerous applications
including sparse DFT [12, 13], sparse Walsh-Hadamard Transform [14, 15, 16],
compressive phase retrieval [17, 18], sparse covariance estimation [19],
sparse polynomial learning and graph sketching [20], and learning mixtures of
sparse linear regressions [21].

The two-stage recovery algorithm we propose is analogous to the peeling
decoder for Low Density Parity Check (LDPC) codes over an erasure channel [22,
23], and the first stage of the algorithm is similar to the one used for
compressed sensing in [2, 1]. However, a key difference from these works is
that our sketching matrix has row weights that scale with the sparsity level
$k$ and therefore, the existing peeling decoder analysis based on density
evolution cannot be applied. We characterize the performance of the algorithm
by obtaining nonasymptotic probability bounds on the number of unknown nonzero
entries after each stage. For the first stage, this is done by establishing
negative association (defined in Section 5.1) between the right node degrees
in the associated bipartite graph, and using concentration inequalities for
negatively associated random variables [24]. For the second stage, we model
the algorithm as a random graph process on another bipartite graph, and we
analyze an alternative random graph process whose evolution is easier to
characterize than the original one.

_Other related work_ : Recovering low-rank matrices from linear measurements
(without sparsity constraints) has been widely studied in the past decade; see
[25] for an overview. A key result in this area is that if the linear
measurement operator satisfies a matrix restricted isometry property, then the
low-rank matrix can be recovered via nuclear-norm minimization [26, 27]. For
matrices that are simultaneously low-rank and sparse (with $k=o(n)$), such
optimization-based approaches are highly sub-optimal with respect to sample
complexity.

Several authors have investigated sketching schemes for sparse matrices [28,
29], and simultaneously sparse and low-rank matrices [7, 30, 31, 32].
Moreover, [33, 34] studied the related problem of recovering a sparse low-rank
covariance matrix from rank-1 sketches of a sample covariance (the sketching
matrices $\boldsymbol{A}_{1},\ldots,\boldsymbol{A}_{m}$ are rank-1). The
recovery algorithms in all these works (usually based on convex programming)
have running time that is polynomial in $n$; the sample complexity also
depends at least logarithmically on $n$.

The problem of estimating sparse eigenvectors has also been widely studied in
the context of sparse PCA [5, 35, 36, 37]. In these works, the principal
eigenvector of the population covariance matrix is assumed to be sparse. The
goal is to recover the principal eigenvector from the sample covariance
matrix. This is distinct from the sketching problem considered here.

Some inverse problems such as phase retrieval [38] and blind deconvolution
[39] can be viewed as recovering a rank-1 signal matrix from linear
measurements. The structure of the linear operator in these problems is
constrained by the applications, whereas the operator in our setting can be
designed flexibly.

###  1.3 Notation

We write $[n]$ for the set of integers $\left\\{1,2,\dots,n\right\\}$,
$\mathbb{N}_{0}$ for $\\{0,1,2,\ldots\\}$, and use i to denote $\sqrt{-1}$.
For a length-$n$ vector $\boldsymbol{u}$, we denote the set of nonzero
locations by ${\sf supp}(\boldsymbol{u}):=\\{\ell\in[n]:u[\ell]\neq 0\\}$. The
notation $x_{n}=o(n^{c})$ is used to denote a positive number $x_{n}$ such
that $\frac{x_{n}}{n^{c}}\to 0$ as $n\to\infty$.

The Bernoulli distribution with parameter $p\in[0,1]$ is denoted by
$\text{Bern}\left(p\right)$, and the Binomial distribution with parameters
$n\in\mathbb{N}$ and $p\in[0,1]$ is denoted by $\text{Bin}\left(n,p\right)$.
$\text{Pois}\left(\lambda\right)$ denotes a Poisson distribution with
parameter $\lambda>0$.

##  2 Baseline sketching scheme and recovery algorithm

We first focus on the noiseless symmetric case in Sections 2.1–2.3, and then
extend the scheme to non-symmetric matrices in Section 2.4. The noisy setting
is discussed in Section 4.

###  2.1 Sketching scheme

Consider the noiseless symmetric sparse, low-rank matrix
$\boldsymbol{X}=\sum_{i=1}^{r}\lambda_{i}\boldsymbol{v}_{i}\boldsymbol{v}_{i}^{T}$
in (2). Let $\tilde{n}=\binom{n}{2}+n$, and let
$\boldsymbol{x}\in{\mathbb{R}}^{\tilde{n}}$ be the vectorized upper-triangular
part of $\boldsymbol{X}$. The sketch is
$\boldsymbol{y}=\boldsymbol{B}\boldsymbol{x}$, with the sketching matrix
$\boldsymbol{B}\in\mathbb{C}^{2R\times\tilde{n}}$ described below. The sample
complexity is $m=2R$ with the exact value of $R=\mathcal{O}(rk^{2}/\ln k)$
given in Theorem 1.

The sketching matrix $\boldsymbol{B}\in\mathbb{C}^{2R\times\tilde{n}}$ is
constructed as in [1] by taking the column-wise Kronecker product of two
matrices: a sparse parity check matrix

$\boldsymbol{H}:=\begin{bmatrix}\boldsymbol{h}_{1}&\boldsymbol{h}_{2}&\dots&\boldsymbol{h}_{\tilde{n}}\end{bmatrix}\in\\{0,1\\}^{R\times\tilde{n}},$
(4)

and a matrix $\boldsymbol{S}$ consisting of the first two rows of an
$\tilde{n}$-point DFT matrix:

$\boldsymbol{S}:=\begin{bmatrix}1&1&1&\ldots&1\\\
1&W&W^{2}&\ldots&W^{\tilde{n}-1}\end{bmatrix}\quad\text{where
}W=\exp\left(\frac{2\pi\textrm{i}}{\tilde{n}}\right).$ (5)

Then,

$\small{\boldsymbol{B}:=\begin{bmatrix}\boldsymbol{h}_{1}\otimes\begin{bmatrix}1\\\
1\end{bmatrix}&\boldsymbol{h}_{2}\otimes\begin{bmatrix}1\\\
W\end{bmatrix}&\dots&\boldsymbol{h}_{\tilde{n}}\otimes\begin{bmatrix}1\\\
W^{\tilde{n}-1}\end{bmatrix}\end{bmatrix}.}$ (6)

For column vectors $\boldsymbol{a},\boldsymbol{b}$ with lengths $n_{a},n_{b}$,
we recall that the Kronecker product $\boldsymbol{a}\otimes\boldsymbol{b}$ is
the length-$(n_{a}+n_{b})$ vector
$[a_{1}\boldsymbol{b}^{T},a_{2}\boldsymbol{b}^{T},\dots
a_{n_{a}}\boldsymbol{b}^{T}]^{T}$. As an example, let $\tilde{n}=6$, $R=4$ and

$\displaystyle\boldsymbol{H}=\begin{bmatrix}1&1&1&0&1&1\\\ 0&0&0&1&1&1\\\
1&0&0&0&0&0\\\ 0&1&1&1&0&0\end{bmatrix},\ \
\boldsymbol{S}=\begin{bmatrix}1&1&\ldots&1\\\ 1&W&\ldots&W^{5}\end{bmatrix}.$
(7)

Then the sketching matrix $\boldsymbol{B}$ is

$\displaystyle\boldsymbol{B}=\begin{bmatrix}1&1&1&0&1&1\\\
1&W&W^{2}&0&W^{4}&W^{5}\\\ 0&0&0&1&1&1\\\ 0&0&0&W^{3}&W^{4}&W^{5}\\\
1&0&0&0&0&0\\\ 1&0&0&0&0&0\\\ 0&1&1&1&0&0\\\
0&W&W^{2}&W^{3}&0&0\end{bmatrix}.$ (8)

We choose $\boldsymbol{H}$ to be column-regular, with each column containing
$d\geq 2$ ones at locations chosen uniformly at random; the example in (7)
uses $d=2$. The sparse matrix $\boldsymbol{H}$ determines which nonzero
entries of $\boldsymbol{x}$ contribute to each entry of the sketch
$\boldsymbol{y}$.

It is convenient to view the sketch
$\boldsymbol{y}=\boldsymbol{B}\boldsymbol{x}\in{\mathbb{C}}^{2R}$ as having
$R$ pairs of entries. We denote the $j$-th pair
$\begin{bmatrix}y_{(2j-1)},y_{(2j)}]\end{bmatrix}^{T}$ by
$\boldsymbol{y}_{j}\in{\mathbb{C}}^{2}$ and call it the $j$-th sketch _bin_.
Let $\boldsymbol{B}_{j}\in{\mathbb{C}}^{2\times\tilde{n}}$ denote the matrix
containing the $j$-th pair of rows in $\boldsymbol{B}$, for $j\in[R]$. We
observe that the $j$-th bin can be expressed as

$\boldsymbol{y}_{j}=\boldsymbol{B}_{j}\boldsymbol{x}=\boldsymbol{S}\boldsymbol{x}^{*}_{j},$
(9)

where the entries of $\boldsymbol{x}^{*}_{j}\in{\mathbb{R}}^{\tilde{n}}$ are

$(x^{*}_{j})_{\ell}=\left\\{\begin{array}[]{ll}x_{\ell}&\ell\in\mathcal{N}(j)\\\
0&\ell\notin\mathcal{N}(j)\end{array}\right.\quad\text{for
}\ell\in[\tilde{n}].$ (10)

Here, $\mathcal{N}(j):=\\{\ell\in[\tilde{n}]:\,H_{j,\ell}=1\\}$ contains the
locations of the ones in the $j$-th row of $\boldsymbol{H}$. In words, (9)
shows that the bin $\boldsymbol{y}_{j}$ is a linear combination of the
$\boldsymbol{x}$ entries that correspond to the ones in the $j$-th row of
$\boldsymbol{H}$. Specifically, $\boldsymbol{y}_{j}$ is the sum of these
$\boldsymbol{x}$ entries, weighted by their respective columns in
$\boldsymbol{S}$. Note that (9) can be recast into two equations in the form
of (1), one for each component in the bin $\boldsymbol{y}_{j}$. Moreover, when
only one nonzero $\boldsymbol{x}$ entry contributes to $\boldsymbol{y}_{j}$,
the DFT structure of $\boldsymbol{S}$ enables exact recovery of that nonzero
from $\boldsymbol{y}_{j}$.

###  2.2 Recovery algorithm for rank-1 symmetric matrices

For ease of exposition, we first describe the recovery algorithm for rank-1
matrices (i.e., $\boldsymbol{X}=\lambda\boldsymbol{v}\boldsymbol{v}^{T}$), and
then explain how to extend the algorithm to tackle rank-$r$ matrices. The
algorithm has two stages, which we refer to as stage A and stage B.

  * •

In stage A, by exploiting the sparsity of $\boldsymbol{X}$, the algorithm
iteratively recovers as many of the nonzero entries in $\boldsymbol{X}$ as
possible from the sketch $\boldsymbol{y}$. This stage is similar to the
compressed sensing recovery algorithm used in [2, 1]. For sufficiently large
$k$, stage A recovers at least one nonzero diagonal entry and a fraction
$k^{-\delta}$ of the nonzero above-diagonal entries in $\boldsymbol{X}$ with
high probability, for a constant $\delta\in(0,1)$.

  * •

In stage B, by exploiting the rank-1 structure of $\boldsymbol{X}$, the
algorithm iteratively recovers the nonzero entries in $\boldsymbol{v}$ from
the partially recovered matrix $\boldsymbol{X}$. Theorem 1 shows that with a
sketch of size $3k^{2}/\ln k$, the algorithm recovers $\boldsymbol{v}$ by the
end of stage B with high probability, for sufficiently large $k$.

We note that the recovery algorithm does not require knowledge of the sparsity
$k$.

####  2.2.1 Stage A of the algorithm

Following the terminology in [1], we classify each of the $R$ bins of
$\boldsymbol{y}$ based on how many nonzero entries of $\boldsymbol{X}$
contribute to the pair of linear combinations defining the bin. Specifically,
a bin $\boldsymbol{y}_{j}$ is referred to as a _zeroton, singleton_ , or
_multiton_ respectively if it involves zero, one, or more than one nonzero
entries of $\boldsymbol{X}$, that is, when $\left\lvert{\sf
supp}(\boldsymbol{x})\cap\mathcal{N}(j)\right\rvert=0,1$ or $>1$ with
$\mathcal{N}(j)$ as defined in (10). This is illustrated by the following
example with the sketching matrix $\boldsymbol{B}$ in (8):

$\displaystyle\boldsymbol{y}\in{\mathbb{C}}^{2R}$
$\displaystyle\quad\quad\quad\quad\boldsymbol{B}\in{\mathbb{C}}^{2R\times\tilde{n}}\quad\quad\quad\boldsymbol{x}\in{\mathbb{R}}^{\tilde{n}}$
$\displaystyle\begin{matrix}\boldsymbol{y}_{1}\\\ \\\ \boldsymbol{y}_{2}\\\
\\\ \boldsymbol{y}_{3}\\\ \\\
\boldsymbol{y}_{4}\end{matrix}\begin{bmatrix}{\color[rgb]{1,0.5,0.5}\definecolor[named]{pgfstrokecolor}{rgb}{1,0.5,0.5}\rule{2.27626pt}{7.11317pt}}{\color[rgb]{0.5,0.5,1}\definecolor[named]{pgfstrokecolor}{rgb}{0.5,0.5,1}\rule{2.27626pt}{7.11317pt}}{\color[rgb]{0.5,1,1}\definecolor[named]{pgfstrokecolor}{rgb}{0.5,1,1}\pgfsys@color@cmyk@stroke{0.5}{0}{0}{0}\pgfsys@color<EMAIL_ADDRESS>{\color[rgb]{1,0.5,0.5}\definecolor[named]{pgfstrokecolor}{rgb}{1,0.5,0.5}\rule{2.27626pt}{7.11317pt}}{\color[rgb]{0.5,0.5,1}\definecolor[named]{pgfstrokecolor}{rgb}{0.5,0.5,1}\rule{2.27626pt}{7.11317pt}}{\color[rgb]{0.5,1,1}\definecolor[named]{pgfstrokecolor}{rgb}{0.5,1,1}\pgfsys@color@cmyk@stroke{0.5}{0}{0}{0}\pgfsys@color<EMAIL_ADDRESS>{\color[rgb]{0.5,1,1}\definecolor[named]{pgfstrokecolor}{rgb}{0.5,1,1}\pgfsys@color@cmyk@stroke{0.5}{0}{0}{0}\pgfsys@color<EMAIL_ADDRESS>{\color[rgb]{0.5,1,1}\definecolor[named]{pgfstrokecolor}{rgb}{0.5,1,1}\pgfsys@color@cmyk@stroke{0.5}{0}{0}{0}\pgfsys@color<EMAIL_ADDRESS>\tiny\square\\\ \tiny\square\\\
{\color[rgb]{1,0.5,0.5}\definecolor[named]{pgfstrokecolor}{rgb}{1,0.5,0.5}\rule{3.55658pt}{7.11317pt}}{\color[rgb]{0.5,0.5,1}\definecolor[named]{pgfstrokecolor}{rgb}{0.5,0.5,1}\rule{3.55658pt}{7.11317pt}}\\\
{\color[rgb]{1,0.5,0.5}\definecolor[named]{pgfstrokecolor}{rgb}{1,0.5,0.5}\rule{3.55658pt}{7.11317pt}}{\color[rgb]{0.5,0.5,1}\definecolor[named]{pgfstrokecolor}{rgb}{0.5,0.5,1}\rule{3.55658pt}{7.11317pt}}\end{bmatrix}=$
$\displaystyle\begin{bmatrix}{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&\tiny\square&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color<EMAIL_ADDRESS>{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&\tiny\square&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color<EMAIL_ADDRESS>\tiny\square&\tiny\square&\tiny\square&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color<EMAIL_ADDRESS>\tiny\square&\tiny\square&\tiny\square&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color<EMAIL_ADDRESS>{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color<EMAIL_ADDRESS>{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color<EMAIL_ADDRESS>\tiny\square&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color<EMAIL_ADDRESS>\tiny\square&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\rule{7.11317pt}{7.11317pt}}&{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color<EMAIL_ADDRESS>{\color[rgb]{1,0.5,0.5}\definecolor[named]{pgfstrokecolor}{rgb}{1,0.5,0.5}\rule{7.11317pt}{7.11317pt}}\\\
{\color[rgb]{0.5,0.5,1}\definecolor[named]{pgfstrokecolor}{rgb}{0.5,0.5,1}\rule{7.11317pt}{7.11317pt}}\\\
\tiny\square\\\ \tiny\square\\\
{\color[rgb]{0.5,1,1}\definecolor[named]{pgfstrokecolor}{rgb}{0.5,1,1}\pgfsys@color@cmyk@stroke{0.5}{0}{0}{0}\pgfsys@color<EMAIL_ADDRESS>X_{11}\\\ x_{2}\equiv X_{12}\\\ x_{3}\equiv X_{13}\\\ x_{4}\equiv X_{22}\\\
x_{5}\equiv X_{23}\\\ x_{6}\equiv X_{33}\end{matrix}\ .$

The white squares ($\square$) represent zeros, and the black squares
($\blacksquare$) represent the DFT coefficients in $\boldsymbol{B}$. The
coloured squares in $\boldsymbol{x}$ represent the nonzero entries. The sketch
$\boldsymbol{y}$ consists of $R=4$ bins, of which $\boldsymbol{y}_{3}$
($\square$) is a zeroton, $\boldsymbol{y}_{2}$ ( ) is a singleton, and
$\boldsymbol{y}_{1}$ ( ) and $\boldsymbol{y}_{4}$ ( ) are multitons.

The first step is to classify each bin $\boldsymbol{y}_{j}$ into one of the
three types. For brevity, we denote the two components of $\boldsymbol{y}_{j}$
as $y_{(2j-1)}$ and $y_{(2j)}$. If both $y_{(2j-1)}$ and $y_{(2j)}$ are zero,
then $\boldsymbol{y}_{j}$ is declared a zeroton. This is because the DFT
coefficients in $\boldsymbol{B}$ ensure that $y_{(2j)}$ can be zero only when
the $\boldsymbol{x}$ entries involved in the linear combination are all zero,
i.e., ${\sf supp}(\boldsymbol{x})\cap\mathcal{N}(j)=\emptyset$. Next, to
identify singletons, (9) and (10) indicate that a singleton takes the form

$\boldsymbol{y}_{j}=\begin{bmatrix}y_{(2j-1)}\\\
y_{(2j)}\end{bmatrix}=x_{\ell}\begin{bmatrix}1\\\ W^{\ell-1}\end{bmatrix},$
(11)

where $x_{\ell}$ is the sole nonzero entry contributing to the bin
$\boldsymbol{y}_{j}$, i.e., ${\sf
supp}(\boldsymbol{x})\cap\mathcal{N}(j)=\ell$. While $y_{(2j-1)}$ directly
captures the value of $x_{\ell}$, the components $y_{(2j-1)}$ and $y_{(2j)}$
together capture the location $\ell$. Therefore, to identify a singleton, the
algorithm computes

$\hat{\ell}=\frac{\arg\left\\{y_{(2j)}/y_{(2j-1)}\right\\}}{2\pi/\tilde{n}}+1,\quad\hat{x}_{\hat{\ell}}=y_{(2j-1)}.$
(12)

If $\left\lvert y_{(2j-1)}\right\rvert=\left\lvert y_{(2j)}\right\rvert$ and
the estimated index $\hat{\ell}$ takes an integer value in
$\left[\tilde{n}\right]$, the algorithm declares $\boldsymbol{y}_{j}$ a
singleton and the $\hat{\ell}$-th entry of $\boldsymbol{x}$ as
$\hat{x}_{\hat{\ell}}$. The algorithm declares $\boldsymbol{y}_{j}$ a multiton
if it is found to be neither a zeroton nor a singleton.

Let $t\in\mathbb{N}_{0}$ denote the iteration number and let $\mathcal{S}(t)$
be the set of singletons at time $t$. At $t=0$, the algorithm identifies the
initial set of singletons $\mathcal{S}(0)$ among the $R$ bins. Then at each
$t\geq 1$, the algorithm picks a singleton uniformly at random from
$\mathcal{S}(t-1)$, and recovers the support of $\boldsymbol{x}$ underlying
the singleton using the index-value pair $(\hat{\ell},\hat{x}_{\hat{\ell}})$
in (12). The algorithm then subtracts (or ‘peels off’) the contribution of the
$\hat{\ell}$-th entry of $\boldsymbol{x}$ from each of the $d$ bins the entry
is involved in, re-categorizes these bins, and updates the set of singletons
to include any new singletons created by the peeling of the $\hat{\ell}$-th
entry. The updated set is $\mathcal{S}(t)$. The algorithm continues until
singletons run out, i.e., when $\mathcal{S}(t)=\emptyset$.

It is useful to visualize the peeling algorithm using a bipartite graph
constructed from the parity check matrix
$\boldsymbol{H}\in\\{0,1\\}^{R\times\tilde{n}}$ (which was used to define the
sketching matrix $\boldsymbol{B}$). As shown in Fig. 1(a), the
$\tilde{k}=\binom{k}{2}+k$ left nodes correspond to the unknown nonzero
entries in $\boldsymbol{X}$ at $t=0$, and the $R$ right nodes correspond to
the bins. Each left node connects to $d$ distinct right nodes. The left nodes
connected to each right node (bin) are those that appear in the linear
constraints defining the bin. We call Fig. 1(a) a ‘pruned’ graph as we have
removed the left nodes detected to be zeros (via zeroton bins). As illustrated
in Figs. 1(b) and 1(c), the algorithm first recovers and peels off ${X}_{24}$
from the singleton $\boldsymbol{y}_{7}$ at $t=1$. This creates a singleton
$\boldsymbol{y}_{6}$ from which ${X}_{44}$ is recovered and peeled off at
$t=2$. The algorithm terminates at $t=2$ as there are no more singletons.

$X_{11}$ $X_{12}$ $X_{13}$ $X_{14}$ $X_{22}$ $X_{23}$ $X_{24}$ $X_{33}$
$X_{34}$ $X_{44}$ $\boldsymbol{y}_{1}$ $\boldsymbol{y}_{2}$
$\boldsymbol{y}_{3}$ $\boldsymbol{y}_{4}$ $\boldsymbol{y}_{5}$
$\boldsymbol{y}_{6}$ $\boldsymbol{y}_{7}$ $\boldsymbol{y}_{8}$
$\binom{k}{2}+k=10$ distinct nonzero entries of $\boldsymbol{X}$ $R=8$ sketch
bins (a) $t=0$.

$X_{11}$ $X_{12}$ $X_{13}$ $X_{14}$ $X_{22}$ $X_{23}$ $X_{24}$ $X_{33}$
$X_{34}$ $X_{44}$ $\boldsymbol{y}_{1}$ $\boldsymbol{y}_{2}$
$\boldsymbol{y}_{3}$ $\boldsymbol{y}_{4}$ $\boldsymbol{y}_{5}$
$\boldsymbol{y}_{6}$ $\boldsymbol{y}_{7}$ $\boldsymbol{y}_{8}$ (b) $t=1$.

$X_{11}$ $X_{12}$ $X_{13}$ $X_{14}$ $X_{22}$ $X_{23}$ $X_{24}$ $X_{33}$
$X_{34}$ $X_{44}$ $\boldsymbol{y}_{1}$ $\boldsymbol{y}_{2}$
$\boldsymbol{y}_{3}$ $\boldsymbol{y}_{4}$ $\boldsymbol{y}_{5}$
$\boldsymbol{y}_{6}$ $\boldsymbol{y}_{7}$ $\boldsymbol{y}_{8}$ (c) $t=2$.

Figure 1: (a): Pruned bipartite graph for stage A corresponding to
$\boldsymbol{H}\in\\{0,1\\}^{R\times\tilde{n}}$. Here,
$\boldsymbol{X}=\lambda\boldsymbol{v}\boldsymbol{v}^{T}$ with $v_{i}\neq 0$
for $i\in[4]$ and $v_{i}=0$ for $5\leq i\leq n$, and $R=8$. (b)–(c): The graph
process that models the recovery of the nonzero entries of $\boldsymbol{X}$.
The faded nodes and edges are those that have been peeled off.

####  2.2.2 Stage B of the algorithm

Recalling that $\boldsymbol{X}=\lambda\boldsymbol{v}\boldsymbol{v}^{T}$, we
write
$\tilde{\boldsymbol{v}}:=\sqrt{\left\lvert\lambda\right\rvert}\boldsymbol{v}$
for the unsigned, unnormalized eigenvector. Stage B of the algorithm uses the
rank-1 structure of $\boldsymbol{X}$ to recover the nonzeros in
$\tilde{\boldsymbol{v}}$ from the nonzero entries of $\boldsymbol{X}$ that
were recovered in stage A. When fully recovered, $\tilde{\boldsymbol{v}}$ is
normalized to give $\boldsymbol{v}$. Since the $j$-th diagonal entry of
$\boldsymbol{X}$ is $\lambda v_{j}^{2}$, the sign of $\lambda$ can be
determined from any nonzero diagonal entry recovered in stage A. We now
describe stage B for $\lambda>0$, and then explain the small adjustment needed
for $\lambda<0$.

Since $\boldsymbol{X}$ is rank-1, its $\binom{k}{2}$ nonzero above-diagonal
entries are pairwise products of the form $\tilde{v}_{i}\tilde{v}_{j}$. The
proof of Theorem 1 shows that in stage A, with high probability, at least a
fraction $k^{-\delta}$ of these pairwise products are recovered, for a
constant $\delta\in(0,1)$. (See Lemma 3.1). Stage B of the algorithm can be
visualized using another bipartite graph shown in Fig. 2(a). The $k$ left
nodes represent the unknown nonzero entries in $\tilde{\boldsymbol{v}}$, and
the right nodes represent the nonzero pairwise products in $\boldsymbol{X}$
recovered in stage A. The right nodes are the product constraints that the
left nodes satisfy. This is a pruned graph because any zero entries in
$\tilde{\boldsymbol{v}}$ or $\boldsymbol{X}$ have been excluded from the
graph.

At $t=0$, the algorithm picks one of the nonzero diagonal entries recovered in
stage A, say $X_{jj}$, and declares $\sqrt{X_{jj}}$ as the value of the left
node $\tilde{v}_{j}$. The algorithm then removes (peels off) the contribution
of $\tilde{v}_{j}$ from the right nodes that $\tilde{v}_{j}$ connects to by
dividing these right nodes by $\sqrt{X_{jj}}$. Fig. 2(b) shows this initial
step with $\tilde{v}_{4}$ being a left node whose squared value was recovered
in stage A and is peeled off at $t=0$. Once the first left node is peeled off,
all of its connecting right nodes reduce from degree-2 nodes to degree-1. For
$t\geq 1$, the algorithm proceeds similarly to stage A, noting that in stage
B, each right node represents the product of its two connecting left nodes
rather than a linear combination. At each step, a degree-1 right node is
picked uniformly at random from the available ones and its connecting left
node recovered and peeled off to update the set of degree-1 right nodes. The
algorithm continues until there are no more degree-1 right nodes. All the
unrecovered entries of $\tilde{\boldsymbol{v}}$ are set to zero. In Fig. 2(c)
and 2(d), the left nodes $\tilde{v}_{1}$ and $\tilde{v}_{5}$ are recovered via
$\tilde{v}_{1}\tilde{v}_{4}$ and $\tilde{v}_{4}\tilde{v}_{5}$, respectively,
and peeled off. Moreover, $\tilde{v}_{2}$ and $\tilde{v}_{3}$ can be recovered
via $\tilde{v}_{2}\tilde{v}_{5}$ and $\tilde{v}_{2}\tilde{v}_{3}$ before
degree-1 right nodes run out.

##### The $\lambda<0$ case

We flip the sign of each right node (pairwise product) before running stage B
of the algorithm as above to recover
$\tilde{\boldsymbol{v}}=\sqrt{\left\lvert\lambda\right\rvert}\boldsymbol{v}$.

$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{v}_{4}$$\tilde{v}_{5}$$\tilde{v}_{1}\tilde{v}_{4}$$\tilde{v}_{2}\tilde{v}_{3}$$\tilde{v}_{2}\tilde{v}_{5}$$\tilde{v}_{4}\tilde{v}_{5}$
$k=5$ nonzero entries of $\tilde{\boldsymbol{v}}$ $4$ nonzero pairwise
products recovered in stage A Degree of each right node = 2 (a) Prior to
$t=0$.

$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{v}_{4}$$\tilde{v}_{5}$$\tilde{v}_{1}\tilde{v}_{4}$$\tilde{v}_{2}\tilde{v}_{3}$$\tilde{v}_{2}\tilde{v}_{5}$$\tilde{v}_{4}\tilde{v}_{5}$
$k=5$ supports of $\tilde{\boldsymbol{v}}$ Right nodes connected to
$\tilde{v}_{4}$ have degree 1; others have degree 2 (b) $t=0$.

$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{v}_{4}$$\tilde{v}_{5}$$\tilde{v}_{1}\tilde{v}_{4}$$\tilde{v}_{2}\tilde{v}_{3}$$\tilde{v}_{2}\tilde{v}_{5}$$\tilde{v}_{4}\tilde{v}_{5}$
$k=5$ supports of $\tilde{\boldsymbol{v}}$ Right nodes connected to
$\tilde{v}_{4}$ have degree 1; others have degree 2 (c) $t=1$.

$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{v}_{4}$$\tilde{v}_{5}$$\tilde{v}_{1}\tilde{v}_{4}$$\tilde{v}_{2}\tilde{v}_{3}$$\tilde{v}_{2}\tilde{v}_{5}$$\tilde{v}_{4}\tilde{v}_{5}$
$k=5$ supports of $\tilde{\boldsymbol{v}}$ Right nodes connected to
$\tilde{v}_{4}$ have degree 1; others have degree 2 (d) $t=2$.

Figure 2: (a): Pruned bipartite graph for stage B with
$\boldsymbol{X}=\tilde{\boldsymbol{v}}\tilde{\boldsymbol{v}}^{T}$, where
$\tilde{v}_{i}\neq 0$ for $i\in[5]$ and $\tilde{v}_{i}=0$ for $6\leq i\leq n$.
(b)–(d): The graph process that models the recovery of the nonzero entries in
$\tilde{\boldsymbol{v}}$ from the nonzero matrix entries recovered in stage A.
The faded nodes and edges are those that have been peeled off. At $t=0$, the
left node $\tilde{v}_{4}$ is peeled off based on the recovered nonzero
diagonal entry $X_{44}$.

###  2.3 Recovery algorithm for rank-$r$ symmetric matrices

##### Disjoint suppports

When the eigenvectors $\\{\boldsymbol{v}_{i}\\}$ are known to have disjoint
supports, the same two-stage algorithm can be used. The only difference is
that the graphical representation for stage B (Fig. 2(a)) will now consist of
$r$ disjoint bipartite graphs. For each of these graphs, stage B of the
algorithm is initialized using a nonzero diagonal entry recovered in stage A.

##### General case

When the supports of the $r$ eigenvectors have overlaps, we cannot use stage B
of the algorithm as it depends on the entries of $\boldsymbol{X}$ being
pairwise products of the eigenvector entries. We therefore take a large enough
sketch to identify all the nonzero entries in $\boldsymbol{X}$ in stage A.
This is equivalent to compressed sensing recovery of the complete vectorized
matrix using the scheme of [1, 2]. We then perform an eigendecomposition on
the recovered submatrix of nonzero entries to obtain the nonzero entries of
$\\{\boldsymbol{v}_{1},\ldots,\boldsymbol{v}_{r}\\}$. This submatrix has size
at most $rk\times rk$, hence the complexity of its eigendecomposition is
$\mathcal{O}((rk)^{3})$, which does not depend on the ambient dimension $n$.
Recovering the entire matrix $\boldsymbol{X}$ in stage A requires a sample
complexity $m=\mathcal{O}(rk^{2})$, which is larger by a factor of $\ln k$
than that needed in the case of disjoint supports.

###  2.4 Recovery of non-symmetric matrices

Consider the non-symmetric case where $\boldsymbol{X}$ has the form
$\boldsymbol{X}=\sum_{i=1}^{r}\sigma_{i}\boldsymbol{u}_{i}\boldsymbol{v}_{i}^{T}$
where $\\{\boldsymbol{u}_{i},\boldsymbol{v}_{i}\\}$ each has at most $k$
nonzero entries. The sketching operator is the same as in the symmetric case,
except that here it acts on the entire matrix $\boldsymbol{X}$. Consider first
the case where the singular vectors $\\{\boldsymbol{u}_{i}\\}$ have disjoint
supports as do $\\{\boldsymbol{v}_{i}\\}$. Based on the sketch, stage A of the
algorithm can be used to recover a small fraction of the nonzero entries in
$\boldsymbol{X}$ like in the symmetric case. Using the recovered nonzero
matrix entries, stage B of the algorithm iteratively recovers the nonzeros in
$\\{\boldsymbol{u}_{i},\boldsymbol{v}_{i}\\}$. While stage B can still be
viewed as a peeling decoder, the associated bipartite graph has a different
structure from the symmetric case. The initialization step is also different.

We take the rank-1 case as an example, where
$\boldsymbol{X}=\sigma\boldsymbol{u}\boldsymbol{v}^{T}$. The peeling algorithm
recovers $\tilde{\boldsymbol{u}}$ and $\tilde{\boldsymbol{v}}$, which are
scaled versions of $\boldsymbol{u}$ and $\boldsymbol{v}$ such that
$\tilde{\boldsymbol{u}}\tilde{\boldsymbol{v}}^{T}=\boldsymbol{X}$; these
vectors are then normalized to obtain $\sigma,\boldsymbol{u}$ and
$\boldsymbol{v}$. Fig. 3(a) illustrates the pruned stage B graph in this
setting. The left nodes represent the nonzero entries in
$\tilde{\boldsymbol{u}}$ and $\tilde{\boldsymbol{v}}$, and the right nodes
represent the nonzero pairwise products recovered in stage A. Each right node
connects to one nonzero in $\tilde{\boldsymbol{u}}$ and one nonzero in
$\tilde{\boldsymbol{v}}$ and equals the product of these two nonzeros.

Since a non-symmetric matrix has no entries in squared form, the peeling
process is initiated by arbitrarily assigning the value 1 to one of the left
nodes who have at least one neighbouring right node. At $t=0$, this left node
is peeled off from the graph, and its connecting right node then have degree
1. In Fig. 3(b), the algorithm recovers $\tilde{u}_{1}$ as 1 and peels it off
at $t=0$. This turns the right nodes $\tilde{u}_{1}\tilde{v}_{1}$ and
$\tilde{u}_{1}\tilde{v}_{3}$ from degree-2 into degree-1.

For $t\geq 1$, like in the symmetric case, the algorithm peels one degree-1
right node at a time along with the left node connecting to the right node,
until no degree-1 right nodes remain in the graph. In Figs. 3(c) and 3(d), at
$t=1$ and $t=2$, $\tilde{v}_{1}$ and $\tilde{v}_{3}$ are recovered and peeled
off sequentially.

$\tilde{u}_{1}$$\tilde{u}_{2}$$\tilde{u}_{3}$$\tilde{u}_{4}$$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{u}_{1}\tilde{v}_{1}$$\tilde{u}_{4}\tilde{v}_{1}$$\tilde{u}_{3}\tilde{v}_{2}$$\tilde{u}_{1}\tilde{v}_{3}$$\tilde{u}_{2}\tilde{v}_{2}$
$4$ nonzero entries of $\tilde{\boldsymbol{u}}$ $3$ nonzero entries of
$\tilde{\boldsymbol{v}}$ $5$ nonzero pairwise products recovered in stage A
Degree of each right node = 2 (a) Prior to $t=0$.

$\tilde{u}_{1}$$\tilde{u}_{2}$$\tilde{u}_{3}$$\tilde{u}_{4}$$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{u}_{1}\tilde{v}_{1}$$\tilde{u}_{4}\tilde{v}_{1}$$\tilde{u}_{3}\tilde{v}_{2}$$\tilde{u}_{1}\tilde{v}_{3}$$\tilde{u}_{2}\tilde{v}_{2}$
Right nodes connected to $\tilde{u}_{1}$ have degree 1; others have degree 2
(b) $t=0$.

$\tilde{u}_{1}$$\tilde{u}_{2}$$\tilde{u}_{3}$$\tilde{u}_{4}$$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{u}_{1}\tilde{v}_{1}$$\tilde{u}_{4}\tilde{v}_{1}$$\tilde{u}_{3}\tilde{v}_{2}$$\tilde{u}_{1}\tilde{v}_{3}$$\tilde{u}_{2}\tilde{v}_{2}$
Right nodes connected to $\tilde{u}_{1}$ have degree 1; others have degree 2
(c) $t=1$.

$\tilde{u}_{1}$$\tilde{u}_{2}$$\tilde{u}_{3}$$\tilde{u}_{4}$$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{u}_{1}\tilde{v}_{1}$$\tilde{u}_{4}\tilde{v}_{1}$$\tilde{u}_{3}\tilde{v}_{2}$$\tilde{u}_{1}\tilde{v}_{3}$$\tilde{u}_{2}\tilde{v}_{2}$
Right nodes connected to $\tilde{u}_{1}$ have degree 1; others have degree 2
(d) $t=2$.

Figure 3: (a): Pruned bipartite graph for stage B when $\boldsymbol{X}$ is
non-symmetric. Here,
$\boldsymbol{X}=\sigma\boldsymbol{u}\boldsymbol{v}^{T}=\tilde{\boldsymbol{u}}\tilde{\boldsymbol{v}}^{T}$
with $\tilde{u}_{j}\neq 0$ for $j\in[4]$, $\tilde{u}_{j}=0$ for $5\leq j\leq
n$, $\tilde{v}_{j}\neq 0$ for $j\in[3]$ and $\tilde{v}_{j}=0$ for $4\leq j\leq
n$. (b)–(d): The graph process that models the recovery of
$\tilde{\boldsymbol{u}}$ and $\tilde{\boldsymbol{v}}$. At $t=0$,
$\tilde{u}_{1}$ is recovered as 1 and peeled off. Following the recovery of
$\tilde{u}_{1}$, the left nodes $\tilde{v}_{1},\tilde{v}_{3}$ (and
$\tilde{u}_{4}$) can be peeled off too.

For rank $r>1$, when both sets of singular vectors $\\{\boldsymbol{u}_{i}\\}$
and $\\{\boldsymbol{v}_{i}\\}$ have disjoint supports, the bipartite graph for
stage B consists of $r$ disjoint subgraphs. The peeling algorithm in this case
is equivalent to the rank-1 algorithm run in parallel on the $r$ subgraphs.

For the general case with $r>1$ and the $\\{\boldsymbol{u}_{i}\\}$ or
$\\{\boldsymbol{v}_{i}\\}$ having overlapping supports, like in the symmetric
case, we use a large enough sketch such that stage A recovers all the nonzeros
in $\boldsymbol{X}$. A singular value decomposition on the recovered nonzero
submatrix gives the nonzeros in $\\{\boldsymbol{u}_{i},\boldsymbol{v}_{i}\\}$.

##  3 Main results

######  Theorem 1 (Noiseless symmetric case).

Consider the matrix
$\boldsymbol{X}=\sum_{i=1}^{r}\lambda_{i}\boldsymbol{v}_{i}\boldsymbol{v}_{i}^{T}$
with each eigenvector $\boldsymbol{v}_{i}$ having $k$ nonzero entries, with
the locations of the nonzeros being uniformly random. For sufficiently large
$k$, the sketching scheme with recovery algorithm described in Sections
2.1–2.3 has the following guarantees.

  1. 1)

For the case $r=1$ or $r>1$ with the supports of $\\{\boldsymbol{v}_{i}\\}$
disjoint, fix $\delta\in(0,1)$ and let $R=dr\tilde{k}/(\delta\ln k)$. Here, we
recall that $\tilde{k}=\binom{k}{2}+k$, and $d\geq 2$ is the number of ones in
each column of $\boldsymbol{H}$ (i.e., the degree of each left node in the
pruned stage A bipartite graph).

Then, with probability at least $1\,-2r\exp(-\frac{1}{30}k^{1-\delta})$, the
two-stage algorithm recovers $\\{\boldsymbol{v}_{i}\\}$ (up to a sign
ambiguity) from the sketch of size $m=2R$ with running time
$\mathcal{O}(rk^{2}/\ln k)$.

  2. 2)

When $r>1$ and the supports of $\\{\boldsymbol{v}_{i}\\}$ may overlap, with
probability at least $1-\mathcal{O}(k^{-2})$, stage A of the algorithm
followed by eigendecomposition of the recovered nonzero submatrix recovers
$\\{\boldsymbol{v}_{i}\\}$ from a sketch of size $m=3r\tilde{k}$ with running
time $\mathcal{O}\left((rk)^{3}\right)$.

Remarks:

  1. 1.

The eigenvectors $\\{\boldsymbol{v}_{i}\\}$ can only be recovered up to a sign
ambiguity since
$\boldsymbol{v}_{i}\boldsymbol{v}_{i}^{T}=(-\boldsymbol{v}_{i})(-\boldsymbol{v}_{i})^{T}$.
Moreover, when there is an eigenvalue with multiplicity $m>1$, the
corresponding eigenvectors can only be recovered up to a rotation within the
$m$-dimensional subspace they span.

  2. 2.

When the eigenvectors have different sparsities $\\{k_{1},\ldots,k_{r}\\}$,
Theorem 1 holds with the sample complexity depending on
$\max\\{k_{1},\ldots,k_{r}\\}$ and the success probability on
$\min\\{k_{1},\ldots,k_{r}\\}$.

  3. 3.

The bound on the success probability in part 1) is not tight. The bound also
does not depend $d$ since this parameter appears in a term in the proof that
is of a smaller order. In the simulations in Sections 3.4 and 4.4, we use
$d=2$.

  4. 4.

In part 1), the first stage of the algorithm is similar to the peeling decoder
for LDPC codes over an erasure channel [22, 23] and the decoder used for
compressed sensing in [2, 1]. However, since our sketching matrix
$\boldsymbol{B}$ (defined via the parity check matrix $\boldsymbol{H}$) has
row weights that scale with the sparsity level $k$, the existing peeling
decoder analysis based on density evolution cannot be applied. See comments
following Lemma 3.1.

######  Theorem 2 (Noiseless non-symmetric case).

Consider the matrix
$\boldsymbol{X}=\sum_{i=1}^{r}\sigma_{i}\boldsymbol{u}_{i}\boldsymbol{v}_{i}^{T}$
with each left singular vector $\boldsymbol{u}_{i}$ having $k$ nonzero entries
and each right singular vector $\boldsymbol{v}_{i}$ having $\beta k$ nonzero
entries, for some constant $\beta\in(0,1]$. Assume that the locations of the
supports are uniformly random. Then, for sufficiently large $k$, the sketching
scheme with recovery algorithm described in Section 2.4 has the following
guarantees.

  1. 1)

For the case $r=1$ or $r>1$ with each set of singular vectors
$\\{\boldsymbol{u}_{i}\\}$ and $\\{\boldsymbol{v}_{i}\\}$ having disjoint
supports, fix $\delta\in(0,\frac{1}{2})$ and let $R=dr\beta k^{2}/(\delta\ln
k)$. Then, with probability at least
$1-2r\exp(-\frac{\beta}{8}k^{1-2\delta})$, the two-stage algorithm recovers
$\\{\boldsymbol{u}_{i},\boldsymbol{v}_{i}\\}$ from the sketch of size $m=2R$
with running time $\mathcal{O}(rk^{2}/\ln k)$.

  2. 2)

When $r>1$ and the $\\{\boldsymbol{u}_{i}\\}$ or $\\{\boldsymbol{v}_{i}\\}$
have overlapping supports, with probability at least $1-\mathcal{O}(k^{-2})$,
stage A of the algorithm followed by singular value decomposition of the
recovered nonzero submatrix recovers
$\\{\boldsymbol{u}_{i},\boldsymbol{v}_{i}\\}$ from a sketch of size
$m=3rk^{2}$ with running time $\mathcal{O}\left((rk)^{3}\right)$.

Remarks similar to those for the symmetric case (below Theorem 1) hold for the
non-symmetric case as well.

###  3.1 Proof of Theorem 1

##### Part 1).

In the setting of part 1), the eigenvectors $\\{\boldsymbol{v}_{i}\\}$ have
disjoint supports with each $\boldsymbol{v}_{i}$ being $k$-sparse. In this
case, the nonzero entries of
$\boldsymbol{X}=\sum_{i=1}^{r}\lambda_{i}\boldsymbol{v}_{i}\boldsymbol{v}_{i}^{T}$
form $r$ disjoint submatrices, each of size $k\times k$. The result of Theorem
1 for this setting is proved via two lemmas, which give high probability
bounds on the number of nonzero matrix entries recovered in stage A and the
number of nonzero eigenvector entries recovered in stage B, respectively.

######  Lemma 3.1.

Consider the setting of part 1) of Theorem 1. Let $\sf{A}$ be the fraction
(out of $r\tilde{k}$) of nonzero entries in the upper triangular part of
$\boldsymbol{X}$ recovered in stage A. Then, with $\delta\in(0,1)$ as defined
in Theorem 1, there exists $\alpha^{*}=k^{-\delta}-o(k^{-\delta})$ such that
for sufficiently large $k$,

$\mathbb{P}\left\\{{\sf{A}}\geq\alpha^{*}\right\\}\geq
1-4\exp\left(-\frac{dr}{24\delta^{2}}\frac{k^{2-\delta}}{\ln^{2}k}\right)\,.$
(13)

Moreover, in the $k\times k$ nonzero submatrix corresponding to
$\lambda_{i}\boldsymbol{v}_{i}\boldsymbol{v}_{i}^{T}$, let ${\sf{A}}_{i}$ be
the fraction (out of $\binom{k}{2}$) of above-diagonal entries recovered and
let $N_{Di}$ be the number of diagonal entries recovered, for $i\in[r]$. Then
there exists $\alpha^{*}_{i}=k^{-\delta}-o(k^{-\delta})$ for $i\in[r]$ such
that for sufficiently large $k$,

$\mathbb{P}\left\\{N_{Di}\geq 1\text{ and
}{\sf{A}}_{i}\geq\alpha^{*}_{i},\;\forall\;i\in[r]\right\\}\geq
1-2r\exp\left(-\frac{1}{2}k^{1-\delta}\right).$ (14)

The proof of the lemma, given in Section 5.2, first establishes that in the
bipartite graph at the start of stage A (Fig. 1(a)), the degrees of the right
nodes are each Binomial with mean $\delta\ln k$ and negatively associated.
(Negative association is defined in Section 5.1). A Chernoff bound for
negatively associated random variables is then used to obtain a high
probability guarantee on the number of left nodes that are connected to
singleton right nodes (and can hence be recovered).

The next lemma shows that the conditional probability of recovering all the
nonzero entries in each eigenvector is close to 1, given the high probability
event in (14).

######  Lemma 3.2.

Consider the setting of part 1) of Theorem 1 and let
$N_{Di},{\sf{A}}_{i},\alpha_{i}^{*}$ be as defined in Lemma 3.1. Let
${\sf{B}}_{i}$ be the fraction (out of $k$) of nonzero entries of
$\boldsymbol{v}_{i}$ that are recovered by the end of stage B. Then, for all
$N_{Di}\geq 1,{\sf{A}}_{i}\geq\alpha^{*}_{i}$ and sufficiently large $k$ we
have:

$\mathbb{P}\left\\{{\sf{B}}_{i}<1\,\middle|\,N_{Di}\,,\,{\sf{A}}_{i}\right\\}<\exp\left(-\frac{1}{30}k^{1-\delta}\right),\quad
i\in[r].$ (15)

The proof of the lemma is given in Section 5.3. Recall that on each subgraph
$i\in[r]$ in stage B, the algorithm sequentially identifies and peels off left
nodes connected to degree-1 right nodes. To successfully recover all the $k$
nonzeros in the corresponding eigenvector $\boldsymbol{v}_{i}$, the residual
graph needs to have at least one degree-1 right node at the end of each
iteration $0\leq t\leq(k-2)$. Conditioned on ${\sf{A}}_{i}$ and the high-
probability event in (14), we show that at the start of each iteration $t$,
the number of degree-1 right nodes connected to each remaining left node is
approximately Binomial with mean ${\sf{A}}_{i}t$ (See Lemma 5.7). This is then
used to show that with high probability, there is at least one degree-1 right
node in each iteration until all the $k$ nonzeros are recovered.

_Proof of part 1) of Theorem 1_: From Lemmas 3.1 and 3.2, for sufficiently
large $k$ we have:

$\begin{split}&\mathbb{P}\left\\{{\sf{B}}_{i}=1,\ \forall i\in[r]\right\\}\\\
&\geq\,\mathbb{E}\left[\,\mathbb{P}\left\\{{\sf{B}}_{i}=1,\ \forall
i\in[r]\mid(N_{Di},\,{\sf{A}}_{i})_{i\in[r]}\right\\}\mathbbm{1}\\{N_{Di}\geq
1\text{ and }{\sf{A}}_{i}\geq\alpha_{i}^{*},\,\forall i\in[r]\\}\right]\\\
&\stackrel{{\scriptstyle(\rm{i})}}{{\geq}}\left(1-\,r\exp\Big{(}-\frac{1}{30}k^{1-\delta}\Big{)}\right)\mathbb{P}\left\\{N_{Di}\geq
1\text{ and }{\sf{A}}_{i}\geq\alpha_{i}^{*},\,\forall i\in[r]\right\\}\\\
&\stackrel{{\scriptstyle(\rm{ii})}}{{\geq}}\left(1-\,r\exp\Big{(}-\frac{1}{30}k^{1-\delta}\Big{)}\right)\left(1\,-\,2r\exp\Big{(}-\frac{1}{2}k^{1-\delta}\Big{)}\right)\\\
&>1-\,2r\exp\Big{(}-\frac{1}{2}k^{1-\delta}\Big{)}\,-\,r\exp\Big{(}-\frac{1}{30}k^{1-\delta}\Big{)}>1-2r\exp\Big{(}-\frac{1}{30}k^{1-\delta}\Big{)},\end{split}$
(16)

where the inequality (i) is obtained from (15) along with a union bound. The
inequality (ii) follows from (14). This completes the proof of part 1) of
Theorem 1. ∎

##### Part 2).

When $\\{\boldsymbol{v}_{i}\\}$ have overlapping supports, recall from Section
2.3 that the algorithm recovers all the nonzero entries of $\boldsymbol{X}$,
and then performs an eigendecomposition on the recovered nonzero submatrix. In
this case, the first stage is equivalent to the compressed sensing recovery of
the vectorized matrix $\boldsymbol{X}$ using the scheme of [1, 2]. Theorem 4
in [1] shows that the compressed sensing scheme can recover a $K$-sparse
vector with probability at least $1-\mathcal{O}(1/K)$ with a sample complexity
of $3K$ and running time $\mathcal{O}(K)$. This result directly yields the
high probability guarantee in part 2) by noting that the number of nonzeros in
the vectorized upper-triangular part of $\boldsymbol{X}$ is bounded below by
$(\binom{k}{2}+k)$ and above by $r({k\choose 2}+k)$.

###  3.2 Proof of Theorem 2

##### Part 1)

Part 1) of the theorem is proved using the following two lemmas which
characterize the high probability performance of stage A and stage B,
respectively.

######  Lemma 3.3.

Consider the setting of part 1) of Theorem 2. Let ${\sf{A}}_{i}$ denote the
fraction of entries recovered in the $k\times\beta k$ nonzero submatrix
corresponding to $\sigma_{i}\boldsymbol{u}_{i}\boldsymbol{v}_{i}^{T}$, for
$i\in[r]$. Then, with $\delta\in(0,\frac{1}{2})$ as defined in Theorem 2,
there exists $\alpha_{i}^{*}=k^{-\delta}-o(k^{-\delta})$ for $i\in[r]$ such
that for sufficiently large $k$,

$\mathbb{P}\left\\{{\sf{A}}_{i}\geq\alpha_{i}^{*}\,,\forall
i\in[r]\right\\}\geq
1-2r\exp\left(-\frac{\beta}{4}k^{\frac{3}{2}-\delta}\right).$ (17)

######  Lemma 3.4.

Consider the setting of part 1) of Theorem 2, and let
${\sf{A}}_{i},\alpha_{i}^{*}$ be as defined in Lemma 3.3. Let
${\sf{B}}^{u}_{i}$ and ${\sf{B}}^{v}_{i}$ be the fraction of nonzeros in
$\boldsymbol{u}_{i}$ and $\boldsymbol{v}_{i}$, respectively, that are
recovered in stage B. Then, for all ${\sf{A}}_{i}\geq\alpha_{i}^{*}$ and
sufficiently large $k$, we have

$\mathbb{P}\left\\{{\sf{B}}^{u}_{i}<1\,\text{ or
}\,{\sf{B}}^{v}_{i}<1\mid{\sf{A}}_{i}\right\\}\leq\exp\left(-\frac{\beta}{8}k^{1-2\delta}\right),\quad
i\in[r].$ (18)

The proofs of Lemmas 3.3 and 3.4 are similar to those of Lemmas 3.1 and 3.2,
respectively. We describe the main steps and highlight the key differences
from the symmetric case in Sections 5.4 and 5.5.

Part 1) of Theorem 2 follows from Lemmas 3.3 and 3.4. Indeed, for sufficiently
large $k$, we have:

$\begin{split}&\mathbb{P}\left\\{{\sf{B}}^{u}_{i}=1\text{ and
}{\sf{B}}^{v}_{i}=1\,,\forall i\in[r]\right\\}\\\
&\geq\mathbb{E}\left[\,\mathbb{P}\left\\{{\sf{B}}^{u}_{i}=1\text{ and
}{\sf{B}}^{v}_{i}=1,\ \forall i\in[r]\mid{\sf{A}}_{i},\
i\in[r]\right\\}\mathbbm{1}\\{{\sf{A}}_{i}\geq\alpha^{*}_{i},\,\forall\,i\in[r]\\}\right]\\\
&\stackrel{{\scriptstyle(\rm{i})}}{{\geq}}\left(1-r\exp\Big{(}-\frac{\beta}{8}k^{1-2\delta}\Big{)}\right)\mathbb{P}\left\\{{\sf{A}}_{i}\geq\alpha^{*}_{i},\,\forall\,i\in[r]\right\\}\\\
&\stackrel{{\scriptstyle(\rm{ii})}}{{\geq}}\left(1-r\exp\Big{(}-\frac{\beta}{8}k^{1-2\delta}\Big{)}\right)\left(1-2r\exp\Big{(}-\frac{\beta}{4}k^{\frac{3}{2}-\delta}\Big{)}\right)\\\
&\geq
1\,-\,r\exp\Big{(}-\frac{\beta}{8}k^{1-2\delta}\Big{)}\,-\,2r\exp\Big{(}-\frac{\beta}{4}k^{\frac{3}{2}-\delta}\Big{)}\geq
1\,-\,2r\exp\Big{(}-\frac{\beta}{8}k^{1-2\delta}\Big{)}.\end{split}$ (19)

Here the inequality (i) is obtained using Lemma 3.4, and (ii) using Lemma 3.3.
∎

##### Part 2).

This proof is similar to the symmetric case (part 2) of Theorem 1). Indeed,
applying [1, Theorem 4] guarantees that with probability
$1-\mathcal{O}(k^{-2})$, stage A of the algorithm recovers all the nonzero
entries in the matrix. The nonzeros in the singular vectors are then obtained
via a singular value decomposition on the recovered nonzero submatrix.

###  3.3 Computational cost of the recovery algorithm

We discuss the running time in the symmetric case, with the non-symmetric case
being analogous.

##### Eigenvectors with disjoint supports.

Consider the setting of part 1) of Theorem 1, where the supports of
$\\{\boldsymbol{v}_{i}\\}$ are disjoint. Lemma 3.1 guarantees that for any
$\delta\in(0,1)$ and sufficiently large $k$, the fraction ${\sf{A}}$ of
nonzero entries of $\boldsymbol{X}$ recovered in stage A is at least
$k^{-\delta}(1-o(1))$ with high probability. We will analyze the complexity of
the two stages assuming that ${\sf{A}}=k^{-\delta}(1-o(1))$. This is without
loss of generality as one can terminate stage A of the algorithm (prematurely)
once a fraction $k^{-\delta}(1-o(1))$ of nonzero entries of $\boldsymbol{X}$
have been recovered. From the proof of Lemma 3.1, this also means that the
fraction of entries recovered in the nonzero submatrix corresponding to
$\lambda_{i}\boldsymbol{v}_{i}\boldsymbol{v}_{i}^{T}$ is ${\sf{A}}_{i}\sim
k^{-\delta}$, for each $i\in[r].$ By Lemma 3.2, this is sufficient for
recovering each $\boldsymbol{v}_{i}$ with high probability.

_Stage A_ : To begin, each of the $R$ sketch bins requires $\mathcal{O}(1)$
numerical operations to be classified via a zeroton test and a singleton test
(specified in (12)). This requires a total of
$\mathcal{O}(R)=\mathcal{O}({rk^{2}}/{\ln k})$ operations. In each peeling
iteration, the contribution of a nonzero matrix entry is subtracted from the
$d$ bins that it is involved in and these bins are re-classified. Since $d$ is
a constant, each peeling iteration requires $\mathcal{O}(1)$ operations. By
the termination assumption above, the number of iterations in stage A is
$\mathcal{O}(rk^{2}{\sf{A}})=\mathcal{O}(rk^{2-\delta})$, corresponding to
$\mathcal{O}(rk^{2-\delta})$ operations. Therefore, the total computational
cost for stage A is $\mathcal{O}({rk^{2}}/{\ln
k})+\mathcal{O}(rk^{2-\delta})=\mathcal{O}({rk^{2}}/{\ln k})$.

_Stage B_ : Recall that the graph at the start of stage B consists of $r$
disjoint subgraphs, each with $k$ left nodes representing the nonzero entries
in the corresponding eigenvector and right nodes representing nonzero pairwise
products recovered in stage A. Lemma 5.6 shows that with high probability, the
degree of each left node in the $i$-th subgraph is
$\mathcal{O}({\sf{A}}_{i}k)$. Since we assumed that ${\sf{A}}_{i}\sim
k^{-\delta}$, the degree of each left node is $\mathcal{O}(k^{1-\delta})$ with
high probability. The algorithm peels off one left node at a time, and the
cost of each peeling iteration is proportional to the number of edges peeled
off during the iteration. Thus, to peel off all the $k$ left nodes from the
$i$-th subgraph, the computational cost is $\mathcal{O}(k^{2-\delta})$. Since
there are $r$ subgraphs, the total computational cost for stage B is
$\mathcal{O}(rk^{2-\delta})$.

Finally, adding the costs for the two stages gives a total computational cost
of $\mathcal{O}(rk^{2}/\ln k)$.

##### Eigenvectors with overlapping supports.

In the setting of part 2) of Theorem 1, the algorithm recovers all the nonzero
matrix entries in stage A. In this case, the number of nonzeros and the number
of bins $R$ are both $\mathcal{O}(rk^{2})$. Therefore, the computational cost
of stage A is $\mathcal{O}(rk^{2})$. The nonzeros in the eigenvectors are then
recovered by an eigendecompostion of the recovered nonzero submatrix. Since
this submatrix has size at most $rk\times rk$, the computational cost of the
eigendecomposition is $\mathcal{O}((rk)^{3})$, which dominates the total cost
of the algorithm.

###  3.4 Numerical results

We investigate the empirical performance of the scheme for both symmetric and
non-symmetric matrices, with $k$-sparse signal vectors
$\\{\boldsymbol{u}_{i},\boldsymbol{v}_{i}\\}$ that have disjoint or
overlapping supports. In the simulations, each signal vector
$\boldsymbol{u}_{i}$ or $\boldsymbol{v}_{i}$ is obtained by first sampling its
$k$ nonzero entries from the mixture of Gaussians $\frac{1}{2}{\sf
N}(-5,1)+\frac{1}{2}{\sf N}(5,1)$ and then normalizing so that the resulting
vector has unit norm. When $\\{\boldsymbol{u}_{i}\\}$ and
$\\{\boldsymbol{v}_{i}\\}$ have overlapping supports, we use the Gram–Schmidt
process to ensure both sets of signal vectors are orthogonal. The number of
ones in each column of the parity check matrix $\boldsymbol{H}$ is chosen to
be $d=2$. We declare exact recovery if: i) the locations of the nonzeros in
the signal vectors are correctly recovered and the values of the nonzeros are
recovered within an absolute deviation of $10^{-7}$ from the ground truth, and
ii) the eigenvalues or singular values are also recovered within an absolute
deviation of $10^{-7}$ from the ground truth.

In each subfigure of Fig. 4, we plot the fraction of trials in which exact
recovery is achieved versus the sketch size $m$, for different sparsity levels
$k$. In Figs. 4(a)–4(b), the dashed lines show the sketch sizes specified in
part 1) of Theorems 1 or 2 for the values of $\delta$ indicated in the
caption. The figures illustrate that for a fixed $\delta$, the success
probability at the sketch size corresponding to the dashed lines increases
with $k$. In Figs. 4(c)–4(d), the dashed lines indicate the sketch size
specified in part 2) of the theorems. The empirical success probability at the
dashed lines similarly increases with $k$.

Each subfigure in Fig. 5 plots the running time of the recovery algorithm
(Matlab implementation) versus the ambient dimension $n$, for three different
sparsity levels $k$ and sketch sizes $m$. The plots confirm that the running
time of the algorithm does not depend on $n$. The running times in different
subfigures are not comparable because the experiments were executed at
different times on a shared machine.

(a) Symmetric, $\\{\boldsymbol{v}_{i}\\}$ with disjoint supports.

(b) Nonsymmetric, $\\{\boldsymbol{u}_{i},\boldsymbol{v}_{i}\\}$ with disjoint
supports.

(c) Symmetric, $\\{\boldsymbol{v}_{i}\\}$ with overlapping supports.

(d) Nonsymmetric, $\\{\boldsymbol{u}_{i},\boldsymbol{v}_{i}\\}$ with
overlapping supports.

Figure 4: Probability of exact recovery ($y$-axis) versus sketch size $m$
($x$-axis), for different sparsity levels $k$. The matrices used in all cases
have rank $r=3$ and size $n\times n$ with $n=10^{4}$. The dashed lines
indicate the sketch sizes stated in Theorems 1 and 2; the dashed lines in (a)
and (b) are determined using $\delta=\frac{5}{7}$ and $\delta=\frac{2}{5}$,
respectively. Results are averaged over 100 trials.

(a)
Symmetric, $r=1$.

(b) Symmetric, $r=3$. $\\{\boldsymbol{v}_{i}\\}$ with disjoint supports.

(c) Symmetric, $r=3$. $\\{\boldsymbol{v}_{i}\\}$ with overlapping supports.

(d) Nonsymmetric, $r=3$. $\\{\boldsymbol{u}_{i},\boldsymbol{v}_{i}\\}$ w.
disjoint supports

Figure 5: Running time (in secs) of the recovery algorithm versus ambient
dimension $n$, for different sparsity levels $k$ and sketch sizes $m$. Exact
recovery is achieved in every trial. Results are averaged over 50 trials and
error bars indicate one standard deviation.

##  4 Recovery in the presence of noise

Consider the general model for the data matrix $\boldsymbol{X}$ in (2):

$\boldsymbol{X}=\boldsymbol{X}_{0}+\boldsymbol{W},\ \text{ where }\
\boldsymbol{X}_{0}=\sum_{i=1}^{r}\lambda_{i}\boldsymbol{v}_{i}\boldsymbol{v}_{i}^{T}\
(\text{symmetric})\ \text{ or }\
\boldsymbol{X}_{0}=\sum_{i=1}^{r}\sigma_{i}\boldsymbol{u}_{i}\boldsymbol{v}_{i}^{T}\
(\text{non-symmetric}).$ (20)

The vectors $\boldsymbol{v}_{i}$ and $\boldsymbol{u}_{i}$ are $k$-sparse, for
$i\in[r]$. In the non-symmetric case, the elements of the noise matrix
$\boldsymbol{W}$ are assumed to be independent and identically distributed
with zero mean and variance $\sigma^{2}$. In the symmetric case, the same
assumption holds for the upper triangular entries, with
$W_{\ell,j}=W_{j,\ell}$ for $j<\ell\in[n]$.

###  4.1 Sketching scheme

We describe the scheme for symmetric matrices, with the non-symmetric case
being analogous. We write $\boldsymbol{x},\boldsymbol{x}_{0},\boldsymbol{w}$
for the vectorized upper-triangular parts of the symmetric matrices
$\boldsymbol{X},\boldsymbol{X}_{0},\boldsymbol{W}$, respectively. As before,
we let $\tilde{n}={n\choose 2}+n$ and $\tilde{k}={k\choose 2}+k$. As in the
noiseless setting (Section 2.1), the sketching matrix $\boldsymbol{B}$ is
formed by combining two matrices: a sparse parity check matrix
$\boldsymbol{H}$, and a matrix $\boldsymbol{S}$ that helps classify each
sketch bin as a zeroton, singleton, or multiton. When there is noise, we
cannot reliably identify zerotons and singletons by simply taking
$\boldsymbol{S}$ to be the first two rows of an $\tilde{n}$-point DFT matrix
as in (5). Therefore, following [1], we choose
$\boldsymbol{S}\in{\mathbb{R}}^{P\times\tilde{n}}$ with

$P=\mathcal{O}\left(\log\left(\frac{n}{k}\right)\right)\ \text{ and }\quad
S_{\kappa,\ell}\stackrel{{\scriptstyle\text{iid}}}{{\sim}}{\sf N}(0,1),\
\kappa\in[P],\,\ell\in[\tilde{n}].$ (21)

The sparse parity check matrix $\boldsymbol{H}\in\\{0,1\\}^{R\times\tilde{n}}$
is chosen in the same way as in the noiseless setting: column-regular with
$d\geq 2$ ones per column in uniformly random locations. The sketching matrix
$\boldsymbol{B}\in{\mathbb{R}}^{PR\times\tilde{n}}$ is then formed as the
column-wise Kronecker product of $\boldsymbol{H}$ and $\boldsymbol{S}$,
computed as illustrated in (7)–(8). The sketch
$\boldsymbol{y}=\boldsymbol{B}\boldsymbol{x}$ has size $m=PR$.

We can view the sketch $\boldsymbol{y}$ as comprising $R$ bins, each with $P$
entries. For $j\in[R]$, we let
$\boldsymbol{y}_{j}:=\begin{bmatrix}y_{(Pj-P+1)},\dots,y_{(Pj)}\end{bmatrix}^{T}\in{\mathbb{R}}^{P}$
be the $j$-th sketch bin. Similarly to (9) and (10), writing
$\boldsymbol{B}_{j}\in{\mathbb{R}}^{P\times\tilde{n}}$ for the $j$-th set of
$P$ rows in $\boldsymbol{B}$, we have

$\displaystyle\boldsymbol{y}_{j}=\boldsymbol{B}_{j}\boldsymbol{x}=\boldsymbol{S}\boldsymbol{x}_{0j}^{*}+\boldsymbol{S}\boldsymbol{w}_{j}^{*},\quad
j\in[R],$ (22)

where the entries of $\boldsymbol{x}_{0j}^{*}\in{\mathbb{R}}^{\tilde{n}}$ and
$\boldsymbol{w}_{j}^{*}\in{\mathbb{R}}^{\tilde{n}}$ are

$(x^{*}_{0j})_{\ell}=\left\\{\begin{array}[]{ll}x_{0\ell}&\ell\in\mathcal{N}(j)\\\
0&\ell\notin\mathcal{N}(j)\end{array}\right.,\quad(w^{*}_{j})_{\ell}=\left\\{\begin{array}[]{ll}w_{\ell}&\ell\in\mathcal{N}(j)\\\
0&\ell\notin\mathcal{N}(j)\end{array}\right.\quad\text{for
}\ell\in\left[\tilde{n}\right].$ (23)

As in (10), $\mathcal{N}(j)=\\{\ell\in[\tilde{n}]:\,H_{j,\ell}=1\\}$ contains
the locations of the ones in the $j$-th row of $\boldsymbol{H}$. The vectors
$\boldsymbol{x}^{*}_{0j}$ and $\boldsymbol{w}^{*}_{j}$ store the signal and
noise entries that contribute to the $j$-th bin. The $j$-th bin
$\boldsymbol{y}_{j}$ is the sum of these entries weighted by their respective
columns in $\boldsymbol{S}$.

###  4.2 Recovery algorithm

We first consider symmetric rank-$r$ matrices where
$\\{\boldsymbol{v}_{1},\ldots,\boldsymbol{v}_{r}\\}$ have disjoint supports.
Fix $\delta\in(0,1)$, let $R=dr\tilde{k}/(\delta\ln k)$, and recall that the
sketch size is $m=PR$.

####  4.2.1 Accumulated noise variance in each bin

The recovery of $\boldsymbol{x}$ crucially depends on the variance of the
noise vector $\boldsymbol{S}\boldsymbol{w}_{j}^{*}$ in each sketch bin
$\boldsymbol{y}_{j}$ (see (22)). Since $\boldsymbol{H}$ has $d$ ones per
column, each entry of $\boldsymbol{x}$ contributes to $d$ bins, chosen
uniformly at random from the $R$ bins in the sketch. Thus we have
$\left\lvert\mathcal{N}(j)\right\rvert\sim\text{Bin}\left(\tilde{n},\frac{d}{R}\right)$.
Noting that
$\mathbb{E}\\{\left\lvert\mathcal{N}(j)\right\rvert\\}=\frac{\tilde{n}d}{R}=({n\choose
2}+n)\frac{\delta\ln k}{r\tilde{k}}$, we apply the tail bound in Lemma 5.1 to
obtain

$\mathbb{P}\left\\{\left\lvert\mathcal{N}(j)\right\rvert\leq\frac{\delta}{4}\frac{n^{2}\ln
k}{rk^{2}}\right\\}\leq\exp\left(-\frac{\delta}{7}\frac{n^{2}\ln
k}{rk^{2}}\right).$ (24)

Writing $\tilde{\boldsymbol{w}}_{j}^{*}:=\boldsymbol{S}\boldsymbol{w}^{*}_{j}$
for the noise vector in (22), the entries in $\tilde{\boldsymbol{w}}_{j}^{*}$
are:

$(\tilde{w}^{*}_{j})_{\kappa}=\sum_{\ell\in\mathcal{N}(j)}{S}_{\kappa,\ell}\,w_{\ell},\quad\text{for
}\kappa\in[P].$ (25)

Recalling that $w_{\ell}$ is zero-mean with variance $\sigma^{2}$ and i.i.d.
across $\ell$, conditioned on $\boldsymbol{S}$ and $\mathcal{N}(j)$, the
vector $\tilde{\boldsymbol{w}}_{j}^{*}$ is zero-mean with covariance matrix
$\boldsymbol{\Sigma}\in\mathbb{R}^{P\times P}$, where

$\displaystyle\Sigma_{\kappa,\kappa}=\sigma^{2}\sum_{\ell\in\mathcal{N}(j)}S_{\kappa,\ell}^{2},\quad\Sigma_{\kappa,\iota}=\sigma^{2}\sum_{\ell\in\mathcal{N}(j)}S_{\kappa,\ell}S_{\iota,\ell}\quad\text{for
}\kappa,\iota\in[P]\text{ and }\kappa\neq\iota.$ (26)

Since $\boldsymbol{S}$ has independent standard Gaussian entries,
$\mathbb{E}[S_{\kappa,\ell}^{2}]=1$ and
$\mathbb{E}[S_{\kappa,\ell}S_{\iota,\ell}]=0$. Moreover, (24) implies that
$\left\lvert\mathcal{N}(j)\right\rvert\to\infty$ almost surely as
$n\to\infty$. Therefore, by the weak law of large numbers,
$\frac{\Sigma_{\kappa,\kappa}}{|\mathcal{N}(j)|}\to\sigma^{2}$ and
$\frac{\Sigma_{\kappa,\iota}}{|\mathcal{N}(j)|}\to 0$ as $n\to\infty$. It
follows that the distribution of the scaled accumulated noise vector
$\frac{\tilde{\boldsymbol{w}}_{j}^{*}}{\sqrt{|\mathcal{N}(j)|}}$ converges to
the Gaussian ${\sf N}(\boldsymbol{0},\sigma^{2}\boldsymbol{I}_{P})$ as
$n\rightarrow\infty$. We will use this limiting distribution to devise zeroton
and singleton tests in stage A of the algorithm.

####  4.2.2 Stage A of the algorithm

##### Zeroton test.

When $\boldsymbol{y}_{j}$ is a zeroton, i.e., when ${\sf
supp}(\boldsymbol{x}_{0})\cap\mathcal{N}(j)=\emptyset$, from (22) and (23) we
have

$\boldsymbol{y}_{j}=\boldsymbol{S}\boldsymbol{w}^{*}_{j}=\tilde{\boldsymbol{w}}_{j}^{*},\quad\text{with
}\frac{\tilde{\boldsymbol{w}}_{j}^{*}}{\sqrt{|\mathcal{N}(j)|}}\ \text{
approximately}\sim{\sf
N}\left(\boldsymbol{0},\sigma^{2}\mathbf{I}_{P}\right).$ (27)

Based on this, a given bin $\boldsymbol{y}_{j}$ is declared a zeroton if

$\frac{1}{P\left\lvert\mathcal{N}(j)\right\rvert}\left\|\boldsymbol{y}_{j}\right\|^{2}\leq\gamma_{0}\sigma^{2},$
(28)

for a suitably chosen constant $\gamma_{0}$. Empirically, $99.7\%$ of the
samples of a Gaussian distribution lie within $3$ standard deviations of the
mean, so we set $\gamma_{0}\geq 3$. In our experiments, $\gamma_{0}=5$ gives
reasonable performance.

##### Singleton test.

When $\boldsymbol{y}_{j}$ is a singleton with ${\sf
supp}(\boldsymbol{x}_{0})\cap\mathcal{N}(j)=\ell$, it takes the form

$\boldsymbol{y}_{j}=x_{0\ell}\boldsymbol{s}_{\ell}+\tilde{\boldsymbol{w}}^{*}_{j},$
(29)

where $\boldsymbol{s}_{\ell}$ is the $\ell$-th column of $\boldsymbol{S}.$
Given any bin $\boldsymbol{y}_{j}$, assuming it is a singleton of the form
(29), the algorithm first estimates the index $\ell$ and the value $x_{0\ell}$
as follows:

$\hat{\ell}=\operatorname*{arg\,min}_{\ell\in\mathcal{N}(j)}\left\|\boldsymbol{y}_{j}-\hat{x}_{0\ell}\boldsymbol{s}_{\ell}\right\|^{2},\
\text{ where }\
\hat{x}_{0\ell}=\frac{\boldsymbol{s}_{\ell}^{T}\boldsymbol{y}_{j}}{\left\|\boldsymbol{s}_{\ell}\right\|^{2}}.$
(30)

The bin $j$ is declared a singleton if

$\frac{1}{P\left\lvert\mathcal{N}(j)\right\rvert}\left\|\boldsymbol{y}_{j}-\hat{x}_{0\hat{\ell}}\,\boldsymbol{s}_{\hat{\ell}}\right\|^{2}\leq\gamma_{1}\sigma^{2},$
(31)

for some constant $\gamma_{1}\leq\gamma_{0}$. The estimated index-value pair
$(\hat{\ell},\hat{x}_{0\hat{\ell}})$ in (30) is the maximum-likelihood
estimate of $(\ell,x_{0\ell})$ corresponding to the singleton bin in (29),
given the limiting distribution of the noise
$\frac{\tilde{\boldsymbol{w}}_{j}^{*}}{\sqrt{|\mathcal{N}(j)|}}$ is ${\sf
N}(\boldsymbol{0},\sigma^{2}\mathbf{I}_{P})$. We note that
$(\hat{\ell},\hat{x}_{0\hat{\ell}})$ will be an accurate estimate of the true
index-value pair $(\ell,x_{0\ell})$ only when the accumulated noise level
$\sqrt{\left\lvert\mathcal{N}(j)\right\rvert}\sigma$ (i.e., the standard
deviation of $\tilde{\boldsymbol{w}}^{*}_{j}$) is much smaller than the signal
magnitude $\left\lvert x_{0\ell}\right\rvert$.

The algorithm declares $\boldsymbol{y}_{j}$ a multiton if it fails both
zeroton and singleton tests.

At $t=0$, the algorithm identifies the initial set of zerotons and singletons
among the $R$ bins using the tests above. Let $\mathcal{S}(t)$ denote the set
of singletons at time $t$. Then at each $t\geq 1$, the algorithm picks a
singleton uniformly at random from $\mathcal{S}(t-1)$, and recovers the
nonzero entry of $\boldsymbol{x}_{0}$ underlying the singleton using the
index-value pair $(\hat{\ell},\hat{x}_{0\hat{\ell}})$ in (30). The algorithm
then peels off the contribution of the $\hat{\ell}$-th entry of
$\boldsymbol{x}_{0}$ from each of the $d$ bins it is involved in, re-
categorizes these bins, and updates the set of singletons to include any new
singletons created by the peeling of the $\hat{\ell}$-th entry. The updated
set of singletons is $\mathcal{S}(t)$. Stage A continues until singletons run
out, i.e., when $\mathcal{S}(t)=\emptyset$. If the supports of
$\\{\boldsymbol{v}_{1},\ldots,\boldsymbol{v}_{r}\\}$ are known to be disjoint,
we proceed to stage B, described in the next subsection.

When the supports of the $r$ eigenvectors have overlaps, we take a large
enough sketch (with $R=\mathcal{O}(rk^{2})$) to identify all the nonzero
entries in $\boldsymbol{X}_{0}$ in stage A, and then perform an
eigendecomposition on the recovered submatrix of nonzero entries to obtain the
nonzero entries of $\\{\boldsymbol{v}_{1},\ldots,\boldsymbol{v}_{r}\\}$.

_Picking the threshold parameters $\gamma_{0},\gamma_{1}$_: Decreasing the
value of $\gamma_{0}$ makes the zeroton test in (28) more stringent. As
$\gamma_{0}$ becomes smaller, zerotons are more likely to be mistaken as
singletons with small signal values. This kind of error results in some zero
entries of $\boldsymbol{X}_{0}$ being mistakenly recovered as small nonzeros.
Since our recovery algorithm crucially relies on the sparsity of
$\boldsymbol{X}_{0}$ to maintain its low complexity, empirically we find that
such errors greatly increase the algorithm’s running time. Moreover, in the
case where $\\{\boldsymbol{v}_{i}\\}$ have disjoint supports, such errors can
lead to stability issues in stage B when peeling steps involve taking ratios
between small signal values. For these reasons, we set $\gamma_{0}=5$ in our
numerical simulations. We set a more stringent threshold for the singleton
test by choosing $\gamma_{1}\leq\gamma_{0}$ to prevent multitons from being
mistaken as singletons.

####  4.2.3 Stage B of the algorithm

In the stage B graph (Fig. 2(a) or Fig. 3(a)), the right nodes are the nonzero
matrix entries (i.e., pairwise products) recovered in stage A. The values of
these recovered entries are no longer exact in the noisy setting. Hence, a
left node connected to multiple degree-1 right nodes will receive different
suggested values from these right nodes. Therefore, instead of peeling off
edges like in Section 2.2.2, we use a message passing algorithm on the stage B
graph. In each iteration, all the left nodes with known values send these
values to their neighboring right nodes. Recall that each right node in the
stage B graph has two edges, and corresponds to the product of the two left
nodes it is connected to. The message a right node sends along each edge is
the ratio of its pairwise product with the incoming message along the other
edge. Each left node then updates its value by taking the average of all
messages it receives. This averaging operation makes the algorithm more robust
to noise.

###  4.3 Sample complexity and computational cost

##### Disjoint supports

We have $R=\mathcal{O}(rk^{2}/\ln k)$ bins and $P=\mathcal{O}(\log(n/k))$
components in each bin, therefore the sample complexity is
$RP=\mathcal{O}(rk^{2}\log(n/k)/\ln k)$, which is higher than the noiseless
case by a factor of $P=\mathcal{O}(\log(n/k))$.

The computational cost of the recovery algorithm is
$\mathcal{O}(n^{2}\log(n/k))$, significantly higher than the
$\mathcal{O}(rk^{2}/\ln k)$ in the noiseless setting. The increase in
complexity is due to the singleton test in (30)–(31), which requires finding a
minimum over the $\left\lvert\mathcal{N}(j)\right\rvert$ left nodes that
connect to bin $\boldsymbol{y}_{j}$, for each $j\in[R]$. From (24) and a
similar upper tail bound, we know that with high probability
$\left\lvert\mathcal{N}(j)\right\rvert$ is of the order $n^{2}\ln k/(rk^{2})$.
Therefore, the total computational cost of the singleton test for the $R$ bins
is $R\left\lvert\mathcal{N}(j)\right\rvert P=\mathcal{O}(n^{2}\log(n/k))$,
which dominates the cost of stage A. When the noise level is small (so that in
stage A, zero matrix entries are not mistakenly recovered as nonzeros), stage
B has the same $\mathcal{O}(rk^{2-\delta})$ complexity as in the noiseless
case. Thus the computational cost of the algorithm is dominated by stage A.

##### General case

The sample complexity is $\mathcal{O}(rk^{2}\log(n/k))$, again higher than the
noiseless setting by a factor of $P=\mathcal{O}(\log(n/k))$. The computational
cost of recovering all the nonzero entries of the signal matrix is
$\mathcal{O}(n^{2}\log(n/k))$ and that of the eigendecomposition of the
recovered nonzero submatrix is $\mathcal{O}((rk)^{3})$. Thus, the total cost
of the algorithm is $\mathcal{O}(\max\\{n^{2}\log(n/k),(rk)^{3}\\})$.

In the noisy setting, the nonzero matrix entries recovered in stage A will not
be exact, which makes it hard to analyze the subsequent step (either stage B
or eigendecomposition) and obtain a rigorous performance guarantee. This can
be addressed by making suitable assumptions on the signal alphabet. The
analysis is conceptually similar to the noiseless setting but with additional
technical challenges, so it is omitted. In the following subsection, we
demonstrate the performance of the scheme in the noisy case via numerical
simulations.

###  4.4 Numerical results

We apply the robust sketching scheme and recovery algorithm to symmetric
matrices $\boldsymbol{X}=\boldsymbol{X}_{0}+\boldsymbol{W}$, with
$\boldsymbol{X}_{0}=\sum_{i=1}^{r}\tilde{\boldsymbol{v}}_{i}\tilde{\boldsymbol{v}}_{i}^{T}$
and $\boldsymbol{W}$ a random symmetric Gaussian matrix. Each unnormalized
eigenvector $\tilde{\boldsymbol{v}}_{i}$ is $k$-sparse with its nonzero
entries drawn uniformly at random from a discrete alphabet $\\{\pm 10,\pm
20,\dots,\pm 50\\}$. When the supports of $\\{\tilde{\boldsymbol{v}}_{i}\\}$
overlap, we adjust in each $\tilde{\boldsymbol{v}}_{i}$ a small number of
nonzero entries to ensure orthogonality among the
$\\{\tilde{\boldsymbol{v}}_{i}\\}$. The resulting eigenvalues of
$\boldsymbol{X}_{0}$ are $\lambda_{i}=\|\tilde{\boldsymbol{v}}_{i}\|^{2}$ for
$i\in[r]$. We measure the signal strength by the expected average eigenvalue
$\lambda:=\mathbb{E}\left[\frac{1}{r}\sum_{i=1}^{r}\lambda_{i}\right]$. The
entries in the noise matrix $\boldsymbol{W}$ are drawn i.i.d. up to symmetry
from ${\sf N}(0,\sigma^{2})$. Let
$\lambda_{W,1}\geq\lambda_{W,2}\geq\dots\geq\lambda_{W,n}$ denote the
eigenvalues of $\boldsymbol{W}$. From Wigner’s semicircle law [40], we have
that $\max_{i\in[r]}\frac{\lambda_{W,i}}{\sqrt{n}}$ concentrates around
$2\sigma$ and that
$\lim_{n\to\infty}\mathbb{E}\left[\frac{1}{n^{2}}\sum_{i=1}^{n}\lambda_{W,i}^{2}\right]=\sigma^{2}.$
Thus, to ensure the right scaling in our signal to noise ratio (SNR), we
define our SNR to be $\lambda/(\sqrt{n}\sigma)$.

Let $\hat{\boldsymbol{X}}_{0}$ denote the estimated signal matrix computed
using the eigenvectors recovered by the algorithm. Fig. 6 plots the normalized
squared error of the recovered matrix (defined as
$\|\hat{\boldsymbol{X}}_{0}-\boldsymbol{X}_{0}\|^{2}_{F}/{\|\boldsymbol{X}_{0}\|^{2}_{F}}$)
against the sketch size $m$, for different SNR values. Recall that the sketch
size is $m=RP$, where $R$ is the number of bins and $P$ is the number of
components in each bin. In Fig. 6(a), the sketch size is varied via $R$, with
$P$ held constant. As $R$ increases, more nonzero matrix entries are recovered
in the first stage of the algorithm. This increases the number of right nodes
in the graph for stage B, allowing more nonzero eigenvector entries to be
recovered which leads to a smaller normalized squared error. In Figs. 6(b) and
6(c), the sketch size is varied via $P$, with the number of bins $R$ held
constant. The $P$ values we vary across are $\lceil c\log(n/k)\rceil$ for
$c=1,2,\dots,6.$ The accuracy of the recovery improves as $P$ increases due to
more reliable zeroton and singleton detection. In all subfigures, we observe
that the normalized squared error decreases with SNR as expected.

(a) $\\{\boldsymbol{v}_{i}\\}$ with disjoint supports. Vary $m$ via $R$;
$P=\lceil 7\log(n/k)\rceil$.

(b) $\\{\boldsymbol{v}_{i}\\}$ with disjoint supports. Vary $m$ via $P$;
$R=5r\tilde{k}/(\ln k)$.

(c) $\\{\boldsymbol{v}_{i}\\}$ with overlapping supports. Vary $m$ via $P$;
$R=10r\tilde{k}$.

Figure 6: Normalized squared error
$\|\hat{\boldsymbol{X}}_{0}-\boldsymbol{X}_{0}\|^{2}_{F}/{\|\boldsymbol{X}_{0}\|^{2}_{F}}$
versus sketch size $m$, for different SNR values computed as
$\lambda/(\sqrt{n}\sigma)$. In all cases, $n=5000$, $k=100$ and $r=3$. The
column weight of the parity check matrix $\boldsymbol{H}$ is $d=2$. The
parameters in the zeroton and singleton tests are $\gamma_{0}=5$ and
$\gamma_{1}=1.5$ in (a)–(b) and $\gamma_{1}=2.5$ in (c). Results are averaged
over 10 sketching matrices and 15 trials each.

##  5 Proofs of main lemmas

In this section, we give the proofs of Lemmas 3.1–3.4.

###  5.1 Preliminaries

We begin with the following tail bound for binomial random variables, which
follows from Sanov’s theorem [41, Theorem 12.4.1].

######  Lemma 5.1 (Tail bound for Binomial).

For $X\sim\text{Bin}\left(n,p\right)$,

$\displaystyle\mathbb{P}\left\\{X\leq x\right\\}$
$\displaystyle\leq\exp\left(-nD_{e}\left(\frac{x}{n}\middle\|p\right)\right)\
\text{ for }0\leq x\leq np,$ (32) $\displaystyle\mathbb{P}\left\\{X\geq
x\right\\}$
$\displaystyle\leq\exp\left(-nD_{e}\left(\frac{x}{n}\middle\|p\right)\right)\
\text{ for }np\leq x\leq n,$ (33)

where $D_{e}(q\|p):=q\ln\frac{q}{p}+(1-q)\ln\frac{1-q}{1-p}$.

We will use the notion of _negative association_ (NA) [42, 24] to handle the
dependence between the degrees of the bipartite graphs we analyze. Roughly
speaking, a collection of random variables $X_{1},X_{2},\ldots,X_{n}$ are
negatively associated if whenever a subset of them is high, then a disjoint
subset of them must be low.

######  Definition 5.1 (Negative Association (NA)).

The random variables $X_{1},X_{2},\dots,X_{n}$ are said to be negatively
associated if for any two disjoint index sets $\mathcal{I,J}\subseteq[n]$ and
two functions
$f:{\mathbb{R}}^{\left\lvert\mathcal{I}\right\rvert}\to{\mathbb{R}}$,
$g:{\mathbb{R}}^{\left\lvert\mathcal{J}\right\rvert}\to{\mathbb{R}}$ both
monotonically increasing or monotonically decreasing,

$\mathbb{E}\left[f(X_{i},i\in\mathcal{I})g(X_{j},j\in\mathcal{J})\right]\leq\mathbb{E}\left[f(X_{i},i\in\mathcal{I})\right]\mathbb{E}\left[g(X_{j},j\in\mathcal{J})\right].$
(34)

######  Lemma 5.2 (Useful properties of NA [42, 24]).

  1. (i)

The union of independent sets of NA random variables is NA.

  2. (ii)

Concordant monotone functions (i.e., all monotonically increasing or all
monotonically decreasing functions) defined on disjoint subsets of a set of NA
random variables are NA.

  3. (iii)

Let $\boldsymbol{x}=\left[x_{1},x_{2},\dots,x_{n}\right]$ be a real-valued
vector and let $\boldsymbol{X}=\left[X_{1},X_{2},\dots,X_{n}\right]$ be a
random vector which takes as values all the $n!$ permutations of
$\boldsymbol{x}$ with equal probabilities. Then $X_{1},X_{2},\dots,X_{n}$ are
NA.

######  Lemma 5.3 (Chernoff bound for NA Bernoulli variables [24]).

Let $X_{1},X_{2},\dots,X_{n}$ be NA random variables with $X_{i}\in\\{0,1\\}$
for $i\in[n]$. Then, $Y=\sum_{i=1}^{n}X_{i}$ satisfies:

$\displaystyle\mathbb{P}\left\\{Y>(1+\epsilon)\mathbb{E}[Y]\right\\}\leq\exp\left(-\frac{\epsilon^{2}\,\mathbb{E}[Y]}{2+\epsilon}\right)\quad\text{
for any }\epsilon\geq 0,$ (35)
$\displaystyle\mathbb{P}\left\\{Y<(1-\epsilon)\mathbb{E}[Y]\right\\}\leq\exp\left(-\frac{\epsilon^{2}\,\mathbb{E}[Y]}{2-\epsilon}\right)\quad\text{
for any }\epsilon\in[0,1].$ (36)

###  5.2 Proof of Lemma 3.1

The pruned graph at the start of stage A has $r\tilde{k}=r(\binom{k}{2}+k)$
left nodes, each with $d$ edges. Fig. 1(a) is an example for $r=1$. The number
of right nodes is $R=dr\tilde{k}/(\delta\ln k)$. Let $Z_{j}$ denote the degree
of the $j$-th right node for $j\in[R]$. Since the total number of edges in the
graph is $dr\tilde{k}$, we have $\sum_{j=1}^{R}Z_{j}=dr\tilde{k}$.

######  Lemma 5.4 (Initial right degrees are Binomial and NA).

For $j\in[R]$, we have
$Z_{j}\sim\text{Bin}\left(r\tilde{k},\frac{d}{R}\right)$. Furthermore,
$Z_{1},Z_{2},\ldots,Z_{R}$ are NA.

###### Proof.

Recall that the bipartite graph at the start of stage A for a rank-$r$ matrix
$\boldsymbol{X}$ has $r\tilde{k}$ left nodes and $R$ right nodes. For
$l\in[r\tilde{k}]$ and $j\in[R]$, let

$Y_{l,j}:=\begin{cases}1&\text{if the $l$-th left node connects to the $j$-th
right node,}\\\ 0&\text{otherwise}.\end{cases}$ (37)

Observe that $Y_{l,j}\sim\text{Bern}\left(\frac{d}{R}\right)$ and that
$Y_{l_{1},j_{1}}\perp Y_{l_{2},j_{2}}$ for $l_{1}\neq l_{2}$. Therefore,

$Z_{j}=\sum_{l=1}^{r\tilde{k}}Y_{l,j}\sim\text{Bin}\left(r\tilde{k},\frac{d}{R}\right)\quad\text{for
}j\in[R].$ (38)

By the construction of the bipartite graph, the vector
$\boldsymbol{Y}_{l}:=\left[Y_{l,1},Y_{l,2},\dots,Y_{l,R}\right]$ contains $d$
ones, distributed uniformly at random among its $R$ entries; the remaining
$(R-d)$ entries are zeros. That is, the joint distribution of the entries of
$\boldsymbol{Y}_{l}$ is a permutation distribution. Hence,
$\left[Y_{l,1},Y_{l,2},\dots,Y_{l,R}\right]$ is NA by Lemma 5.2(iii). Since
the vectors $\boldsymbol{Y}_{l}$ are mutually independent for
$l\in[r\tilde{k}]$, the concatenated vector $[Y_{l,j}\ \forall
l\in[r\tilde{k}],\ \forall j\in[R]]$ is NA by Lemma 5.2(i). Furthermore, from
(38), we note that $Z_{1},\dots,Z_{R}$ are increasing functions defined on
disjoint subsets of the $Y_{l,j}$’s. Thus, by Lemma 5.2(ii),
$Z_{1},\dots,Z_{R}$ are NA. ∎

The number of left nodes (i.e. nonzero matrix entries) recovered in stage A is
at least as large as the number of singletons at the start of stage A, which
equals $\sum_{j=1}^{R}\mathbbm{1}{\left\\{Z_{j}=1\right\\}}$. We now use Lemma
5.4 to obtain a high probability bound on this sum.

######  Lemma 5.5 (Bound on number of singletons).

Let $Z\sim\text{Bin}\left(r\tilde{k},\frac{d}{R}\right)$. Then, for
sufficiently large $k$, the number of singletons
$\sum_{j=1}^{R}\mathbbm{1}{\left\\{Z_{j}=1\right\\}}$ satisfies:

$\mathbb{P}\left\\{\left\lvert\frac{1}{R}\sum_{j=1}^{R}\mathbbm{1}{\left\\{Z_{j}=1\right\\}}-\mathbb{P}\left\\{Z=1\right\\}\right\rvert\geq
k^{-\delta}\right\\}<4\exp\left(-\frac{dr}{24\delta^{2}}\frac{k^{2-\delta}}{\ln^{2}k}\right).$
(39)

###### Proof.

We first note

$\frac{1}{R}\sum_{j=1}^{R}\mathbbm{1}{\left\\{Z_{j}=1\right\\}}-\mathbb{P}\left\\{Z=1\right\\}=\frac{1}{R}\sum_{j=1}^{R}\mathbbm{1}{\left\\{Z_{j}\leq
1\right\\}}-\mathbb{P}\left\\{Z\leq
1\right\\}\,+\,\mathbb{P}\left\\{Z=0\right\\}-\frac{1}{R}\sum_{j=1}^{R}\mathbbm{1}{\left\\{Z_{j}=0\right\\}}\,.$

Applying the triangle inequality and a union bound, we have

$\displaystyle\mathbb{P}\left\\{\left\lvert\frac{1}{R}\sum_{j=1}^{R}\mathbbm{1}{\left\\{Z_{j}=1\right\\}}-\mathbb{P}\left\\{Z=1\right\\}\right\rvert\geq
k^{-\delta}\right\\}$
$\displaystyle\leq\mathbb{P}\left\\{\left\lvert\frac{1}{R}\sum_{j=1}^{R}\mathbbm{1}{\left\\{Z_{j}\leq
1\right\\}}-\mathbb{P}\left\\{Z\leq
1\right\\}\right\rvert\geq\frac{k^{-\delta}}{2}\right\\}+\mathbb{P}\left\\{\left\lvert\frac{1}{R}\sum_{j=1}^{R}\mathbbm{1}{\left\\{Z_{j}=0\right\\}}-\mathbb{P}\left\\{Z=0\right\\}\right\rvert\geq\frac{k^{-\delta}}{2}\right\\}.$
(40)

For brevity, denote the two terms in (40) by $T_{1}$ and $T_{2}$,
respectively. Since $Z_{1},\dots,Z_{R}$ are NA (by Lemma 5.4) and
$\mathbbm{1}{\left\\{Z_{j}\leq 1\right\\}}$ is monotonic in $Z_{j}$, by Lemma
5.2(ii), the random variables $\mathbbm{1}{\left\\{Z_{j}\leq 1\right\\}}$ for
$j\in[R]$ are NA. Thus, applying the Chernoff bound in Lemma 5.3 and recalling
that $R=dr\tilde{k}/(\delta\ln k)$, we obtain that for sufficiently large $k$

$\displaystyle T_{1}\leq
2\exp\left(-\frac{R}{3(2k^{\delta})^{2}\mathbb{P}\left\\{Z\leq
1\right\\}}\right)<2\exp\left(-\frac{dr}{24\delta^{2}}\frac{k^{2-\delta}}{\ln^{2}k}\right),$
(41)

where we have used the fact that
$\mathbb{P}\left\\{Z=0\right\\}=k^{-\delta}(1-o(1))$ and
$\mathbb{P}\left\\{Z=1\right\\}=(\delta\ln k)k^{-\delta}(1-o(1))$. Similarly,
using that $\mathbbm{1}{\left\\{Z_{j}=0\right\\}}$ is monotonic in $Z_{j}$, we
deduce that

$\displaystyle T_{2}\leq
2\exp\left(-\frac{R}{3(2k^{\delta})^{2}\mathbb{P}\left\\{Z=0\right\\}}\right)<2\exp\left(-\frac{dr}{24\delta}\frac{k^{2-\delta}}{\ln
k}\right).$ (42)

Combining the upper bounds on $T_{1}$ and $T_{2}$ yields (39). ∎

We now prove (13) from (39). Using $R=dr\tilde{k}/(\delta\ln k)$ and
$\mathbb{P}\left\\{Z=1\right\\}=(\delta\ln k)k^{-\delta}(1-o(1))$, we deduce
from (39) that there exists a non-negative sequence $c_{k,\delta}=o(1)$ such
that for sufficiently large $k$,

$\mathbb{P}\left\\{\sum_{j=1}^{R}\mathbbm{1}{\left\\{Z_{j}=1\right\\}}\leq
dr\tilde{k}\cdot
k^{-\delta}(1-c_{k,\delta})\right\\}<4\exp\left(-\frac{dr}{24\delta^{2}}\frac{k^{2-\delta}}{\ln^{2}k}\right).$
(43)

Since each left node connects to at most $d$ singletons, the number of left
nodes recoverable from just the singletons in the initial graph is at least
$\frac{1}{d}\sum_{j=1}^{R}\mathbbm{1}{\\{Z_{j}=1\\}}$. This number is smaller
than the total number of left nodes recoverable in stage A because additional
singletons may be created during the peeling process. Thus, with ${\sf{A}}$
denoting the total fraction of left nodes (out of $r\tilde{k}$) recovered in
stage A, using (43) we have:

$\displaystyle\mathbb{P}\left\\{{\sf{A}}r\tilde{k}\leq
r\tilde{k}k^{-\delta}(1-c_{k,\delta})\right\\}$
$\displaystyle\leq\mathbb{P}\left\\{\frac{1}{d}\sum_{j=1}^{R}\mathbbm{1}{\left\\{Z_{j}=1\right\\}}\leq
r\tilde{k}k^{-\delta}(1-c_{k,\delta})\right\\}$
$\displaystyle<4\exp\left(-\frac{dr}{24\delta^{2}}\frac{k^{2-\delta}}{\ln^{2}k}\right).$
(44)

Rearranging (5.2) yields (13).

We now prove (14). For $l\in[r\tilde{k}]$, let

$\bar{V}_{l}:=\begin{cases}1&\text{if the $l$-th left node is recovered by the
end of stage A},\\\ 0&\text{otherwise}.\end{cases}$ (45)

Since the recovered left nodes are distributed uniformly at random among all
left nodes, each $\bar{V}_{l}\sim\text{Bern}\left({\sf{A}}\right)$. Moreover,
conditioned on $\sum_{l=1}^{r\tilde{k}}\bar{V}_{l}={\sf{A}}r\tilde{k}$, the
$\\{\bar{V}_{l}\\}$ are NA by Lemma 5.2(iii). Let $\mathcal{S}_{Di}$ and
$\mathcal{S}_{i}$ be the set of indices of the left nodes that represent the
diagonal and above-diagonal entries in the $i$-th nonzero submatrix
(corresponding to $\lambda_{i}\boldsymbol{v}_{i}\boldsymbol{v}_{i}^{T}$),
respectively, with $\left\lvert\mathcal{S}_{Di}\right\rvert=k$ and
$\left\lvert\mathcal{S}_{i}\right\rvert=\binom{k}{2}$.

The number of diagonal entries (out of $k$) recovered in the $i$-th nonzero
submatrix is $N_{Di}=\sum_{l\in\mathcal{S}_{Di}}\bar{V}_{l}$, which can be
bounded as follows using the Chernoff bound in Lemma 5.3. Recalling that
$\alpha^{*}=k^{-\delta}-o(k^{-\delta})$, we have for any fixed
$\epsilon\in(0,\frac{1}{5})$ and sufficiently large $k$,

$\displaystyle\mathbb{P}\left\\{N_{Di}<1\,\middle|\,{\sf{A}}=\alpha^{*}\right\\}$
$\displaystyle\leq\mathbb{P}\left\\{N_{Di}<k\alpha^{*}\epsilon\,\middle|\,{\sf{A}}=\alpha^{*}\right\\}$
$\displaystyle\leq\mathbb{P}\left\\{\sum_{l\in\mathcal{S}_{Di}}\bar{V}_{l}<k\mathbb{E}[\bar{V}_{l}\mid{\sf{A}}]\epsilon\,\middle|\,{\sf{A}}=\alpha^{*}\right\\}$
$\displaystyle\leq\exp\left(-\frac{(1-\epsilon)^{2}k\alpha^{*}}{1+\epsilon}\right)\leq\exp\left(-\frac{1}{2}k^{1-\delta}\right).$
(46)

The fraction of above-diagonal entries (out of $\binom{k}{2}$) recovered in
the $i$-th nonzero submatrix is
${\sf{A}}_{i}=\frac{1}{\binom{k}{2}}\sum_{l\in\mathcal{S}_{i}}\bar{V}_{l}$.
Let $\alpha^{*}_{i}=\alpha^{*}(1-k^{-\frac{1}{4}})$. Then using the Chernoff
bound in Lemma 5.3, we have for sufficiently large $k$,

$\displaystyle\mathbb{P}\left\\{{\sf{A}}_{i}<\alpha^{*}_{i}\,\middle|\,{\sf{A}}=\alpha^{*}\right\\}$
$\displaystyle=\mathbb{P}\left\\{\frac{1}{\binom{k}{2}}\sum_{l\in\mathcal{S}_{i}}\bar{V}_{l}<\alpha^{*}\left(1-k^{-\frac{1}{4}}\right)\,\middle|\,{\sf{A}}=\alpha^{*}\right\\}$
$\displaystyle=\mathbb{P}\left\\{\frac{1}{\binom{k}{2}}\sum_{l\in\mathcal{S}_{i}}\bar{V}_{l}<\mathbb{E}[\bar{V}_{l}\mid{\sf{A}}]\left(1-k^{-\frac{1}{4}}\right)\,\middle|\,{\sf{A}}=\alpha^{*}\right\\}$
$\displaystyle\leq\exp\left(-\frac{\left(k^{-\frac{1}{4}}\right)^{2}\binom{k}{2}\alpha^{*}}{2}\right)\leq\exp\left(-\frac{1}{8}k^{\frac{3}{2}-\delta}\right).$
(47)

Finally, we obtain (14) using (46) and (47) as follows. For sufficiently large
$k$, we have

$\displaystyle\mathbb{P}\left\\{N_{Di}\geq 1\ \text{and}\
{\sf{A}}_{i}\geq\alpha^{*}_{i},\ \forall\ i\in[r]\right\\}$
$\displaystyle\geq\mathbb{P}\left\\{{\sf{A}}\geq\alpha^{*}\right\\}\mathbb{P}\left\\{N_{Di}\geq
1\text{ and }{\sf{A}}_{i}\geq\alpha^{*}_{i},\ \forall\
i\in[r]\mid{\sf{A}}=\alpha^{*}\right\\}$
$\displaystyle\stackrel{{\scriptstyle(\rm{i})}}{{\geq}}\mathbb{P}\left\\{{\sf{A}}\geq\alpha^{*}\right\\}\left[1-\sum_{i=1}^{r}\mathbb{P}\left\\{N_{Di}<1\,\middle|\,{\sf{A}}=\alpha^{*}\right\\}-\sum_{i=1}^{r}\mathbb{P}\left\\{{\sf{A}}_{i}<\alpha_{i}^{*}\,\middle|\,{\sf{A}}=\alpha^{*}\right\\}\right]$
$\displaystyle\stackrel{{\scriptstyle(\rm{ii})}}{{\geq}}\left[1-4\exp\left(-\frac{dr}{24\delta^{2}}\frac{k^{2-\delta}}{\ln^{2}k}\right)\right]\left[1-r\exp\left(-\frac{1}{2}k^{1-\delta}\right)-r\exp\left(-\frac{1}{8}k^{\frac{3}{2}-\delta}\right)\right]$
$\displaystyle\geq 1-2r\exp\left(-\frac{1}{2}k^{1-\delta}\right),$ (48)

where (i) is obtained using a union bound, and (ii) by applying (13), (46) and
(47). ∎

###  5.3 Proof of Lemma 3.2

Since the eigenvectors $\\{\boldsymbol{v}_{1},\ldots,\boldsymbol{v}_{r}\\}$
have non-overlapping supports, the nonzero entries of
$\boldsymbol{X}=\sum_{i=1}^{r}\lambda_{i}\boldsymbol{v}_{i}\boldsymbol{v}_{i}^{T}$
form $r$ disjoint submatrices, each of size $k\times k$. Therefore, the graph
for stage B consists of $r$ disjoint bipartite subgraphs, one corresponding to
each submatrix. Specifically, the $i$-th subgraph has $k$ left nodes
representing the unknown nonzeros in the $i$-th eigenvector and
${\sf{A}}_{i}\binom{k}{2}$ right nodes representing the nonzero pairwise
products recovered in the $i$-th submatrix in stage A. (See Fig. 2(a) for an
example subgraph.) The peeling algorithm described in Section 2.2.2 is applied
to each of these subgraphs.

$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{v}_{4}$$\tilde{v}_{5}$$\tilde{v}_{1}\tilde{v}_{4}$$\tilde{v}_{2}\tilde{v}_{3}$$\tilde{v}_{2}\tilde{v}_{3}$$\tilde{v}_{4}\tilde{v}_{5}$
$k=5$ nonzero entries of $\tilde{\boldsymbol{v}}_{i}$
${\sf{A}}_{i}\binom{k}{2}=4$ right nodes Degree of each right node = 2 (a)
Prior to $t=0$.

$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{v}_{4}$$\tilde{v}_{5}$$\tilde{v}_{1}\tilde{v}_{4}$$\tilde{v}_{2}\tilde{v}_{3}$$\tilde{v}_{2}\tilde{v}_{3}$$\tilde{v}_{4}\tilde{v}_{5}$
$k=5$ supports of $\tilde{\boldsymbol{v}}$ Right nodes connected to
$\tilde{v}_{4}$ have degree 1; others have degree 2 (b) $t=0$.

$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{v}_{4}$$\tilde{v}_{5}$$\tilde{v}_{1}\tilde{v}_{4}$$\tilde{v}_{2}\tilde{v}_{3}$$\tilde{v}_{2}\tilde{v}_{3}$$\tilde{v}_{4}\tilde{v}_{5}$
$k=5$ supports of $\tilde{\boldsymbol{v}}$ Right nodes connected to
$\tilde{v}_{4}$ have degree 1; others have degree 2 (c) $t=1$.

$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{v}_{4}$$\tilde{v}_{5}$$\tilde{v}_{1}\tilde{v}_{4}$$\tilde{v}_{2}\tilde{v}_{3}$$\tilde{v}_{2}\tilde{v}_{3}$$\tilde{v}_{4}\tilde{v}_{5}$
$k=5$ supports of $\tilde{\boldsymbol{v}}$ Right nodes connected to
$\tilde{v}_{4}$ have degree 1; others have degree 2 (d) $t=2$.

Figure 7: Alternative graph process showing the recovery of
$\tilde{\boldsymbol{v}}_{i}$ for
$\boldsymbol{X}=\sum_{i=1}^{r}\tilde{\boldsymbol{v}}_{i}\tilde{\boldsymbol{v}}_{i}^{T}$
given ${\sf{A}}_{i}$ and ${N_{D}}_{i}\geq 1$ (as the counterpart of Fig. 2).
The double occurrence of $\tilde{v}_{2}\tilde{v}_{3}$ as right nodes
highlights that the right nodes are sampled _with replacement_ from the
$\binom{k}{2}$ nonzero pairwise products in the nonzero submatrix
corresponding to $\tilde{\boldsymbol{v}}_{i}\tilde{\boldsymbol{v}}_{i}^{T}$.
The faded nodes and edges are those that have been peeled off.

Consider the $i$-th subgraph, given ${\sf{A}}_{i}$ and the event that
${N_{D}}_{i}\geq 1$. The ${\sf{A}}_{i}\binom{k}{2}$ right nodes can be seen as
drawn uniformly at random without replacement from the $\binom{k}{2}$ nonzero
pairwise products in the $i$-th submatrix. This creates dependence between the
degrees of the left nodes, which makes the evolution of the random graph
process in Fig. 2 difficult to characterize. We therefore consider an
alternative graph process, in which the ${\sf{A}}_{i}\binom{k}{2}$ right nodes
at the start of stage B are drawn uniformly at random _with replacement_ from
the $\binom{k}{2}$ nonzero pairwise products in the $i$-th submatrix; see Fig.
7(a) for an illustration. Stage B of the algorithm proceeds in the same way on
the alternative graph as the original. At $t=0$, a left node is recovered
based on one of the nonzero diagonal entries recovered in stage A (Fig. 7(b)).
Then, for each $t\geq 1$, a degree-1 right node is picked uniformly at random
from the available ones and its connecting left node recovered and peeled off
to update the set of degree-1 right nodes. This continues till there are no
more degree-1 right nodes (Figs. 7(c) and 7(d)).

Since the right nodes in the alternative graph are sampled with replacement,
the number of distinct right nodes in the alternative graph at the start of
stage B can be no larger than ${\sf{A}}_{i}\binom{k}{2}$. Any repeated right
nodes do not play a role in the recovery algorithm. Let
${\sf{B}}_{\text{alt},i}$ denote the fraction of left nodes recovered when
there remain no degree-1 right nodes in the alternative graph process. Since
the number of distinct right nodes in the initial alternative graph is no
larger than that in the original one, we have:

$\mathbb{P}\left\\{{\sf{B}}_{\text{alt},i}<1\,\middle|\,N_{Di},\,{\sf{A}}_{i}\right\\}\geq\mathbb{P}\left\\{{\sf{B}}_{i}<1\,\middle|\,N_{Di},\,{\sf{A}}_{i}\right\\},\quad
i\in[r],\qquad\forall\,N_{Di}\geq 1,\ {\sf{A}}_{i}\geq\alpha^{*}_{i}.$ (49)

We will prove (15) by showing that

$\mathbb{P}\left\\{{\sf{B}}_{\text{alt},i}<1\,\middle|\,N_{Di},\,{\sf{A}}_{i}\right\\}<\exp\left(-\frac{1}{30}k^{1-\delta}\right),\quad
i\in[r],\qquad\forall\,N_{Di}\geq 1,\ {\sf{A}}_{i}\geq\alpha^{*}_{i}.$ (50)

For brevity, we shall refer to ‘prior to $t=0$’ as ‘at $t=-1$’. Given
${\sf{A}}_{i}$, let $\tilde{Z}_{l}$ be the degree of the $l$-th left node in
the alternative graph at $t=-1$ (i.e. Fig. 7(a)), for $l\in[k]$. Since each
right node has degree-2 at $t=-1$, the total number of edges at $t=-1$ is
$\sum_{l=1}^{k}\tilde{Z}_{l}=2{\sf{A}}_{i}\binom{k}{2}$.

######  Lemma 5.6 (Initial left degrees are Binomial).

For $l\in[k]$,

$\tilde{Z}_{l}\sim\text{Bin}\left({\sf{A}}_{i}\binom{k}{2},\frac{2}{k}\right).$
(51)

Let

$\mu:={\sf{A}}_{i}\binom{k}{2}\frac{2}{k}={\sf{A}}_{i}(k-1).$ (52)

Then, for ${\sf{A}}_{i}\geq\alpha^{*}_{i}$, we have for sufficiently large
$k$,

$\displaystyle\mathbb{P}\left\\{\tilde{Z}_{l}\leq\frac{\mu}{2}\right\\}<\exp\left(-\frac{1}{10}k^{1-\delta}\right),\qquad\mathbb{P}\left\\{\tilde{Z}_{l}\geq\frac{3\mu}{2}\right\\}<\exp\left(-\frac{1}{10}k^{1-\delta}\right).$
(53)

###### Proof.

The result (51) follows from a standard counting argument, noting that the
right nodes are picked uniformly at random with replacement and that each
right node has degree 2 at $t=-1$. Result (53) is obtained using the Binomial
tail bound in Lemma 5.1 and simplifying the resulting expression. ∎

Recall that one left node is recovered in each iteration $t\geq 0$ of the
peeling algorithm, until there are no more degree-1 right nodes left. Let
$\pi(t)$ denote the index of the left node recovered in iteration $t$. Thus,
$\pi:\left\\{0,1,\dots,k-1\right\\}\rightarrow\left\\{1,2,\dots,k\right\\}$.
If the decoding process has not terminated after iteration $(t-1)$, for $1\leq
t\leq(k-1)$, then $t$ left nodes out of $k$ have been peeled off. For the
algorithm to proceed, we need at least one of the remaining $(k-t)$ left nodes
to be connected to at least one degree-1 right node. The following lemma
obtains the distribution of the number of degree-1 right nodes connected to
each remaining left node, after each iteration of the peeling algorithm.

######  Lemma 5.7 (Number of degree-1 right nodes connected to each unpeeled
left node).

Consider the alternative graph process on subgraph $i\in[r]$, given
${\sf{A}}_{i}$ and ${N_{D}}_{i}$. Suppose that the peeling algorithm on this
subgraph has not terminated after iteration $(t-1)$, for $1\leq t\leq(k-1)$.
Let $S_{l}^{(t)}$ be the number of degree-1 right nodes connected to the
$l$-th remaining left node at the start of iteration $t$, for
$l\in[k]\setminus\left\\{\pi(0),\ldots,\pi(t-1)\right\\}$. Then

$S_{l}^{(t)}\sim\text{Bin}\left({\sf{A}}_{i}\binom{k}{2},\frac{t}{\binom{k}{2}}\right).$
(54)

Moreover, for ${\sf{A}}_{i}\geq\alpha^{*}_{i}$ and sufficiently large $k$, the
following holds for $1\leq t\leq\frac{k-1}{20}$:

$\mathbb{P}\left\\{S_{l}^{(t)}\geq\frac{\mu}{5}\mid S_{l}^{(t)}\geq
1\right\\}<\exp\left(-\frac{1}{16}k^{1-\delta}\right),$ (55)

where $\mu={\sf{A}}_{i}(k-1)$ as defined in (52).

###### Proof.

For $l\in[k]$ and $j\in[{\sf{A}}_{i}\binom{k}{2}]$, let

$Y_{l,j}:=\begin{cases}1&\text{if the $l$-th left node connects to the $j$-th
right node at $t=-1$},\\\ 0&\text{otherwise}.\end{cases}$ (56)

At $t=-1$, since each right node has degree $2$ (see Fig. 7(a)), the vector
$\left[Y_{1,j},Y_{2,j},\dots,Y_{k,j}\right]$ contains exactly two ones
distributed uniformly at random, and the remaining entries are zero. It
follows that $\sum_{l=1}^{k}Y_{l,j}=2$,
$Y_{l,j}\sim\text{Bern}\left(\frac{2}{k}\right)$ and $Y_{l_{1},j}\cdot
Y_{l_{2},j}\sim\text{Bern}(1/\binom{k}{2})$ for $l_{1}\neq l_{2}$. Moreover,
by construction of the alternative graph, the edges of each right node are
independent from those of the others, so $Y_{l_{1},j_{1}}\perp
Y_{l_{2},j_{2}}$ for $j_{1}\neq j_{2}$.

In iteration $t=0$, the first left node $v_{\pi(0)}$ is peeled off. Consider
one of the $(k-1)$ remaining left nodes after iteration $t=0$, say $v_{l}$.
The number of degree-1 right nodes connected to $v_{l}$ is equal to the number
of right nodes which were connected to $v_{l}$ and $v_{\pi(0)}$ at $t=-1$.
That is, at the start of iteration $t=1$, we have

$\displaystyle S_{l}^{(1)}=\sum_{j=1}^{{\sf{A}}_{i}\binom{k}{2}}Y_{l,j}\cdot
Y_{\pi(0),j}\stackrel{{\scriptstyle\rm{(i)}}}{{\sim}}\text{Bin}\left({\sf{A}}_{i}\binom{k}{2},\frac{1}{\binom{k}{2}}\right),\qquad
l\in[k]\setminus\pi(0),$ (57)

where (i) holds because $Y_{l,j}\cdot
Y_{\pi(0),j}\stackrel{{\scriptstyle\text{iid}}}{{\sim}}\text{Bern}(1/\binom{k}{2})$
for $l\neq\pi(0)$ and $j\in[{\sf{A}}_{i}\binom{k}{2}]$.

More generally, when $t$ left nodes have been peeled off by the end of
iteration $(t-1)$, consider one of the remaining left nodes, say $v_{l}$. Then
the number of degree-1 right nodes connected to $v_{l}$ is the number of right
nodes which were originally connected to $v_{l}$ and some recovered
$v_{\iota}$, i.e., $v_{\iota}$ for
$\iota\in\left\\{\pi(0),\ldots,\pi(t-1)\right\\}$. Hence, we deduce that

$S_{l}^{(t)}=\sum_{j=1}^{{\sf{A}}_{i}\binom{k}{2}}\sum_{\kappa=0}^{t-1}Y_{l,j}\cdot
Y_{\pi(\kappa),j}\stackrel{{\scriptstyle(\rm{i})}}{{\sim}}\text{Bin}\left({\sf{A}}_{i}\binom{k}{2},\frac{t}{\binom{k}{2}}\right),\qquad
l\in[k]\setminus\\{\pi(0),\ldots,\pi(t-1)\\},$ (58)

where (i) is obtained by noting that for
$l\notin\left\\{\pi(0),\ldots,\pi(t-1)\right\\}$ and each $j$, the sum
$\sum_{\kappa=0}^{t-1}Y_{l,j}Y_{\pi(\kappa),j}\sim\text{Bern}(t/\binom{k}{2})$,
and that these sums are mutually independent across
$j\in[{\sf{A}}_{i}\binom{k}{2}]$. This gives the result in (54).

We next prove (55) for ${\sf{A}}_{i}\geq\alpha^{*}_{i}$ using the distribution
in (54). For $1\leq t\leq\frac{k-1}{20}$ and sufficiently large $k$, applying
the Binomial tail bound in Lemma 5.1, we obtain

$\displaystyle\mathbb{P}\left\\{S_{l}^{(t)}\geq\frac{\mu}{5}\mid
S_{l}^{(t)}\geq 1\right\\}$
$\displaystyle\leq\frac{\exp\left(-{\sf{A}}_{i}\binom{k}{2}D_{e}\Big{(}\frac{\mu/5}{{\sf{A}}_{i}\binom{k}{2}}\middle\|\frac{t}{\binom{k}{2}}\Big{)}\right)}{1-\left(1-t/\binom{k}{2}\right)^{{\sf{A}}_{i}\binom{k}{2}}}$
$\displaystyle\leq\frac{\exp\left(-{\sf{A}}_{i}\binom{k}{2}D_{e}\Big{(}\frac{\mu/5}{{\sf{A}}_{i}\binom{k}{2}}\middle\|\frac{(k-1)/20}{\binom{k}{2}}\Big{)}\right)}{1-\left(1-1/\binom{k}{2}\right)^{{\sf{A}}_{i}\binom{k}{2}}}$
$\displaystyle\stackrel{{\scriptstyle(\rm{i})}}{{<}}\frac{\exp\left(-k^{1-\delta}/8\right)}{\alpha^{*}_{i}/2}<\exp\left(-\frac{1}{16}k^{1-\delta}\right).$
(59)

For the inequality (i), the numerator is obtained by simplifying the relative
entropy term and noting that ${\sf{A}}_{i}\geq\alpha^{*}_{i}$; the denominator
is obtained using
$(1-1/\binom{k}{2})^{{\sf{A}}_{i}\binom{k}{2}}\leq(1-1/\binom{k}{2})^{\alpha^{*}_{i}\binom{k}{2}}\leq
e^{-\alpha^{*}_{i}}<1-\alpha^{*}_{i}/2$ (where the last step holds since
$\alpha^{*}_{i}\in(0,1)$). ∎

######  Lemma 5.8 (Degree-1 right nodes don’t run out in initial iterations).

Consider the alternative graph process on subgraph $i\in[r]$, given
${\sf{A}}_{i}\geq\alpha^{*}_{i}$ and ${N_{D}}_{i}\geq 1$. Let
$C_{1}(t)\in\mathbb{N}_{0}$ denote the number of degree-1 right nodes in the
residual subgraph after iteration $t$ of the peeling algorithm. Let
$t_{0}:=\frac{k-1}{20}$, and define the event

$\mathcal{E}:=\left\\{C_{1}(t)>0,\ \text{ for }\ 0\leq t\leq t_{0}\right\\}.$
(60)

Then, for sufficiently large $k$,

$\displaystyle\mathbb{P}\left\\{\mathcal{E}\right\\}>1-\exp\left(-\frac{1}{20}k^{1-\delta}\right).$
(61)

###### Proof.

Recall that the initial left degrees (at $t=-1$), denoted by
$\\{\tilde{Z}_{l}\\}_{l\in[k]}$, each has a Binomial distribution given in
(51). Recall also that $\pi(t)$ is the index of the left node recovered in
iteration $t$. At $t=0$, the first left node $v_{\pi(0)}$ is peeled off,
creating $\tilde{Z}_{\pi(0)}$ degree-1 right nodes. Thus, we have

$C_{1}(0)=\tilde{Z}_{\pi(0)}.$ (62)

More generally, after iteration $(t-1)$, when $t$ left nodes have been peeled
off, consider one of the $(k-t)$ remaining left nodes, say $v_{l}$, for some
$l\in[k]\setminus\left\\{\pi(0),\dots,\pi(t-1)\right\\}$. At the start of
iteration $t$, among the $\tilde{Z}_{l}$ edges of $v_{l}$, the number of edges
connected to degree-1 right nodes is $S_{l}^{(t)}$ whose distribution is given
in (54). When one of the $(k-t)$ remaining left nodes, denoted by
$v_{\pi(t)}$, is peeled off in iteration $t$, the number of degree-1 right
nodes in the residual graph can be expressed as:

$\displaystyle C_{1}(t)$ $\displaystyle=C_{1}(t-1)+\\#\text{ degree-$2$ right
nodes reduced to degree-$1$ due to the removal of }v_{\pi(t)}$
$\displaystyle\quad-\\#\text{ degree-$1$ right nodes reduced to degree-$0$ due
to the removal of }v_{\pi(t)}$
$\displaystyle=C_{1}(t-1)+\left(\tilde{Z}_{\pi(t)}-S_{\pi(t)}^{(t)}\right)-S_{\pi(t)}^{(t)}$
$\displaystyle=C_{1}(t-1)+\left(\tilde{Z}_{\pi(t)}-2S_{\pi(t)}^{(t)}\right).$
(63)

The marginal distributions of $\tilde{Z}_{\pi(t)}$ and $S_{\pi(t)}^{(t)}$ are
given by Lemmas 5.6 and 5.7, which also guarantee that with high probability,
$\tilde{Z}_{\pi(t)}>\frac{\mu}{2}$ and $S_{\pi(t)}^{(t)}<\frac{\mu}{5}$ for
$t\leq t_{0}$.

By the same reasoning as above, we can write $C_{1}(t-1)$ in terms of
$C_{1}(t-2)$ and so on to obtain

$\displaystyle C_{1}(t)$
$\displaystyle=C_{1}(t-1)+\left(\tilde{Z}_{\pi(t)}-2S_{\pi(t)}^{(t)}\right)$
$\displaystyle=C_{1}(0)+\sum_{\kappa=1}^{t}\left(\tilde{Z}_{\pi(\kappa)}-2S_{\pi(\kappa)}^{(\kappa)}\right)$
$\displaystyle=\tilde{Z}_{\pi(0)}+\sum_{\kappa=1}^{t}\left(\tilde{Z}_{\pi(\kappa)}-2S_{\pi(\kappa)}^{(\kappa)}\right),$
(64)

where the last step follows from (62). We next show that $C_{1}(t)>0$ with
high probability for $0\leq t\leq t_{0}$ by bounding each term on the RHS of
(64). From (64), we note that

$\displaystyle\left\\{\tilde{Z}_{\pi(t)}>\frac{\mu}{2}\;\text{and}\;S_{\pi(t)}^{(t)}<\frac{\mu}{5},\
\forall\ t\leq t_{0}\right\\}$
$\displaystyle\Rightarrow\left\\{C_{1}(t)>\frac{\mu}{2}+t\left(\frac{\mu}{2}-\frac{2\mu}{5}\right),\
\forall\ t\in[0,t_{0}]\right\\}$ $\displaystyle\Rightarrow\left\\{C_{1}(t)>0,\
\forall\ t\in[0,t_{0}]\right\\}.$ (65)

It follows that

$\displaystyle\mathbb{P}\left\\{\mathcal{E}\right\\}$
$\displaystyle=\mathbb{P}\left\\{C_{1}(t)>0,\ \forall\ t\in[0,t_{0}]\right\\}$
$\displaystyle\geq\mathbb{P}\left\\{\tilde{Z}_{\pi(t)}>\frac{\mu}{2}\;\text{and}\;S_{\pi(t)}^{(t)}<\frac{\mu}{5},\
\forall\ t\leq t_{0}\,\middle|\,S_{\pi(t)}^{(t)}\geq 1,\ \forall\
t\in[1,t_{0}]\right\\},$ (66)

where the condition $S_{\pi(t)}^{(t)}\geq 1,\ \forall\ t\in[1,t_{0}]$ captures
the fact that the left node recovered in each iteration $t$ was connected to
at least one degree-1 right node at $t$. Applying a union bound gives

$\mathbb{P}\left\\{\mathcal{E}\right\\}\geq
1-\sum_{t=0}^{t_{0}}\mathbb{P}\left\\{\tilde{Z}_{\pi(t)}\leq\frac{\mu}{2}\right\\}-\sum_{t=1}^{t_{0}}\mathbb{P}\left\\{S_{\pi(t)}^{(t)}\geq\frac{\mu}{5}\mid
S_{\pi(t)}^{(t)}\geq 1\right\\}.$ (67)

Each summand in the second term can be bounded using (53) and each summand in
the third term can be bounded using (55) for $t\leq t_{0}=\frac{k-1}{20}$. We
therefore obtain

$\displaystyle\mathbb{P}\left\\{\mathcal{E}\right\\}$
$\displaystyle>1-(t_{0}+1)\exp\left(-\frac{1}{10}k^{1-\delta}\right)-t_{0}\exp\left(-\frac{1}{16}k^{1-\delta}\right)$
$\displaystyle>1-\left(\frac{k-1}{20}+1\right)\exp\left(-\frac{1}{10}k^{1-\delta}\right)-\frac{k-1}{20}\exp\left(-\frac{1}{16}k^{1-\delta}\right)$
$\displaystyle>1-\exp\left(-\frac{1}{20}k^{1-\delta}\right),$ (68)

where the last inequality holds for sufficiently large $k$. ∎

We now prove Lemma 3.2 by proving (50) based on Lemmas 5.7 and 5.8. In
particular, we show that after $t_{0}=\frac{k-1}{20}$ iterations, every
remaining left node in the residual subgraph is connected to at least one
degree-1 right node with high probability and so can be recovered by the
peeling algorithm.

###### Proof of Lemma 3.2.

Consider the alternative graph process on subgraph $i\in[r]$. Conditioned on
the event

$\mathcal{E}=\left\\{C_{1}(t)>0,\ \forall\ t\in[0,t_{0}]\right\\},$

$(t_{0}+1)$ left nodes (out of $k$) have been recovered and peeled off by
iteration $t_{0}$. Consider one of the remaining left nodes after iteration
$t_{0}$, i.e., $v_{l}$ with
$l\in[k]\setminus\left\\{\pi(0),\dots,\pi(t_{0})\right\\}$. We show that with
high probability, $v_{l}$ is connected to at least one degree-1 right node and
is therefore recoverable in the forthcoming steps. From Lemma 5.7, the number
of degree-1 right nodes connected to $v_{l}$ after iteration $t_{0}$, denoted
by $S_{l}^{(t_{0}+1)}$, is distributed as
$S_{l}^{(t_{0}+1)}\sim\text{Bin}\left({\sf{A}}_{i}\binom{k}{2},(t_{0}+1)/\binom{k}{2}\right)$.
Thus, using $t_{0}=\frac{k-1}{20}$ and ${\sf{A}}_{i}\geq\alpha_{i}^{*}$, we
have the bound

$\displaystyle\mathbb{P}\left\\{S_{l}^{(t_{0}+1)}=0\mid\mathcal{E}\right\\}=\left(1-\frac{t_{0}+1}{\binom{k}{2}}\right)^{{\sf{A}}_{i}\binom{k}{2}}\leq\exp\left(-(t_{0}+1)\alpha^{*}_{i}\right)<\exp\left(-\frac{1}{25}k^{1-\delta}\right).$
(69)

Applying a union bound gives

$\displaystyle\mathbb{P}\left\\{S_{l}^{(t_{0}+1)}\geq 1,\ \forall\
l\in[k]\setminus\left\\{\pi(0),\ldots,\pi(t_{0})\right\\}\mid\mathcal{E}\right\\}$
$\displaystyle>1-\left(k-(t_{0}+1)\right)\exp\left(-\frac{1}{25}k^{1-\delta}\right)$
$\displaystyle>1-k\exp\left(-\frac{1}{25}k^{1-\delta}\right).$ (70)

We recall that ${\sf{B}}_{\text{alt},i}$ is the fraction of left nodes
recovered when degree-1 right nodes run out in the alternative graph process
on subgraph $i\in[r]$. Then, using (61) and (5.3), for all $N_{Di}\geq 1$ and
${\sf{A}}_{i}\geq\alpha^{*}_{i}$, we have for sufficiently large $k$:

$\displaystyle\mathbb{P}\left\\{{\sf{B}}_{\text{alt},i}=1\,\middle|\,N_{Di},\,{\sf{A}}_{i}\right\\}$
$\displaystyle\geq\mathbb{P}\left\\{\mathcal{E}\right\\}\mathbb{P}\left\\{S_{l}^{(t_{0}+1)}\geq
1,\ \forall\
l\in[k]\setminus\left\\{\pi(0),\ldots,\pi(t_{0})\right\\}\mid\mathcal{E}\right\\}$
$\displaystyle>\left[1-\exp\left(-\frac{1}{20}k^{1-\delta}\right)\right]\left[1-k\exp\left(-\frac{1}{25}k^{1-\delta}\right)\right]$
$\displaystyle>1-\exp\left(-\frac{1}{30}k^{1-\delta}\right).$ (71)

This proves (50), and the result of Lemma 3.2 follows via (49). ∎

###  5.4 Proof of Lemma 3.3

At the start of stage A, the pruned bipartite graph for a non-symmetric matrix
$\boldsymbol{X}$ as defined in Lemma 3.3 has $r\beta k^{2}$ left nodes and
$R=dr\beta k^{2}/(\delta\ln k)$ right nodes. The proof of (17) is along the
same lines as that of (14). Specifically, we first obtain a concentration
inequality similar to (39) on the number of singleton right nodes in the
initial graph. The subsequent steps are very similar to (43)–(5.2) and
(47)–(5.2) and are omitted for brevity.

###  5.5 Proof of Lemma 3.4

Recall that the vectors in each of the sets
$\\{\boldsymbol{u}_{1},\ldots,\boldsymbol{u}_{r}\\}$ and
$\\{\boldsymbol{v}_{1},\ldots,\boldsymbol{v}_{r}\\}$ have $k$ and $\beta k$
disjoint supports, respectively, for some constant $\beta\in(0,1]$. Thus, the
nonzero entries of
$\boldsymbol{X}=\sum_{i=1}^{r}\sigma_{i}\boldsymbol{u}_{i}\boldsymbol{v}_{i}^{T}=\sum_{i=1}^{r}\tilde{\boldsymbol{u}}_{i}\tilde{\boldsymbol{v}}_{i}^{T}$
form $r$ disjoint submatrices, each of size $k\times\beta k$. The initial
graph for stage B consists of $r$ disjoint bipartite subgraphs. The peeling
algorithm in stage B is applied to each of these subgraphs. In the $i$-th
subgraph, the $(k+\beta k)$ left nodes represent the unknown nonzeros in
$\tilde{\boldsymbol{u}}_{i}$ and $\tilde{\boldsymbol{v}}_{i}$. The
${\sf{A}}_{i}\beta k^{2}$ right nodes represent the nonzero pairwise products
recovered in the $i$-th nonzero submatrix in stage A.

As in the proof of the symmetric case, we will consider an alternative initial
graph in which for each subgraph $i$, the ${\sf{A}}_{i}\beta k^{2}$ right
nodes are drawn uniformly at random _with replacement_ from the set of $\beta
k^{2}$ nonzero pairwise products in the $i$-th nonzero submatrix. See Fig.
8(a) for an illustration. The number of distinct right nodes for the
alternative initial graph can be no larger than that for the original one.
Moreover, as illustrated in Fig. 8, following the recovery of the first
nonzero entry of $\tilde{\boldsymbol{u}}_{i}$, the process on the alternative
graph peels off all the recoverable entries of $\tilde{\boldsymbol{v}}_{i}$
before peeling any nonzero entries of $\tilde{\boldsymbol{u}}_{i}$ that have
become recoverable along the way. The order of peeling does not affect the
total number of nonzeros recoverable in $\tilde{\boldsymbol{u}}_{i}$ and
$\tilde{\boldsymbol{v}}_{i}$ before degree-1 right nodes run out.

Writing ${\sf{B}}^{u}_{\text{alt},i}$ and ${\sf{B}}^{v}_{\text{alt},i}$ for
the fractions of nonzeros recoverable in $\tilde{\boldsymbol{u}}_{i}$ and
$\tilde{\boldsymbol{v}}_{i}$ in the alternative process, we therefore have

$\mathbb{P}\left\\{\\{{\sf{B}}^{u}_{\text{alt},i}<1\\}\,\cup\,\\{{\sf{B}}^{v}_{\text{alt},i}<1\\}\mid{\sf{A}}_{i}\right\\}\geq\mathbb{P}\left\\{\\{{\sf{B}}^{u}_{i}<1\\}\,\cup\,\\{{\sf{B}}^{v}_{i}<1\\}\mid{\sf{A}}_{i}\right\\}.$
(72)

To prove (18), we will show that for ${\sf{A}}_{i}\geq\alpha_{i}^{*}$ and
sufficiently large $k$:

$\mathbb{P}\left\\{\\{{\sf{B}}^{u}_{\text{alt},i}<1\\}\,\cup\,\\{{\sf{B}}^{v}_{\text{alt},i}<1\\}\mid{\sf{A}}_{i}\right\\}\leq\exp\left(-\frac{\beta}{8}k^{1-2\delta}\right),\quad
i\in[r].$ (73)

$\tilde{u}_{1}$$\tilde{u}_{2}$$\tilde{u}_{3}$$\tilde{u}_{4}$$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{u}_{1}\tilde{v}_{1}$$\tilde{u}_{4}\tilde{v}_{1}$$\tilde{u}_{3}\tilde{v}_{2}$$\tilde{u}_{1}\tilde{v}_{3}$$\tilde{u}_{3}\tilde{v}_{2}$
$k=4$ nonzero entries of $\tilde{\boldsymbol{u}}_{i}$ $\beta k=3$ nonzero
entries of $\tilde{\boldsymbol{v}}_{i}$ ${\sf{A}}_{i}\beta k^{2}$ right nodes
Degree of each right node = 2 (a) Prior to $t=0$.

$\tilde{u}_{1}$$\tilde{u}_{2}$$\tilde{u}_{3}$$\tilde{u}_{4}$$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{u}_{1}\tilde{v}_{1}$$\tilde{u}_{4}\tilde{v}_{1}$$\tilde{u}_{3}\tilde{v}_{2}$$\tilde{u}_{1}\tilde{v}_{3}$$\tilde{u}_{3}\tilde{v}_{2}$
Right nodes connected to $\tilde{u}_{1}$ have degree 1; others have degree 2
(b) $t=0$.

$\tilde{u}_{1}$$\tilde{u}_{2}$$\tilde{u}_{3}$$\tilde{u}_{4}$$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{u}_{1}\tilde{v}_{1}$$\tilde{u}_{4}\tilde{v}_{1}$$\tilde{u}_{3}\tilde{v}_{2}$$\tilde{u}_{1}\tilde{v}_{3}$$\tilde{u}_{3}\tilde{v}_{2}$
Right nodes connected to $\tilde{u}_{1}$ have degree 1; others have degree 2
(c) $t=1$.

$\tilde{u}_{1}$$\tilde{u}_{2}$$\tilde{u}_{3}$$\tilde{u}_{4}$$\tilde{v}_{1}$$\tilde{v}_{2}$$\tilde{v}_{3}$$\tilde{u}_{1}\tilde{v}_{1}$$\tilde{u}_{4}\tilde{v}_{1}$$\tilde{u}_{3}\tilde{v}_{2}$$\tilde{u}_{1}\tilde{v}_{3}$$\tilde{u}_{3}\tilde{v}_{2}$
Right nodes connected to $\tilde{u}_{1}$ have degree 1; others have degree 2
(d) $t=2$.

Figure 8: Alternative graph process showing the recovery of
$\tilde{\boldsymbol{u}}_{i}$ and $\tilde{\boldsymbol{v}}_{i}$ given
${\sf{A}}_{i}$ (as the counterpart of Fig. 3). The double occurrence of
$\tilde{u}_{3}\tilde{v}_{2}$ as right nodes highlights that the right nodes
are sampled _with replacement_ from the $\beta k^{2}$ nonzero pairwise
products in the nonzero submatrix corresponding to
$\tilde{\boldsymbol{u}}_{i}\tilde{\boldsymbol{v}}_{i}^{T}$. (b): At $t=0$,
$\tilde{u}_{1}$ is recovered as 1 and peeled off. (c)–(d): Following the
recovery of $\tilde{u}_{1}$, we peel off all the recoverable $\tilde{v}_{j}$’s
(i.e., $\tilde{v}_{1},\tilde{v}_{3}$) before peeling any $\tilde{u}_{l}$’s
that have become recoverable along the way (i.e., $\tilde{u}_{4}$).

The following lemma specifies the degree distribution of each left node at the
beginning of stage B. The proof is along the same lines as that of Lemma 5.6
and is omitted.

######  Lemma 5.9 (Initial left degrees are Binomial).

Let $Z^{u}_{l}$ be the degree of the left node representing the $l$-th nonzero
of $\tilde{\boldsymbol{u}}_{i}$, and $Z^{v}_{j}$ be the degree of the left
node representing the $j$-th nonzero of $\tilde{\boldsymbol{v}}_{i}$, for
$l\in[k]$ and $j\in[\beta k]$. Then,

$\displaystyle Z^{u}_{l}\sim\text{Bin}\left({\sf{A}}_{i}\beta
k^{2},\frac{1}{k}\right),\quad Z^{v}_{j}\sim\text{Bin}\left({\sf{A}}_{i}\beta
k^{2},\frac{1}{\beta k}\right).$ (74)

Then, for ${\sf{A}}_{i}\geq\alpha^{*}_{i}$ and sufficiently large $k$:

$\displaystyle\mathbb{P}\left\\{Z^{u}_{l}\notin\left[\frac{{\sf{A}}_{i}\beta
k}{2},\,\frac{3{\sf{A}}_{i}\beta
k}{2}\right]\,\middle|\,{\sf{A}}_{i}\right\\}<2\exp\left(-\frac{\beta}{10}k^{1-\delta}\right),\quad$
(75)
$\displaystyle\mathbb{P}\left\\{Z^{v}_{j}\notin\left[\frac{{\sf{A}}_{i}k}{2},\,\frac{3{\sf{A}}_{i}k}{2}\right]\,\middle|\,{\sf{A}}_{i}\right\\}<2\exp\left(-\frac{1}{10}k^{1-\delta}\right).$
(76)

Define the event

$\mathcal{T}_{0}:=\left\\{Z^{u}_{l}\in\left[\frac{{\sf{A}}_{i}\beta
k}{2},\,\frac{3{\sf{A}}_{i}\beta k}{2}\right],\
\forall\,l\in[k]\right\\}\cap\left\\{Z^{v}_{j}\in\left[\frac{{\sf{A}}_{i}k}{2},\,\frac{3{\sf{A}}_{i}k}{2}\right],\
\forall\,j\in[\beta k]\right\\}.$ (77)

Then for ${\sf{A}}_{i}\geq\alpha^{*}_{i}$ and sufficiently large $k$:

$\mathbb{P}\left\\{\mathcal{T}_{0}\mid{\sf{A}}_{i}\right\\}>1-2k\exp\left(-\frac{\beta}{10}k^{1-\delta}\right)-2\beta
k\exp\left(-\frac{1}{10}k^{1-\delta}\right).$ (78)

The recovery of each nonzero entry of $\tilde{\boldsymbol{u}}_{i}$ or
$\tilde{\boldsymbol{v}}_{i}$ is counted as one iteration of the algorithm. The
recovery of the first nonzero entry of $\tilde{\boldsymbol{u}}_{i}$ or
$\tilde{\boldsymbol{v}}_{i}$ corresponds to iteration $t=0$. The next lemma
specifies the distribution of degree-1 right nodes connected to each remaining
left node after each iteration of the algorithm. The proof is similar to that
of Lemma 5.7 and is omitted for brevity.

######  Lemma 5.10 (Number of degree-1 right nodes connected to each unpeeled
left node).

Suppose that the peeling algorithm on the $i$-th subgraph has not terminated
after iteration $(t-1)$, for $1\leq t\leq(k+\beta k-1)$. Let
$\mathcal{U}(t-1)$ and $\mathcal{V}(t-1)$ be the set of indices of the
nonzeros recovered in $\tilde{\boldsymbol{u}}_{i}$ and
$\tilde{\boldsymbol{v}}_{i}$ by the end of iteration $(t-1)$. At the start of
iteration $t$, let $S^{(u,t)}_{l}$ be the number of degree-1 right nodes
connected to the remaining left node $\tilde{u}_{l}$ and let $S^{(v,t)}_{j}$
be that connected to the remaining left node $\tilde{v}_{j}$, where
$l\in[k]\setminus\mathcal{U}(t-1)$ and $j\in[\beta
k]\setminus\mathcal{V}(t-1)$. Then,

$S^{(u,t)}_{l}\sim\text{Bin}\left({\sf{A}}_{i}\beta
k^{2},\frac{\left\lvert\mathcal{V}(t-1)\right\rvert}{\beta k^{2}}\right),\quad
S^{(v,t)}_{j}\sim\text{Bin}\left({\sf{A}}_{i}\beta
k^{2},\frac{\left\lvert\mathcal{U}(t-1)\right\rvert}{\beta k^{2}}\right).$
(79)

In the following, we refer to the left nodes representing the nonzero entries
of $\tilde{\boldsymbol{u}}_{i}$ as $\tilde{u}$-left nodes, and those
representing nonzero entries of $\tilde{\boldsymbol{v}}_{i}$ as
$\tilde{v}$-left nodes. We assume without loss of generality that the peeling
algorithm is initialized by setting a $\tilde{u}$-left node to 1. As
illustrated in Fig. 8(b), following the recovery of the first $\tilde{u}$-left
node, all the right nodes connected to this node reduce to degree-1. The next
lemma gives a high probability lower bound on the number of degree-1 right
nodes created.

######  Lemma 5.11 (Recovery of the first left node).

Denote by $N$ the number of distinct $\tilde{v}$-left nodes connected to the
degree-1 right nodes created by the recovery of the first $\tilde{u}$-left
node. Then, there exists $b_{k,\delta}=o(k^{1-\delta})$ such that for
${\sf{A}}_{i}\geq\alpha^{*}_{i}$ and sufficiently large $k$, we have

$\mathbb{P}\left\\{N\geq\frac{\beta}{2}k^{1-\delta}-b_{k,\delta}\,\middle|\,\mathcal{T}_{0},{\sf{A}}_{i}\right\\}\geq
1-2\exp\left(-\frac{\beta}{2}k^{1-\frac{3}{2}\delta}\right).$ (80)

Here $\mathcal{T}_{0}$ is the event defined in (77).

###### Proof.

Let $Z^{u}$ be the degree of the first $\tilde{u}$-left node recovered. Then
the recovery of this $\tilde{u}$-left node creates $Z^{u}$ degree-1 right
nodes. Consider a $\tilde{v}$-left node, say $\tilde{v}_{j}$, and let
$S^{v}_{j}$ be the number of degree-1 right nodes connected to $\tilde{v}_{j}$
after the first $\tilde{u}$-left node is peeled off. Then, conditioned on
$Z^{u}$, we have

$S^{v}_{j}:=\sum_{s=1}^{Z^{u}}\mathbbm{1}{\\{\text{The $s$-th degree-1 right
node is connected to
}\tilde{v}_{j}\\}}\sim\text{Bin}\left(Z^{u},\frac{1}{\beta
k}\right),\quad\text{for }j\in[\beta k],$ (81)

and we can express $N$ (defined in the lemma statement) as

$N:=\sum_{j=1}^{\beta k}\mathbbm{1}{\\{S^{v}_{j}\geq 1\\}}.$ (82)

Let $V_{s}$ be the index of the $\tilde{v}$-left node that the $s$-th degree-1
right node is connected to, for $s\in[Z^{u}]$. Let $N_{0}=\mathbb{E}[N\mid
Z^{u},{\sf{A}}_{i}]$ and

$N_{s}=\mathbb{E}[N\mid V_{1},\dots,V_{s},\,Z^{u},{\sf{A}}_{i}],\qquad
s\in[Z^{u}].$ (83)

Conditioned on $Z^{u}$, the sequence $\\{N_{s}\\}$ is a Doob Martingale, with
$N_{Z^{u}}=N$. Writing $N=f(V_{1},\dots,V_{Z^{u}})$, we note that $f$ is
1-Lipschitz since changing the connection of a single right node cannot change
$N$ by more than one. Moreover, by the construction of the alternative graph,
$V_{1},\ldots,V_{Z^{u}}$ are independent given ($Z^{u},{\sf{A}}_{i}$).
Therefore, by McDiarmid’s inequality [43, Sec. 13.5], we have

$\mathbb{P}\left\\{\left\lvert N-\mathbb{E}[N\mid
Z^{u},{\sf{A}}_{i}]\right\rvert\geq
k^{-\frac{\delta}{4}}Z^{u}\,\,\middle|\,\,Z^{u},\,{\sf{A}}_{i}\right\\}\leq
2\exp\left(-2k^{-\frac{\delta}{2}}Z^{u}\right).$ (84)

Hence

$\mathbb{P}\left\\{N\,\geq\,\mathbb{E}[N\mid
Z^{u},{\sf{A}}_{i}]-k^{-\frac{\delta}{4}}Z^{u}\,\,\middle|\,\,Z^{u},\,{\sf{A}}_{i}\right\\}\geq
1-2\exp\left(-2k^{-\frac{\delta}{2}}Z^{u}\right).$ (85)

From (81) and (82), we obtain that

$\mathbb{E}[N\mid Z^{u},{\sf{A}}_{i}]=\beta k\,\mathbb{P}\left\\{S^{v}_{j}\geq
1\right\\}=\beta k\left[1-\left(1-\frac{1}{\beta
k}\right)^{Z^{u}}\right]\geq\beta k\left(1-\exp\Big{(}-\frac{Z^{u}}{\beta
k}\Big{)}\right).$ (86)

Furthermore, given ${\sf{A}}_{i}\geq\alpha_{i}^{*}=k^{-\delta}-o(k^{-\delta})$
and conditioning on $\mathcal{T}_{0}$, we have

$Z^{u}\geq\frac{\beta}{2}k^{1-\delta}-o(k^{1-\delta}),\qquad\mathbb{E}[N\mid
Z^{u},{\sf{A}}_{i}]-k^{-\frac{\delta}{4}}Z^{u}\geq\frac{\beta}{2}k^{1-\delta}-o(k^{1-\delta}).$
(87)

Using (87) in (85) gives the result in (80). ∎

###### Proof of Lemma 3.4.

With $N$ as defined in Lemma 5.11, at the end of iteration $N$ of the peeling
algorithm, one $\tilde{u}$-left node and $N$ $\tilde{v}$-left nodes have been
recovered. Adopting the notation in Lemma 5.10, this means
$\left\lvert\mathcal{U}(N)\right\rvert=1$ and
$\left\lvert\mathcal{V}(N)\right\rvert=N$. Lemma 5.10 implies that at the
start of iteration $(N+1)$, the number of degree-1 right nodes connected to
each remaining $\tilde{u}$-left node is

$S^{(u,N+1)}_{l}\sim\text{Bin}\left({\sf{A}}_{i}\beta k^{2},\frac{N}{\beta
k^{2}}\right),\quad\text{for }l\in[k]\setminus\mathcal{U}(N).$ (88)

Recall the definition of $\mathcal{T}_{0}$ from (77) and let

$\mathcal{T}_{1}:=\left\\{N\geq\frac{\beta}{2}k^{1-\delta}-b_{k,\delta}\right\\},\qquad\mathcal{T}_{2}:=\left\\{S^{(u,N+1)}_{l}\geq
1,\ \forall\,l\in[k]\setminus\mathcal{U}(N)\right\\}.$ (89)

Here $b_{k,\delta}=o(k^{1-\delta})$ is the same as in Lemma 5.11. Note that
$\mathcal{T}_{2}$ is the event that at the start of iteration $(N+1)$, each of
the $(k-1)$ remaining $\tilde{u}$-left nodes is connected to at least one
degree-1 right node and can therefore be recovered. Analogously to (69) and
(5.3), we can show that for ${\sf{A}}_{i}\geq\alpha^{*}_{i}$ and sufficiently
large $k$,

$\mathbb{P}\left\\{\mathcal{T}_{2}\mid\mathcal{T}_{0},\,\mathcal{T}_{1},\,{\sf{A}}_{i}\right\\}>1-(k-1)\exp\left(-\frac{\beta}{4}k^{1-2\delta}\right).$
(90)

Since the recovery of each left node corresponds to one iteration, conditioned
on the event $\mathcal{T}_{2}$, we have
$\left\lvert\mathcal{U}(N+k-1)\right\rvert=k$ and
$\left\lvert\mathcal{V}(N+k-1)\right\rvert=N$.

At the start of iteration $(N+k)$, all the left nodes that remain are
$\tilde{v}$-left nodes, so all the right nodes that remain are degree-1.
Moreover, for ${\sf{A}}_{i}\geq\alpha^{*}_{i}$ and sufficiently large $k$,
conditioning on event $\mathcal{T}_{0}$ ensures that each $\tilde{v}$-left
node is connected to at least one right node. Therefore, all the remaining
$\tilde{v}$-left nodes can be recovered. If the events
$\mathcal{T}_{0},\mathcal{T}_{1},\mathcal{T}_{2}$ all hold, then the algorithm
successfully recovers all the nonzeros of $\tilde{\boldsymbol{u}}_{i}$ and
$\tilde{\boldsymbol{v}}_{i}$. Indeed, we have for
${\sf{A}}_{i}\geq\alpha^{*}_{i}$ and sufficiently large $k$:

$\displaystyle\mathbb{P}\left\\{\\{{\sf{B}}^{u}_{\text{alt},i}=1\\}\,\cap\,\\{{\sf{B}}^{v}_{\text{alt},i}=1\\}\,\mid\,{\sf{A}}_{i}\right\\}$
$\displaystyle\geq\mathbb{P}\left\\{\mathcal{T}_{0},\mathcal{T}_{1},\mathcal{T}_{2}\mid{\sf{A}}_{i}\right\\}$
$\displaystyle=\mathbb{P}\left\\{\mathcal{T}_{0}\mid{\sf{A}}_{i}\right\\}\mathbb{P}\left\\{\mathcal{T}_{1}\mid\mathcal{T}_{0},{\sf{A}}_{i}\right\\}\mathbb{P}\left\\{\mathcal{T}_{2}\mid\mathcal{T}_{0},\mathcal{T}_{1},{\sf{A}}_{i}\right\\}$