Joint Sparse Recovery Using Signal Space Matching Pursuit

*Junhan Kim, ^†Jian Wang, *Luong Trung Nguyen, and *Byonghyo Shim, ^∗Department of Electrical and Computer Engineering, Seoul National University, Seoul, Korea ^†School of Data Science, Fudan University, Shanghai, China Email: ^∗{junhankim, ltnguyen, bshim}@islab.snu.ac.kr, ^†[email protected] This work was supported in part by the MSIT (Ministry of Science and ICT), Korea, under the ITRC (Information Technology Research Center) support program (IITP-2020-2017-0-01637) supervised by the IITP (Institute for Information & communications Technology Promotion) and in part by the Samsung Research Funding & Incubation Center for Future Technology of Samsung Electronics under Grant SRFC-IT1901-17.This paper was presented in part at the IEEE International Symposium on Information Theory, Colorado, USA, June, 2018 [1]. (Junhan Kim and Jian Wang equally contributed to this work.)(Corresponding author: Byonghyo Shim.)

Abstract

In this paper, we put forth a new joint sparse recovery algorithm called signal space matching pursuit (SSMP). The key idea of the proposed SSMP algorithm is to sequentially investigate the support of jointly sparse vectors to minimize the subspace distance to the residual space. Our performance guarantee analysis indicates that SSMP accurately reconstructs any row $K$ -sparse matrix of rank $r$ in the full row rank scenario if the sampling matrix $\mathbf{A}$ satisfies $\operatornamewithlimits{krank}(\mathbf{A})\geq K+1$ , which meets the fundamental minimum requirement on $\mathbf{A}$ to ensure exact recovery. We also show that SSMP guarantees exact reconstruction in at most $K-r+\lceil\frac{r}{L}\rceil$ iterations, provided that $\mathbf{A}$ satisfies the restricted isometry property (RIP) of order $L(K-r)+r+1$ with

\delta_{L(K-r)+r+1}<\max\left\{\frac{\sqrt{r}}{\sqrt{K+\frac{r}{4}}+\sqrt{\frac{r}{4}}},\frac{\sqrt{L}}{\sqrt{K}+1.15\sqrt{L}}\right\},

where $L$ is the number of indices chosen in each iteration. This implies that the requirement on the RIP constant becomes less restrictive when $r$ increases. Such behavior seems to be natural but has not been reported for most of conventional methods. We also show that if $r=1$ , then by running more than $K$ iterations, the performance guarantee of SSMP can be improved to $\delta_{\lfloor 7.8K\rfloor}\leq 0.155$ . Furthermore, we show that under a suitable RIP condition, the reconstruction error of SSMP is upper bounded by a constant multiple of the noise power, which demonstrates the robustness of SSMP to measurement noise. Finally, from extensive numerical experiments, we show that SSMP outperforms conventional joint sparse recovery algorithms both in noiseless and noisy scenarios.

Index Terms:

Joint sparse recovery, multiple measurement vectors (MMV), subspace distance, rank aware order recursive matching pursuit (RA-ORMP), orthogonal least squares (OLS), restricted isometry property (RIP)

I Introduction

In recent years, sparse signal recovery has received considerable attention in image processing, seismology, data compression, source localization, wireless communication, machine learning, to name just a few [3, 2, 4, 5]. The main goal of sparse signal recovery is to reconstruct a high dimensional $K$ -sparse vector $\mathbf{x}\in\mathbb{R}^{n}$ ( $\|\mathbf{x}\|_{0}\leq K\ll n$ where $\|\mathbf{x}\|_{0}$ denotes the number of nonzero elements in $\mathbf{x}$ ) from its compressed linear measurements

\mathbf{y}=\mathbf{A}\mathbf{x},

(1)

where $\mathbf{A}\in\mathbb{R}^{m\times n}$ ( $m\ll n$ ) is the sampling (sensing) matrix. In various applications, such as wireless channel estimation [6, 5], sub-Nyquist sampling of multiband signals [7, 8], angles of departure and arrival (AoD and AoA) estimation in mmWave communication systems [9], and brain imaging [10], we encounter a situation where multiple measurement vectors (MMV) of a group of jointly sparse vectors are available. By the jointly sparse vectors, we mean multiple sparse vectors having a common support (the index set of nonzero entries). In this situation, one can dramatically improve reconstruction accuracy by recovering all the desired sparse vectors simultaneously [11, 12]. The problem to reconstruct a group $\{\mathbf{x}_{i}\}_{i=1}^{r}$ of jointly $K$ -sparse vectors¹¹1In the sequel, we assume that $\mathbf{x}_{1},{\color[rgb]{0,0,1}{\ldots}},\mathbf{x}_{r}$ are linearly independent. is often referred to as the joint sparse recovery problem [17]. Let $\mathbf{y}_{i}=\mathbf{A}\mathbf{x}_{i}$ be the measurement vector of $\mathbf{x}_{i}$ acquired through the sampling matrix $\mathbf{A}$ . Then the system model describing the MMV can be expressed as

\mathbf{Y}=\mathbf{A}\mathbf{X},

(2)

where $\mathbf{Y}=[\mathbf{y}_{1},\ldots,\mathbf{y}_{r}]$ and $\mathbf{X}=[\mathbf{x}_{1},\ldots,\mathbf{x}_{r}]$ .

The main task of joint sparse recovery problems is to identify the common support shared by the unknown sparse signals. Once the support is determined accurately, then the system model in (2) can be reduced to $r$ independent overdetermined systems and thus the solutions can be found via the conventional least squares (LS) approach. The problem to identify the support can be formulated as

\begin{split}&\underset{S\subset\{1,\ldots,n\}}{\min}\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ |S|\\ &\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \text{s.t.}\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \mathbf{y}_{1},\ldots,\mathbf{y}_{r}\in\mathcal{R}(\mathbf{A}_{S}),\end{split}

(3)

where $\mathbf{A}_{S}$ is the submatrix of $\mathbf{A}$ that contains the columns indexed by $S$ and $\mathcal{R}(\mathbf{A}_{S})$ is the subspace spanned by the columns of $\mathbf{A}_{S}$ . A naive way to solve (3) is to search over all possible subspaces spanned by $\mathbf{A}$ . Such a combinatorial search, however, is exhaustive and thus infeasible for most practical scenarios. Over the years, various approaches to address this problem have been proposed [11, 12, 13, 14, 15, 17, 16]. Roughly speaking, these approaches can be classified into two categories: 1) those based on greedy search principles (e.g., simultaneous orthogonal matching pursuit (SOMP) [13] and rank aware order recursive matching pursuit (RA-ORMP) [16]) and 2) those relying on convex optimization (e.g., mixed norm minimization [14]). Hybrid approaches combining greedy search techniques and conventional methods such as multiple signal classification (MUSIC) have also been proposed (e.g., compressive MUSIC (CS-MUSIC) [15] and subspace-augmented MUSIC (SA-MUSIC) [17]).

In this paper, we put forth a new algorithm, called signal space matching pursuit (SSMP), to improve the quality of joint sparse recovery. Basically, the SSMP algorithm identifies columns of $\mathbf{A}$ that span the measurement space. By the measurement space, we mean the subspace $\mathcal{R}(\mathbf{Y})$ spanned by the measurement vectors $\mathbf{y}_{1},\ldots,\mathbf{y}_{r}$ . Towards this end, we recast (3) as

\begin{split}&\underset{S\subset\{1,\ldots,n\}}{\min}\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ |S|\\ &\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \text{s.t.}\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \mathcal{R}(\mathbf{Y})\subseteq\mathcal{R}(\mathbf{A}_{S}).\end{split}

(4)

In solving this problem, SSMP exploits the notion of subspace distance which measures the closeness between two vector spaces (see Definition 1 in Section II-A) [18]. Specifically, SSMP chooses multiple, say $L$ , columns of $\mathbf{A}$ that minimize the subspace distance to the measurement space in each iteration.

The contributions of this paper are as follows:

1)

We propose a new joint sparse recovery algorithm called SSMP (Section II). From the simulation results, we show that SSMP outperforms conventional techniques both in noiseless and noisy scenarios (Section VI). Specifically, in the noiseless scenario, the critical sparsity (the maximum sparsity level at which exact reconstruction is ensured [22]) of SSMP is about 1.5 times higher than those obtained by conventional techniques. In the noisy scenario, SSMP performs close to the oracle least squares (Oracle-LS) estimator²²2Oracle-LS is the estimator that provides the best achievable bound using prior knowledge on the common support. in terms of the mean square error (MSE) when the signal-to-noise ratio (SNR) is high.

We analyze a condition under which the SSMP algorithm exactly reconstructs any group of jointly $K$ -sparse vectors in at most $K$ iterations in the noiseless scenario (Section III).

–

In the full row rank scenario ( $r=K$ ),³³3Note that $\mathbf{X}$ has at most $K$ nonzero rows so that $r=\operatornamewithlimits{rank}(\mathbf{X})\leq K$ . In this sense, we refer to the case where $r=K$ as the full row rank scenario. we show that SSMP accurately recovers any group $\{\mathbf{x}_{i}\}_{i=1}^{r}$ of jointly $K$ -sparse vectors in $\lceil\frac{K}{L}\rceil$ iterations if the sampling matrix $\mathbf{A}$ satisfies (Theorem 1(i))

$\displaystyle\operatornamewithlimits{krank}(\mathbf{A})$ $\displaystyle\geq K+1,$

where $\operatornamewithlimits{krank}(\mathbf{A})$ is the maximum number of columns in $\mathbf{A}$ that are linearly independent. This implies that under a mild condition on $\mathbf{A}$ (any $m$ columns of $\mathbf{A}$ are linearly independent), SSMP guarantees exact reconstruction with $m=K+1$ measurements, which meets the fundamental minimum requirement on the number of measurements to ensure perfect recovery of $\{\mathbf{x}_{i}\}_{i=1}^{r}$ [16].

–

In the rank deficient scenario ( $r<K$ ), we show that SSMP reconstructs $\{\mathbf{x}_{i}\}_{i=1}^{r}$ accurately in at most $K-r+\lceil\frac{r}{L}\rceil$ iterations if $\mathbf{A}$ satisfies the restricted isometry property (RIP [27]) of order $L(K-r)+r+1$ with (Theorem 1(ii))

\displaystyle\delta_{L(K-r)+r+1}

\displaystyle<\max\left\{\frac{\sqrt{r}}{\sqrt{K+\frac{r}{4}}+\sqrt{\frac{r}{4}}},\frac{\sqrt{L}}{\sqrt{K}+1.15\sqrt{L}}\right\}.

(5)

Using the monotonicity of the RIP constant (see Lemma 3), one can notice that the requirement on the RIP constant becomes less restrictive when the number $r$ of (linearly independent) measurement vectors increases. This behavior seems to be natural but has not been reported for conventional methods such as SOMP [13] and mixed norm minimization [14]. In particular, if $r$ is on the order of $K$ , e.g., $r=\lceil\frac{K}{2}\rceil$ , then (5) is satisfied under

\delta_{L(K-\lceil\frac{K}{2}\rceil)+\lceil\frac{K}{2}\rceil+1}<\frac{1}{2},

which implies that SSMP ensures exact recovery with overwhelming probability as long as the number of random measurements scales linearly with $K\log\frac{n}{K}$ [27, 31].

3)

We analyze the performance of SSMP in the scenario where the observation matrix $\mathbf{Y}$ is contaminated by noise (Section IV). Specifically, we show that under a suitable RIP condition, the reconstruction error of SSMP is upper bounded by a constant multiple of the noise power (Theorems 2 and 3), which demonstrates the robustness of SSMP to measurement noise.
4)

As a special case, when $r=1$ , we establish the performance guarantee of SSMP running more than $K$ iterations (Section V). Specifically, we show that SSMP exactly recovers any $K$ -sparse vector in $\max\{K,\lfloor\frac{8K}{L}\rfloor\}$ iterations under (Theorem 5)

$\delta_{\lfloor 7.8K\rfloor}\leq 0.155.$

In contrast to (5), this bound is a constant and unrelated to the sparsity $K$ . This implies that even when $r$ is not on the order of $K$ , SSMP guarantees exact reconstruction with $\mathcal{O}(K\log\frac{n}{K})$ random measurements by running slightly more than $K$ iterations.

We briefly summarize the notations used in this paper.

•

Let $\Omega=\{1,{\color[rgb]{0,0,1}{\ldots}},n\}$ ;
•

For $J\subset\Omega$ , $|J|$ is the cardinality of $J$ and $\Omega\setminus J$ is the set of all indices in $\Omega$ but not in $J$ ;
•

For a vector $\mathbf{x}$ , $\mathbf{x}_{J}\in\mathbb{R}^{|J|}$ is the restriction of $\mathbf{x}$ to the elements indexed by $J$ ;
•

We refer to $\mathbf{X}$ having at most $K$ nonzero rows as a row $K$ -sparse signal and define its support $S$ as the index set of its nonzero rows;⁴⁴4For example, if $\mathbf{x}_{1}=[1\ 0\ 0\ 2]^{\prime}$ and $\mathbf{x}_{2}=[2\ 0\ 0\ 1]^{\prime}$ , then $\mathbf{X}=[\mathbf{x}_{1}\ \mathbf{x}_{2}]$ is a row $2$ -sparse matrix with support $\{1,4\}$ .
•

The submatrix of $\mathbf{X}$ containing the rows indexed by $J$ is denoted by $\mathbf{X}^{J}$ ;

We use the following notations for a general matrix $\mathbf{H}\in\mathbb{R}^{m\times n}$ .

•

The $i$ -th column of $\mathbf{H}$ is denoted by $\mathbf{h}_{i}\in\mathbb{R}^{m}$ ;
•

The $j$ -th row of $\mathbf{H}$ is denoted by $\mathbf{h}^{j}\in\mathbb{R}^{n}$ ;
•

The submatrix of $\mathbf{H}$ containing the columns indexed by $J$ is denoted by $\mathbf{H}_{J}$ ;
•

If $\mathbf{H}_{J}$ has full column rank, $\mathbf{H}_{J}^{\dagger}=(\mathbf{H}_{J}^{\prime}\mathbf{H}_{J})^{-1}\mathbf{H}_{J}^{\prime}$ is the pseudoinverse of $\mathbf{H}_{J}$ where $\mathbf{H}_{J}^{\prime}$ is the transpose of $\mathbf{H}_{J}$ ;
•

We define $\mathcal{R}(\mathbf{H})$ as the column space of $\mathbf{H}$ ;
•

The Frobenius norm and the spectral norm of $\mathbf{H}$ are denoted by $\|\mathbf{H}\|_{F}$ and $\|\mathbf{H}\|_{2}$ , respectively;
•

The mixed $\ell_{1,2}$ -norm of $\mathbf{H}$ is denoted by $\|\mathbf{H}\|_{1,2}$ , i.e., $\|\mathbf{H}\|_{1,2}=\sum_{j=1}^{m}\|\mathbf{h}^{j}\|_{2}$ ;
•

The minimum, maximum, and $i$ -th largest singular values of $\mathbf{H}$ are denoted by $\sigma_{\min}(\mathbf{H})$ , $\sigma_{\max}(\mathbf{H})$ , and $\sigma_{i}(\mathbf{H})$ , respectively;
•

The orthogonal projections onto a subspace $\mathcal{V}\subset\mathbb{R}^{m}$ and its orthogonal complement $\mathcal{V}^{\perp}$ are denoted by $\mathbf{P}_{\mathcal{V}}$ and $\mathbf{P}_{\mathcal{V}}^{\perp}$ , respectively. For simplicity, we write $\mathbf{P}_{J}$ and $\mathbf{P}_{J}^{\perp}$ instead of $\mathbf{P}_{\mathcal{R}(\mathbf{A}_{J})}$ and $\mathbf{P}_{\mathcal{R}(\mathbf{A}_{J})}^{\perp}$ , respectively.

II The Proposed SSMP Algorithm

As mentioned, we use the notion of subspace distance in solving our main problem (4). In this section, we briefly introduce the definition of subspace distance and its properties, and then describe the proposed SSMP algorithm.

II-A Preliminaries

We begin with the definition of subspace distance. In a nutshell, the subspace distance is minimized if two subspaces coincide with each other and maximized if two subspaces are perpendicular to each other.

Definition 1 (Subspace distance [18]).

Let $\mathcal{V}$ and $\mathcal{W}$ be subspaces in $\mathbb{R}^{m}$ , and let $\{\mathbf{v}_{1},\ldots,\mathbf{v}_{p}\}$ and $\{\mathbf{w}_{1},\ldots,\mathbf{w}_{q}\}$ be orthonormal bases of $\mathcal{V}$ and $\mathcal{W}$ , respectively. Then the subspace distance $\operatornamewithlimits{dist}(\mathcal{V},\mathcal{W})$ between $\mathcal{V}$ and $\mathcal{W}$ is

\displaystyle\operatornamewithlimits{dist}(\mathcal{V},\mathcal{W})

\displaystyle=\sqrt{\max\{p,q\}-\underset{i=1}{\overset{p}{\sum}}\underset{j=1}{\overset{q}{\sum}}|\langle\mathbf{v}_{i},\mathbf{w}_{j}\rangle|^{2}}.

(6)

Refer to caption — Figure 1: Subspace distance between one-dimensional subspaces.

As a special case, suppose $\mathcal{V}$ and $\mathcal{W}$ are one-dimensional subspaces in $\mathbb{R}^{m}$ (i.e., $\mathcal{V}$ and $\mathcal{W}$ are two lines in $\mathbb{R}^{m}$ ). Further, let $\{\mathbf{v}\}$ and $\{\mathbf{w}\}$ be orthonormal bases of $\mathcal{V}$ and $\mathcal{W}$ , respectively, and $\theta$ be the angle between $\mathbf{v}$ and $\mathbf{w}$ . Then the subspace distance $\operatornamewithlimits{dist}(\mathcal{V},\mathcal{W})$ between $\mathcal{V}$ and $\mathcal{W}$ is (see Fig. 1)

$\displaystyle\operatornamewithlimits{dist}(\mathcal{V},\mathcal{W})$	$\displaystyle\overset{(a)}{=}\sqrt{1-\|\langle\mathbf{v},\mathbf{w}\rangle\|^{2}}$
	$\displaystyle=\sqrt{1-(\\|\mathbf{v}\\|_{2}\\|\mathbf{w}\\|_{2}\cos\theta)^{2}}$
	$\displaystyle=\sin\theta,$	(7)

where (a) is because $\mathcal{V}$ and $\mathcal{W}$ are one-dimensional (i.e., $p=q=1$ in (6)). One can easily see that $\operatornamewithlimits{dist}(\mathcal{V},\mathcal{W})$ is maximal when $\theta=\frac{\pi}{2}$ (i.e., $\mathcal{V}\perp\mathcal{W}$ ) and $\operatornamewithlimits{dist}(\mathcal{V},\mathcal{W})=0$ if and only if $\theta=0$ (i.e., $\mathcal{V}=\mathcal{W}$ ). In fact, $\operatornamewithlimits{dist}(\mathcal{V},\mathcal{W})$ is maximized when $\mathcal{V}$ and $\mathcal{W}$ are orthogonal and $\operatornamewithlimits{dist}(\mathcal{V},\mathcal{W})=0$ if and only if $\mathcal{V}$ and $\mathcal{W}$ coincide with each other [18]. Also, note that the subspace distance is a proper metric between subspaces since it satisfies the following three properties of a metric [18, 19]:

(i)

$\operatornamewithlimits{dist}(\mathcal{V},\mathcal{W})\geq 0$ for any subspaces $\mathcal{V},\mathcal{W}\subset\mathbb{R}^{m}$ , and the equality holds if and only if $\mathcal{V}=\mathcal{W}$ .
(ii)

$\operatornamewithlimits{dist}(\mathcal{V},\mathcal{W})=\operatornamewithlimits{dist}(\mathcal{W},\mathcal{V})$ for any subspaces $\mathcal{V},\mathcal{W}\subset\mathbb{R}^{m}$ .
(iii)

$\operatornamewithlimits{dist}(\mathcal{U},\mathcal{W})\leq\operatornamewithlimits{dist}(\mathcal{U},\mathcal{V})+\operatornamewithlimits{dist}(\mathcal{V},\mathcal{W})$ for any subspaces $\mathcal{U},\mathcal{V},\mathcal{W}\subset\mathbb{R}^{m}$ .

Exploiting the subspace distance, (4) can be reformulated as

\begin{split}&\underset{S\subset\Omega}{\min}\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ |S|\\ &\leavevmode\nobreak\ \text{s.t.}\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \operatornamewithlimits{dist}\left(\mathcal{R}(\mathbf{Y}),\mathbf{P}_{S}\mathcal{R}(\mathbf{Y})\right)=0,\end{split}

(8)

where $\mathbf{P}_{S}\mathcal{R}(\mathbf{Y})=\{\mathbf{P}_{S}\mathbf{z}:\mathbf{z}\in\mathcal{R}(\mathbf{Y})\}$ .

Lemma 1.

Problems (4) and (8) are equivalent.

Proof.

It suffices to show that two constraints $\mathcal{R}(\mathbf{Y})\subseteq\mathcal{R}(\mathbf{A}_{S})$ and $\operatornamewithlimits{dist}\left(\mathcal{R}(\mathbf{Y}),\mathbf{P}_{S}\mathcal{R}(\mathbf{Y})\right)=0$ are equivalent. If $\operatornamewithlimits{dist}\left(\mathcal{R}(\mathbf{Y}),\mathbf{P}_{S}\mathcal{R}(\mathbf{Y})\right)=0$ , then $\mathcal{R}(\mathbf{Y})=\mathbf{P}_{S}\mathcal{R}(\mathbf{Y})\subseteq\mathcal{R}(\mathbf{A}_{S})$ by the property (i) of the subspace distance. Conversely, if $\mathcal{R}(\mathbf{Y})\subseteq\mathcal{R}(\mathbf{A}_{S})$ , then $\mathbf{P}_{S}\mathcal{R}(\mathbf{Y})=\mathcal{R}(\mathbf{Y})$ and thus $\operatornamewithlimits{dist}\left(\mathcal{R}(\mathbf{Y}),\mathbf{P}_{S}\mathcal{R}(\mathbf{Y})\right)=0$ . ∎

It is worth mentioning that problem (8) can be extended to the noisy scenario by relaxing the constraint $\operatornamewithlimits{dist}\left(\mathcal{R}(\mathbf{Y}),\mathbf{P}_{S}\left(\mathcal{R}(\mathbf{Y})\right)\right)=0$ to $\operatornamewithlimits{dist}\left(\mathcal{R}(\mathbf{Y}),\mathbf{P}_{S}\left(\mathcal{R}(\mathbf{Y})\right)\right)\leq\epsilon$ for some properly chosen threshold $\epsilon>0$ .

II-B Algorithm Description

The proposed SSMP algorithm solves problem (8) using the greedy principle. Note that the greedy principle has been popularly used in sparse signal recovery for its computational simplicity and competitive performance [20, 21, 22, 23]. In a nutshell, SSMP sequentially investigates the support to minimize the subspace distance to the residual space.

In the first iteration, SSMP chooses multiple, say $L$ , columns of the sampling matrix $\mathbf{A}$ that minimize the subspace distance to the measurement space $\mathcal{R}(\mathbf{Y})$ . Towards this end, SSMP computes the subspace distance

d_{i}=\operatornamewithlimits{dist}\left(\mathcal{R}(\mathbf{Y}),\mathbf{P}_{\{i\}}\mathcal{R}(\mathbf{Y})\right)

between $\mathcal{R}(\mathbf{Y})$ and its orthogonal projection $\mathbf{P}_{\{i\}}\mathcal{R}(\mathbf{Y})$ onto the subspace spanned by each column $\mathbf{a}_{i}$ of $\mathbf{A}$ . Let $0\leq d_{i_{1}}\leq d_{i_{2}}\leq\ldots\leq d_{i_{n}}$ , then SSMP chooses $i_{1},\ldots,i_{L}$ (the indices corresponding to the $L$ smallest subspace distances). In other words, the estimated support $S^{1}=\{i_{1},\ldots,i_{L}\}$ is given by

	$\displaystyle S^{1}$	$\displaystyle=\underset{I:I\subset\Omega,\|I\|=L}{\arg\min}\sum_{i\in I}d_{i}$
		$\displaystyle=\underset{I:I\subset\Omega,\|I\|=L}{\arg\min}\sum_{i\in I}\operatornamewithlimits{dist}\left(\mathcal{R}(\mathbf{Y}),\mathbf{P}_{\{i\}}\mathcal{R}(\mathbf{Y})\right).$		(9)

For example, let $\mathbf{y}\in\mathbb{R}^{2}$ , $\mathbf{A}\in\mathbb{R}^{2\times 3}$ , and $\theta_{i}$ be the angle between $\mathcal{R}(\mathbf{y})$ and $\mathbf{P}_{\{i\}}\mathcal{R}(\mathbf{y})$ (see Fig. 2). Also, suppose SSMP picks up two indices in each iteration (i.e., $L=2$ ). Then, by (7), $d_{i}=\operatornamewithlimits{dist}\left(\mathcal{R}(\mathbf{y}),\mathbf{P}_{\{i\}}\mathcal{R}(\mathbf{y})\right)=\sin\theta_{i}$ so that $d_{2}<d_{3}<d_{1}$ (see $\theta_{2}<\theta_{3}<\theta_{1}$ in Fig. 2) and $S^{1}=\{2,3\}$ . After updating the support set $S^{1}$ , SSMP computes the estimate $\mathbf{X}^{1}$ of the desired signal $\mathbf{X}$ by solving the LS problem:

\displaystyle\mathbf{X}^{1}

\displaystyle=\underset{\mathbf{U}:\operatornamewithlimits{supp}(\mathbf{U})\subset S^{1}}{\arg\min}\|\mathbf{Y}-\mathbf{A}\mathbf{U}\|_{F}.

Note that $(\mathbf{X}^{1})^{S^{1}}=\mathbf{A}_{S^{1}}^{\dagger}\mathbf{Y}$ and $(\mathbf{X}^{1})^{\Omega\setminus S^{1}}=\mathbf{0}_{|\Omega\setminus S^{1}|\times r}$ , where $\mathbf{0}_{d_{1}\times d_{2}}$ is the $(d_{1}\times d_{2})$ -dimensional zero matrix. Finally, SSMP generates the residual matrix

\mathbf{R}^{1}=\mathbf{Y}-\mathbf{A}\mathbf{X}^{1}=\mathbf{Y}-\mathbf{A}_{S^{1}}(\mathbf{X}^{1})^{S^{1}}=\mathbf{Y}-\mathbf{A}_{S^{1}}\mathbf{A}_{S^{1}}^{\dagger}\mathbf{Y}=\mathbf{P}_{S^{1}}^{\perp}\mathbf{Y},

which will be used as an observation matrix for the next iteration.

The general iteration of SSMP is similar to the first iteration except for the fact that the residual space is used instead of the measurement space $\mathcal{R}(\mathbf{Y})$ in the support identification. In other words, SSMP identifies the columns of $\mathbf{A}$ that minimize the subspace distance to the residual space $\mathcal{R}(\mathbf{R}^{k-1})$ in the $k$ -th iteration. Let $I^{k}$ be the set of $L$ indices newly chosen in the $k$ -th iteration, then

\displaystyle I^{k}

\displaystyle=\underset{I:I\subset\Omega\setminus S^{k-1},|I|=L}{\arg\min}\sum_{i\in I}\operatornamewithlimits{dist}\left(\mathcal{R}(\mathbf{R}^{k-1}),\mathbf{P}_{S^{k-1}\cup\{i\}}\mathcal{R}(\mathbf{R}^{k-1})\right).

(10)

These set of operations are repeated until the iteration number reaches the maximum value $k_{\max}=\min\{K,\lfloor\frac{m}{L}\rfloor\}$ ⁵⁵5In the estimation step, SSMP computes $(\mathbf{X}^{k})^{S^{k}}=\mathbf{A}_{S^{k}}^{\dagger}\mathbf{Y}=(\mathbf{A}_{S^{k}}^{\prime}\mathbf{A}_{S^{k}})^{-1}\mathbf{A}_{S^{k}}^{\prime}\mathbf{Y}$ . In order to perform this task, $\mathbf{A}_{S^{k}}$ should have full column rank and thus $|S^{k}|=kL\leq m$ . or the pre-defined stopping criterion

\operatornamewithlimits{dist}(\mathcal{R}(\mathbf{Y}),\mathbf{P}_{S^{k}}\mathcal{R}(\mathbf{Y}))\leq\epsilon

is satisfied ( $\epsilon$ is a stopping threshold).

Suppose SSMP runs $k$ iterations in total and more than $K$ indices are chosen (i.e., $|S^{k}|>K$ ). Then, by identifying the $K$ rows of $\mathbf{X}^{k}$ with largest $\ell_{2}$ -norms, the estimated support $S^{k}$ is pruned to its subset $\widehat{S}$ consisting of $K$ elements, i.e.,

\displaystyle\widehat{S}

\displaystyle=\underset{J:|J|=K}{\arg\min}\leavevmode\nobreak\ \|\mathbf{X}^{k}-(\mathbf{X}^{k})^{J}\|_{F}.

Finally, SSMP computes the row $K$ -sparse estimate $\widehat{\mathbf{X}}$ of $\mathbf{X}$ by solving the LS problem on $\widehat{S}$ . This pruning step, also known as debiasing, helps to reduce the reconstruction error $\|\mathbf{X}-\widehat{\mathbf{X}}\|_{F}$ in the noisy scenario [34].

II-C Refined Identification Rule

As shown in (10), in the $k$ -th iteration, SSMP computes the subspace distance between the residual space $\mathcal{R}(\mathbf{R}^{k-1})$ and its projection space $\mathbf{P}_{S^{k-1}\cup\{i\}}\mathcal{R}(\mathbf{R}^{k-1})$ for all remaining indices $i$ . Since this operation requires $n-(k-1)L$ projection operators $\mathbf{P}_{S^{k-1}\cup\{i\}}$ , it would clearly be computationally prohibitive, especially for large $n$ . In view of this, it is of importance to come up with a computationally efficient way to perform the support identification. In the following proposition, we introduce a simplified selection rule under which the number of projections required in each identification step is independent of the signal dimension $n$ .

Proposition 1.

Consider the system model in (2). Suppose the SSMP algorithm chooses $L$ indices in each iteration. Then the set $I^{k+1}$ of $L$ indices chosen in the $(k+1)$ -th iteration satisfies

\displaystyle I^{k+1}

\displaystyle=\underset{I:I\subset\Omega\setminus S^{k},|I|=L}{\arg\max}\leavevmode\nobreak\ \underset{i\in I}{\sum}\left\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}}{\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\|_{2}}\right\|_{2}.

(11)

Proof.

See Appendix A. ∎

Note that the selection rule in (11) requires only two projection operators $\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}$ and $\mathbf{P}_{S^{k}}^{\perp}$ in each identification step. One can also observe from (11) that SSMP performs the identification step by one simple matrix multiplication followed by the selection of $L$ columns with largest $\ell_{2}$ -norms. Specifically, if $\mathbf{B}\in\mathbb{R}^{m\times n}$ is the $\ell_{2}$ -normalized counterpart of $\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}$ , i.e.,

\mathbf{b}_{i}=\begin{cases}\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}}{\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\|_{2}},&i\in\Omega\setminus S^{k},\\ \mathbf{0}_{m\times 1},&i\in S^{k},\end{cases}

then SSMP picks $L$ column indices of $\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{B}$ with largest $\ell_{2}$ -norms. In Algorithm 1, we summarize the proposed SSMP algorithm.

Algorithm 1 The SSMP algorithm

0: sampling matrix

\mathbf{A}\in\mathbb{R}^{m\times n}

, observation matrix

\mathbf{Y}\in\mathbb{R}^{m\times r}

, sparsity

K

, number

L

of selected indices per iteration, and stopping threshold

\epsilon

0: iteration counter

k=0

, estimated support

S^{0}=\emptyset

, and residual matrix

\mathbf{R}^{0}=\mathbf{Y}

1: While

k<\min\{K,\lfloor\frac{m}{L}\rfloor\}

and

\operatornamewithlimits{dist}(\mathcal{R}(\mathbf{Y}),\mathbf{P}_{S^{k}}\mathcal{R}(\mathbf{Y}))>\epsilon

k=k+1

;

3: Construct the

\ell_{2}

-normalized counterpart

\mathbf{B}

\mathbf{P}_{S^{k-1}}^{\perp}\mathbf{A}

, i.e.,

\mathbf{b}_{i}=\mathbf{P}_{S^{k-1}}^{\perp}\mathbf{a}_{i}/\|\mathbf{P}_{S^{k-1}}^{\perp}\mathbf{a}_{i}\|_{2}

for each of

i\in\Omega\setminus S^{k-1}

and

\mathbf{b}_{i}=\mathbf{0}_{m\times 1}

for each of

i\in S^{k-1}

;

4: Identify the indices

\phi(1),\ldots,\phi(L)

L

columns in

\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k-1})}\mathbf{B}

with largest

\ell_{2}

-norms;

S^{k}=S^{k-1}\cup\{\phi(1),\ldots,\phi(L)\}

;

\mathbf{X}^{k}=\underset{\mathbf{U}:\operatornamewithlimits{supp}(\mathbf{U})\subset S^{k}}{\arg\min}\|\mathbf{Y}-\mathbf{A}\mathbf{U}\|_{F}

;

\mathbf{R}^{k}=\mathbf{Y}-\mathbf{A}\mathbf{X}^{k}

;

8: end while

\widehat{S}=\underset{J:|J|=K}{\arg\min}\|\mathbf{X}^{k}-(\mathbf{X}^{k})^{J}\|_{F}

and

\widehat{\mathbf{X}}

satisfying

\widehat{\mathbf{X}}^{\widehat{S}}=\mathbf{A}_{\widehat{S}}^{\dagger}\mathbf{Y}

and

\widehat{\mathbf{X}}^{\Omega\setminus\widehat{S}}=\mathbf{0}_{|\Omega\setminus\widehat{S}|\times r}

III Exact Joint Sparse Recovery via SSMP

In this section, we analyze a sufficient condition under which SSMP recovers any row $K$ -sparse matrix accurately in the noiseless scenario. In this case, we set the stopping threshold $\epsilon$ to zero. Also, we assume that the number $L$ of indices chosen in each iteration satisfies $L\leq\min\{K,\frac{m}{K}\}$ . Then SSMP performs at most $K$ iterations before stopping (see Algorithm 1). Finally, we assume that the sampling matrix $\mathbf{A}$ has unit $\ell_{2}$ -norm columns since the performance of SSMP is not affected by the $\ell_{2}$ -normalization of columns in $\mathbf{A}$ .⁶⁶6Note that the term $\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}/\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\|_{2}$ in (11) is not affected by the $\ell_{2}$ -normalization of columns in $\mathbf{A}$ .

III-A Definition and Lemmas

In our analysis, we employ the RIP framework, a widely used tool to analyze sparse recovery algorithms [25, 21, 22, 26, 23, 24].

Definition 2 (RIP [27]).

A matrix $\mathbf{A}\in\mathbb{R}^{m\times n}$ is said to satisfy the RIP of order $K$ if there exists a constant $\delta\in(0,1)$ such that

\displaystyle(1-\delta)\|\mathbf{x}\|_{2}^{2}\leq\|\mathbf{A}\mathbf{x}\|_{2}^{2}\leq(1+\delta)\|\mathbf{x}\|_{2}^{2}

(12)

for any $K$ -sparse vector $\mathbf{x}\in\mathbb{R}^{n}$ . In particular, the minimum value of $\delta$ satisfying (12) is called the RIP constant and denoted by $\delta_{K}$ .

We next introduce several lemmas useful in our analysis.

Lemma 2 ([17, Lemma A.2]).

Let $\mathbf{A}\in\mathbb{R}^{m\times n}$ and $S,J\subset\Omega$ with $S\setminus J\neq\emptyset$ , then

\sigma_{\min}(\mathbf{P}_{J}^{\perp}\mathbf{A}_{S\setminus J})\geq\sigma_{\min}(\mathbf{A}_{S\cup J}).

Lemma 2 implies that if $\mathbf{A}_{S\cup J}$ has full column rank, so does the projected matrix $\mathbf{P}_{J}^{\perp}\mathbf{A}_{S\setminus J}$ . The next lemma describes the monotonicity of the RIP constant.

Lemma 3 ([22, Lemma 1]).

If a matrix $\mathbf{A}$ satisfies the RIP of orders $K_{1}$ and $K_{2}$ $(K_{1}\leq K_{2})$ , then $\delta_{K_{1}}\leq\delta_{K_{2}}$ .

Lemma 4 ([28, Lemma 1]).

Let $\mathbf{A}\in\mathbb{R}^{m\times n}$ and $S,J\subset\Omega$ . If $\mathbf{A}$ satisfies the RIP of order $|S\cup J|$ , then for any $\mathbf{z}\in\mathbb{R}^{|S\setminus J|}$ ,

\displaystyle(1-\delta_{|S\cup J|})\|\mathbf{z}\|_{2}^{2}\leq\|\mathbf{P}_{J}^{\perp}\mathbf{A}_{S\setminus J}\mathbf{z}\|_{2}^{2}\leq(1+\delta_{|S\cup J|})\|\mathbf{z}\|_{2}^{2}.

Lemma 4 implies that if $\mathbf{A}$ satisfies the RIP of order $|S\cup J|$ , then the projected matrix $\mathbf{P}_{J}^{\perp}\mathbf{A}_{S\setminus J}$ obeys the RIP of order $|S\setminus J|$ and the corresponding RIP constant $\delta_{|S\setminus J|}(\mathbf{P}_{J}^{\perp}\mathbf{A}_{S\setminus J})$ satisfies $\delta_{|S\setminus J|}(\mathbf{P}_{J}^{\perp}\mathbf{A}_{S\setminus J})\leq\delta_{|S\cup J|}$ .

Recall that the SSMP algorithm picks the indices of the $L$ largest elements in $\{\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}:i\in\Omega\setminus S^{k}\}$ in the $(k+1)$ -th iteration, where $\mathbf{B}$ is the $\ell_{2}$ -normalized counterpart of $\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}$ (see Algorithm 1). This implies that SSMP chooses at least one support element in the $(k+1)$ -th iteration if and only if the largest element in $\{\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}:i\in S\setminus S^{k}\}$ is larger than the $L$ -th largest element in $\{\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}:i\in\Omega\setminus(S\cup S^{k})\}$ . The following lemma provides a lower bound of $\max_{i\in S\setminus S^{k}}\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}$ .

Lemma 5.

Consider the system model in (2). Let $S^{k}$ be the estimated support and $\mathbf{R}^{k}$ be the residual generated in the k-th iteration of the SSMP algorithm. Also, let $\mathbf{B}$ be the $\ell_{2}$ -normalized counterpart of $\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}$ . If $\mathbf{A}_{S\cup S^{k}}$ has full column rank and $|S\setminus S^{k}|>0$ , then

	$\displaystyle\max_{i\in S\setminus S^{k}}\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\\|_{2}^{2}$	$\displaystyle\geq\frac{1}{\|S\setminus S^{k}\|}\sum_{i\in S\setminus S^{k}}\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\\|_{2}^{2}$
		$\displaystyle\geq\frac{1}{\|S\setminus S^{k}\|}\sum_{i=1}^{d}\sigma_{\|S\setminus S^{k}\|+1-i}^{2}(\mathbf{B}_{S\setminus S^{k}}),$		(13)

where $d=\operatornamewithlimits{rank}(\mathbf{X}^{S\setminus S^{k}})$ .

Proof.

See Appendix C. ∎

The following lemmas play an important role in bounding the $L$ -th largest element in $\{\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}:i\in\Omega\setminus(S\cup S^{k})\}$ .

Lemma 6 ([29, Lemma 3]).

Let $\mathbf{A}\in\mathbb{R}^{m\times n}$ and $S,J\subset\Omega$ . If $\mathbf{A}$ satisfies the RIP of order $|S\cup J|+1$ , then for any $i\in\Omega\setminus(S\cup J)$ ,

\displaystyle\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{J}^{\perp}\mathbf{A}_{S\setminus J})}^{\perp}\mathbf{P}_{J}^{\perp}\mathbf{a}_{i}\|_{2}^{2}

\displaystyle\geq(1-\delta_{|S\cup J|+1}^{2})\|\mathbf{P}_{J}^{\perp}\mathbf{a}_{i}\|_{2}^{2}.

Lemma 7 ([28, Lemma 2]).

Let $\mathbf{A}\in\mathbb{R}^{m\times n}$ and $S,J,\Lambda\subset\Omega$ with $(S\cup J)\cap\Lambda=\emptyset$ . If $\mathbf{A}$ satisfies the RIP of order $|S\cup J|+|\Lambda|$ , then for any $\mathbf{z}\in\mathbb{R}^{|S\setminus J|}$ ,

\|\mathbf{A}_{\Lambda}^{\prime}\mathbf{P}_{J}^{\perp}\mathbf{A}_{S\setminus J}\mathbf{z}\|_{2}\leq\delta_{|S\cup J|+|\Lambda|}\|\mathbf{z}\|_{2}.

III-B Performance Guarantee of SSMP

We now analyze a condition of SSMP to guarantee exact reconstruction of any row $K$ -sparse signal in $K$ iterations. In the sequel, we say that SSMP is successful in the $k$ -th iteration if at least one support element is chosen (i.e., $I^{k}\cap S\neq\emptyset$ ). First, we present a condition under which SSMP chooses at least $K-r$ support elements in the first $K-r$ iterations (i.e., $|S\cap S^{K-r}|\geq K-r$ ).

Proposition 2.

Consider the system model in (2), where $\mathbf{A}$ has unit $\ell_{2}$ -norm columns and any $r$ nonzero rows of $\mathbf{X}$ are linearly independent. Let $L$ be the number of indices chosen in each iteration of the SSMP algorithm. If $\mathbf{A}$ satisfies the RIP of order $L(K-r)+r+1$ with

\displaystyle\delta_{L(K-r)+r+1}

\displaystyle<\max\left\{\frac{\sqrt{r}}{\sqrt{K+\frac{r}{4}}+\sqrt{\frac{r}{4}}},\frac{\sqrt{L}}{\sqrt{K}+1.15\sqrt{L}}\right\},

(14)

then SSMP picks at least $K-r$ support elements in the first $K-r$ iterations.

Proof.

We show that $|S\cap S^{k}|\geq k$ for each of $k\in\{0,\ldots,K-r\}$ . First, we consider the case where $k=0$ . This case is trivial since $S^{0}=\emptyset$ and thus

|S\cap S^{0}|=0.

Next, we assume that $|S\cap S^{k}|\geq k$ for some integer $k\leavevmode\nobreak\ (0\leq k<K-r)$ . In other words, we assume that the SSMP algorithm chooses at least $k$ support elements in the first $k$ iterations. In particular, we consider the case where $|S\cap S^{k}|=k$ , since otherwise $|S\cap S^{k+1}|\geq|S\cap S^{k}|\geq k+1$ . Under this assumption, we show that SSMP picks at least one support element in the $(k+1)$ -th iteration. As mentioned, SSMP is successful in the $(k+1)$ -th iteration if and only if the largest element $p_{1}$ in $\{\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}\}_{i\in S\setminus S^{k}}$ is larger than the $L$ -th largest element $q_{L}$ in $\{\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}\}_{i\in\Omega\setminus(S\cup S^{k})}$ . In our proof, we build a lower bound of $p_{1}$ and an upper bound of $q_{L}$ and then show that the former is larger than the latter under (14).

$\bullet$ Lower bound of $p_{1}$ :

Note that $|S\setminus S^{k}|=|S|-|S\cap S^{k}|=K-k>r$ . Then, $\operatornamewithlimits{rank}(\mathbf{X}^{S\setminus S^{k}})\geq r$ since any $r$ nonzero rows of $\mathbf{X}$ are linearly independent. Also, note that $\operatornamewithlimits{rank}(\mathbf{X}^{S\setminus S^{k}})\leq r$ since $\mathbf{X}^{S\setminus S^{k}}$ consists of $r$ columns. As a result, we have

\displaystyle\operatornamewithlimits{rank}(\mathbf{X}^{S\setminus S^{k}})

\displaystyle=r.

(15)

In addition, since

	$\displaystyle\|S\cup S^{k}\|$	$\displaystyle=\|S^{k}\|+\|S\setminus S^{k}\|=Lk+K-k$
		$\displaystyle\leq(L-1)(K-r-1)+K=L(K-r-1)+r+1,$		(16)

$\mathbf{A}$ satisfies the RIP of order $|S\cup S^{k}|$ . Then, by Lemma 5, we have

\displaystyle p_{1}^{2}

\displaystyle=\max_{i\in S\setminus S^{k}}\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}^{2}\geq\frac{r\sigma_{\min}^{2}(\mathbf{B}_{S\setminus S^{k}})}{K-k}.

(17)

Let $\mathbf{D}=\text{diag}\{\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\|_{2}:i\in S\setminus S^{k}\}$ , then

$\displaystyle\sigma_{\min}^{2}(\mathbf{B}_{S\setminus S^{k}})$	$\displaystyle\hskip 1.05273pt=\sigma_{\min}^{2}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}\mathbf{D}^{-1})$
	$\displaystyle\hskip 1.05273pt\geq\sigma_{\min}^{2}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}})\sigma_{\min}^{2}(\mathbf{D}^{-1})$
	$\displaystyle\hskip 1.05273pt=\frac{\sigma_{\min}^{2}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}})}{\max_{i\in S\setminus S^{k}}\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\\|_{2}^{2}}$
	$\displaystyle\overset{(a)}{\geq}\sigma_{\min}^{2}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}),$	(18)

where (a) is because $\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\|_{2}^{2}\leq\|\mathbf{a}_{i}\|_{2}^{2}=1$ for each of $i\in S\setminus S^{k}$ . Note that $\mathbf{A}$ satisfies the RIP of order $|S\cup S^{k}|$ . Then, by Lemma 4, the projected matrix $\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}$ obeys the RIP of order $|S\setminus S^{k}|$ and the corresponding RIP constant $\delta_{|S\setminus S^{k}|}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}})$ satisfies

\displaystyle\delta_{|S\setminus S^{k}|}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}})

\displaystyle\leq\delta_{|S\cup S^{k}|}.

(19)

Also, by the definition of the RIP (see Definition 2), we have

\displaystyle\sigma_{\min}^{2}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}})

\displaystyle\geq 1-\delta_{|S\setminus S^{k}|}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}).

(20)

Finally, by combining (17)-(20), we obtain

\displaystyle p_{1}^{2}

\displaystyle\geq\frac{r(1-\delta_{|S\cup S^{k}|})}{K-k}\geq\frac{r(1-\delta_{L(K-r)+r+1})}{K},

(21)

where the last inequality follows from Lemma 3.

$\bullet$ Upper bound of $q_{L}$ :

We build an upper bound of $q_{L}$ in two different ways and combine the results.

First, let $\psi_{l}$ be the index corresponding to the $l$ -th largest element $q_{l}$ in $\left\{\left\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\right\|_{2}\right\}_{i\in\Omega\setminus(S\cup S^{k})}$ and $\Lambda=\{\psi_{1},\ldots,\psi_{L}\}$ . Since $\mathcal{R}(\mathbf{R}^{k})=\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}\mathbf{X}^{S\setminus S^{k}})\subseteq\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}})$ , we have

$\displaystyle q_{L}^{2}$	$\displaystyle=\left\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}}{\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}\\|_{2}}\right\\|_{2}^{2}$
	$\displaystyle\leq\left\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}})}\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}}{\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}\\|_{2}}\right\\|_{2}^{2}$
	$\displaystyle=1-\left\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}})}^{\perp}\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}}{\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}\\|_{2}}\right\\|_{2}^{2}.$	(22)

Also, since $|S\cup S^{k}|+1\leq|S\cup S^{k}|+L\leq L(K-r)+r+1$ by (16), $\mathbf{A}$ satisfies the RIP of order $|S\cup S^{k}|+1$ and thus

	$\displaystyle\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}})}^{\perp}\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}\\|_{2}^{2}$	$\displaystyle\overset{(a)}{\geq}(1-\delta_{\|S\cup S^{k}\|+1}^{2})\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}\\|_{2}^{2}$		(23)
		$\displaystyle\overset{(b)}{\geq}(1-\delta_{L(K-r)+r+1}^{2})\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}\\|_{2}^{2},$		(24)

where (a) and (b) follow from Lemmas 6 and 3, respectively. By combining (22) and (24), we obtain

\displaystyle q_{L}^{2}

\displaystyle\leq\delta_{L(K-r)+r+1}^{2}.

(25)

We next obtain an upper bound of $q_{L}$ in a different way. Since $q_{L}$ is the $L$ -th largest element, we have

$\displaystyle q_{L}^{2}$	$\displaystyle\hskip 1.05273pt\leq\frac{1}{L}(q_{1}^{2}+\ldots+q_{L}^{2})$
	$\displaystyle\hskip 1.05273pt=\frac{1}{L}\underset{l=1}{\overset{L}{\sum}}\left\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{l}}}{\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{l}}\\|_{2}}\right\\|_{2}^{2}$
	$\displaystyle\overset{(a)}{\leq}\frac{1}{L(1-\delta_{\|S^{k}\|+1}^{2})}\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{\Lambda}\\|_{F}^{2}$
	$\displaystyle\hskip 1.05273pt=\frac{1}{L(1-\delta_{\|S^{k}\|+1}^{2})}\\|\mathbf{A}_{\Lambda}^{\prime}\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}\mathbf{U}\\|_{F}^{2}$
	$\displaystyle\overset{(b)}{\leq}\frac{\delta_{\|S\cup S^{k}\|+L}^{2}\\|\mathbf{U}\\|_{F}^{2}}{L(1-\delta_{\|S^{k}\|+1}^{2})},$	(26)

where (a) is because $\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{l}}\|_{2}^{2}\geq 1-\delta_{|S^{k}|+1}^{2}$ by Lemma 6, $\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}\mathbf{U}$ is an orthonormal basis of $\mathcal{R}(\mathbf{R}^{k})\leavevmode\nobreak\ (\subseteq\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}))$ , and (b) follows from Lemma 7. Also, since $\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}\mathbf{U}$ is an orthonormal basis of $\mathcal{R}(\mathbf{R}^{k})$ and $\operatornamewithlimits{rank}(\mathbf{R}^{k})=r$ ,⁷⁷7Note that $\operatornamewithlimits{rank}(\mathbf{R}^{k})=\operatornamewithlimits{rank}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}\mathbf{X}^{S\setminus S^{k}})=\operatornamewithlimits{rank}(\mathbf{X}^{S\setminus S^{k}})=r$ , where the last equality follows from (15). we have

\displaystyle r=\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}\mathbf{U}\|_{F}^{2}\overset{(a)}{\geq}(1-\delta_{|S\cup S^{k}|})\|\mathbf{U}\|_{F}^{2},

(27)

where (a) follows from Lemma 4. Using this together with (26), we have

	$\displaystyle q_{L}^{2}$	$\displaystyle\leq\frac{r\delta_{\|S\cup S^{k}\|+L}^{2}}{L(1-\delta_{\|S^{k}\|+1}^{2})(1-\delta_{\|S\cup S^{k}\|})}$
		$\displaystyle\leq\frac{r\delta_{L(K-r)+r+1}^{2}}{L(1-\delta_{L(K-r)+r+1}^{2})(1-\delta_{L(K-r)+r+1})},$		(28)

where the last inequality follows from Lemma 3.

Finally, from (25) and (28), we obtain the following upper bound of $q_{L}$ :

\displaystyle q_{L}^{2}

\displaystyle\leq\min\left\{\delta_{L(K-r)+r+1}^{2},\frac{r\delta_{L(K-r)+r+1}^{2}}{L(1-\delta_{L(K-r)+r+1}^{2})(1-\delta_{L(K-r)+r+1})}\right\}.

(29)

$\bullet$ When is $p_{1}>q_{L}$ ?

From (21) and (29), we have

	$\displaystyle p_{1}^{2}-q_{L}^{2}$
		$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \geq\frac{r(1-\delta_{L(K-r)+r+1})}{K}-\min\left\{\delta_{L(K-r)+r+1}^{2},\frac{r\delta_{L(K-r)+r+1}^{2}}{L(1-\delta_{L(K-r)+r+1}^{2})(1-\delta_{L(K-r)+r+1})}\right\}.$		(30)

One can easily check that under (14), the right-hand side of (30) is strictly larger than zero. As a result, $p_{1}>q_{L}$ , and hence SSMP is successful in the $(k+1)$ -th iteration. ∎

Thus far, we have shown that SSMP picks at least $K-r$ support elements in the first $K-r$ iterations under (14). We next analyze the performance of SSMP when at least $K-r$ support elements are chosen.

Proposition 3.

Consider the system model in (2), where $\mathbf{A}$ has unit $\ell_{2}$ -norm columns and any $r$ nonzero rows of $\mathbf{X}$ are linearly independent. Let $L$ be the number of indices chosen in each iteration of the SSMP algorithm. Suppose SSMP picks at least $K-r$ support elements in the first $k$ iterations (i.e., $|S\cap S^{k}|\geq K-r$ ). If $\mathbf{A}$ satisfies

\displaystyle\operatornamewithlimits{krank}(\mathbf{A})

\displaystyle\geq|S\cup S^{k}|+1,

(31)

then SSMP chooses $\min\{L,|S\setminus S^{k}|\}$ support elements in the $(k+1)$ -th iteration.

Proof.

In a nutshell, we will show that


	$\displaystyle\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\\|_{2}=1,\leavevmode\nobreak\ \forall i\in S\setminus S^{k},$		(32a)
	$\displaystyle\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\\|_{2}<1,\leavevmode\nobreak\ \forall i\in\Omega\setminus(S\cup S^{k}).$		(32b)

If this argument holds, then $\min\{L,|S\setminus S^{k}|\}$ support elements are chosen since SSMP picks the indices of the $L$ largest elements in $\{\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}\}_{i\in\Omega\setminus S^{k}}$ in the $(k+1)$ -th iteration (see Algorithm 1).

$\bullet$ Proof of (32a):

Since $|S\setminus S^{k}|=K-|S\cap S^{k}|\leq r$ and any $r$ nonzero rows of $\mathbf{X}$ are linearly independent, we have

\displaystyle\operatornamewithlimits{rank}(\mathbf{X}^{S\setminus S^{k}})

\displaystyle=|S\setminus S^{k}|.

Then, by Lemma 5, we have

	$\displaystyle\frac{1}{\|S\setminus S^{k}\|}\sum_{i\in S\setminus S^{k}}\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\\|_{2}^{2}$	$\displaystyle\geq\frac{1}{\|S\setminus S^{k}\|}\sum_{i=1}^{\|S\setminus S^{k}\|}\sigma_{\|S\setminus S^{k}\|+1-i}^{2}(\mathbf{B}_{S\setminus S^{k}})$
		$\displaystyle=\frac{\\|\mathbf{B}_{S\setminus S^{k}}\\|_{F}^{2}}{\|S\setminus S^{k}\|}=1.$		(33)

Also, since $\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}\leq\|\mathbf{b}_{i}\|_{2}=1$ for each of $i\in S\setminus S^{k}$ , we have

\displaystyle\frac{1}{|S\setminus S^{k}|}\sum_{i\in S\setminus S^{k}}\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}^{2}

\displaystyle\leq 1.

(34)

By combining (33) and (34), we obtain

\frac{1}{|S\setminus S^{k}|}\sum_{i\in S\setminus S^{k}}\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}^{2}=1,

which in turn implies (32a).

$\bullet$ Proof of (32b):

If $\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}=1\leavevmode\nobreak\ (\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\|_{2}=\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\|_{2})$ for some incorrect index $i\in\Omega\setminus(S\cup S^{k})$ , then

\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\in\mathcal{R}(\mathbf{R}^{k})\subset\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}),

which implies that the matrix $\mathbf{P}_{S^{k}}^{\perp}[\mathbf{A}_{S\setminus S^{k}}\ \mathbf{a}_{i}]$ does not have full column rank. This is a contradiction since

\sigma_{\min}(\mathbf{P}_{S^{k}}^{\perp}[\mathbf{A}_{S\setminus S^{k}}\ \mathbf{a}_{i}])\overset{(a)}{\geq}\sigma_{\min}(\mathbf{A}_{(S\cup S^{k})\cup\{i\}})\overset{(b)}{>}0,

where (a) and (b) follow from Lemma 2 and (31), respectively.

Therefore, $\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}<1$ for all of $i\in\Omega\setminus(S\cup S^{k})$ . ∎

We are now ready to establish a sufficient condition for SSMP to guarantee exact reconstruction of any row $K$ -sparse matrix.

Theorem 1.

Consider the system model in (2), where $\mathbf{A}$ has unit $\ell_{2}$ -norm columns and any $r$ nonzero rows of $\mathbf{X}$ are linearly independent. Let $L\leavevmode\nobreak\ (L\leq\min\{K,\frac{m}{K}\})$ be the number of indices chosen in each iteration of the SSMP algorithm. Then, SSMP exactly reconstructs $\mathbf{X}$ from $\mathbf{Y}=\mathbf{A}\mathbf{X}$ in at most $K-r+\lceil\frac{r}{L}\rceil$ iterations if one of the following conditions is satisfied:

(i)

$r=K$ and $\mathbf{A}$ satisfies

$\displaystyle\operatornamewithlimits{krank}(\mathbf{A})$ $\displaystyle\geq K+1.$ (35)
(ii)

$r<K$ and $\mathbf{A}$ satisfies the RIP of order $L(K-r)+r+1$ with (14).

Proof.

We show that SSMP picks all support elements in at most $K-r+\lceil\frac{r}{L}\rceil$ iterations under (i) or (ii). It is worth pointing out that even if several incorrect indices are added, SSMP still reconstructs $\mathbf{X}$ accurately as long as all the support elements are chosen [23, eq. (11)].

$\bullet$ Case 1: $r=K$ and $\mathbf{A}$ satisfies (35).

In this case, it suffices to show that $\min\{L,|S\setminus S^{k}|\}$ support elements are chosen in the $(k+1)$ -th iteration for each of $k\in\{0,\ldots,\lceil\frac{K}{L}\rceil-1\}$ .

First, if $k=0$ , then $S^{k}=\emptyset$ and thus

	$\displaystyle\|S\cap S^{k}\|=0\geq K-r,$
	$\displaystyle\operatornamewithlimits{krank}(\mathbf{A})\geq K+1=\|S\cup S^{0}\|+1.$

Therefore, SSMP picks $\min\{L,|S\setminus S^{0}|\}$ support elements in the first iteration by Proposition 3.

Next, we assume that the argument holds up to $k=\alpha$ ( $0\leq\alpha<\lceil\frac{K}{L}\rceil-1$ ). Then, since

\displaystyle|S\setminus S^{k}|\geq K-\alpha L\geq K-\left(\left\lceil\frac{K}{L}\right\rceil-2\right)L>L

for each of $k\in\{0,\ldots,\alpha\}$ , $L\leavevmode\nobreak\ (=\min\{L,|S\setminus S^{k}|\})$ support elements are chosen in each of the first $\alpha+1$ iterations. In other words, SSMP does not choose an incorrect index until the $(\alpha+1)$ -th iteration (i.e., $S^{\alpha+1}\subset S$ ). Thus, we have

	$\displaystyle\|S\cap S^{\alpha+1}\|=L(\alpha+1)\geq K-r,$
	$\displaystyle\operatornamewithlimits{krank}(\mathbf{A})\geq K+1=\|S\cup S^{\alpha+1}\|+1,$

and then $\min\{L,|S\setminus S^{\alpha+1}|\}$ support elements are chosen in the ( $\alpha+2$ )-th iteration by Proposition 3.

$\bullet$ Case 2: $r<K$ and $\mathbf{A}$ satisfies the RIP with (14).

In this case, we have $|S\cap S^{K-r}|\geq K-r$ by Proposition 2. In other words, SSMP picks at least $K-r$ support elements during the first $K-r$ iterations. Then, since

	$\displaystyle\|S\cup S^{K-r}\|+1$	$\displaystyle=\|S\|+\|S^{K-r}\|-\|S\cap S^{K-r}\|+1$
		$\displaystyle\leq K+L(K-r)-(K-r)+1$
		$\displaystyle=L(K-r)+r+1,$

$\operatornamewithlimits{krank}(\mathbf{A})\geq|S\cup S^{K-r}|+1$ ,⁸⁸8Note that since $\mathbf{A}$ satisfies the RIP of order $L(K-r)+r+1$ , any $p\leavevmode\nobreak\ (p\leq L(K-r)+r+1)$ columns of $\mathbf{A}$ are linearly independent. and one can deduce (in a similar way to the case 1) that SSMP picks the rest of $|S\setminus S^{K-r}|$ support elements by running $\lceil\frac{|S\setminus S^{K-r}|}{L}\rceil\leavevmode\nobreak\ (\leq\lceil\frac{r}{L}\rceil)$ additional iterations. As a result, SSMP picks all support elements and reconstructs $\mathbf{X}$ accurately in at most $K-r+\lceil\frac{r}{L}\rceil$ iterations. ∎

Remark 1.

The assumption that any $r$ nonzero rows of $\mathbf{X}$ are linearly independent is fairly mild since it applies to many naturally acquired signals. For example, any random matrix whose entries are drawn i.i.d. from a continuous probability distribution (e.g., Gaussian, uniform, exponential, and chi-square) obeys this assumption [17, 30].

Theorem 1 indicates that SSMP does not require any RIP condition to guarantee exact reconstruction in the full row rank scenario $(r=K)$ . In [16, Theorem 2], it has been shown that

\displaystyle\operatornamewithlimits{krank}(\mathbf{A})

\displaystyle\geq 2K-r+1

(36)

is the fundamental minimum requirement on $\mathbf{A}$ to ensure exact joint sparse recovery. Combining this with Theorem 1, one can see that SSMP guarantees exact reconstruction with the minimum requirement on $\mathbf{A}$ in the full row rank scenario. This is in contrast to conventional joint sparse recovery algorithms such as SOMP [13], M-ORMP [11], and mixed norm minimization [14], which require additional conditions on $\mathbf{A}$ (e.g., null space property) to guarantee exact reconstruction (see [16, Theorems 4 and 5]). In addition, if the sampling matrix $\mathbf{A}\in\mathbb{R}^{m\times n}$ satisfies $\operatornamewithlimits{krank}(\mathbf{A})=m$ , then SSMP recovers any row $K$ -sparse signal accurately with $m=K+1$ measurements in the full row rank scenario, which meets the fundamental minimum number of measurements to ensure exact joint sparse recovery [16].

Furthermore, one can see that the requirement (14) on the RIP constant becomes less restrictive when the number $r$ of (linearly independent) measurement vectors increases, since $\delta_{L(K-r)+r+1}$ decreases with $r$ and the upper bound in (14) increases with $r$ . Such behavior seems to be natural but has not been reported for conventional methods such as SOMP and mixed norm minimization.

In Table I, we summarize performance guarantees of various joint sparse recovery algorithms including SSMP. One can observe that SSMP is very competitive both in full row rank and rank deficient scenarios.

Table I: Performance Guarantees of SSMP and Conventional Techniques

	Full Row Rank ( $r=K$ )	Rank Deficient ( $r<K$ )
	Is $\operatornamewithlimits{krank}(\mathbf{A})\geq K+1$ sufficient?	Is there any known guarantee that improves with $r$ ?
SOMP [13]	No	No
M-ORMP [11]	No	No
$\ell_{1}/\ell_{2}$ -norm minimization [14]	No	No
CS-MUSIC [15]	Yes	No
SA-MUSIC [17]	Yes	Yes ( $\delta_{K+1}<\frac{r}{K+r}$ [17])
SSMP	Yes	Yes (see (14))

It is well-known that a random matrix $\mathbf{A}\in\mathbb{R}^{m\times n}$ whose entries are drawn i.i.d. from a Gaussian distribution $\mathcal{N}(0,\frac{1}{m})$ satisfies the RIP of order $K$ with $\delta_{K}\leq\epsilon\in(0,1)$ with overwhelming probability $\chi$ , provided that

\displaystyle m

\displaystyle\geq\frac{C_{\chi}K\log\frac{n}{K}}{\epsilon^{2}},

(37)

where $C_{\chi}$ is the constant depending on $\chi$ [27, 31]. When combined with Theorem 1, one can notice that SSMP requires a smaller number of (random Gaussian) measurements for exact joint sparse recovery as $r$ increases. In particular, if $r$ is on the order of $K$ , e.g., $r=\lceil\frac{K}{2}\rceil$ , then (14) is satisfied under

\displaystyle\delta_{L(K-\lceil\frac{K}{2}\rceil)+\lceil\frac{K}{2}\rceil+1}

\displaystyle<\frac{1}{2}.

This implies that SSMP accurately recovers any row $K$ -sparse matrix in at most $K$ iterations with overwhelming probability as long as the number of random measurements scales linearly with $K\log\frac{n}{K}$ .

We would like to mention that when analyzing the number of measurements ensuring exact joint sparse recovery, probabilistic approaches have been popularly used. For example, it has been shown in [30, Theorem 9], [32], [33, Table I] that if $r=\mathcal{O}(\lceil\log n\rceil)$ measurement vectors are available, then $m=\mathcal{O}(K)$ measurements are sufficient for exact joint sparse recovery. Main benefit of our result, when compared to the previous results, is that it holds uniformly for all sampling matrices and row sparse signals. For example, the result in [30] holds only for a Gaussian or Bernoulli sampling matrix and a fixed row sparse signal that is independent of the sampling matrix. Also, results in [32, 33] are based on an asymptotic analysis where the dimension $n$ and sparsity level $K$ of a desired sparse signal go to infinity. In contrast, we put our emphasis on the finite-size problem model ( $n,K<\infty$ ) so that our result is more realistic and comprehensive.

III-C Connection With Previous Efforts

First, we consider the SSMP algorithm in the single measurement vector (SMV) scenario (i.e., $r=1$ ). By Proposition 1, the set $I^{k+1}$ of $L$ indices chosen in the $(k+1)$ -th iteration is

	$\displaystyle I^{k+1}$	$\displaystyle=\underset{I:I\subset\Omega\setminus S^{k},\|I\|=L}{\arg\max}\leavevmode\nobreak\ \underset{i\in I}{\sum}\left\|\left\langle\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}}{\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\\|_{2}},\frac{\mathbf{r}^{k}}{\\|\mathbf{r}^{k}\\|_{2}}\right\rangle\right\|$
		$\displaystyle=\underset{I:I\subset\Omega\setminus S^{k},\|I\|=L}{\arg\max}\leavevmode\nobreak\ \underset{i\in I}{\sum}\left\|\left\langle\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}}{\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\\|_{2}},\mathbf{r}^{k}\right\rangle\right\|,$		(38)

where $\mathbf{r}^{k}$ is the residual defined as $\mathbf{r}^{k}=\mathbf{P}_{S^{k}}^{\perp}\mathbf{y}$ . One can see that the selection rule of SSMP simplifies to the support identification rule of the multiple orthogonal least squares (MOLS) algorithm when $r=1$ [34, Proposition 1]. In view of this, SSMP can also be considered as an extension of MOLS to the MMV scenario. Moreover, since MOLS reduces to the conventional OLS algorithm when it chooses one index in each iteration [35, 36, 34], SSMP includes OLS as a special case when $r=L=1$ . Using these connections with Theorem 1, one can establish the performance guarantees of MOLS and OLS, respectively, as follows:


	$\displaystyle\delta_{LK-L+2}<\frac{\sqrt{L}}{\sqrt{K}+1.15\sqrt{L}},\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ L>1,$		(39a)
	$\displaystyle\delta_{K+1}<\frac{1}{\sqrt{K+\frac{1}{4}}+\frac{1}{2}},\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ L=1.$		(39b)

In [34], it has been shown that MOLS accurately recovers any $K$ -sparse vector in at most $K$ iterations under


	$\displaystyle\delta_{LK}<\frac{\sqrt{L}}{\sqrt{K}+2\sqrt{L}},\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ L>1,$		(40a)
	$\displaystyle\delta_{K+1}<\frac{1}{\sqrt{K}+2},\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ L=1.$		(40b)

Clearly, the proposed guarantees (39a) and (39b) are less restrictive than (40a) and (40b), respectively. Furthermore, we would like to mention that there exists a $K$ -sparse vector that cannot be recovered by OLS running $K$ iterations under $\delta_{K+1}=\frac{1}{\sqrt{K+\frac{1}{4}}}$ [37, Example 2], which implies that a sufficient condition of OLS running $K$ iterations cannot be less restrictive than

\displaystyle\delta_{K+1}

\displaystyle<\frac{1}{\sqrt{K+\frac{1}{4}}}.

(41)

One can see that the gap between (39b) and (41) is very small and vanishes for large $K$ , which demonstrates the near-optimality of (39b) for OLS running $K$ iterations.

Next, we consider the case where the SSMP algorithm picks one index in each iteration (i.e., $L=1$ ). Then the selection rule in (11) simplifies to

\displaystyle s^{k+1}

\displaystyle=\underset{i\in\Omega\setminus S^{k}}{\arg\max}\left\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}}{\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\|_{2}}\right\|_{2}.

(42)

In this case, SSMP reduces to the RA-ORMP algorithm [16]. Exploiting the relationship between SSMP and RA-ORMP, one can deduce from Theorem 1 that RA-ORMP exactly recovers any row $K$ -sparse matrix of rank $r$ in $K$ iterations under

\displaystyle\delta_{K+1}

\displaystyle<\frac{\sqrt{r}}{\sqrt{K+\frac{r}{4}}+\sqrt{\frac{r}{4}}},

(43)

which is consistent with the best known guarantee for RA-ORMP [29]. It is worth mentioning that (43) is a near-optimal recovery condition of RA-ORMP running $K$ iterations, since there exists a row $K$ -sparse matrix of rank $r$ that cannot be recovered by RA-ORMP running $K$ iterations under $\delta_{K+1}\geq\sqrt{\frac{r}{K}}$ [29, Theorem 2]. In Table II, we summarize the relationship between the proposed SSMP, MOLS, OLS, and RA-ORMP algorithms and their performance guarantees.

Table II: Relationship between the proposed SSMP, MOLS, OLS, and RA-ORMP algorithms

	Connection With SSMP	Performance Guarantee
MOLS [34]	SSMP when $r=1$	$\delta_{LK-L+2}<\frac{\sqrt{L}}{\sqrt{K}+1.15\sqrt{L}}$
OLS [35]	SSMP when $r=L=1$	$\delta_{K+1}<\frac{1}{\sqrt{K+\frac{1}{4}}+\frac{1}{2}}$
RA-ORMP [16]	SSMP when $L=1$	$\delta_{K+1}<\frac{\sqrt{r}}{\sqrt{K+\frac{r}{4}}+\sqrt{\frac{r}{4}}}$

IV Robustness of SSMP to Measurement Noise

Thus far, we have focused on the performance guarantee of SSMP in the noiseless scenario. In this section, we analyze the performance of SSMP in the more realistic scenario where the observation matrix $\mathbf{Y}$ is contaminated by noise $\mathbf{W}\in\mathbb{R}^{m\times r}$ :

\displaystyle\mathbf{Y}

\displaystyle=\mathbf{A}\mathbf{X}+\mathbf{W}.

(44)

Here, the noise matrix $\mathbf{W}$ is assumed to be bounded (i.e., $\|\mathbf{W}\|_{F}\leq\epsilon$ for some $\epsilon>0$ ) or Gaussian. In this paper, we exclusively consider the bounded scenario, but our analysis can be easily extended to the Gaussian noise scenario after small modifications (see [38, Lemma 3]).

In our analysis, we employ the Frobenius norm $\|\mathbf{X}-\widehat{\mathbf{X}}\|_{F}$ of the reconstruction error as a performance measure since exact recovery of $\mathbf{X}$ is not possible in the noisy scenario. Also, we assume that the number $L$ of indices chosen in each iteration satisfies $L\leq\min\{K,\frac{m}{K}\}$ . Then SSMP continues to perform an iteration until the iteration number $k$ reaches $k_{\max}=K$ or $\operatornamewithlimits{dist}(\mathcal{R}(\mathbf{Y}),\mathbf{P}_{S^{k}}\mathcal{R}(\mathbf{Y}))\leq\epsilon$ for some $k<K$ (see Algorithm 1). The following theorem presents an upper bound of $\|\mathbf{X}-\widehat{\mathbf{X}}\|_{F}$ when SSMP is terminated by the condition $\operatornamewithlimits{dist}(\mathcal{R}(\mathbf{Y}),\mathbf{P}_{S^{k}}\mathcal{R}(\mathbf{Y}))\leq\epsilon$ .

Theorem 2.

Consider the system model in (44) where $\|\mathbf{W}\|_{F}$ is bounded. Suppose SSMP picks $L\leavevmode\nobreak\ (L\leq\min\{K,\frac{m}{K}\})$ indices in each iteration and $\operatornamewithlimits{dist}(\mathcal{R}(\mathbf{Y}),\mathbf{P}_{S^{k}}\mathcal{R}(\mathbf{Y}))\leq\epsilon$ for some $k<K$ . Also, suppose $\mathbf{A}$ satisfies the RIP of order $\max\{Lk+K,2K\}$ . Then the output $\widehat{\mathbf{X}}$ of SSMP satisfies

\displaystyle\|\mathbf{X}-\widehat{\mathbf{X}}\|_{F}

\displaystyle\leq\frac{2\sigma_{\max}(\mathbf{Y})\epsilon\sqrt{1+\delta_{2K}}+2(\sqrt{1+\delta_{2K}}+\sqrt{1-\delta_{Lk+K}})\|\mathbf{W}\|_{F}}{\sqrt{(1-\delta_{Lk+K})(1-\delta_{2K})}}.

In particular, when $\epsilon=\|\mathbf{W}\|_{F}/\sigma_{\max}(\mathbf{Y})$ , $\widehat{\mathbf{X}}$ satisfies

\displaystyle\|\mathbf{X}-\widehat{\mathbf{X}}\|_{F}

\displaystyle\leq\frac{(4\sqrt{1+\delta_{2K}}+2\sqrt{1-\delta_{Lk+K}})\|\mathbf{W}\|_{F}}{\sqrt{(1-\delta_{Lk+K})(1-\delta_{2K})}}.

(45)

Proof.

Recall that if SSMP chooses more than $K$ indices (i.e., $|S^{k}|>K$ ), then by identifying the $K$ rows of $\mathbf{X}^{k}$ with largest $\ell_{2}$ -norms, $S^{k}$ is pruned to its subset $\widehat{S}$ consisting of $K$ elements (see Algorithm 1). Let $\mathbf{Z}^{k}$ be the row $K$ -sparse matrix defined as $(\mathbf{Z}^{k})^{\widehat{S}}=(\mathbf{X}^{k})^{\widehat{S}}$ and $(\mathbf{Z}^{k})^{\Omega\setminus\widehat{S}}=\mathbf{0}_{|\Omega\setminus\widehat{S}|\times r}$ . Then, one can show that (see Appendix D)

\displaystyle\|\mathbf{Z}^{k}-\mathbf{X}\|_{F}

\displaystyle\leq\frac{2(\|\mathbf{R}^{k}\|_{F}+\|\mathbf{W}\|_{F})}{\sqrt{1-\delta_{Lk+K}}}

(46)

and

\displaystyle\|\mathbf{Z}^{k}-\mathbf{X}\|_{F}

\displaystyle\geq\frac{\sqrt{1-\delta_{2K}}\|\mathbf{X}-\widehat{\mathbf{X}}\|_{F}-2\|\mathbf{W}\|_{F}}{\sqrt{1+\delta_{2K}}}.

(47)

Also, since $\operatornamewithlimits{dist}(\mathcal{R}(\mathbf{Y}),\mathbf{P}_{S^{k}}\mathcal{R}(\mathbf{Y}))\leq\epsilon$ , we have

$\displaystyle\epsilon$	$\displaystyle\overset{(a)}{\geq}\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{Y}(\mathbf{Y}^{\prime}\mathbf{Y})^{-1/2}\\|_{F}$
	$\displaystyle\hskip 1.05273pt=\\|\mathbf{R}^{k}(\mathbf{Y}^{\prime}\mathbf{Y})^{-1/2}\\|_{F}$
	$\displaystyle\hskip 1.05273pt\geq\sigma_{\min}((\mathbf{Y}^{\prime}\mathbf{Y})^{-1/2})\\|\mathbf{R}^{k}\\|_{F}$
	$\displaystyle\hskip 1.05273pt=\frac{\\|\mathbf{R}^{k}\\|_{F}}{\sigma_{\max}(\mathbf{Y})},$	(48)

where (a) follows from (112) in Appendix A. By combining (46)-(48), we obtain the desired result. ∎

Theorem 2 implies that if the SSMP algorithm is terminated by the condition

\operatornamewithlimits{dist}(\mathcal{R}(\mathbf{Y}),\mathbf{P}_{S^{k}}\mathcal{R}(\mathbf{Y}))\leq\frac{\|\mathbf{W}\|_{F}}{\sigma_{\max}(\mathbf{Y})},

then $\|\mathbf{X}-\widehat{\mathbf{X}}\|_{F}$ is upper bounded by a constant multiple of the noise power $\|\mathbf{W}\|_{F}$ , which demonstrates the robustness of SSMP to the measurement noise. One can also deduce from (45) that $\mathbf{X}$ is recovered accurately (i.e., $\widehat{\mathbf{X}}=\mathbf{X}$ ) if SSMP finishes before running $K$ iterations in the noiseless scenario.

We next consider the case where SSMP finishes after running $K$ iterations. In our analysis, we first establish a condition of SSMP choosing all support elements (i.e., $S\subseteq S^{K}$ ) and then derive an upper bound of the reconstruction error $\|\mathbf{X}-\widehat{\mathbf{X}}\|_{F}$ under the obtained condition. The following proposition presents a condition under which SSMP picks at least one support element in the $(k+1)$ -th iteration.

Proposition 4.

Consider the system model in (44), where $\mathbf{A}$ has unit $\ell_{2}$ -norm columns, any $r$ nonzero rows of $\mathbf{X}$ are linearly independent, and $\|\mathbf{W}\|_{F}$ is bounded. Let $S^{k}$ be the estimated support generated in the $k$ -th iteration of the SSMP algorithm and $L$ be the number of indices chosen in each iteration. Suppose there exists at least one remaining support element after the $k$ -th iteration (i.e., $|S\setminus S^{k}|>0$ ). Also, suppose $\mathbf{A}$ satisfies the RIP of order $|S\cup S^{k}|+L$ and $\bar{\eta}=\|\mathbf{P}_{\mathcal{R}(\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{Y})}\|_{2}$ obeys

\displaystyle\bar{\eta}

\displaystyle<\sqrt{\frac{1-\delta_{|S\cup S^{k}|}}{1+\delta_{|S\cup S^{k}|}}}.

(49)

Then, the following statements hold:

(i)

If $|S\setminus S^{k}|>r$ and

	$\displaystyle\frac{2\bar{\eta}\sqrt{1+\delta_{\|S\cup S^{k}\|}}}{\sqrt{1-\delta_{\|S\cup S^{k}\|}}-\bar{\eta}\sqrt{1+\delta_{\|S\cup S^{k}\|}}}$
		$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ <\sqrt{\frac{r(1-\delta_{\|S\cup S^{k}\|})}{K}}-\min\left\{\delta_{\|S\cup S^{k}\|+1},\sqrt{\frac{r\delta^{2}_{\|S\cup S^{k}\|+L}}{L(1-\delta^{2}_{\|S^{k}\|+1})(1-\delta_{\|S\cup S^{k}\|})}}\right\},$		(50)

then SSMP chooses at least one support element in the $(k+1)$ -th iteration.

(ii)

If $|S\setminus S^{k}|\leq r$ and

\displaystyle\frac{2\bar{\eta}\sqrt{1+\delta_{|S\cup S^{k}|}}}{\sqrt{1-\delta_{|S\cup S^{k}|}}-\bar{\eta}\sqrt{1+\delta_{|S\cup S^{k}|}}}

\displaystyle<1-\delta_{|S\cup S^{k}|+1},

(51)

then SSMP picks $\min\{L,|S\setminus S^{k}|\}$ support elements in the $(k+1)$ -th iteration.

Proof.

We consider the following two cases: 1) $|S\setminus S^{k}|>r$ and 2) $|S\setminus S^{k}|\leq r$ .

1)

$|S\setminus S^{k}|>r$ :

Recall that SSMP chooses at least one support element in the $(k+1)$ -th iteration if the largest element $p_{1}$ in $\{\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}\}_{i\in S\setminus S^{k}}$ is larger than the $L$ -th largest element $q_{L}$ in $\{\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}\}_{i\in\Omega\setminus(S\cup S^{k})}$ , where $\mathbf{B}$ is the $\ell_{2}$ -normalized counterpart of $\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}$ (see Algorithm 1). In our proof, we construct a lower bound of $p_{1}$ and an upper bound of $q_{L}$ and then show that the former is larger than the latter under (49) and (50).

$\bullet$ Lower bound of $p_{1}$ :

Note that for each of $i\in S\setminus S^{k}$ ,

$\displaystyle\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\\|_{2}$	$\displaystyle\overset{}{=}\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}\mathbf{b}_{i}-(\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})})\mathbf{b}_{i}\\|_{2}$
	$\displaystyle\overset{(a)}{\geq}\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}\mathbf{b}_{i}\\|_{2}-\\|(\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})})\mathbf{b}_{i}\\|_{2}$
	$\displaystyle\hskip 1.05273pt\geq\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}\mathbf{b}_{i}\\|_{2}-\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\\|_{2},$	(52)

where (a) is from the triangle inequality. Thus, $p_{1}=\max_{i\in S\setminus S^{k}}\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}$ satisfies

	$\displaystyle p_{1}$	$\displaystyle\overset{}{\geq}\underset{i\in S\setminus S^{k}}{\max}\hskip 1.42262pt\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}\mathbf{b}_{i}\\|_{2}-\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\\|_{2}$
		$\displaystyle\overset{}{\geq}\sqrt{\frac{r(1-\delta_{\|S\cup S^{k}\|})}{K}}-\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\\|_{2},$		(53)

where the last inequality follows from (21).

$\bullet$ Upper bound of $q_{L}$ :

We construct an upper bound of $q_{L}$ in two different ways and then combine the results.

First, we note that for each of $i\in\Omega\setminus(S\cup S^{k})$ , $\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}$ satisfies

$\displaystyle\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\\|_{2}$	$\displaystyle\overset{}{=}\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}\mathbf{b}_{i}-(\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})})\mathbf{b}_{i}\\|_{2}$
	$\displaystyle\overset{(a)}{\leq}\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}\mathbf{b}_{i}\\|_{2}+\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\\|_{2}$	(54)
	$\displaystyle\overset{(b)}{\leq}\delta_{\|S\cup S^{k}\|+1}+\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\\|_{2},$	(55)

where (a) follows from the triangle inequality and (b) is from (22) and (23). Therefore, the $L$ -th largest element $q_{L}$ in $\{\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}\}_{i\in\Omega\setminus(S\cup S^{k})}$ also satisfies

\displaystyle q_{L}

\displaystyle\leq\delta_{|S\cup S^{k}|+1}+\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\|_{2}.

(56)

We next derive an upper bound of $q_{L}$ in a different way. Let $\psi_{l}$ be the index corresponding to the $l$ -th largest element $q_{l}$ in $\{\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}\}_{i\in\Omega\setminus(S\cup S^{k})}$ , then $q_{L}$ satisfies

$\displaystyle q_{L}$	$\displaystyle\overset{}{\leq}\frac{1}{L}(q_{1}+\ldots+q_{L})$
	$\displaystyle\overset{(a)}{\leq}\frac{1}{L}\underset{l=1}{\overset{L}{\sum}}\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}\mathbf{b}_{\psi_{l}}\\|_{2}+\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\\|_{2}$
	$\displaystyle\overset{(b)}{\leq}\sqrt{\frac{1}{L}\underset{l=1}{\overset{L}{\sum}}\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}\mathbf{b}_{\psi_{l}}\\|_{2}^{2}}+\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\\|_{2}$
	$\displaystyle\overset{(c)}{\leq}\sqrt{\frac{r\delta_{\|S\cup S^{k}\|+L}^{2}}{L(1-\delta_{\|S^{k}\|+1}^{2})(1-\delta_{\|S\cup S^{k}\|})}}+\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\\|_{2},$	(57)

where (a) is due to (54), (b) follows from the Cauchy-Schwarz inequality, and (c) is from (26) and (27).

Finally, by combining (56) and (57), we have the following upper bound of $q_{L}$ :

\displaystyle q_{L}

\displaystyle\leq\min\left\{\delta_{|S\cup S^{k}|+1},\sqrt{\frac{r\delta_{|S\cup S^{k}|+L}^{2}}{L(1-\delta_{|S^{k}|+1}^{2})(1-\delta_{|S\cup S^{k}|})}}\right\}+\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\|_{2}.

(58)

$\bullet$ When is $p_{1}>q_{L}$ ?

From (53) and (58), we have

\displaystyle p_{1}-q_{L}

\displaystyle\geq\sqrt{\frac{r(1-\delta_{|S\cup S^{k}|})}{K}}-\min\left\{\delta_{|S\cup S^{k}|+1},\sqrt{\frac{r\delta_{|S\cup S^{k}|+L}^{2}}{L(1-\delta_{|S^{k}|+1}^{2})(1-\delta_{|S\cup S^{k}|})}}\right\}-2\eta_{k},

(59)

where $\eta_{k}=\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\|_{2}$ . Also, it has been shown in [17, Proposition 7.6] that if the condition number $\kappa(\mathbf{A}_{S\cup S^{k}})$ of the matrix $\mathbf{A}_{S\cup S^{k}}$ obeys⁹⁹9In our case, (60) is satisfied since $\kappa(\mathbf{A}_{S\cup S^{k}})\bar{\eta}\overset{(a)}{<}\frac{\sigma_{\max}(\mathbf{A}_{S\cup S^{k}})}{\sigma_{\min}(\mathbf{A}_{S\cup S^{k}})}\sqrt{\frac{1-\delta_{|S\cup S^{k}|}}{1+\delta_{|S\cup S^{k}|}}}\overset{(b)}{\leq}1$ , where (a) is due to (49) and (b) follows from the definition of the RIP.

\displaystyle\kappa(\mathbf{A}_{S\cup S^{k}})\bar{\eta}

\displaystyle<1,

(60)

then $\eta_{k}=\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{Y})}\|_{2}$ satisfies¹⁰¹⁰10In [17, Proposition 7.6], it has been shown that (61) holds when $S^{k}\subset S$ . After small modifications, the proof can readily be extended to the case where $S^{k}\not\subset S$ .

$\displaystyle\eta_{k}$	$\displaystyle\hskip 1.05273pt\leq\frac{\bar{\eta}\cdot\kappa(\mathbf{A}_{S\cup S^{k}})}{1-\bar{\eta}\cdot\kappa(\mathbf{A}_{S\cup S^{k}})}$	(61)
	$\displaystyle\hskip 1.05273pt=\frac{\bar{\eta}\sigma_{\max}(\mathbf{A}_{S\cup S^{k}})}{\sigma_{\min}(\mathbf{A}_{S\cup S^{k}})-\bar{\eta}\sigma_{\max}(\mathbf{A}_{S\cup S^{k}})}$
	$\displaystyle\overset{(a)}{\leq}\frac{\bar{\eta}\sqrt{1+\delta_{\|S\cup S^{k}\|}}}{\sqrt{1-\delta_{\|S\cup S^{k}\|}}-\bar{\eta}\sqrt{1+\delta_{\|S\cup S^{k}\|}}},$	(62)

where (a) follows from the definition of the RIP. Then, by combining (59) and (62), we have

	$\displaystyle p_{1}-q_{L}$	$\displaystyle\geq\sqrt{\frac{r(1-\delta_{\|S\cup S^{k}\|})}{K}}-\min\left\{\delta_{\|S\cup S^{k}\|+1},\sqrt{\frac{r\delta_{\|S\cup S^{k}\|+L}^{2}}{L(1-\delta_{\|S^{k}\|+1}^{2})(1-\delta_{\|S\cup S^{k}\|})}}\right\}$
		$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ -\frac{2\bar{\eta}\sqrt{1+\delta_{\|S\cup S^{k}\|}}}{\sqrt{1-\delta_{\|S\cup S^{k}\|}}-\bar{\eta}\sqrt{1+\delta_{\|S\cup S^{k}\|}}}.$		(63)

One can easily check that under (50), the right-hand side of (63) is strictly larger than zero. Therefore, $p_{1}>q_{L}$ , and hence SSMP picks at least one support element (the index corresponding to $p_{1}$ ) in the $(k+1)$ -th iteration.

2)

$|S\setminus S^{k}|\leq r$ :

By combining (32a), (52), and (62), we have

\displaystyle\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}

\displaystyle\geq 1-\frac{\bar{\eta}\sqrt{1+\delta_{|S\cup S^{k}|}}}{\sqrt{1-\delta_{|S\cup S^{k}|}}-\bar{\eta}\sqrt{1+\delta_{|S\cup S^{k}|}}},\leavevmode\nobreak\ \forall i\in S\setminus S^{k}.

Also, by combining (55) and (62), we have

\displaystyle\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{j}\|_{2}

\displaystyle\leq\delta_{|S\cup S^{k}|+1}+\frac{\bar{\eta}\sqrt{1+\delta_{|S\cup S^{k}|}}}{\sqrt{1-\delta_{|S\cup S^{k}|}}-\bar{\eta}\sqrt{1+\delta_{|S\cup S^{k}|}}},\leavevmode\nobreak\ \forall j\in\Omega\setminus(S\cup S^{k}).

Then, for any support element $i\in S\setminus S^{k}$ and any incorrect index $j\in\Omega\setminus(S\cup S^{k})$ , we obtain

\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}>\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{j}\|_{2}

by (51). Therefore, the SSMP algorithm picks $\min\{L,|S\setminus S^{k}|\}$ support elements in the $(k+1)$ -th iteration. ∎

Proposition 4 indicates that if more than $r$ support elements remain after the $k$ -th iteration, then SSMP picks at least one support element in the $(k+1)$ -th iteration under (50). This in turn implies that SSMP chooses at least $K-r$ support elements in the first $K-r$ iterations, provided that the sampling matrix $\mathbf{A}$ obeys the RIP of order $L(K-r)+r+1$ and the corresponding RIP constant $\delta$ satisfies

\displaystyle\sqrt{\frac{r(1-\delta)}{K}}-\min\left\{\delta,\sqrt{\frac{r\delta^{2}}{L(1-\delta^{2})(1-\delta)}}\right\}

\displaystyle>\frac{2\bar{\eta}\sqrt{1+\delta}}{\sqrt{1-\delta}-\bar{\eta}\sqrt{1+\delta}}.

(64)

In particular, in the noiseless case ( $\|\mathbf{W}\|_{F}=0$ ), $\bar{\eta}=0$ so that (64) is satisfied under (14). Thus, SSMP picks at least $K-r$ support elements in the first $K-r$ iterations under (14), which coincides with the result in Proposition 2.

The next proposition presents a relationship between $\bar{\eta}$ and noise.

Proposition 5.

Consider the system model in (44) where $\|\mathbf{W}\|_{F}$ is bounded. If $\sigma_{\min}(\mathbf{A}\mathbf{X})>\sigma_{\max}\|\mathbf{W}\|_{F}$ , then $\bar{\eta}=\|\mathbf{P}_{\mathcal{R}(\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{Y})}\|_{2}$ satisfies

\displaystyle\bar{\eta}

\displaystyle\leq\left(\frac{\sigma_{\min}(\mathbf{A}\mathbf{X})}{\sigma_{\max}(\mathbf{W})}-1\right)^{-1}.

(65)

Proof.

Let $\mathcal{U}=\mathcal{R}(\mathbf{A}\mathbf{X})$ and $\mathcal{V}=\mathcal{R}(\mathbf{A}\mathbf{X}+\mathbf{W})$ , then it is well-known that [39, p. 275]

\displaystyle\|\mathbf{P}_{\mathcal{U}}-\mathbf{P}_{\mathcal{V}}\|_{2}

\displaystyle=\max\big{\{}\|\mathbf{P}_{\mathcal{U}}^{\perp}\mathbf{P}_{\mathcal{V}}\|_{2},\|\mathbf{P}_{\mathcal{V}}^{\perp}\mathbf{P}_{\mathcal{U}}\|_{2}\big{\}}.

Now, what remains is to show that $\max\{\|\mathbf{P}_{\mathcal{U}}^{\perp}\mathbf{P}_{\mathcal{V}}\|_{2},\|\mathbf{P}_{\mathcal{V}}^{\perp}\mathbf{P}_{\mathcal{U}}\|_{2}\}\leq\eta$ , where

\eta=\left(\frac{\sigma_{\min}(\mathbf{A}\mathbf{X})}{\sigma_{\max}(\mathbf{W})}-1\right)^{-1}=\frac{\sigma_{\max}(\mathbf{W})}{\sigma_{\min}(\mathbf{A}\mathbf{X})-\sigma_{\max}(\mathbf{W})}.

$\bullet$ $\|\mathbf{P}_{\mathcal{U}}^{\perp}\mathbf{P}_{\mathcal{V}}\|_{2}\leq\eta$ ?

Since

\|\mathbf{P}_{\mathcal{U}}^{\perp}\mathbf{P}_{\mathcal{V}}\|_{2}=\sup_{\mathbf{v}\in\mathcal{V},\|\mathbf{v}\|_{2}=1}\inf_{\mathbf{u}\in\mathcal{U}}\|\mathbf{v}-\mathbf{u}\|_{2},

it suffices to show that $\inf_{\mathbf{u}\in\mathcal{U}}\|\mathbf{v}-\mathbf{u}\|_{2}\leq\eta$ for any unit vector $\mathbf{v}$ in $\mathcal{V}$ . Let $\mathbf{v}=(\mathbf{A}\mathbf{X}+\mathbf{W})\bar{\mathbf{v}}$ be an arbitrary unit vector in $\mathcal{V}=\mathcal{R}(\mathbf{A}\mathbf{X}+\mathbf{W})$ , then

$\displaystyle 1$	$\displaystyle\hskip 1.05273pt=\\|(\mathbf{A}\mathbf{X}+\mathbf{W})\bar{\mathbf{v}}\\|_{2}$
	$\displaystyle\overset{(a)}{\geq}\\|\mathbf{A}\mathbf{X}\bar{\mathbf{v}}\\|_{2}-\\|\mathbf{W}\bar{\mathbf{v}}\\|_{2}$
	$\displaystyle\hskip 1.05273pt\geq(\sigma_{\min}(\mathbf{A}\mathbf{X})-\sigma_{\max}(\mathbf{W}))\\|\bar{\mathbf{v}}\\|_{2},$	(66)

where (a) follows from the triangle inequality. Note that $\inf_{\mathbf{u}\in\mathcal{U}}\|\mathbf{v}-\mathbf{u}\|_{2}\leq\|\mathbf{v}-\bar{\mathbf{u}}\|_{2}$ for any $\bar{\mathbf{u}}\in\mathcal{U}=\mathcal{R}(\mathbf{A}\mathbf{X})$ . In particular, when $\bar{\mathbf{u}}=\mathbf{A}\mathbf{X}\bar{\mathbf{v}}$ , we have

	$\displaystyle\inf_{\mathbf{u}\in\mathcal{U}}\\|\mathbf{v}-\mathbf{u}\\|_{2}$	$\displaystyle\leq\\|\mathbf{v}-\bar{\mathbf{u}}\\|_{2}=\\|\mathbf{W}\bar{\mathbf{v}}\\|_{2}$
		$\displaystyle\leq\sigma_{\max}(\mathbf{W})\\|\bar{\mathbf{v}}\\|_{2}$
		$\displaystyle\leq\frac{\sigma_{\max}(\mathbf{W})}{\sigma_{\min}(\mathbf{A}\mathbf{X})-\sigma_{\max}(\mathbf{W})}=\eta,$

where the last inequality follows from (66).

$\bullet$ $\|\mathbf{P}_{\mathcal{V}}^{\perp}\mathbf{P}_{\mathcal{U}}\|_{2}\leq\eta$ ?

Let $\mathbf{u}=\mathbf{A}\mathbf{X}\bar{\mathbf{u}}$ be an arbitrary unit vector in $\mathcal{U}=\mathcal{R}(\mathbf{A}\mathbf{X})$ , then $\|\bar{\mathbf{u}}\|_{2}\leq 1/\sigma_{\min}(\mathbf{A}\mathbf{X})$ . Also, let $\bar{\mathbf{v}}=(\mathbf{A}\mathbf{X}+\mathbf{W})\bar{\mathbf{u}}\in\mathcal{V}$ , then

	$\displaystyle\underset{\mathbf{v}\in\mathcal{V}}{\inf}\hskip 0.85358pt\\|\mathbf{u}-\mathbf{v}\\|_{2}$	$\displaystyle\leq\\|\mathbf{u}-\bar{\mathbf{v}}\\|_{2}$
		$\displaystyle\leq\sigma_{\max}(\mathbf{W})\\|\bar{\mathbf{u}}\\|_{2}$
		$\displaystyle\leq\frac{\sigma_{\max}(\mathbf{W})}{\sigma_{\min}(\mathbf{A}\mathbf{X})}$
		$\displaystyle\leq\eta.$

Since $\mathbf{u}$ is an arbitrary unit vector in $\mathcal{U}$ , we have $\|\mathbf{P}_{\mathcal{V}}^{\perp}\mathbf{P}_{\mathcal{U}}\|_{2}=\sup_{\mathbf{u}\in\mathcal{U},\|\mathbf{u}\|_{2}=1}\inf_{\mathbf{v}\in\mathcal{V}}\|\mathbf{u}-\mathbf{v}\|_{2}\leq\eta$ , which is the desired result. ∎

Since $\sigma_{\max}(\mathbf{W})\leq\|\mathbf{W}\|_{F}$ , it is clear from Proposition 5 that

\displaystyle\bar{\eta}

\displaystyle\leq\left(\frac{\sigma_{\min}(\mathbf{A}\mathbf{X})}{\|\mathbf{W}\|_{F}}-1\right)^{-1}.

One can observe that the upper bound increases with the noise power $\|\mathbf{W}\|_{F}$ . In particular, if $\|\mathbf{W}\|_{F}=0$ , then $\bar{\eta}=0$ , which in turn implies that the measurement space $\mathcal{R}(\mathbf{Y})$ coincides with the signal space $\mathcal{R}(\mathbf{A}\mathbf{X})$ .

Having the results of Propositions 4 and 5 in hand, we are now ready to establish a condition under which SSMP picks all support elements.

Theorem 3.

Consider the system model in (44), where $\mathbf{A}$ has unit $\ell_{2}$ -norm columns, any $r$ nonzero rows of $\mathbf{X}$ are linearly independent, and $\|\mathbf{W}\|_{F}$ is bounded. Suppose SSMP chooses $L$ $(L\leq\min\{K,\frac{m}{K}\})$ indices in each iteration. Also, suppose $\mathbf{A}$ obeys the RIP of order $L(K-r)+r+1$ and the corresponding RIP constant $\delta$ satisfies

\displaystyle 0\leq\eta=\left(\frac{\sigma_{\min}(\mathbf{A}\mathbf{X})}{\sigma_{\max}(\mathbf{W})}-1\right)^{-1}<\sqrt{\frac{1-\delta}{1+\delta}}.

(67)

Then SSMP picks all support elements in at most $K-r+\lceil\frac{r}{L}\rceil$ iterations if one of the following conditions is satisfied:

(i)

$r=K$ and $\delta$ satisfies

$\displaystyle 1-\delta$ $\displaystyle>\frac{2\eta\sqrt{1+\delta}}{\sqrt{1-\delta}-\eta\sqrt{1+\delta}}.$ (68)

(ii)

$r<K$ and $\delta$ satisfies

\displaystyle\sqrt{\frac{r(1-\delta)}{K}}-\min\left\{\delta,\sqrt{\frac{r\delta^{2}}{L(1-\delta^{2})(1-\delta)}}\right\}

\displaystyle>\frac{2\eta\sqrt{1+\delta}}{\sqrt{1-\delta}-\eta\sqrt{1+\delta}}.

(69)

Proof.

By Propositions 4 and 5, the SSMP algorithm chooses at least $K-r$ support elements in the first $K-r$ iterations under (69). Furthermore, similar to the proof of Theorem 1, one can show that if SSMP picks at least $K-r$ support elements in the first $K-r$ iterations, then SSMP chooses the remaining support elements by running $\lceil\frac{|S\setminus S^{K-r}|}{L}\rceil\leavevmode\nobreak\ (\leq\lceil\frac{r}{L}\rceil)$ additional iterations under (68). Also, using $r,L\leq K$ , one can easily show that

\sqrt{\frac{r(1-\delta)}{K}}-\min\left\{\delta,\sqrt{\frac{r\delta^{2}}{L(1-\delta^{2})(1-\delta)}}\right\}<1-\delta,

and thus (68) is satisfied under (69). By combining these results, we can conclude that SSMP picks all support elements in at most $K-r+\lceil\frac{r}{L}\rceil$ iterations if (i) or (ii) holds. ∎

In the noiseless scenario ( $\|\mathbf{W}\|_{F}=0$ ), $\eta=0$ so that (69) is satisfied under (14). Combining this with Theorem 3, one can see that SSMP chooses all support elements and recovers $\mathbf{X}$ accurately in at most $K-r+\lceil\frac{r}{L}\rceil$ iterations under (14), which is consistent with the result in Theorem 1. One can also infer from Theorem 3 that all support elements are chosen if

\displaystyle\eta

\displaystyle<\sqrt{\frac{1-\delta}{1+\delta}}\frac{f(\delta,r)}{2+f(\delta,r)},

(70)

where

f(\delta,r)=\sqrt{\frac{r(1-\delta)}{K}}-\min\left\{\delta,\sqrt{\frac{r\delta^{2}}{L(1-\delta^{2})(1-\delta)}}\right\}.

Note that $f(\delta,r)$ is a decreasing function of $\delta$ and thus the upper bound in (70) also decreases with $\delta$ . Then, since $\eta$ decreases with $\sigma_{\min}(\mathbf{A}\mathbf{X})/\sigma_{\max}(\mathbf{W})$ , the RIP condition in (70) becomes less restrictive when $\sigma_{\min}(\mathbf{A}\mathbf{X})/\sigma_{\max}(\mathbf{W})$ increases. Furthermore, note that

f(\delta,r)=\max\left\{\sqrt{\frac{r(1-\delta)}{K}}-\delta,\sqrt{r}\left(\sqrt{\frac{1-\delta}{K}}-\sqrt{\frac{\delta^{2}}{L(1-\delta^{2})(1-\delta)}}\right)\right\},

increases with the number $r$ of (linearly independent) measurement vectors so that the upper bound in (70) also increases with $r$ . Thus, when $r$ increases, the requirement on the RIP constant becomes less restrictive (see Fig. 3). This behavior seems to be natural but has not been reported for conventional joint sparse recovery algorithms such as SOMP, M-ORMP, and mixed norm minimization techniques [13, 11, 14, 16]. Moreover, we note that if all support elements are chosen, then the output $\widehat{\mathbf{X}}$ of SSMP satisfies¹¹¹¹11We note that the result in (71) is the extension of the result in [34, eq. (19)], which is related to the SMV version of SSMP, to the MMV scenario. This extension can be easily done by taking similar steps to the proofs of (46) and (47) in Appendix D.

\displaystyle\begin{cases}\|\mathbf{X}-\widehat{\mathbf{X}}\|_{F}\leq\frac{\|\mathbf{W}\|_{F}}{\sqrt{1-\delta_{K}}},&L=1,\\ \|\mathbf{X}-\widehat{\mathbf{X}}\|_{F}\leq\left(1+\sqrt{\frac{1+\delta_{2K}}{1-\delta_{LK}}}\right)\frac{2\|\mathbf{W}\|_{F}}{\sqrt{1-\delta_{2K}}},&L>1.\end{cases}

(71)

This means that the reconstruction error $\|\mathbf{X}-\widehat{\mathbf{X}}\|_{F}$ of SSMP is upper bounded by a constant multiple of the noise power $\|\mathbf{W}\|_{F}$ , which clearly demonstrates the robustness of the SSMP algorithm to measurement noise.

Finally, it is worth mentioning that our analysis can readily be extended to the scenario where the input signal $\mathbf{X}$ is approximately row $K$ -sparse (a.k.a. row compressible), meaning that $0<\|\mathbf{X}-\mathbf{X}_{K}\|_{F}\leq\rho\|\mathbf{X}\|_{F}$ for some small $\rho>0$ .¹²¹²12 $\mathbf{X}_{K}$ is the matrix obtained from $\mathbf{X}$ by maintaining the $K$ rows with largest $\ell_{2}$ -norms and setting all the other rows to the zero vector. In this case, one can establish a condition ensuring the robustness of the SSMP algorithm by 1) partitioning the observation matrix $\mathbf{Y}$ as

\displaystyle\mathbf{Y}

\displaystyle=\mathbf{A}\mathbf{X}_{K}+(\mathbf{A}(\mathbf{X}-\mathbf{X}_{K})+\mathbf{W}),

(72)

2) considering $\widehat{\mathbf{W}}=\mathbf{A}(\mathbf{X}-\mathbf{X}_{K})+\mathbf{W}$ as modified measurement noise satisfying (see Appendix E)

\displaystyle\|\widehat{\mathbf{W}}\|_{F}

\displaystyle\leq\sqrt{1+\delta_{K}}\left(\|\mathbf{X}-\mathbf{X}_{K}\|_{F}+\frac{1}{\sqrt{K}}\|\mathbf{X}-\mathbf{X}_{K}\|_{1,2}\right)+\|\mathbf{W}\|_{F},

(73)

and then 3) applying Theorems 2 and 3 to the system model in (72).

V SSMP Running More Than $K$ Iterations

Thus far, we have analyzed the performance of SSMP running at most $K$ iterations. Our result in Theorem 1 implies that if the number $r$ of measurement vectors is on the order of $K$ (e.g., $r=\lceil\frac{K}{2}\rceil$ ), then SSMP recovers any row $K$ -sparse matrix accurately with overwhelming probability as long as the number $m$ of random measurements scales linearly with $K\log\frac{n}{K}$ . However, if $r$ is not on the order of $K$ (e.g., $r=1$ ), then the upper bound of the proposed guarantee (14) is inversely proportional to $\sqrt{K}$ , which requires that $m$ should scale with $K^{2}\log\frac{n}{K}$ (see (37)).

In the compressed sensing literature, there have been some efforts to improve the performance guarantee by running an algorithm more than $K$ iterations [41, 42, 45, 43, 44]. For example, it has been shown in [42, Theorem 6.25] that the optimal performance guarantee $\delta_{K+1}<\frac{1}{\sqrt{K+1}}$ of the conventional OMP algorithm running $K$ iterations [40] can be relaxed to $\delta_{13K}<\frac{1}{6}$ if it runs $12K$ iterations. In this section, we show that if $r=1$ , then by running more than $K$ iterations, SSMP ensures exact reconstruction with $\mathcal{O}(K\log\frac{n}{K})$ random measurements.

V-A Main Results

In this subsection, we provide our main results. In the next theorem, we demonstrate that if $\gamma$ support elements remain after some iterations, then under a suitable RIP condition, SSMP picks the remaining support elements by running a specified number of additional iterations.

Theorem 4.

Consider the SSMP algorithm and the system model in (1) where $\mathbf{A}$ has unit $\ell_{2}$ -norm columns. Let $L$ be the number of indices chosen in each iteration and $\gamma$ be the number of remaining support elements after the $k$ -th iteration. Let $c$ be an integer such that $c\geq 2$ . If $\mathbf{A}$ obeys the RIP of order $Lk+\lfloor\gamma(1+4c-\frac{4c}{e^{c}-1})\rfloor$ and the corresponding RIP constant $\delta$ satisfies


	$\displaystyle c\geq-\frac{2}{(1-\delta)^{2}}\log\frac{1}{2},$		(74a)
	$\displaystyle c\geq-\frac{1}{(1-\delta)^{2}}\log\left(\frac{1}{2}-\sqrt{\frac{\delta}{2(1+\delta)}}\right),$		(74b)
	$\displaystyle c>-\frac{1}{(1-\delta)^{2}}\log\left(\frac{1}{2}-\frac{\delta}{(1+\delta)(1-\delta)^{2}}\right),$		(74c)

then

\displaystyle S\subset S^{k+\max\{\gamma,\lfloor\frac{4c\gamma}{L}\rfloor\}}.

(75)

Proof.

See Section V-B. ∎

Theorem 4 implies that if $\gamma$ support elements remain, then under (74a)-(74c), SSMP chooses all these elements by running $\max\{\gamma,\lfloor\frac{4c\gamma}{L}\rfloor\}$ additional iterations. In particular, when $\gamma=K$ , we obtain the following result.

Theorem 5.

Consider the system model in (1) where $\mathbf{A}$ has unit $\ell_{2}$ -norm columns. Let $L$ be the number of indices chosen in each iteration of the SSMP algorithm and $c$ be an integer such that $c\geq 2$ . Suppose $\mathbf{A}$ obeys the RIP of order $\lfloor K(1+4c-\frac{4c}{e^{c}-1})\rfloor$ and the corresponding RIP constant $\delta$ satisfies (74a)-(74c). Then the SSMP algorithm accurately recovers $\mathbf{x}$ from $\mathbf{y}=\mathbf{A}\mathbf{x}$ in $\max\{K,\lfloor\frac{4cK}{L}\rfloor\}$ iterations.

Note that if $c=2$ , then (74a)-(74c) are satisfied under $\delta\leq 0.167$ , $\delta\leq 0.155$ , and $\delta<0.185$ , respectively. Combining this with Theorem 5, one can see that SSMP ensures exact recovery of any $K$ -sparse vector in $\max\{K,\lfloor\frac{8K}{L}\rfloor\}$ iterations under

\displaystyle\delta_{\lfloor 7.8K\rfloor}

\displaystyle\leq 0.155.

(76)

In particular, our result indicates that the OLS algorithm, a special case of SSMP for $r=L=1$ (see Table II), recovers any $K$ -sparse vector accurately in $8K$ iterations under (76). In Table III, we summarize the performance guarantee of SSMP for different choices of $c$ .

The beneficial point of (76) is that the upper bound is a constant and unrelated to the sparsity $K$ . This implies that by running slightly more than $K$ iterations, SSMP accurately recovers any $K$ -sparse vector with overwhelming probability with $\mathcal{O}(K\log\frac{n}{K})$ random measurements (see (37)). This is in contrast to the guarantees (39a) and (39b) of SSMP running $K$ iterations, which require that the number of random measurements should scale with $K^{2}\log\frac{n}{K}$ .

Table III: The Result in Theorem 5 With Different

c

	Number of Iterations	Performance Guarantee
$c=2$	$\max\left\{K,\left\lfloor\frac{8K}{L}\right\rfloor\right\}$	$\delta_{\lfloor 7.8K\rfloor}\leq 0.155$
$c=3$	$\max\left\{K,\left\lfloor\frac{12K}{L}\right\rfloor\right\}$	$\delta_{\lfloor 12.4K\rfloor}<0.235$
$c=4$	$\max\left\{K,\left\lfloor\frac{16K}{L}\right\rfloor\right\}$	$\delta_{\lfloor 16.8K\rfloor}<0.263$
$c=5$	$\max\left\{K,\left\lfloor\frac{20K}{L}\right\rfloor\right\}$	$\delta_{\lfloor 20.9K\rfloor}<0.281$
$c=6$	$\max\left\{K,\left\lfloor\frac{24K}{L}\right\rfloor\right\}$	$\delta_{25K}<0.291$

V-B Proof of Theorem 4

The proof of Theorem 4 is based on induction on the number $\gamma$ of remaining support elements. We note that this proof technique is similar in spirit to the works of Zhang [41] and Foucart and Rauhut [42].

First, we consider the case where $\gamma=0$ . In this case, $S\setminus S^{k}=\emptyset$ and thus

\displaystyle S\subset S^{k}=S^{k+\max\{\gamma,\lfloor\frac{4c\gamma}{L}\rceil\}}.

Next, we assume that the argument holds up to an integer $\gamma-1\leavevmode\nobreak\ (\gamma\geq 1)$ . In other words, we assume that if the number $t$ of remaining support elements is less than $\gamma-1$ , then SSMP chooses all these elements by running $\max\{t,\lfloor\frac{4ct}{L}\rfloor\}$ additional iterations. Under this assumption, we show that if $\gamma$ support elements remain after the $k$ -th iteration, then

\displaystyle S\subset S^{k+\max\{\gamma,\lfloor\frac{4c\gamma}{L}\rfloor\}}.

(77)

1)

Preliminaries

Before we proceed, we define notations used in our analysis. Without loss of generality, let $\Gamma^{k}=S\setminus S^{k}=\{1,\ldots,\gamma\}$ and $|x_{1}|\geq\ldots\geq|x_{\gamma}|$ . We define the subset $\Gamma_{\tau}^{k}$ of $\Gamma^{k}$ as

\displaystyle\Gamma^{k}_{\tau}

\displaystyle=\begin{cases}\emptyset,&\tau=0,\\ \{1,\ldots,2^{\tau-1}L\},&\tau=1,\ldots,\max\{0,\lceil\log_{2}\frac{\gamma}{L}\rceil\},\\ \Gamma^{k},&\tau=\max\{0,\lceil\log_{2}\frac{\gamma}{L}\rceil\}+1.\end{cases}

(78)

Let

\displaystyle\sigma

\displaystyle=\frac{1}{2}\exp(c(1-\delta)^{2}),

(79)

then $\sigma\geq 2$ by (74a). Also, let $N$ be the integer such that


$\displaystyle\\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma^{k}_{0}}\\|_{2}^{2}$	$\displaystyle<\sigma\\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma^{k}_{1}}\\|_{2}^{2},$	(80a)
$\displaystyle\\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma^{k}_{1}}\\|_{2}^{2}$	$\displaystyle<\sigma\\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma^{k}_{2}}\\|_{2}^{2},$	(80b)
	$\displaystyle\leavevmode\nobreak\ \hskip 1.99168pt\vdots$
$\displaystyle\\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma^{k}_{N-2}}\\|_{2}^{2}$	$\displaystyle<\sigma\\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma^{k}_{N-1}}\\|_{2}^{2},$	(80c)
$\displaystyle\\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma^{k}_{N-1}}\\|_{2}^{2}$	$\displaystyle\geq\sigma\\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma^{k}_{N}}\\|_{2}^{2}.$	(80d)

If (80d) holds for $N=1$ , then we take $N=1$ . We note that $N$ always exists, since $\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma^{k}_{\tau}}\|_{2}=0$ when $\tau=\max\{0,\lceil\log_{2}\frac{\gamma}{L}\rceil\}+1$ so that (80d) holds at least for $N=\max\{0,\lceil\log_{2}\frac{\gamma}{L}\rceil\}+1$ . From (80a)-(80d), one can easily see that

\displaystyle\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\|_{2}^{2}

\displaystyle\leq\sigma^{N-1-\tau}\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{N-1}^{k}}\|_{2}^{2}

(81)

for each of $\tau\in\{0,\ldots,N\}$ . Also, if $N\geq 2$ , then we have

	$\displaystyle\gamma$	$\displaystyle\overset{(a)}{>}\frac{2\sigma-1}{2\sigma-2}2^{N-2}L$
		$\displaystyle\overset{(b)}{>}\frac{e^{c}-1}{e^{c}-2}2^{N-2}L,$		(82)

where (a) follows from [45, eq. (21)] and (b) is because $\sigma<\frac{1}{2}e^{c}$ by (79). We next provide lemmas useful in our proof.

Lemma 8.

For any integer $l$ satisfying $l\geq k$ and $\tau\in\{1,\ldots,\max\{0,\lceil\log_{2}\frac{\gamma}{L}\rceil\}+1\}$ , the residual of the SSMP algorithm satisfies

\displaystyle\|\mathbf{r}^{l}\|_{2}^{2}-\|\mathbf{r}^{l+1}\|_{2}^{2}

\displaystyle\geq\frac{(1-\delta_{|S^{l}|+1}^{2})(1-\delta_{|\Gamma_{\tau}^{k}\cup S^{l}|})}{\lceil\frac{|\Gamma_{\tau}^{k}|}{L}\rceil(1+\delta_{L})}\left(\|\mathbf{r}^{l}\|_{2}^{2}-\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\|_{2}^{2}\right).

(83)

Proof.

See Appendix F. ∎

Lemma 9.

For any integer $l\geq k$ , $\Delta l>0$ , and $\tau\in\{1,\ldots,\max\{0,\lceil\log_{2}\frac{\gamma}{L}\rceil\}+1\}$ , the residual $r^{l+\Delta l}$ generated in the $(l+\Delta l)$ -th iteration of the SSMP algorithm satisfies

\displaystyle\|\mathbf{r}^{l+\Delta l}\|_{2}^{2}-\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\|_{2}^{2}

\displaystyle\leq C_{\tau,l,\Delta l}\left(\|\mathbf{r}^{l}\|_{2}^{2}-\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\|_{2}^{2}\right),

(84)

where

\displaystyle C_{\tau,l,\Delta l}

\displaystyle=\exp\left(-\frac{\Delta l(1-\delta_{|S^{l+\Delta l-1}|+1}^{2})(1-\delta_{|\Gamma_{\tau}^{k}\cup S^{l+\Delta l-1}|})}{\lceil\frac{|\Gamma_{\tau}^{k}|}{L}\rceil(1+\delta_{L})}\right).

(85)

Proof.

See Appendix G. ∎

2)

Sketch of Proof

We now proceed to the proof of (77). In our proof, we consider the following two cases: i) $N\geq 2$ and ii) $N=1$ .

i) $N\geq 2$ : In this case, one can show that (see justifications in Section V-C)

\displaystyle\|\mathbf{x}_{\Gamma^{k_{N}}}\|_{2}^{2}

\displaystyle\leq\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{N-1}^{k}}\|_{2}^{2},

(86)

where $\Gamma^{k_{N}}=S\setminus S^{k_{N}}$ and

\displaystyle k_{N}

\displaystyle=k+c\sum_{\tau=1}^{N}\left\lceil\frac{|\Gamma^{k}_{\tau}|}{L}\right\rceil.

(87)

This implies that $|\Gamma^{k_{N}}|\leq|\Gamma^{k}\setminus\Gamma_{N-1}^{k}|=\gamma-2^{N-2}L$ , since otherwise

\|\mathbf{x}_{\Gamma^{k_{N}}}\|_{2}^{2}>|x_{2^{N-2}L+1}|^{2}+\ldots+|x_{\gamma}|^{2}=\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{N-1}^{k}}\|_{2}^{2},

which is a contradiction to (86). In other words, at most $\gamma-2^{N-2}L$ support elements remain after the $k_{N}$ -th iteration. Then, by the induction hypothesis, SSMP picks the remaining support elements by running $\max\{\gamma-2^{N-2}L,\lfloor\frac{4c(\gamma-2^{N-2}L)}{L}\rfloor\}$ additional iterations, i.e.,

\displaystyle S\subset S^{k_{N}+\max\{\gamma-2^{N-2}L,\lfloor\frac{4c(\gamma-2^{N-2}L)}{L}\rfloor\}}.

(88)

Note that

$\displaystyle k_{N}+\max\left\{\gamma-2^{N-2}L,\left\lfloor\frac{4c(\gamma-2^{N-2}L)}{L}\right\rfloor\right\}$
	$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \overset{(a)}{<}k+c2^{N}+\max\left\{\gamma-2^{N-2}L,\left\lfloor\frac{4c\gamma}{L}\right\rfloor-c2^{N}\right\}$
	$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \hskip 1.05273pt=k+\max\left\{\gamma-2^{N-2}(L-4c),\left\lfloor\frac{4c\gamma}{L}\right\rfloor\right\},$	(89)

where (a) is because

k_{N}=k+c\sum_{\tau=1}^{N}\left\lceil\frac{|\Gamma^{k}_{\tau}|}{L}\right\rceil\leq k+c\sum_{\tau=1}^{N}2^{\tau-1}<k+c2^{N}.

Then, by combining (88) and (89), we have

\displaystyle S

\displaystyle\subset S^{k+\max\left\{\gamma-2^{N-2}(L-4c),\left\lfloor\frac{4c\gamma}{L}\right\rfloor\right\}}.

(90)

Also, note that¹³¹³13If $L\geq 4c$ , then it is trivial. If $L<4c$ , then $\gamma+2^{N-2}(4c-L)\leq\gamma+2^{\lceil\log_{2}\frac{\gamma}{L}\rceil-1}(4c-L)<\frac{4c\gamma}{L}$ so that $\max\left\{\gamma-2^{N-2}(L-4c),\left\lfloor\frac{4c\gamma}{L}\right\rfloor\right\}=\lfloor\frac{4c\gamma}{L}\rfloor\leq\max\left\{\gamma,\left\lfloor\frac{4c\gamma}{L}\right\rfloor\right\}$ .

\displaystyle\max\left\{\gamma-2^{N-2}(L-4c),\left\lfloor\frac{4c\gamma}{L}\right\rfloor\right\}

\displaystyle\leq\max\left\{\gamma,\left\lfloor\frac{4c\gamma}{L}\right\rfloor\right\}.

Using this together with (90), we obtain (77), which completes the proof.

ii) $N=1$ : In this case, one can show that (see justifications in Section V-D)

\displaystyle\|\mathbf{x}_{\Gamma^{k+1}}\|_{2}^{2}

\displaystyle<\|\mathbf{x}_{\Gamma^{k}}\|_{2}^{2},

(91)

which in turn implies that $|\Gamma^{k+1}|<|\Gamma^{k}|=\gamma$ . In other words, at most $\gamma-1$ support elements remain after the $(k+1)$ -th iteration. Then, by the induction hypothesis, SSMP picks the remaining support elements by running $\max\{\gamma-1,\lfloor\frac{4c(\gamma-1)}{L}\rfloor\}$ iterations, i.e.,

\displaystyle S\subset S^{k+1+\max\left\{\gamma-1,\left\lfloor\frac{4c(\gamma-1)}{L}\right\rfloor\right\}}.

(92)

Also, one can show that¹⁴¹⁴14If $L\leq 4c$ , then it is trivial. If $L>4c$ , then $\frac{4c\gamma}{L}+\frac{L-4c}{L}<\gamma+(1-\frac{4c}{L})$ so that $\lfloor\frac{4c\gamma}{L}+\frac{L-4c}{L}\rfloor\leq\gamma$ . As a result, $\max\{\gamma,\lfloor\frac{4c\gamma}{L}+\frac{L-4c}{L}\rfloor\}=\gamma\leq\max\{\gamma,\lfloor\frac{4c\gamma}{L}\rfloor\}$ .

	$\displaystyle k+1+\max\left\{\gamma-1,\left\lfloor\frac{4c(\gamma-1)}{L}\right\rfloor\right\}$	$\displaystyle=k+\max\left\{\gamma,\left\lfloor\frac{4c\gamma}{L}+\frac{L-4c}{L}\right\rfloor\right\}$
		$\displaystyle\leq k+\max\left\{\gamma,\left\lfloor\frac{4c\gamma}{L}\right\rfloor\right\}.$		(93)

By combining (92) and (93), we obtain (77), which is the desired result.

V-C Proof of (86)

In our proof, we build an upper bound of $\|\mathbf{x}_{\Gamma^{k_{N}}}\|_{2}^{2}$ and a lower bound of $\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{N-1}^{k}}\|_{2}^{2}$ and then show that the former is smaller than the latter under (74b).

$\bullet$ Upper bound of $\|\mathbf{x}_{\Gamma^{k_{N}}}\|_{2}^{2}$ :

Note that

$\displaystyle\\|\mathbf{r}^{k_{N}}\\|_{2}^{2}$	$\displaystyle\hskip 1.05273pt=\\|\mathbf{P}_{S^{k_{N}}}^{\perp}\mathbf{A}_{S\setminus S^{k_{N}}}\mathbf{x}_{S\setminus S^{k_{N}}}\\|_{2}^{2}$
	$\displaystyle\overset{(a)}{\geq}(1-\delta_{\|S\cup S^{k_{N}}\|})\\|\mathbf{x}_{S\setminus S^{k_{N}}}\\|_{2}^{2}$
	$\displaystyle\hskip 1.05273pt=(1-\delta_{\|S\cup S^{k_{N}}\|})\\|\mathbf{x}_{\Gamma^{k_{N}}}\\|_{2}^{2},$	(94)

where (a) is from Lemma 4. Also, one can show that (see justifications in Appendix H)

\displaystyle|S\cup S^{k_{N}}|

\displaystyle\leq Lk+\left\lfloor\gamma\left(1+4c-\frac{4c}{e^{c}-1}\right)\right\rfloor,

(95)

and thus the RIP constant $\delta$ of order $Lk+\lfloor\gamma(1+4c-\frac{4c}{e^{c}-1})\rfloor$ satisfies

\displaystyle\delta_{|S\cup S^{k_{N}}|}

\displaystyle\leq\delta

(96)

by Lemma 3. Using this together with (94), we obtain

\displaystyle\|\mathbf{x}_{\Gamma^{k_{N}}}\|_{2}^{2}

\displaystyle\leq\frac{\|\mathbf{r}^{k_{N}}\|_{2}^{2}}{1-\delta}.

(97)

$\bullet$ Lower bound of $\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma^{k}_{N-1}}\|_{2}^{2}$ :

Let $k_{0}=k$ and $k_{i}=k+c\sum_{\tau=1}^{i}\lceil\frac{|\Gamma_{\tau}^{k}|}{L}\rceil$ for each of $i\in\{1,\ldots,N\}$ . Then, by applying Lemma 9 with $\tau=i$ , $l=k_{i-1}$ , and $\Delta l=k_{i}-k_{i-1}$ , we have

\displaystyle\|\mathbf{r}^{k_{i}}\|_{2}^{2}-\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{i}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{i}^{k}}\|_{2}^{2}

\displaystyle\leq C_{i,k_{i-1},k_{i}-k_{i-1}}\left(\|\mathbf{r}^{k_{i-1}}\|_{2}^{2}-\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{i}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{i}^{k}}\|_{2}^{2}\right)

(98)

for each of $i\in\{1,\ldots,N\}$ . Note that

$\displaystyle C_{i,k_{i-1},k_{i}-k_{i-1}}$	$\displaystyle\overset{(a)}{=}\exp\left(-\frac{c(1-\delta_{\|S^{k_{i}-1}\|+1}^{2})(1-\delta_{\|\Gamma_{i}^{k}\cup S^{k_{i}-1}\|})}{1+\delta_{L}}\right)$
	$\displaystyle\overset{(b)}{\leq}\exp\left(-c(1-\delta_{\|S\cup S^{k_{N}}\|})^{2}\right)$
	$\displaystyle\overset{(c)}{\leq}\exp\left(-c(1-\delta)^{2}\right),$	(99)

where (a), (b), and (c) follow from (85), Lemma 3, and (96) respectively. By combining (98) and (99), we have

\displaystyle\|\mathbf{r}^{k_{i}}\|_{2}^{2}

\displaystyle\hskip 1.05273pt\leq\exp\left(-c(1-\delta)^{2}\right)\|\mathbf{r}^{k_{i-1}}\|_{2}^{2}+\left(1-\exp\left(-c(1-\delta)^{2}\right)\right)\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{i}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{i}^{k}}\|_{2}^{2}

(100)

for each of $i\in\{1,\ldots,N\}$ .¹⁵¹⁵15If $\|\mathbf{r}^{k_{i-1}}\|_{2}^{2}-\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{i}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{i}^{k}}\|_{2}^{2}<0$ , then (100) holds since $\|\mathbf{r}^{k_{i}}\|_{2}\leq\|\mathbf{r}^{k_{i-1}}\|_{2}^{2}$ due to the orthogonal projection at each iteration of SSMP and $\|\mathbf{r}^{k_{i-1}}\|_{2}^{2}<\exp\left(-c(1-\delta)^{2}\right)\|\mathbf{r}^{k_{i-1}}\|_{2}^{2}+\left(1-\exp\left(-c(1-\delta)^{2}\right)\right)\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{i}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{i}^{k}}\|_{2}^{2}$ . Now, after taking similar steps to [43, p. 4201], one can show that

\displaystyle\|\mathbf{r}^{k_{N}}\|_{2}^{2}

\displaystyle\leq 4(1+\delta)\exp\left(-c(1-\delta)^{2}\right)\left(1-\exp\left(-c(1-\delta)^{2}\right)\right)\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{N-1}^{k}}\|_{2}^{2},

which is equivalent to

\displaystyle\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{N-1}^{k}}\|_{2}^{2}

\displaystyle\geq\frac{\|\mathbf{r}^{k_{N}}\|_{2}^{2}}{4(1+\delta)\exp\left(-c(1-\delta)^{2}\right)\left(1-\exp\left(-c(1-\delta)^{2}\right)\right)}.

(101)

$\bullet$ When is $\|\mathbf{x}_{\Gamma^{k_{N}}}\|_{2}^{2}\leq\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma^{k}_{N-1}}\|_{2}^{2}$ ?

From (97) and (101), we have

\displaystyle\frac{\|\mathbf{x}_{\Gamma^{k_{N}}}\|_{2}^{2}}{\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{N-1}^{k}}\|_{2}^{2}}

\displaystyle\leq\frac{4(1+\delta)\exp\left(-c(1-\delta)^{2}\right)\left(1-\exp\left(-c(1-\delta)^{2}\right)\right)}{1-\delta}.

(102)

One can easily check that under (74b), the right-hand side of (102) is smaller than one, which completes the proof.

V-D Proof of (91)

In our proof, we build an upper bound of $\|\mathbf{x}_{\Gamma^{k+1}}\|_{2}^{2}$ and a lower bound of $\|\mathbf{x}_{\Gamma^{k}}\|_{2}^{2}$ and then show that the former is smaller than the latter under (74c).

$\bullet$ Upper bound of $\|\mathbf{x}_{\Gamma^{k+1}}\|_{2}^{2}$ :

By taking similar steps to the proof of (97), we can show that

\displaystyle\|\mathbf{x}_{\Gamma^{k+1}}\|_{2}^{2}

\displaystyle\leq\frac{\|\mathbf{r}^{k+1}\|_{2}^{2}}{1-\delta}.

(103)

$\bullet$ Lower bound of $\|\mathbf{x}_{\Gamma^{k}}\|_{2}^{2}$ :

From Lemma 8, we have

$\displaystyle\\|\mathbf{r}^{k}\\|_{2}^{2}-\\|\mathbf{r}^{k+1}\\|_{2}^{2}$	$\displaystyle\hskip 1.05273pt\geq\frac{(1-\delta_{\|S^{k}\|+1}^{2})(1-\delta_{\|\Gamma_{1}^{k}\cup S^{k}\|})}{\lceil\frac{\|\Gamma_{1}^{k}\|}{L}\rceil(1+\delta_{L})}\left(\\|\mathbf{r}^{k}\\|_{2}^{2}-\\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{1}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{1}^{k}}\\|_{2}^{2}\right)$
	$\displaystyle\overset{(a)}{=}\frac{(1-\delta_{\|S^{k}\|+1}^{2})(1-\delta_{\|\Gamma_{1}^{k}\cup S^{k}\|})}{1+\delta_{L}}\left(\\|\mathbf{r}^{k}\\|_{2}^{2}-\\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{1}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{1}^{k}}\\|_{2}^{2}\right)$
	$\displaystyle\overset{(b)}{\geq}(1-\delta)^{2}\left(\\|\mathbf{r}^{k}\\|_{2}^{2}-\\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{1}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{1}^{k}}\\|_{2}^{2}\right),$	(104)

where (a) is because $\frac{|\Gamma_{1}^{k}|}{L}\leq 1$ and (b) follows from Lemma 3.¹⁶¹⁶16Again, if $\|\mathbf{r}^{k}\|_{2}^{2}-\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{1}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{1}^{k}}\|_{2}^{2}<0$ , then (104) holds since $\|\mathbf{r}^{k}\|_{2}^{2}-\|\mathbf{r}^{k+1}\|_{2}^{2}\geq 0$ by the orthogonal projection at each iteration. After re-arranging terms, we have

\displaystyle\|\mathbf{r}^{k+1}\|_{2}^{2}

\displaystyle\hskip 1.05273pt\leq\left(1-(1-\delta)^{2}\right)\|\mathbf{r}^{k}\|_{2}^{2}+(1-\delta)^{2}\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{1}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{1}^{k}}\|_{2}^{2}.

(105)

Note that

	$\displaystyle\\|\mathbf{r}^{k}\\|_{2}^{2}$	$\displaystyle\hskip 1.05273pt=\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{\Gamma^{k}}\mathbf{x}_{\Gamma^{k}}\\|_{2}^{2}\leq\\|\mathbf{A}_{\Gamma^{k}}\mathbf{x}_{\Gamma^{k}}\\|_{2}^{2}$
		$\displaystyle\overset{(a)}{\leq}(1+\delta_{\|\Gamma^{k}\|})\\|\mathbf{x}_{\Gamma^{k}}\\|_{2}^{2}\overset{(b)}{\leq}(1+\delta)\\|\mathbf{x}_{\Gamma^{k}}\\|_{2}^{2},$		(106)

where (a) and (b) follow from the RIP and Lemma 3, respectively. Also, note that

$\displaystyle\\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{1}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{1}^{k}}\\|_{2}^{2}$	$\displaystyle\hskip 1.05273pt\leq(1+\delta_{\|\Gamma^{k}\setminus\Gamma_{1}^{k}\|})\\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{1}^{k}}\\|_{2}^{2}$
	$\displaystyle\overset{(a)}{\leq}(1+\delta)\\|\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{1}^{k}}\\|_{2}^{2}$
	$\displaystyle\overset{(b)}{\leq}\frac{1+\delta}{\sigma}\\|\mathbf{x}_{\Gamma^{k}}\\|_{2}^{2}$
	$\displaystyle\overset{(c)}{=}2(1+\delta)\exp\left(-c(1-\delta)^{2}\right)\\|\mathbf{x}_{\Gamma^{k}}\\|_{2}^{2},$	(107)

where (a), (b), and (c) follow from Lemma 3, (81), and (79), respectively. By combining (105)-(107), we have

\displaystyle\|\mathbf{x}_{\Gamma^{k}}\|_{2}^{2}

\displaystyle\geq\frac{\|\mathbf{r}^{k+1}\|_{2}^{2}}{(1+\delta)\left(1-(1-\delta)^{2}+2(1-\delta)^{2}\exp(-c(1-\delta)^{2})\right)}.

(108)

$\bullet$ When is $\|\mathbf{x}_{\Gamma^{k+1}}\|_{2}^{2}<\|\mathbf{x}_{\Gamma^{k}}\|_{2}^{2}$ ?

From (103) and (108), we have

\displaystyle\frac{\|\mathbf{x}_{\Gamma^{k+1}}\|_{2}^{2}}{\|\mathbf{x}_{\Gamma^{k}}\|_{2}^{2}}

\displaystyle\leq\frac{(1+\delta)\left(1-(1-\delta)^{2}+2(1-\delta)^{2}\exp(-c(1-\delta)^{2})\right)}{1-\delta}.

(109)

One can easily show that under (74c), the right-hand side of (109) is smaller than one, which completes the proof.

VI Simulation Results

In this section, we study the performance of the proposed SSMP algorithm through empirical simulations. In our simulations, we use an $m\times n$ sampling matrix $\mathbf{A}$ ( $m=64,\ n=512$ ) whose entries are drawn i.i.d. from a Gaussian distribution $\mathcal{N}\left(0,\frac{1}{m}\right)$ . For each $K$ , we generate a row $K$ -sparse matrix $\mathbf{X}\in\mathbb{R}^{n\times r}$ whose support is uniformly chosen at random. Nonzero entries of $\mathbf{X}$ are drawn i.i.d. from a standard Gaussian distribution or binary ( $\pm 1$ ) random distribution. We refer to these two types of signals as the Gaussian signal and the 2-ary pulse amplitude modulation (2-PAM) signal, respectively. In our simulations, the following joint sparse recovery algorithms are considered:

1.

SOMP [13]
2.

M-ORMP [11]
3.

$\ell_{1}/\ell_{2}$ -norm minimization technique [14]
4.

CS-MUSIC [15]
5.

RA-OMP [16]
6.

RA-ORMP [16]
7.

SSMP ( $L=2$ )
8.

SSMP ( $L=3$ )

VI-A Noiseless Scenario

In this subsection, we study the empirical performance of SSMP in the noiseless scenario. In this case, the observation matrix $\mathbf{Y}\in\mathbb{R}^{m\times r}$ follows the system model in (2). We perform $5,000$ independent trials for each point of the algorithm and compute the exact reconstruction ratio (ERR) defined as [22, 23]

\operatornamewithlimits{ERR}=\frac{\text{number of exact reconstructions}}{\text{number of total trials}}.

By comparing the critical sparsity (the maximum sparsity level at which exact reconstruction is ensured [22]), recovery accuracy of different algorithms can be evaluated.

First, we study the recovery performance in the full row rank case. In Fig. 4, we plot the ERR performance as a function of sparsity $K$ . We can observe that the proposed SSMP algorithm shows a perfect performance (i.e., ERR=1) regardless of the sparsity and the type of a row sparse signal. We also observe that RA-ORMP, which can be viewed as a special case of SSMP for $L=1$ (see Table II), achieves an excellent performance. This is because the simulations are performed in the scenario where $\operatornamewithlimits{krank}(\mathbf{A})\geq K+1$ , and SSMP guarantees exact reconstruction in this scenario if $r=K$ (see Theorem 1). On the other hand, conventional algorithms such as SOMP, M-ORMP, $\ell_{1}/\ell_{2}$ -norm minimization, and RA-OMP are imperfect when $K$ is large, which agrees with the result that these algorithms do not uniformly guarantee exact recovery under $\operatornamewithlimits{krank}(\mathbf{A})\geq K+1$ (see Table I).

Next, we investigate the recovery performance of SSMP in the rank deficient case ( $r<K$ ). In Fig. 5, we plot the ERR performance as a function of $K$ . In general, we observe that the ERR performance improves with the number $r$ of measurement vectors. We also observe that SSMP outperforms other joint sparse recovery algorithms in terms of the critical sparsity for both Gaussian and 2-PAM signals. For example, when $r=5$ and the desired signal is Gaussian, the critical sparsity of SSMP is $1.5$ times higher than those obtained by conventional recovery algorithms (see Fig. 5(a)). In Fig. 6, we plot the ERR performance as a function of the number $m$ of measurements. In these simulations, we set the sparsity level to $K=50$ , for which none of recovery algorithms uniformly guarantees exact recovery. Overall, we observe that ERR improves with $m$ and the number of measurements required for exact reconstruction decreases with $r$ . From Fig. 6(a), we observe that SSMP recovers any row $K$ -sparse signal accurately when $m\geq 105$ , while other algorithms require more than $145$ measurements. Interestingly, from Fig. 6(b), we observe that SSMP ensures perfect recovery with $95$ measurements, which meets the fundamental minimum number of measurements ( $m=2K-r+1=95$ ) required for exact joint sparse recovery [16].

Finally, we study the empirical performance of SSMP in the scenario where the desired signal is approximately row sparse. Recall that $\mathbf{X}$ is approximately row $K$ -sparse if $\|\mathbf{X}-\mathbf{X}_{K}\|_{F}\leq\rho\|\mathbf{X}\|_{F}$ for some small $\rho$ . For an approximately row $K$ -sparse signal, we define the support as the index set of the $K$ rows with largest $\ell_{2}$ -norms. For each $K$ , we generate an approximately row $K$ -sparse matrix $\mathbf{X}\in\mathbb{R}^{n\times r}$ whose support $S$ is uniformly chosen at random. The elements of $\mathbf{X}^{S}$ and $\mathbf{X}^{\Omega\setminus S}$ are drawn i.i.d. from Gaussian distributions $\mathcal{N}(0,1)$ and $\mathcal{N}(0,\sigma^{2})$ , respectively. In our simulations, $\sigma^{2}$ is set to $\sigma^{2}=\frac{\rho^{2}K}{(n-K)(1-\rho^{2})}$ so that

\frac{\mathbb{E}[\|\mathbf{X}-\mathbf{X}_{K}\|_{F}^{2}]}{\mathbb{E}[\|\mathbf{X}\|_{F}^{2}]}=\frac{(n-K)r\sigma^{2}}{(n-K)r\sigma^{2}+Kr}=\rho^{2}.

As a performance measure for this scenario, we employ the exact support recovery ratio (ESRR):

\operatornamewithlimits{ESRR}=\frac{\text{number of exact support recovery}}{\text{number of total trials}}.

In Fig. 7, we plot the ESRR performance as a function of $K$ when $r=7$ . In general, we observe that the ESRR performance degrades with $\rho$ . In particular, one can see that the ESRR performance of CS-MUSIC degrades severely with $\rho$ . For example, if $\rho=0$ , then the critical sparsity of CS-MUSIC is 20 (see Fig. 5(c)). However, if $\rho=0.01$ , then the ESRR of CS-MUSIC is less than one even when $K=10$ (see Fig. 7(a)). On the other hand, critical sparsities of other algorithms remain the same when $\rho$ increases from $0$ to $0.05$ . We also observe that the proposed SSMP algorithm is very competitive in recovering approximately row sparse signals. Specifically, the critical sparsity of SSMP is $38$ , which is about $1.5$ times higher than critical sparsities of other approaches.

VI-B Noisy Scenario

In this subsection, we study the empirical performance of SSMP in the noisy scenario. In this case, the observation matrix $\mathbf{Y}$ follows the system model in (44), and we employ the MSE as a performance measure where

\operatornamewithlimits{MSE}=\frac{1}{nr}\|\mathbf{X}-\widehat{\mathbf{X}}\|_{F}^{2}.

For each simulation point of the algorithm, we perform 10,000 independent trials and average the MSE. In our simulations, we set the stopping threshold $\epsilon$ of SSMP to $\epsilon=\|\mathbf{W}\|_{F}/\sigma_{\max}(\mathbf{Y})$ as in Theorem 2.

In Fig. 8, we plot the MSE performance as a function of SNR (in dB) which is defined as $\operatornamewithlimits{SNR}=10\log_{10}\frac{\|\mathbf{A}\mathbf{X}\|_{F}^{2}}{\|\mathbf{W}\|_{F}^{2}}$ . In these simulations, the entries of $\mathbf{W}$ are generated at random from Gaussian distribution $\mathcal{N}(0,\frac{K}{m}10^{-\frac{\text{SNR}}{10}})$ .¹⁷¹⁷17One can easily check that $\mathbb{E}[(\mathbf{A}\mathbf{X})_{ij}^{2}]=\frac{K}{m}$ , since each component of $\mathbf{A}$ is generated independently from $\mathcal{N}(0,\frac{1}{m})$ and $\mathbf{X}$ is a row $K$ -sparse matrix whose nonzero entries are drawn independently from $\mathcal{N}(0,1)$ . Then, from the definition of SNR, we have $\mathbb{E}[\mathbf{W}_{ij}^{2}]=\frac{K}{m}10^{-\frac{\operatornamewithlimits{SNR}}{10}}$ . As a benchmark, we also plot the MSE performance of the Oracle-LS estimator, the best possible estimator using prior knowledge on the support. From the figures, we observe that the MSE performance of SSMP improves linearly with SNR, but the performance improvement of conventional algorithms diminishes with SNR. In particular, SSMP performs close to the Oracle-LS estimator in the high SNR regime ( $\operatornamewithlimits{SNR}\geq 20$ dB).

VII Conclusion and Discussion

In this paper, we proposed a new joint sparse recovery algorithm called signal space matching pursuit (SSMP) that greatly improves the reconstruction accuracy over conventional techniques. The key idea behind SSMP is to sequentially investigate the support of a row sparse matrix to minimize the subspace distance to the residual space. Our theoretical analysis indicates that under a mild condition on the sampling matrix, SSMP exactly reconstructs any row $K$ -sparse matrix $\mathbf{X}$ of rank $K$ using $m=K+1$ measurements, which meets the fundamental minimum number of measurements to ensure perfect recovery of $\mathbf{X}$ . We also showed that SSMP guarantees exact reconstruction in at most $K-r+\lceil\frac{r}{L}\rceil$ iterations, provided that $\mathbf{A}$ satisfies the RIP of order $L(K-r)+r+1$ with

\displaystyle\delta_{L(K-r)+r+1}

\displaystyle<\max\left\{\frac{\sqrt{r}}{\sqrt{K+\frac{r}{4}}+\sqrt{\frac{r}{4}}},\frac{\sqrt{L}}{\sqrt{K}+1.15\sqrt{L}}\right\}.

(110)

This implies that the requirement on the RIP constant becomes less restrictive when the number $r$ of measurement vectors increases. Such behavior seems to be natural but has not been reported for most of conventional methods. The result in (110) also implies that if $r$ is on the order of $K$ , then SSMP ensures perfect recovery with overwhelming probability as long as the number of random measurements scales linearly with $K\log\frac{n}{K}$ . We further showed that if $r=1$ , then by running $\max\{K,\lfloor\frac{8K}{L}\rfloor\}$ iterations, the guarantee (110) can be significantly improved to $\delta_{\lfloor 7.8K\rfloor}\leq 0.155$ . This implies that although $r$ is not on the order of $K$ , SSMP guarantees exact reconstruction with $\mathcal{O}(K\log\frac{n}{K})$ random measurements by running slightly more than $K$ iterations. Moreover, we showed that under a proper RIP condition, the reconstruction error of SSMP is upper bounded by a constant multiple of the noise power, which demonstrates the robustness of SSMP to measurement noise. Finally, from numerical experiments, we confirmed that SSMP outperforms conventional joint sparse recovery algorithms both in noiseless and noisy scenarios.

Finally, we would like to mention some future directions. Firstly, in this work, the number $L$ of indices chosen in each iteration is fixed. It would be more flexible and useful if this constraint is relaxed to the variable length. To achieve this goal, a deliberately designed thresholding operation is needed. Secondly, in analyzing a condition under which SSMP chooses $K-r$ support elements in the first $K-r$ iterations, we only considered the scenario where SSMP picks at least one support element in each iteration. One can see that although SSMP fails to choose a support element in some iterations, it can still choose $K-r$ support elements by identifying multiple support elements at a time. It would be an interesting future work to improve the proposed guarantee (110) for this more realistic scenario. Lastly, our result in Theorem 1 implies that if $r$ is not on the order of $K$ , then SSMP running $K$ iterations requires that the number $m$ of random measurements should scale with $K^{2}\log\frac{n}{K}$ , which is worse than the conventional linear scaling of $m$ ( $m=\mathcal{O}(K\log\frac{n}{K})$ ). When $r=1$ , we can overcome this limitation by running more than $K$ iterations. Extension of our analysis for general $r$ to obtain an improved performance guarantee of SSMP is also an interesting research direction.

Appendix A Proof of Proposition 1

Proof.

From (10), the set $I^{k+1}$ of $L$ indices chosen in the ( $k+1$ )-th iteration of SSMP satisfies

\displaystyle I^{k+1}

\displaystyle=\underset{I:I\subset\Omega\setminus S^{k},|I|=L}{\arg\min}\sum_{i\in I}\operatornamewithlimits{dist}\left(\mathcal{R}(\mathbf{R}^{k}),\mathbf{P}_{S^{k}\cup\{i\}}\mathcal{R}(\mathbf{R}^{k})\right).

(111)

From the definition of subspace distance, one can show that (see justifications in Appendix B)

\displaystyle\operatornamewithlimits{dist}\left(\mathcal{R}(\mathbf{R}^{k}),\mathbf{P}_{S^{k}\cup\{i\}}\mathcal{R}(\mathbf{R}^{k})\right)

\displaystyle=\|\mathbf{P}_{S^{k}\cup\{i\}}^{\perp}\mathbf{U}\|_{F},

(112)

where $\mathbf{U}=[\mathbf{u}_{1},\ldots,\mathbf{u}_{d}]$ is an orthonormal basis of $\mathcal{R}(\mathbf{R}^{k})$ . By combining (111) and (112), we have

	$\displaystyle I^{k+1}$	$\displaystyle=\underset{I:I\subset\Omega\setminus S^{k},\|I\|=L}{\arg\min}\sum_{i\in I}\\|\mathbf{P}_{S^{k}\cup\{i\}}^{\perp}\mathbf{U}\\|_{F}$
		$\displaystyle=\underset{I:I\subset\Omega\setminus S^{k},\|I\|=L}{\arg\max}\sum_{i\in I}\\|\mathbf{P}_{S^{k}\cup\{i\}}\mathbf{U}\\|_{F}.$		(113)

Also, note that $\|\mathbf{P}_{S^{k}\cup\{i\}}\mathbf{u}_{l}\|_{2}^{2}$ can be decomposed as [36, eq. (12)]

\displaystyle\|\mathbf{P}_{S^{k}\cup\{i\}}\mathbf{u}_{l}\|_{2}^{2}

\displaystyle=\|\mathbf{P}_{S^{k}}\mathbf{u}_{l}\|_{2}^{2}+\left|\left\langle\mathbf{u}_{l},\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}}{\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\|_{2}}\right\rangle\right|^{2}.

Using this together with (113), we obtain

	$\displaystyle I^{k+1}$	$\displaystyle=\underset{I:I\subset\Omega\setminus S^{k},\|I\|=L}{\arg\max}\sum_{i\in I}\sqrt{\sum_{l=1}^{d}\left\|\left\langle\mathbf{u}_{l},\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}}{\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\\|_{2}}\right\rangle\right\|^{2}}$
		$\displaystyle=\underset{I:I\subset\Omega\setminus S^{k},\|I\|=L}{\arg\max}\sum_{i\in I}\left\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}}{\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\\|_{2}}\right\\|_{2},$

which is the desired result. ∎

Appendix B Proof of (112)

Proof.

For notational convenience, let $\mathcal{V}=\mathcal{R}(\mathbf{A}_{S^{k}\cup\{i\}})$ . Then, our job is to show that

\displaystyle\operatornamewithlimits{dist}\left(\mathcal{R}(\mathbf{R}^{k}),\mathbf{P}_{\mathcal{V}}\mathcal{R}(\mathbf{R}^{k})\right)

\displaystyle=\|\mathbf{P}_{\mathcal{V}}^{\perp}\mathbf{U}\|_{F},

(114)

where $\mathbf{U}=[\mathbf{u}_{1},\ldots,\mathbf{u}_{d}]$ is an orthonormal basis of $\mathcal{R}(\mathbf{R}^{k})$ . Let $\{\mathbf{v}_{1},\ldots,\mathbf{v}_{e}\}$ be an orthonormal basis of $\mathbf{P}_{\mathcal{V}}\mathcal{R}(\mathbf{R}^{k})$ ( $e\leq d$ ). Then we have

	$\displaystyle\operatornamewithlimits{dist}\left(\mathcal{R}(\mathbf{R}^{k}),\mathbf{P}_{\mathcal{V}}\mathcal{R}(\mathbf{R}^{k})\right)$	$\displaystyle\overset{(a)}{=}\sqrt{d-\underset{i=1}{\overset{d}{\sum}}\underset{j=1}{\overset{e}{\sum}}\|\langle\mathbf{u}_{i},\mathbf{v}_{j}\rangle\|^{2}}$
		$\displaystyle=\sqrt{d-\underset{i=1}{\overset{d}{\sum}}\\|\mathbf{P}_{\mathbf{P}_{\mathcal{V}}\mathcal{R}(\mathbf{R}^{k})}\mathbf{u}_{i}\\|_{2}^{2}},$		(115)

where (a) is from the definition of subspace distance (see (6)). Also, note that

\mathbf{P}_{\mathbf{P}_{\mathcal{V}}\mathcal{R}(\mathbf{R}^{k})}\mathbf{u}_{i}=\mathbf{P}_{\mathbf{P}_{\mathcal{V}}\mathcal{R}(\mathbf{R}^{k})}(\mathbf{P}_{\mathcal{V}}\mathbf{u}_{i}+\mathbf{P}_{\mathcal{V}}^{\perp}\mathbf{u}_{i})=\mathbf{P}_{\mathcal{V}}\mathbf{u}_{i}.

Using this together with (115), we have (114). ∎

Appendix C Proof of Lemma 5

Proof.

Note that

	$\displaystyle\max_{i\in S\setminus S^{k}}\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\\|_{2}^{2}$	$\displaystyle\geq\frac{1}{\|S\setminus S^{k}\|}\sum_{i\in S\setminus S^{k}}\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\\|_{2}^{2}$
		$\displaystyle=\frac{1}{\|S\setminus S^{k}\|}\\|\mathbf{B}_{S\setminus S^{k}}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}^{\perp}\mathbf{B}_{S\setminus S^{k}}\\|_{F}^{2}.$		(116)

Let $t=\operatornamewithlimits{rank}(\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}^{\perp}\mathbf{B}_{S\setminus S^{k}})$ . Then, by Eckart-Young theorem [46], we have

	$\displaystyle\\|\mathbf{B}_{S\setminus S^{k}}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}^{\perp}\mathbf{B}_{S\setminus S^{k}}\\|_{F}^{2}$	$\displaystyle\geq\sum_{i=t+1}^{\|S\setminus S^{k}\|}\sigma_{i}^{2}(\mathbf{B}_{S\setminus S^{k}})$
		$\displaystyle=\sum_{i=1}^{\|S\setminus S^{k}\|-t}\sigma_{\|S\setminus S^{k}\|+1-i}^{2}(\mathbf{B}_{S\setminus S^{k}}).$		(117)

By combining (116) and (117), we obtain

\displaystyle\max_{i\in S\setminus S^{k}}\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\|_{2}^{2}

\displaystyle\geq\frac{1}{|S\setminus S^{k}|}\sum_{i=1}^{|S\setminus S^{k}|-t}\sigma_{|S\setminus S^{k}|+1-i}^{2}(\mathbf{B}_{S\setminus S^{k}}).

(118)

We now take a look at $t=\operatornamewithlimits{rank}(\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}^{\perp}\mathbf{B}_{S\setminus S^{k}})$ . Since $\mathbf{R}^{k}=\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}\mathbf{X}^{S\setminus S^{k}}=\mathbf{B}_{S\setminus S^{k}}\mathbf{D}\mathbf{X}^{S\setminus S^{k}}$ where $\mathbf{D}=\text{diag}\{\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{i}\|_{2}:i\in S\setminus S^{k}\}$ , we have $\mathcal{R}(\mathbf{R}^{k})+\mathcal{R}(\mathbf{B}_{S\setminus S^{k}})=\mathcal{R}(\mathbf{B}_{S\setminus S^{k}})$ and thus

\displaystyle\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}+\mathbf{P}_{\mathcal{R}(\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}^{\perp}\mathbf{B}_{S\setminus S^{k}})}

\displaystyle=\mathbf{P}_{\mathcal{R}(\mathbf{B}_{S\setminus S^{k}})}.

As a result, we have

\displaystyle t

\displaystyle=\operatornamewithlimits{rank}(\mathbf{B}_{S\setminus S^{k}})-\operatornamewithlimits{rank}(\mathbf{R}^{k}).

(119)

Note that since $\mathbf{A}_{S\cup S^{k}}$ has full column rank, the projected matrix $\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}$ also has full column rank by Lemma 4, and so does its $\ell_{2}$ -normalized counterpart $\mathbf{B}_{S\setminus S^{k}}$ (i.e. $\operatornamewithlimits{rank}(\mathbf{B}_{S\setminus S^{k}})=|S\setminus S^{k}|$ ). Also, since $\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}$ has full column rank, $\operatornamewithlimits{rank}(\mathbf{R}^{k})=\operatornamewithlimits{rank}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}\mathbf{X}^{S\setminus S^{k}})=\operatornamewithlimits{rank}(\mathbf{X}^{S\setminus S^{k}})=d$ . Combining these together with (118) and (119), we obtain the desired result (13). ∎

Appendix D Proofs of (46) and (47)

D-A Proof of (46)

From the triangle inequality, we have

\displaystyle\|\mathbf{Z}^{k}-\mathbf{X}\|_{F}

\displaystyle\leq\|\mathbf{Z}^{k}-\mathbf{X}^{k}\|_{F}+\|\mathbf{X}^{k}-\mathbf{X}\|_{F}.

(120)

Note that since $\widehat{S}$ is the set of indices corresponding to the $K$ rows of $\mathbf{X}^{k}$ with largest $\ell_{2}$ -norms (see Algorithm 1) and $\mathbf{Z}^{k}$ is the row $K$ -sparse matrix satisfying $(\mathbf{Z}^{k})^{\widehat{S}}=(\mathbf{X}^{k})^{\widehat{S}}$ and $(\mathbf{Z}^{k})^{\Omega\setminus\widehat{S}}=\mathbf{0}_{|\Omega\setminus\widehat{S}|\times r}$ , $\mathbf{Z}^{k}$ is the best row $K$ -sparse approximation of $\mathbf{X}^{k}$ , i.e.,

\displaystyle\mathbf{Z}^{k}

\displaystyle=\underset{\mathbf{U}:|\text{supp}(\mathbf{U})|\leq K}{\arg\min}\|\mathbf{X}^{k}-\mathbf{U}\|_{F}.

(121)

Then, we have $\|\mathbf{X}^{k}-\mathbf{Z}^{k}\|_{F}\leq\|\mathbf{X}^{k}-\mathbf{X}\|_{F}$ . Using this together with (120), we obtain

\displaystyle\|\mathbf{Z}^{k}-\mathbf{X}\|_{F}

\displaystyle\leq 2\|\mathbf{X}^{k}-\mathbf{X}\|_{F}\leq\frac{2\|\mathbf{A}(\mathbf{X}^{k}-\mathbf{X})\|_{F}}{\sqrt{1-\delta_{Lk+K}}}

(122)

where the last inequality is because $\mathbf{A}$ obeys the RIP of order $\max\{Lk+K,2K\}$ and $\mathbf{X}^{k}$ and $\mathbf{X}$ are row $Lk$ - and $K$ -sparse, respectively. Also, since $\mathbf{R}^{k}=\mathbf{Y}-\mathbf{A}\mathbf{X}^{k}=\mathbf{A}(\mathbf{X}-\mathbf{X}^{k})+\mathbf{W}$ , we have

\displaystyle\|\mathbf{A}(\mathbf{X}^{k}-\mathbf{X})\|_{F}

\displaystyle=\|\mathbf{R}^{k}-\mathbf{W}\|_{F}\leq\|\mathbf{R}^{k}\|_{F}+\|\mathbf{W}\|_{F},

(123)

where the last inequality follows from the triangle inequality. By combining (122) and (123), we obtain the desired result in (46).

D-B Proof of (47)

Note that $\mathbf{Z}^{k}-\mathbf{X}$ is row $2K$ -sparse since $\mathbf{Z}^{k}$ and $\mathbf{X}$ are row $K$ -sparse matrices. Then, we have

$\displaystyle\\|\mathbf{Z}^{k}-\mathbf{X}\\|_{F}$	$\displaystyle\geq\frac{\\|\mathbf{A}(\mathbf{Z}^{k}-\mathbf{X})\\|_{F}}{\sqrt{1+\delta_{2K}}}$
	$\displaystyle=\frac{\\|\mathbf{A}\mathbf{Z}^{k}-(\mathbf{Y}-\mathbf{W})\\|_{F}}{\sqrt{1+\delta_{2K}}}$
	$\displaystyle\geq\frac{\\|\mathbf{A}\mathbf{Z}^{k}-\mathbf{Y}\\|_{F}-\\|\mathbf{W}\\|_{F}}{\sqrt{1+\delta_{2K}}},$	(124)

where the last inequality follows from the triangle inequality. Also, by noting that $\widehat{\mathbf{X}}=\arg\min_{\mathbf{U}:\operatornamewithlimits{supp}(\mathbf{U})\subset\widehat{S}}\|\mathbf{Y}-\mathbf{A}\mathbf{U}\|_{F}$ (see Algorithm 1) and $\operatornamewithlimits{supp}(\mathbf{Z}^{k})\subset\widehat{S}$ , we have

$\displaystyle\\|\mathbf{A}\mathbf{Z}^{k}-\mathbf{Y}\\|_{F}$	$\displaystyle\geq\\|\mathbf{A}\widehat{\mathbf{X}}-\mathbf{Y}\\|_{F}$
	$\displaystyle=\\|\mathbf{A}(\widehat{\mathbf{X}}-\mathbf{X})-\mathbf{W}\\|_{F}$
	$\displaystyle\overset{(a)}{\geq}\\|\mathbf{A}(\widehat{\mathbf{X}}-\mathbf{X})\\|_{F}-\\|\mathbf{W}\\|_{F}$
	$\displaystyle\overset{(b)}{\geq}\sqrt{1-\delta_{2K}}\\|\widehat{\mathbf{X}}-\mathbf{X}\\|_{F}-\\|\mathbf{W}\\|_{F},$	(125)

where (a) follows from the triangle inequality and (b) is because $\widehat{\mathbf{X}}$ and $\mathbf{X}$ are row $K$ -sparse and thus $\widehat{\mathbf{X}}-\mathbf{X}$ is a row $2K$ -sparse matrix. Finally, by combining (124) and (125), we obtain the desired result in (47).

Appendix E Proof of (73)

Proof.

Since $\|\widehat{\mathbf{W}}\|_{F}=\|\mathbf{A}(\mathbf{X}-\mathbf{X}_{K})+\mathbf{W}\|_{F}\leq\|\mathbf{A}(\mathbf{X}-\mathbf{X}_{K})\|_{F}+\|\mathbf{W}\|_{F}$ by the triangle inequality, it suffices to show that

\displaystyle\|\mathbf{A}(\mathbf{X}-\mathbf{X}_{K})\|_{F}

\displaystyle\leq\sqrt{1+\delta_{K}}\left(\|\mathbf{X}-\mathbf{X}_{K}\|_{F}+\frac{\|\mathbf{X}-\mathbf{X}_{K}\|_{1,2}}{\sqrt{K}}\right).

(126)

Note that if $\mathbf{A}$ satisfies the RIP of order $K$ , then for any (not necessarily sparse) vector $\mathbf{z}$ , we have [21, Proposition 3.5]

\displaystyle\|\mathbf{A}\mathbf{z}\|_{2}

\displaystyle\leq\sqrt{1+\delta_{K}}\left(\|\mathbf{z}\|_{2}+\frac{\|\mathbf{z}\|_{1}}{\sqrt{K}}\right).

(127)

Then, we have

$\displaystyle\\|\mathbf{A}(\mathbf{X}-\mathbf{X}_{K})\\|_{F}^{2}$	$\displaystyle=\sum_{i=1}^{r}\\|\mathbf{A}(\mathbf{X}-\mathbf{X}_{K})_{i}\\|_{2}^{2}$
	$\displaystyle\overset{(a)}{\leq}(1+\delta_{K})\sum_{i=1}^{r}\left(\\|(\mathbf{X}-\mathbf{X}_{K})_{i}\\|_{2}+\frac{\\|(\mathbf{X}-\mathbf{X}_{K})_{i}\\|_{1}}{\sqrt{K}}\right)^{2}$
	$\displaystyle\overset{(b)}{\leq}(1+\delta_{K})\left(\sqrt{\sum_{i=1}^{r}\\|(\mathbf{X}-\mathbf{X}_{K})_{i}\\|_{2}^{2}}+\sqrt{\sum_{i=1}^{r}\frac{\\|(\mathbf{X}-\mathbf{X}_{K})_{i}\\|_{1}^{2}}{K}}\right)^{2}$
	$\displaystyle=(1+\delta_{K})\left(\\|\mathbf{X}-\mathbf{X}_{K}\\|_{F}+\sqrt{\sum_{i=1}^{r}\frac{\\|(\mathbf{X}-\mathbf{X}_{K})_{i}\\|_{1}^{2}}{K}}\right)^{2},$	(128)

where (a) and (b) are from (127) and Minkowski inequality, respectively. Let $\mathbf{\Phi}=\mathbf{X}-\mathbf{X}_{K}$ , then we have

$\displaystyle\\|\mathbf{X}-\mathbf{X}_{K}\\|_{1,2}^{2}$	$\displaystyle=\\|\mathbf{\Phi}\\|_{1,2}^{2}=\left(\sum_{j=1}^{n}\\|\bm{\phi}^{j}\\|_{2}\right)^{2}$
	$\displaystyle=\sum_{j=1}^{n}\\|\bm{\phi}^{j}\\|_{2}^{2}+\sum_{k\neq l}\\|\bm{\phi}^{k}\\|_{2}\\|\bm{\phi}^{l}\\|_{2}$
	$\displaystyle\overset{(a)}{\geq}\sum_{j=1}^{n}\\|\bm{\phi}^{j}\\|_{2}^{2}+\sum_{k\neq l}\langle\|\bm{\phi}^{k}\|,\|\bm{\phi}^{l}\|\rangle=\sum_{i=1}^{r}\left(\sum_{j=1}^{n}\|\phi_{ji}\|^{2}+\sum_{k\neq l}\|\phi_{ki}\|\|\phi_{li}\|\right)$
	$\displaystyle=\sum_{i=1}^{r}\\|\bm{\phi}_{i}\\|_{1}^{2}=\sum_{i=1}^{r}\\|(\mathbf{X}-\mathbf{X}_{K})_{i}\\|_{1}^{2},$	(129)

where (a) follows from Cauchy-Schwarz inequality. By combining (128) and (129), we obtain (126), which completes the proof. ∎

Appendix F Proof of Lemma 8

Proof.

For any integer $l$ such that $l\geq k$ , we have [45, eq. (C.10)]

	$\displaystyle\\|\mathbf{r}^{l}\\|_{2}^{2}-\\|\mathbf{r}^{l+1}\\|_{2}^{2}$	$\displaystyle\hskip 1.05273pt\geq\frac{1}{1+\delta_{L}}\\|\mathbf{A}_{I^{l+1}}^{\prime}\mathbf{r}^{l}\\|_{2}^{2}$
		$\displaystyle\overset{(a)}{=}\frac{1}{1+\delta_{L}}\sum_{i\in I^{l+1}}\|\langle\mathbf{P}_{S^{l}}^{\perp}\mathbf{a}_{i},\mathbf{r}^{l}\rangle\|^{2},$		(130)

where $I^{l+1}$ is the set of indices chosen in the $(l+1)$ -th iteration of the SSMP algorithm and (a) is because $\mathbf{P}_{S^{l}}^{\perp}=(\mathbf{P}_{S^{l}}^{\perp})^{\prime}=(\mathbf{P}_{S^{l}}^{\perp})^{2}$ . Also, note that

$\displaystyle\sum_{i\in I^{l+1}}\|\langle\mathbf{P}_{S^{l}}^{\perp}\mathbf{a}_{i},\mathbf{r}^{l}\rangle\|^{2}$	$\displaystyle\hskip 1.05273pt=\sum_{i\in I^{l+1}}\left\|\left\langle\frac{\mathbf{P}_{S^{l}}^{\perp}\mathbf{a}_{i}}{\\|\mathbf{P}_{S^{l}}^{\perp}\mathbf{a}_{i}\\|_{2}},\mathbf{r}^{l}\right\rangle\right\|^{2}\\|\mathbf{P}_{S^{l}}^{\perp}\mathbf{a}_{i}\\|_{2}^{2}$
	$\displaystyle\overset{(a)}{\geq}(1-\delta_{\|S^{l}\|+1}^{2})\sum_{i\in I^{l+1}}\left\|\left\langle\frac{\mathbf{P}_{S^{l}}^{\perp}\mathbf{a}_{i}}{\\|\mathbf{P}_{S^{l}}^{\perp}\mathbf{a}_{i}\\|_{2}},\mathbf{r}^{l}\right\rangle\right\|^{2}$
	$\displaystyle\overset{(b)}{=}(1-\delta_{\|S^{l}\|+1}^{2})\max_{I:I\subset\Omega\setminus S^{l},\|I\|=L}\sum_{i\in I}\left\|\left\langle\frac{\mathbf{P}_{S^{l}}^{\perp}\mathbf{a}_{i}}{\\|\mathbf{P}_{S^{l}}^{\perp}\mathbf{a}_{i}\\|_{2}},\mathbf{r}^{l}\right\rangle\right\|^{2}$
	$\displaystyle\overset{(c)}{\geq}(1-\delta_{\|S^{l}\|+1}^{2})\max_{I:I\subset\Omega\setminus S^{l},\|I\|=L}\\|\mathbf{A}_{I}^{\prime}\mathbf{r}^{l}\\|_{2}^{2},$	(131)

where (a) and (b) follow from Lemma 6 and (38), respectively, and (c) is because $\|\mathbf{P}_{S^{l}}^{\perp}\mathbf{a}_{i}\|_{2}\leq\|\mathbf{a}_{i}\|_{2}=1$ for each of $i\in\Omega\setminus S^{l}$ . Furthermore, it has been shown in [45, eq. (C.5)] that

\displaystyle\max_{I:I\subset\Omega\setminus S^{l},|I|=L}\|\mathbf{A}_{I}^{\prime}\mathbf{r}^{l}\|_{2}^{2}

\displaystyle\geq\frac{1-\delta_{|\Gamma_{\tau}^{k}\cup S^{l}|}}{\lceil\frac{|\Gamma_{\tau}^{k}|}{L}\rceil}\left(\|\mathbf{r}^{l}\|_{2}^{2}-\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\|_{2}^{2}\right).

(132)

By combining (130)-(132), we obtain (83), which is the desired result. ∎

Appendix G Proof of Lemma 9

Proof.

For each of $l^{\prime}\in\{l,\ldots,l+\Delta l-1\}$ , we have

$\displaystyle\\|\mathbf{r}^{l^{\prime}}\\|_{2}^{2}-\\|\mathbf{r}^{l^{\prime}+1}\\|_{2}^{2}$
	$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \overset{(a)}{\geq}\frac{(1-\delta_{\|S^{l^{\prime}}\|+1}^{2})(1-\delta_{\|\Gamma_{\tau}^{k}\cup S^{l^{\prime}}\|})}{\lceil\frac{\|\Gamma_{\tau}^{k}\|}{L}\rceil(1+\delta_{L})}\left(\\|\mathbf{r}^{l^{\prime}}\\|_{2}^{2}-\\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\\|_{2}^{2}\right)$
	$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \overset{(b)}{\geq}\left(1-\exp\left(-\frac{(1-\delta_{\|S^{l^{\prime}}\|+1}^{2})(1-\delta_{\|\Gamma_{\tau}^{k}\cup S^{l^{\prime}}\|})}{\lceil\frac{\|\Gamma_{\tau}^{k}\|}{L}\rceil(1+\delta_{L})}\right)\right)\left(\\|\mathbf{r}^{l^{\prime}}\\|_{2}^{2}-\\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\\|_{2}^{2}\right)$
	$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \overset{(c)}{\geq}\left(1-\exp\left(-\frac{(1-\delta_{\|S^{l+\Delta l-1}\|+1}^{2})(1-\delta_{\|\Gamma_{\tau}^{k}\cup S^{l+\Delta l-1}\|})}{\lceil\frac{\|\Gamma_{\tau}^{k}\|}{L}\rceil(1+\delta_{L})}\right)\right)\left(\\|\mathbf{r}^{l^{\prime}}\\|_{2}^{2}-\\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\\|_{2}^{2}\right),$	(133)

where (a) is from Lemma 8, (b) is because $t>1-e^{-t}$ for $t>0$ , and (c) is from Lemma 3.¹⁸¹⁸18If $\|\mathbf{r}^{l^{\prime}}\|_{2}^{2}-\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma^{k}_{\tau}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma^{k}_{\tau}}\|_{2}^{2}<0$ , then (133) clearly holds true since $\|\mathbf{r}^{l^{\prime}}\|_{2}^{2}-\|\mathbf{r}^{l^{\prime}+1}\|_{2}^{2}\geq 0$ due to the orthogonal projection at each iteration of SSMP. By subtracting both sides of (133) by $\|\mathbf{r}^{l^{\prime}}\|_{2}^{2}-\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\|_{2}^{2}$ , we obtain

	$\displaystyle\\|\mathbf{r}^{l^{\prime}+1}\\|_{2}^{2}-\\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\\|_{2}^{2}$
		$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \hskip 1.05273pt\leq\exp\left(-\frac{(1-\delta_{\|S^{l+\Delta l-1}\|+1}^{2})(1-\delta_{\|\Gamma_{\tau}^{k}\cup S^{l+\Delta l-1}\|})}{\lceil\frac{\|\Gamma_{\tau}^{k}\|}{L}\rceil(1+\delta_{L})}\right)\left(\\|\mathbf{r}^{l^{\prime}}\\|_{2}^{2}-\\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\\|_{2}^{2}\right).$		(134)

By plugging $l^{\prime}=l,\ldots,l+\Delta l-1$ into (134), we have


	$\displaystyle\\|\mathbf{r}^{l+1}\\|_{2}^{2}-\\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\\|_{2}^{2}$
	$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leq\exp\left(-\frac{(1-\delta_{\|S^{l+\Delta l-1}\|+1}^{2})(1-\delta_{\|\Gamma_{\tau}^{k}\cup S^{l+\Delta l-1}\|})}{\lceil\frac{\|\Gamma_{\tau}^{k}\|}{L}\rceil(1+\delta_{L})}\right)\left(\\|\mathbf{r}^{l}\\|_{2}^{2}-\\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\\|_{2}^{2}\right),$		(135a)
	$\displaystyle\\|\mathbf{r}^{l+2}\\|_{2}^{2}-\\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\\|_{2}^{2}$
	$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leq\exp\left(-\frac{(1-\delta_{\|S^{l+\Delta l-1}\|+1}^{2})(1-\delta_{\|\Gamma_{\tau}^{k}\cup S^{l+\Delta l-1}\|})}{\lceil\frac{\|\Gamma_{\tau}^{k}\|}{L}\rceil(1+\delta_{L})}\right)\left(\\|\mathbf{r}^{l+1}\\|_{2}^{2}-\\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\\|_{2}^{2}\right),$		(135b)
	$\displaystyle\hskip 227.62204pt\vdots$
	$\displaystyle\\|\mathbf{r}^{l+\Delta l}\\|_{2}^{2}-\\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\\|_{2}^{2}$
	$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leq\exp\left(-\frac{(1-\delta_{\|S^{l+\Delta l-1}\|+1}^{2})(1-\delta_{\|\Gamma_{\tau}^{k}\cup S^{l+\Delta l-1}\|})}{\lceil\frac{\|\Gamma_{\tau}^{k}\|}{L}\rceil(1+\delta_{L})}\right)\left(\\|\mathbf{r}^{l+\Delta l-1}\\|_{2}^{2}-\\|\mathbf{A}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\mathbf{x}_{\Gamma^{k}\setminus\Gamma_{\tau}^{k}}\\|_{2}^{2}\right).$		(135c)

Finally, by combining (135a)-(135c), we obtain (84), which is the desired result. ∎

Appendix H Proof of (95)

Proof.

Note that

\displaystyle|S\cup S^{k_{N}}|

\displaystyle=|S^{k_{N}}|+|\Gamma^{k_{N}}|\leq|S^{k_{N}}|+|\Gamma^{k}|=Lk_{N}+\gamma.

(136)

Also, from (87), we have

	$\displaystyle k_{N}$	$\displaystyle\hskip 1.05273pt\leq k+c\sum_{\tau=1}^{N}2^{\tau-1}$
		$\displaystyle\hskip 1.05273pt<k+c2^{N}$
		$\displaystyle\overset{(a)}{<}k+\frac{4c\gamma}{L}\left(1-\frac{1}{e^{c}-1}\right)$

where (a) is from (82). By combining this together with (136) and by noting that $|S\cup S^{k_{N}}|$ is an integer, we obtain

\displaystyle|S\cup S^{k_{N}}|

\displaystyle\leq Lk+\left\lfloor\gamma\left(1+4c-\frac{4c}{e^{c}-1}\right)\right\rfloor,

which is the desired result. ∎

References

[1] J. Kim and B. Shim, “Multiple orthogonal least squares for joint sparse recovery,” in Proc. IEEE Int. Symp. Inform. Theory, Vail, CO, USA, Jun. 2018, pp. 61–65.
[2] E. J. Candès and M. B. Wakin, “An introduction to compressive sampling,” IEEE Signal Process. Mag., vol. 25, no. 2, pp. 21–30, Mar. 2008.
[3] D. L. Donoho, “Compressed sensing,” IEEE Trans. Inform. Theory, vol. 52, no. 4, pp. 1289–1306, Apr. 2006.
[4] V. Papyan, Y. Romano, J. Sulam, and M. Elad, “Theoretical foundations of deep learning via sparse representations: A multilayer sparse model and its connection to convolutional neural networks,” IEEE Signal Process. Mag., vol. 35, no. 4, pp. 72–89, Jul. 2018.
[5] J. W. Choi, B. Shim, Y. Ding, B. Rao, and D. I. Kim, “Compressed sensing for wireless communications: Useful tips and tricks,” IEEE Commun. Surveys Tuts., vol. 19, no. 3, pp. 1527–1550, 3rd Quart., 2017.
[6] J. W. Choi and B. Shim, “Statistical recovery of simultaneously sparse time-varying signals from multiple measurement vectors,” IEEE Trans. Signal Process., vol. 63, no. 22, pp. 6136–6148, Nov. 2015.
[7] P. Feng and Y. Bresler, “Spectrum-blind minimum-rate sampling and reconstruction of multiband signals,” in Proc. IEEE Int. Conf. Acoust., Speech, Signal Process., Atlanta, GA, USA, May 1996, pp. 1668–1691.
[8] M. Mishali and Y. C. Eldar, “From theory to practice: Sub-Nyquist sampling of sparse wideband analog signals,” IEEE J. Sel. Topics Signal. Process., vol. 4, no. 2, pp. 375–391, Apr. 2010.
[9] Y. Y. Lee, C. H. Wang, and Y. H. Huang, “A hybrid RF/baseband precoding processor based on parallel-index-selection matrix-inversion-bypass simultaneous orthogonal matching pursuit for millimeter wave MIMO systems,” IEEE Trans. Signal Process., vol. 63, no. 2, pp. 305–317, Jan. 2015.
[10] O. Lee, J. M. Kim, Y. Bresler, and J. C. Ye, “Compressive diffuse optical tomography: Noniterative exact reconstruction using joint sparsity,” IEEE Trans. Med. Imag., vol. 30, no. 5, pp. 1129–1142, May 2011.
[11] S. F. Cotter, B. D. Rao, K. Engan, and K. Kreutz-Delgado, “Sparse solutions to linear inverse problems with multiple measurement vectors,” IEEE Trans. Signal Process., vol. 53, no. 7, pp. 2477–2488, Jul. 2005.
[12] J. Chen and X. Huo, “Theoretical results on sparse representations of multiple-measurement vectors,” IEEE Trans. Signal Process., vol. 54, no. 12, pp. 4634–4643, Dec. 2006.
[13] J. A. Tropp, A. C. Gilbert, and M. J. Strauss, “Algorithms for simultaneous sparse approximation. part I: Greedy pursuit,” Signal Process., vol. 86, no. 3, pp. 572–588, Mar. 2006.
[14] J. A. Tropp, “Algorithms for simultaneous sparse approximation. part II: Convex relaxation,” Signal Process., vol. 86, no. 3, pp. 589–602, Mar. 2006.
[15] J. M. Kim, O. K. Lee, and J. C. Ye, “Compressive MUSIC: Revisiting the link between compressive sensing and array signal processing,” IEEE Trans. Inform. Theory, vol. 58, no. 1, pp. 278–301, Jan. 2012.
[16] M. E. Davies and Y. C. Eldar, “Rank awareness in joint sparse recovery,” IEEE Trans. Inform. Theory, vol. 58, no. 2, pp. 1135–1146, Feb. 2012.
[17] K. Lee, Y. Bresler, and M. Junge, “Subspace methods for joint sparse recovery,” IEEE Trans. Inform. Theory, vol. 58, no. 6, pp. 3613–3641, Jun. 2012.
[18] L. Wang, X. Wang, and J. Feng, “Subspace distance analysis with application to adaptive Bayesian algorithm for face recognition,” Pattern Recognition, vol. 39, no. 3, pp. 456–464, Mar. 2006.
[19] X. Sun, L. Wang, and J. Feng, “Further results on the subspace distance,” Pattern Recognition, vol. 40, no. 1, pp. 328–329, Jan. 2007.
[20] Y. C. Pati, R. Rezaiifar, and P. S. Krishnaprasad, “Orthogonal matching pursuit: Recursive function approximation with applications to wavelet decomposition,” in Proc. 27th Annu. Asilomar Conf. Signals, Syst., and Comput., Pacific Grove, CA, USA, Nov. 1993, pp. 40–44.
[21] D. Needell and J. A. Tropp, “CoSaMP: Iterative signal recovery from incomplete and inaccurate samples,” Appl. Comput. Harmon. Anal., vol. 26, no. 3, pp. 301–321, May 2009.
[22] W. Dai and O. Milenkovic, “Subspace pursuit for compressive sensing signal reconstruction,” IEEE Trans. Inform. Theory, vol. 55, no. 5, pp. 2230–2249, May 2009.
[23] J. Wang, S. Kwon, and B. Shim, “Generalized orthogonal matching pursuit,” IEEE Trans. Signal Process., vol. 60, no. 12, pp. 6202–6216, Dec. 2012.
[24] S. Kwon, J. Wang, and B. Shim, “Multipath matching pursuit,” IEEE Trans. Inform. Theory, vol. 60, no. 5, pp. 2986–3001, May 2014.
[25] E. J. Candès, “The restricted isometry property and its implications for compressed sensing,” Comptes Rendus Mathematique, vol. 346, no. 9-10, pp. 589–592, May 2008.
[26] M. A. Davenport and M. B. Wakin, “Analysis of orthogonal matching pursuit using the restricted isometry property,” IEEE Trans. Inform. Theory, vol. 56, no. 9, pp. 4395–4401, Sep. 2010.
[27] E. J. Candès and T. Tao, “Decoding by linear programming,” IEEE Trans. Inform. Theory, vol. 51, no. 12, pp. 4203–4215, Dec. 2005.
[28] B. Li, Y. Shen, Z. Wu, and J. Li, “Sufficient conditions for generalized orthogonal matching pursuit in noisy case,” Signal Process., vol. 108, pp. 111–123, Mar. 2015.
[29] J. Kim and B. Shim, “Nearly sharp restricted isometry condition of rank aware order recursive matching pursuit,” to appear in Proc. IEEE Int. Symp. Inform. Theory, Paris, France, Jul. 2019.
[30] J. D. Blanchard and M. E. Davies, “Recovery guarantees for rank aware pursuits,” IEEE Signal Process. Lett., vol. 19, no. 7, pp. 427–430, Jul. 2012.
[31] R. Baraniuk, M. Davenport, R. DeVore, and M. Wakin, “A simple proof of the restricted isometry property for random matrices,” Construct. Approx., vol. 28, no. 3, pp. 253–263, Dec. 2008.
[32] Y. Jin and B. D. Rao, “Support recovery of sparse signals in the presence of multiple measurement vectors,” IEEE Trans. Inform. Theory, vol. 59, no. 5, pp. 3139–3157, May 2013.
[33] S. Park, N. Y. Yu, and H. N. Lee, “An information-theoretic study for joint sparsity pattern recovery with different sensing matrices,” IEEE Trans. Inform. Theory, vol. 63, no. 9, pp. 5559–5571, Sept. 2017.
[34] J. Wang and P. Li, “Recovery of sparse signals using multiple orthogonal least squares,” IEEE Trans. Signal Process., vol. 65, no. 8, pp. 2049–2062, Apr. 2017.
[35] S. Chen, S. A. Billings, and W. Luo, “Orthogonal least squares methods and their application to non-linear system identification,” Int. J. Control, vol. 50, no. 5, pp. 1873–1896, 1989.
[36] L. Rebollo-Neira and D. Lowe, “Optimized orthogonal matching pursuit approach,” IEEE Signal Process. Lett., vol. 9, no. 4, pp. 137–140, Apr. 2002.
[37] J. Wen, J. Wang, and Q. Zhang, “Nearly optimal bounds for orthogonal least squares,” IEEE Trans. Signal Process., vol. 65, no. 20, pp. 5347–5356, Oct. 2017.
[38] T. T. Cai and L. Wang, “Orthogonal matching pursuit for sparse signal recovery with noise,” IEEE Trans. Inform. Theory, vol. 57, no. 7, pp. 4680–4688, Jul. 2011.
[39] P. A. Wedin, “On angles between subspaces of a finite dimensional inner product space,” Matrix Pencils, pp. 263–285, 1983.
[40] Q. Mo, “A sharp restricted isometry constant bound of orthogonal matching pursuit,” arXiv:1501.01708, 2015.
[41] T. Zhang, “Sparse recovery with orthogonal matching pursuit under RIP,” IEEE Trans. Inform. Theory, vol. 57, no. 9, pp. 6215–6221, Sep. 2011.
[42] S. Foucart and H. Rauhut, A mathematical introduction to compressive sensing, in Applied and Numerical Harmonic Analysis, Birkhauser, 2013.
[43] J. Wang and B. Shim, “Exact recovery of sparse signals using orthogonal matching pursuit: How many iterations do we need?,” IEEE Trans. Signal Process., vol. 64, no. 16, pp. 4194–4202, Aug. 2016.
[44] A. Cohen, W. Dahmen, and R. DeVore, “Orthogonal matching pursuit under the restricted isometry property,” Construct. Approx., vol. 45, no. 1, pp. 113–127, Feb. 2017.
[45] J. Wang, S. Kwon, P. Li, and B. Shim, “Recovery of sparse signals via generalized orthogonal matching pursuit: A new analysis,” IEEE Trans. Signal Process., vol. 64, no. 4, pp. 1076–1089, Feb. 2016.
[46] C. Eckart and G. Young, “The approximation of one matrix by another of lower rank,” Psychometrika, vol. 1, no. 3, pp. 211–218, 1936.

$\displaystyle q_{L}^{2}$	$\displaystyle=\left\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}}{\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}\\|_{2}}\right\\|_{2}^{2}$
	$\displaystyle\leq\left\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}})}\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}}{\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}\\|_{2}}\right\\|_{2}^{2}$
	$\displaystyle=1-\left\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}})}^{\perp}\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}}{\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}\\|_{2}}\right\\|_{2}^{2}.$	(22)

	$\displaystyle\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}})}^{\perp}\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}\\|_{2}^{2}$	$\displaystyle\overset{(a)}{\geq}(1-\delta_{\|S\cup S^{k}\|+1}^{2})\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}\\|_{2}^{2}$		(23)
		$\displaystyle\overset{(b)}{\geq}(1-\delta_{L(K-r)+r+1}^{2})\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{L}}\\|_{2}^{2},$		(24)

$\displaystyle q_{L}^{2}$	$\displaystyle\hskip 1.05273pt\leq\frac{1}{L}(q_{1}^{2}+\ldots+q_{L}^{2})$
	$\displaystyle\hskip 1.05273pt=\frac{1}{L}\underset{l=1}{\overset{L}{\sum}}\left\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\frac{\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{l}}}{\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{a}_{\psi_{l}}\\|_{2}}\right\\|_{2}^{2}$
	$\displaystyle\overset{(a)}{\leq}\frac{1}{L(1-\delta_{\|S^{k}\|+1}^{2})}\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{\Lambda}\\|_{F}^{2}$
	$\displaystyle\hskip 1.05273pt=\frac{1}{L(1-\delta_{\|S^{k}\|+1}^{2})}\\|\mathbf{A}_{\Lambda}^{\prime}\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}_{S\setminus S^{k}}\mathbf{U}\\|_{F}^{2}$
	$\displaystyle\overset{(b)}{\leq}\frac{\delta_{\|S\cup S^{k}\|+L}^{2}\\|\mathbf{U}\\|_{F}^{2}}{L(1-\delta_{\|S^{k}\|+1}^{2})},$	(26)

$\displaystyle\epsilon$	$\displaystyle\overset{(a)}{\geq}\\|\mathbf{P}_{S^{k}}^{\perp}\mathbf{Y}(\mathbf{Y}^{\prime}\mathbf{Y})^{-1/2}\\|_{F}$
	$\displaystyle\hskip 1.05273pt=\\|\mathbf{R}^{k}(\mathbf{Y}^{\prime}\mathbf{Y})^{-1/2}\\|_{F}$
	$\displaystyle\hskip 1.05273pt\geq\sigma_{\min}((\mathbf{Y}^{\prime}\mathbf{Y})^{-1/2})\\|\mathbf{R}^{k}\\|_{F}$
	$\displaystyle\hskip 1.05273pt=\frac{\\|\mathbf{R}^{k}\\|_{F}}{\sigma_{\max}(\mathbf{Y})},$	(48)

$\displaystyle\\|\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\mathbf{b}_{i}\\|_{2}$	$\displaystyle\overset{}{=}\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}\mathbf{b}_{i}-(\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})})\mathbf{b}_{i}\\|_{2}$
	$\displaystyle\overset{(a)}{\geq}\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}\mathbf{b}_{i}\\|_{2}-\\|(\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})})\mathbf{b}_{i}\\|_{2}$
	$\displaystyle\hskip 1.05273pt\geq\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}\mathbf{b}_{i}\\|_{2}-\\|\mathbf{P}_{\mathcal{R}(\mathbf{P}_{S^{k}}^{\perp}\mathbf{A}\mathbf{X})}-\mathbf{P}_{\mathcal{R}(\mathbf{R}^{k})}\\|_{2},$	(52)

Joint Sparse Recovery Using Signal Space Matching Pursuit

Abstract

Index Terms:

I Introduction

II The Proposed SSMP Algorithm

II-A Preliminaries

Definition 1 (Subspace distance [18]).

Lemma 1.

Proof.

II-B Algorithm Description

II-C Refined Identification Rule

Proposition 1.

Proof.

III Exact Joint Sparse Recovery via SSMP

III-A Definition and Lemmas

Definition 2 (RIP [27]).

Lemma 2 ([17, Lemma A.2]).

Lemma 3 ([22, Lemma 1]).

Lemma 4 ([28, Lemma 1]).

Lemma 5.

Proof.

Lemma 6 ([29, Lemma 3]).

Lemma 7 ([28, Lemma 2]).

III-B Performance Guarantee of SSMP

Proposition 2.

Proof.

Proposition 3.

Proof.

Theorem 1.

Proof.

Remark 1.

III-C Connection With Previous Efforts

IV Robustness of SSMP to Measurement Noise

Theorem 2.

Proof.

Proposition 4.

Proof.

Proposition 5.

Proof.

Theorem 3.

Proof.

V SSMP Running More Than KK Iterations

V-A Main Results

Theorem 4.

Proof.

Theorem 5.

V-B Proof of Theorem 4

Lemma 8.

Proof.

Lemma 9.

Proof.

V-C Proof of (86)

V-D Proof of (91)

VI Simulation Results

VI-A Noiseless Scenario

VI-B Noisy Scenario

VII Conclusion and Discussion

Appendix A Proof of Proposition 1

Proof.

Appendix B Proof of (112)

Proof.

Appendix C Proof of Lemma 5

Proof.

Appendix D Proofs of (46) and (47)

D-A Proof of (46)

D-B Proof of (47)

Appendix E Proof of (73)

Proof.

Appendix F Proof of Lemma 8

Proof.

Appendix G Proof of Lemma 9

Proof.

Appendix H Proof of (95)

Proof.

References

V SSMP Running More Than $K$ Iterations