Signal Recovery from Unlabeled Samples

TL;DR

This paper investigates recovering signals from unlabeled, noisy linear measurements with preserved order, establishing a duality with compressed sensing and proposing an RIP-based recovery method with theoretical and empirical validation.

Contribution

It introduces a duality between unlabeled ordered sampling and compressed sensing, and develops a low-complexity recovery algorithm with theoretical guarantees.

Findings

Successful recovery requires more samples than signal dimension.

The proposed algorithm achieves stable recovery under RIP conditions.

Empirical results show phase-transition behavior similar to compressed sensing.

Abstract

In this paper, we study the recovery of a signal from a set of noisy linear projections (measurements), when such projections are unlabeled, that is, the correspondence between the measurements and the set of projection vectors (i.e., the rows of the measurement matrix) is not known a priori. We consider a special case of unlabeled sensing referred to as Unlabeled Ordered Sampling (UOS) where the ordering of the measurements is preserved. We identify a natural duality between this problem and classical Compressed Sensing (CS), where we show that the unknown support (location of nonzero elements) of a sparse signal in CS corresponds to the unknown indices of the measurements in UOS. While in CS it is possible to recover a sparse signal from an under-determined set of linear equations (less equations than the signal dimension), successful recovery in UOS requires taking more samples than the dimension of the signal. Motivated by this duality, we develop a Restricted Isometry Property (RIP) similar to that in CS. We also design a low-complexity Alternating Minimization algorithm that achieves a stable signal recovery under the established RIP. We analyze our proposed algorithm for different signal dimensions and number of measurements theoretically and investigate its performance empirically via simulations. The results are reminiscent of phase-transition similar to that occurring in CS.

Figures

Equations

$$\mathbf{x}=\mathbf{S}\mathbf{B}\mathbf{y}+\mathbf{w},$$

$$\displaystyle{\bf L}_{p,q}=\left\{\begin{array}[]{ll}1&(p,q)=(i_{l},l)\text{ for }l\in[m],\\ 0&\text{otherwise.}\end{array}\right.$$

$$\mathbf{x}=\mathrm{vec}(\mathbf{S}\mathbf{B}\mathbf{y})+\mathbf{w}=({\mathbf{y}}^{\mathsf{T}}\otimes \mathbf{S})b+\mathbf{w},$$

$$(1-ϵ)\Vert \mathbf{H}{\Vert }_{\mathsf{F}}^2\le \Vert \mathbf{H}b{\Vert }^2\le (1+ϵ)\Vert \mathbf{H}{\Vert }_{\mathsf{F}}^2$$

$$\displaystyle{\cal E}=\bigcup_{{\bf H}\in{\cal H}}\Big{\{}{\mathbb{b}}:\big{|}\|{\bf H}{\mathbb{b}}\|^{2}-\|{\bf H}\|_{\sf F}^{2}\big{|}>{\epsilon}\|{\bf H}\|_{\sf F}^{2}\Big{\}}.$$

$$\displaystyle\Big{\{}{\bf B}:(1-{\epsilon})m\|{\bf y}\|^{2}\leq\|{\bf S}{\bf B}{\bf y}\|^{2}\leq(1+{\epsilon})m\|{\bf y}\|^{2}\Big{\}},$$

$$\displaystyle{\mathbb{P}}\Big{[}\Big{\{}{\bf B}:\big{|}\|{\bf S}{\bf B}{\bf y}\|^{2}-m\|{\bf y}\|^{2}\big{|}>\frac{{\epsilon}}{2}m\|{\bf y}\|^{2}\Big{\}}\Big{]}\leq 2e^{-cm{\epsilon}^{2}}.$$

$$\underset{\mathbf{g}\in {\mathcal{N}}_{ϵ}}{\max }\Vert \mathbf{S}\mathbf{B}\mathbf{g}\Vert \le \underset{\mathbf{y}\in {\mathcal{B}}_2^k}{\max }\Vert \mathbf{S}\mathbf{B}\mathbf{y}\Vert \le (1+\frac{ϵ}{2})\underset{\mathbf{g}\in {\mathcal{N}}_{ϵ}}{\max }\Vert \mathbf{S}\mathbf{B}\mathbf{g}\Vert ,$$

$$\displaystyle{\mathbb{P}}\Big{[}$$

$$\displaystyle\leq{\mathbb{P}}\Big{[}\bigcup_{{\bf g}\in{\cal N}_{{\epsilon}}}\Big{\{}{\bf B}:\big{|}\|{\bf S}{\bf B}{\bf g}\|^{2}-m\|{\bf g}\|^{2}\big{|}>m\|{\bf g}\|^{2}\frac{{\epsilon}}{2}\Big{\}}\Big{]}$$

$$\stackrel{(i)}{\le }|{\mathcal{N}}_{ϵ}|e^{-cm{ϵ}^2}=2(1+\frac{2}{ϵ}{)}^k{e}^{-cm{ϵ}^2}.$$

$$\displaystyle{\mathbb{P}}\Big{[}$$

$$\le 2|\mathcal{S}|(1+\frac{2}{ϵ}{)}^k\left(\genfrac{}{}{0ex}{}{n}{m}\right)e^{-cm{ϵ}^2}=2(1+\frac{2}{ϵ}{)}^k(\genfrac{}{}{0ex}{}{n}{m})e^{-cm{ϵ}^2},$$

$$(1-\delta )\Vert \mathbf{H}-{\mathbf{H}}^{\prime }{\Vert }_{\mathsf{F}}^2\le \Vert (\mathbf{H}-{\mathbf{H}}^{\prime })b{\Vert }^2\le (1+\delta )\Vert \mathbf{H}-{\mathbf{H}}^{\prime }{\Vert }_{\mathsf{F}}^2,$$

$$(\mathbf{H}-{\mathbf{H}}^{\prime })b=\mathbf{S}\mathbf{B}\mathbf{y}-{\mathbf{S}}^{\prime }\mathbf{B}{\mathbf{y}}^{\prime }=(0,0,\dots ,0,g{)}^{\mathsf{T}},$$

$$d_{\mathbf{H},{\mathbf{H}}^{\prime }}$$

$$=\sqrt{m\Vert \mathbf{y}{\Vert }^2+m\Vert {\mathbf{y}}^{\prime }{\Vert }^2-2m{\nu }_{\mathbf{S},{\mathbf{S}}^{\prime }}⟨\mathbf{y},{\mathbf{y}}^{\prime }⟩},$$

$$⟨\mathbf{H},{\mathbf{H}}^{\prime }⟩$$

$$=⟨\mathbf{y},{\mathbf{y}}^{\prime }⟩⟨\mathbf{S},{\mathbf{S}}^{\prime }⟩,$$

$${\nu }_{\mathbf{S},{\mathbf{S}}^{\prime }}=\frac{⟨\mathbf{S},{\mathbf{S}}^{\prime }⟩}{m}=\frac{\mathsf{tr}({\mathbf{S}}^{\mathsf{T}}{\mathbf{S}}^{\prime })}{m}\in [0,1],$$

$$(1-\delta )\Vert \mathbf{H}-{\mathbf{H}}^{\prime }{\Vert }_{\mathsf{F}}^2\le \Vert (\mathbf{H}-{\mathbf{H}}^{\prime })b{\Vert }^2\le (1+\delta )\Vert \mathbf{H}-{\mathbf{H}}^{\prime }{\Vert }_{\mathsf{F}}^2,$$

$$\frac{k}{n}\log (1+\frac{2}{\delta })+\frac{1}{n}\log \left(\genfrac{}{}{0ex}{}{n}{m}\right)-\frac{cm{\delta }^2{\mu }^2}{2n}<0.$$

$$\kappa \log (1+\frac{2}{\delta })+h(\theta )-\frac{c{\delta }^2{\mu }^2}{2}(1-\theta )<0,$$

$$\frac{n}{k}\approx \frac{2\log (1+\frac{2}{\delta })}{c{\delta }^2{\mu }^2}$$

$$\hat{\mathbf{H}}={\mathrm{argmin}}_{\mathbf{H}\in \mathcal{H}}\Vert \mathbf{H}b-\check{\mathbf{x}}\Vert .$$

$$\Vert \hat{\mathbf{H}}b-\check{\mathbf{x}}\Vert \le \Vert \check{\mathbf{H}}b-\check{\mathbf{x}}\Vert =\Vert \mathbf{w}\Vert .$$

$$\Vert \hat{\mathbf{H}}{\Vert }_{\mathsf{F}}$$

$$\stackrel{(ii)}{\le }\Vert \check{\mathbf{H}}{\Vert }_{\mathsf{F}}\frac{\sqrt{1+ϵ}}{\sqrt{1-ϵ}}(1+\frac{2}{\sqrt{\mathsf{snr}}}),$$

$$\max \{\Vert \hat{\mathbf{H}}{\Vert }_{\mathsf{F}},\Vert \check{\mathbf{H}}{\Vert }_{\mathsf{F}}\}\le \Vert \check{\mathbf{H}}{\Vert }_{\mathsf{F}}\frac{\sqrt{1+ϵ}}{\sqrt{1-ϵ}}(1+\frac{2}{\sqrt{\mathsf{snr}}}).$$

$$\frac{\Vert \hat{\mathbf{H}}-\check{\mathbf{H}}{\Vert }_{\mathsf{F}}}{\Vert \check{\mathbf{H}}{\Vert }_{\mathsf{F}}}=\frac{d_{\hat{\mathbf{H}},\check{\mathbf{H}}}}{\Vert \check{\mathbf{H}}{\Vert }_{\mathsf{F}}}\le \mu \frac{\sqrt{1+ϵ}}{\sqrt{1-ϵ}}(1+\frac{2}{\sqrt{\mathsf{snr}}}).$$

$$\sqrt{1-\delta }\frac{\Vert \check{\mathbf{H}}-\hat{\mathbf{H}}{\Vert }_{\mathsf{F}}}{\Vert \check{\mathbf{H}}{\Vert }_{\mathsf{F}}}$$

$$\stackrel{(b)}{\le }\frac{2\Vert \mathbf{w}\Vert }{\Vert \check{\mathbf{H}}b\Vert \sqrt{1-ϵ}}=\frac{2}{\sqrt{\mathsf{snr}(1-ϵ)}},$$

$$\frac{\Vert \check{\mathbf{H}}-\hat{\mathbf{H}}{\Vert }_{\mathsf{F}}}{\Vert \check{\mathbf{H}}{\Vert }_{\mathsf{F}}}\le \frac{2}{\sqrt{\mathsf{snr}(1-\delta )(1-ϵ)}}\stackrel{(i)}{\le }\frac{2}{\sqrt{\mathsf{snr}(1-2\delta )}},$$

$$\hat{\mathbf{H}}\in \{\mathbf{H}:\Vert \mathbf{H}b-\check{\mathbf{x}}\Vert \le \eta \Vert \mathbf{w}\Vert \},$$

Full text

I Introduction

The recovery of a vector-valued signal ${\bf y}\in{\mathbb{R}}^{k}$ from a set of linear and possibly noisy measurements ${\bf x}={\bf B}{\bf y}+{\bf w}$, with an $n\times k$ measurement matrix ${\bf B}$, is the classical problem of linear regression in statistical inference and is arguably the most widely-studied problem in statistics, mathematics, and computer science. For $n\geq k$, one has an over-determined set of noisy linear equations, and the Maximum Likelihood (ML) estimate of ${\bf y}$ under the additive Gaussian noise ${\bf w}$ is given by the well-known method of least squares. For $n<k$, in contrast, one deals with an under-determined set of linear equations, which is only solvable if some additional a priori information about the signal ${\bf y}$ is available. For example, the whole field of Compressed Sensing (CS) deals with the recovery of the signal ${\bf y}$ when it is sparse or more generally compressible, i.e., it has only a few significantly large coefficients when represented in a suitable basis [1, 2, 3].

Almost all the past research in linear regression mainly deals with exploiting the underlying signal structure, whereas it is generally assumed that the regression matrix ${\bf B}$ is fully known. In practice, the matrix ${\bf B}$ is implemented through a measurement device, where due to physical limitations, there might be some uncertainty or mismatch between the intended matrix ${\bf B}$ and the one realized via measurement devices. This has resulted in a surge of interest in generalized linear regression problems in which the matrix ${\bf B}$ is mismatched [4, 5] or is known only up to some unknown transformation. In this paper, we are interested in the unlabeled sensing case, where the observation model is given by

[TABLE]

where ${\bf B}$ is a completely known $n\times k$ matrix, and where ${\bf S}$ is an unknown $m\times n$ matrix with 0-1 components that has only a single $1$ in each row and samples (selects) $m$ out of $n$ elements of ${\bf B}{\bf y}$. Although ${\bf S}$ is not known, it belongs to an a priori known set of selection matrices ${\cal S}$. Identifying ${\bf S}$ in (1), therefore, corresponds to associating the measurements ${\bf x}$ to their corresponding linear projections in ${\bf B}$. Once ${\bf S}$ is identified, (1) reduces to a linear regression problem whose solution can be obtained via standard techniques.

In [6, 7], a variant of this problem was studied when ${\cal S}$ is the set of all permutations of $n$ measurements taken by ${\bf B}$. It was shown that if the measurement matrix ${\bf B}$ is generated randomly, any arbitrary $k$-dim signal can be recovered from a set of $n=2k$ noiseless unlabeled measurements, where this bound was shown to be tight. In [8], the authors studied a similar problem but rather than recovering the unknown signal ${\bf y}$, they obtained a scaling law of the Signal-to-Noise Ratio (SNR) required for identifying the unknown permutation matrix in ${\cal S}$, where they showed that the required SNR increases logarithmically with the signal dimension. A Multiple Measurement version of the problem defined by ${\bf y}^{i}={\bf S}{\bf B}{\bf x}^{i}+{\bf w}^{i}$, $i=1,\dots,l$, where the signals ${\bf x}^{i}$ might vary but ${\bf S}$ remains the same across all $l$ measurements, was recently studied in [9, 10] and was proved to yield a signal recovery at a finite SNR provided that $l$ is sufficiently large.

A practical scenario well-modeled by (1) is sampling in the presence of jitter [11], in which ${\cal S}$ consists of 0-1 matrices with 1s as their diagonal elements and with some off-diagonal 1s representing the location of the jittered samples. A similar problem arises in molecular channels [12] and in the reconstruction of phase-space dynamics of linear/nonlinear dynamical systems [13] because of synchronization issues. Unlabeled regression (1) is also encountered in noncooperative multi-target tracking, e.g., in radar, where the receiver only observes the unlabeled data associated to all the targets, thus, ${\cal S}$ consists of the set of all permutations corresponding to all possible data-target associations [14]. A quite similar scenario arises in robotics and a well-known classical problem is Simultaneous Localization and Mapping (SLAM) where robots measure their relative coordinates and recovery of the underlying geometry requires suitably permuting the data.

A different line of work well-modeled by (1) is the genome assembly problem from shotgun reads [15, 16] in which a sequence $y\in\{{\mathtt{A}},{\mathtt{T}},{\mathtt{G}},{\mathtt{C}}\}^{k}$ of length $k$ is assembled from an unknown permutation of its sub-sequences called “reads”. Designing efficient recovery (assembly) algorithms is still an active research area (see [16] and the refs. therein.). Communication over the classical noisy deletion channel is another example of (1), where ${\bf B}$ represents the linear encoding matrix with elements belonging to a finite field ${\mathbb{F}}_{q}$, for some prime number $q$, where ${\bf y}\in{\mathbb{F}}_{q}^{k}$ denotes the $k$-dim vector containing $k$ information symbols, where ${\bf w}\in{\mathbb{F}}_{q}^{m}$ is the additive noise of the channel, and where ${\cal S}$ consists of all selection operators that keep only $m\leq n$ out of $n$ encoded symbols in ${\bf B}{\bf y}$ while preserving their order. In particular, in contrast with the erasure channel, where the location of erased symbols is known, in a deletion channel the location of deleted symbols is missing. This makes designing good encoding and decoding techniques as well as identifying the capacity quite challenging [17].

I-A Contribution

Since satisfactory efficient algorithms for solving the unlabeled sampling problem in (1) are yet unknown, we make progress towards this goal by addressing a relevant subproblem that we refer to as Unlabeled Ordered Sampling (UOS), where ${\cal S}$ is the set of all 0-1 ordered sampling matrices that select only $m$ out of $n$ components of ${\bf B}{\bf y}$ for some $m\leq n$ while preserving their relative order. This occurs in many practical scenarios (e.g., in a deletion channel) and is illustrated in Fig. 1. We discover a duality between this problem and the CS problem [1, 2, 3], where the unknown location of samples in the former corresponds in a natural way to the unknown location of nonzero coefficients of the signal in the latter. To the best of our knowledge, this is the first paper addressing the underlying duality connection between the unlabeled sensing in (1) and classical CS. Designing a low-complexity algorithm for recovering the desired signal from its unlabeled samples is generally considered to be a challenging problem [6, 7, 8, 11]. In particular, in most practically relevant situations, ${\cal S}$ is a very large set such that a naive exhaustive search over ${\cal S}$ would be unfeasible. In this paper, however, we are able to exploit the underlying ordering in the UOS case to design an efficient Alternating Minimization recovery algorithm ($\mathtt{AltMin}$). We analyze the performance of our proposed algorithm theoretically and investigate its performance empirically via simulations.

I-B Notation

We denote vectors by boldface small letters (e.g., ${\bf x}$), matrices by boldface capital letters (e.g., ${\bf X}$), scalars by non-boldface letters (e.g., $x$), and sets by capital calligraphic letters (e.g., ${\cal X}$). For integers $k_{1},k_{2}\in\mathbb{Z}$, we use the shorthand notation [{\text{k_{1}$\,:\,$k_{2}}}]=\big{\{}k_{1},k_{1}+1,\dots,k_{2}\}, where the set is empty when $k_{2}<k_{1}$, and use the simplified notation $[k_{1}]$ for [{\text{1$\,:\,$k_{1}}}]. We denote the $k$-th component of a vector ${\bf x}$ by ${\bf x}_{k}$ and a subvector of ${\bf x}$ with indices in the range $k_{1}$ : $k_{2}$ by {\bf x}_{{\text{k_{1}$\,:\,$k_{2}}}}. We denote the element of a matrix ${\bf X}$ at row $r$ and column $c$ by ${\bf X}_{r,c}$ and use an indexing notation similar to that for vectors for submatrices of ${\bf X}$, namely, ${\bf X}_{r,c}$, {\bf X}_{r,{\text{c_{1}$\,:\,$c_{2}}}}, {\bf X}_{{\text{r_{1}$\,:\,$r_{2}}},c}, and {\bf X}_{{\text{r_{1}$\,:\,$r_{2}}},{\text{c_{1}$\,:\,$c_{2}}}}. We denote the Kronecker product of a $p\times q$ matrix ${\bf X}$ and an $r\times s$ matrix ${\bf Y}$ by a $pr\times qs$ matrix ${\bf X}\otimes{\bf Y}$ that has $p\times q$ blocks with the $ij$-th block, $i\in[p],j\in[q]$, given by the $r\times s$ matrix ${\bf X}_{i,j}{\bf Y}$. We use $\mathsf{tr}(.)$ for the trace and $(.)^{\sf T}$ for the transpose operator. We represent the inner product between two vectors ${\bf x}$ and ${\bf y}$ and two matrices ${\bf X}$ and ${\bf Y}$, with ${\langle{\bf x},{\bf y}\rangle}={\bf x}^{\sf T}{\bf y}$ and ${\langle{\bf X},{\bf Y}\rangle}=\mathsf{tr}({\bf X}^{\sf T}{\bf Y})$ respectively, and use $\|{\bf x}\|=\sqrt{{\langle{\bf x},{\bf x}\rangle}}$ for the $l_{2}$ norm of a vector ${\bf x}$, and $\|{\bf X}\|_{\sf F}=\sqrt{{\langle{\bf X},{\bf X}\rangle}}$ for the Frobenius norm of a matrix ${\bf X}$. We denote a diagonal matrix with the elements $(x_{1},\dots,x_{k})$ by ${\text{diag}}(x_{1},\dots,x_{k})$ and the identity matrix of the order $k$ with ${\bf I}_{k}$. We represent a sequence of vectors and sequence of matrices by upper indices, e.g., ${\bf x}^{1},{\bf x}^{2},\cdots$ and ${\bf X}^{1},{\bf X}^{2},\cdots$. We denote a Gaussian distribution with a mean $\mu$ and a variance $\sigma^{2}$ by ${\sf N}(\mu,\sigma^{2})$. We use $O(.)$ and $o(.)$ for the big-O and the small-O respectively.

II Statement of the Problem

In this section, we start from the more familiar CS problem. Then, we introduce the UOS and make a duality connection between the two.

II-A Basic Setup

Let $n$ and $m$ be positive integers with $m\leq n$, and let ${\cal I}=\{i_{1},\dots,i_{m}\}\subseteq[n]$ be a subset of $[n]$ consisting of ordered elements of $[n]$ satisfying $i_{l}<i_{l+1}$. We define the lift-up operator associated to ${\cal I}$ as a linear map from ${\mathbb{R}}^{m}$ to ${\mathbb{R}}^{n}$ given by the $\{0,1\}$-valued $n\times m$ tall matrix ${\bf L}$ with

[TABLE]

The operator ${\bf L}$ embeds $m$ components of ${\bf x}=(x_{1},\dots,x_{m})^{\sf T}$ into the index set ${\cal I}$ in the $n$-dim vector ${\bf L}{\bf x}$, while keeping their relative order, i.e., $({\bf L}{\bf x})_{i_{l}}=x_{l}$ for $l\in[m]$, and fills the rest with [math]. For example, for $n=4$, $m=3$, ${\bf x}=(x_{1},x_{2},x_{3})$, and ${\cal I}=\{1,3,4\}$, we have ${\bf L}{\bf x}=(x_{1},0,x_{2},x_{3})$. We define the collection of all ${n\choose m}$ lift-up operators by ${\cal L}_{n,m}$.

II-B Compressed Sensing

In CS [1, 2, 3], the goal is to recover a sparse or compressible signal by taking less measurements than its dimension. We call a signal ${\bf z}\in{\mathbb{R}}^{n}$ $m$-sparse ($m$-compressible) if it has only $m$ nonzero (significantly large) components. For simplicity, we focus on sparse rather than compressible signals. We fix $m,n$ with $m\leq n$ as before. We define an instance of the CS problem for an $m$-sparse signal ${\bf z}\in{\mathbb{R}}^{n}$ by the triple $({\bf x},{\bf L},{\bf A})$, where ${\bf x}\in{\mathbb{R}}^{m}$ denotes the nonzero elements of ${\bf z}$, where ${\bf L}\in{\cal L}_{n,m}$ encodes the location of these nonzero elements, and where ${\bf A}$ is the $k\times n$ matrix whose rows correspond to $k$ linear measurements. The $m$-sparse signal ${\bf z}$ is generated by embedding the $m$ components of ${\bf x}$ via the lift-up matrix ${\bf L}$ as ${\bf z}={\bf L}{\bf x}$, where it is seen that ${\bf z}$, albeit being $n$-dim, has at most $m$ nonzero components. This is illustrated in Fig. 2. In CS, the sparse signal ${\bf z}$ is sampled via the measurement matrix ${\bf A}$ producing $k$ measurements ${\bf y}={\bf A}{\bf z}$. The goal is to recover the unknown $({\bf x},{\bf L})$ or equivalently ${\bf z}={\bf L}{\bf x}$, from the known $({\bf y},{\bf A})$. The crucial idea in CS is that by exploiting the underlying sparsity, ${\bf z}$ can be recovered by taking less samples than its embedding dimension $n$. The practically interesting regime of parameters in CS is given by $m\leq k\leq n$, where the number of measurements $k$ is more than the number of nonzero coefficients of the signal $m$ but much less than the embedding dimension $n$.

II-C Unlabeled Ordered Sampling

By changing the role of ${\bf y}$ and ${\bf x}$ in the CS problem, we obtain an instance of UOS in (1) as the dual problem with ${\bf B}={\bf A}^{\sf T}$ and ${\bf S}={\bf L}^{\sf T}$. This is also illustrated in Fig. 2. In UOS, a $k$-dim signal ${\bf y}$ is oversampled via the tall $n\times k$ matrix ${\bf B}={\bf A}^{\sf T}$, which gives $n$ measurements ${\bf z}={\bf A}^{\sf T}{\bf y}$. The resulting over-complete measurements ($n\geq k$) in ${\bf z}$ are subsampled by the $m\times n$ matrix ${\bf S}={\bf L}^{\sf T}$, which selects only $m$ out of $n$ measurements in ${\bf z}$ and yields the unlabeled samples ${\bf x}$. The goal is to recover the unknown signal ${\bf y}$ from the known $({\bf x},{\bf B})$, without knowing ${\bf S}$ explicitly. Compared with CS, where the support (location of the nonzero values) of ${\bf z}$ is unknown, in UOS  the indices/labels of the measurements are missing. However, it is not difficult to check that, due to the special structure of ${\bf S}$, the relative order of the measurements is still preserved. We define the set of all such selection matrices by ${\cal S}_{m,n}:=\{{\bf L}^{\sf T}:{\bf L}\in{\cal L}_{n,m}\}$. In contrast with the lift-up matrices in ${\cal L}_{n,m}$, which embed signals with a lower dimension $m$ in a higher dimension $n$, the selection matrices in ${\cal S}_{m,n}$ reduce the dimensionality by sampling only $m$ out of $n$ components of their input (while keeping the order). The practically relevant regime of parameters for UOS is given by $k\leq m\leq n$, where one needs to take more measurements $m$ than the signal dimension $k$ to overcome the uncertainty caused by the unlabeled sensing. Motivated by the duality between CS and UOS, in the next section, we develop a Restricted Isometry Property (RIP), which resembles its counterpart in CS [3] and guarantees a stable signal recovery in UOS. We also design an efficient low-complexity algorithm able to recover the target signal under the RIP.

III Restricted Isometry Property

III-A Basic Setup

For the rest of the paper, we focus on UOS  illustrated in Fig. 1. We consider a $k$-dim signal ${\bf y}$ and an $n$-dim vector of measurements ${\bf z}={\bf B}{\bf y}$ taken via an $n\times k$ matrix ${\bf B}$, where ${\bf B}={\bf A}^{\sf T}$ with ${\bf A}$ being the $k\times n$ matrix in the CS variant (see Fig. 2). An instance of UOS  is given by the triple $({\bf y},{\bf S},{\bf B})$, where the goal is to recover the unknown signal ${\bf y}$ from the noisy unlabeled measurements ${\bf x}={\bf S}{\bf B}{\bf y}+{\bf w}$ without having any explicit knowledge about ${\bf S}$ except that ${\bf S}\in{\cal S}_{m,n}$. This corresponds to a variant of unlabeled sensing problem in (1) with the set of possible transformations ${\cal S}$ given by ${\cal S}_{m,n}$. In the rest of the paper, for simplicity, we drop the index $m,n$ and denote ${\cal S}_{m,n}$ by ${\cal S}$. Using the ${\rm vec}$ notation, we can write

[TABLE]

where ${\mathbb{b}}={\rm vec}({\bf B})$ denotes the vector obtained by stacking the columns of ${\bf B}$ and where we used the well-known identity ${{\rm vec}({\bf X}{\bf Y}{\bf Z})=({\bf Z}^{\sf T}\otimes{\bf X}){\rm vec}({\bf Y})}$ for matrices ${\bf X},{\bf Y},{\bf Z}$ of appropriate dimensions. We will use the notations ${\bf S}{\bf B}{\bf y}$ and $({\bf y}^{\sf T}\otimes{\bf S}){\mathbb{b}}$ interchangeably across the paper. Note that ${\mathbb{b}}$ induces a linear map $({\bf y}^{\sf T}\otimes{\bf S})\mapsto({\bf y}^{\sf T}\otimes{\bf S}){\mathbb{b}}$ from the signal set ${{\cal H}:=\{{\bf y}^{\sf T}\otimes{\bf S}:{\bf y}\in{\mathbb{R}}^{k},{\bf S}\in{\cal S}\}}$ into ${\mathbb{R}}^{m}$. We define the Signal-to-Noise Ratio (SNR) in (5) by ${\mathsf{snr}}=\frac{\|({\bf y}^{\sf T}\otimes{\bf S}){\mathbb{b}}\|^{2}}{\|{\bf w}\|^{2}}=\frac{\|{\bf S}{\bf B}{\bf y}\|^{2}}{\|{\bf w}\|^{2}}$.

III-B Restricted Isometry Property on ${\cal H}$

Our goal is to recover the desired signal ${\bf y}$. Since the signal set ${\cal H}$ is unbounded, a requirement is that at least $\|{\bf y}\|$ should be feasibly recovered from the measurements $({\bf y}^{\sf T}\otimes{\bf S}){\mathbb{b}}$. A sufficient condition for this is the Restricted Isometry Property (RIP) over ${\cal H}$, which resembles a similar property in CS [3].

Definition 1** (RIP over ${\cal H}$)**

Let ${\bf B}$ be a matrix and let ${\mathbb{b}}={\rm vec}({\bf B})$. The linear map ${\bf H}\mapsto{\bf H}{\mathbb{b}}$ induced by ${\mathbb{b}}$ satisfies the RIP over ${\cal H}$ with a constant ${\epsilon}\in(0,1)$ if

[TABLE]

holds for all ${\bf H}\in{\cal H}$. $\lozenge$

We prove that under suitable conditions on $m,n,k$, we can obtain a matrix ${\bf B}$ satisfying the RIP in (6) by sampling components of ${\bf B}$ i.i.d. from the Gaussian distribution. This is summarized in the following proposition.

Proposition 1

Let ${\bf B}$ be a random matrix with i.i.d. ${\sf N}(0,1)$ components and let ${\mathbb{b}}={\rm vec}({\bf B})$. Then, there is a constant $c>0$ such that ${\bf B}$ satisfies RIP over ${\cal H}$ with a probability larger than $1-2(1+\frac{2}{{\epsilon}})^{k}{n\choose m}e^{-cm{\epsilon}^{2}}$. $\square$

Proof:

Let us define the following probability event (for the random realization of ${\mathbb{b}}$):

[TABLE]

We need to prove that ${\mathbb{P}}[{\cal E}]\leq 2(1+\frac{2}{{\epsilon}})^{k}{n\choose m}e^{-cm{\epsilon}^{2}}$. Since ${\cal E}$ consists of the union of a continuum of events labeled with ${\bf H}\in{\cal H}$, the conventional union bound is not immediately applicable. So, we need to first quantize the labeling set ${\cal H}$ appropriately and then apply the union bound. We do this by using the net argument (see [18, 19] for further discussion on using the net technique). We start by first deriving a concentration bound for a fixed ${\bf H}={\bf y}^{\sf T}\otimes{\bf S}\in{\cal H}$. From ${\bf H}{\mathbb{b}}={\bf S}{\bf B}{\bf y}$ and ${\mathbb{E}}[\|{\bf H}{\mathbb{b}}\|^{2}]=\|{\bf H}\|_{{\sf F}}^{2}$, the subevent of (7) corresponding to the given ${\bf H}$ can be written as

[TABLE]

where we used the fact $\|{\bf H}\|_{\sf F}^{2}=m\|{\bf y}\|^{2}$. Note that for a fixed ${\bf S}$ and ${\bf y}$, the vector ${\bf S}{\bf B}{\bf y}$ is an $m$-dim vector with i.i.d. ${\sf N}(0,1)$ components. From the Gaussian concentration result [20], there is a $c>0$ such that

[TABLE]

We first generalize this concentration result to all the vectors ${\bf y}$ inside the unit $l_{2}$ ball ${{\cal B}_{2}^{k}:=\{{\bf y}\in{\mathbb{R}}^{k}:\|{\bf y}\|\leq 1\}}$. We consider a discrete $\frac{{\epsilon}}{2}$-net (grid) of minimal size $N$ over ${\cal B}_{2}^{k}$ denoted by ${\cal N}_{{\epsilon}}=\{{\bf g}^{1},\dots,{\bf g}^{N}\}$ that satisfies $\max_{{\bf y}\in{\cal B}_{2}^{k}}\min_{{\bf g}\in{\cal N}_{{\epsilon}}}\|{\bf g}-{\bf y}\|\leq\frac{{\epsilon}}{2}$ (see [19] for further details). Consider the set of $N$ spheres with centers belonging to the net ${\cal N}_{{\epsilon}}$ each having a radius $\frac{{\epsilon}}{2}$. All these spheres lie inside a sphere of radius $1+\frac{{\epsilon}}{2}$ centered at the origin. Thus, using the volume inequality $N(\frac{{\epsilon}}{2})^{k}{\text{vol}}({\cal B}_{2}^{k})\leq(1+\frac{{\epsilon}}{2})^{k}{\text{vol}}({\cal B}_{2}^{k})$ for the $k$-dim unit $l_{2}$ ball ${\cal B}_{2}^{k}$, we obtain that such a minimal net consists of at most $N\leq(1+\frac{2}{{\epsilon}})^{k}$ points. We also have (see Lemma 2.19. in [19])

[TABLE]

which implies that the operator norm of ${\bf S}{\bf B}$ can be well approximated via the points in ${\cal N}_{{\epsilon}}$. From (10), we obtain

[TABLE]

where in $(i)$ we applied the union bound over ${\cal N}_{{\epsilon}}$ and used the bound (9), which holds for any ${\bf y}$ and in particular for any ${\bf g}\in{\cal N}_{\epsilon}$. Finally applying the union bound over all ${n\choose m}$ possible selection matrices ${\bf S}\in{\cal S}$, we have

[TABLE]

which from (7) implies that ${\mathbb{P}}[{\cal E}]\leq 2(1+\frac{2}{{\epsilon}})^{k}{n\choose m}e^{-cm{\epsilon}^{2}}$. This completes the proof.     \sqcap$$\sqcup

III-C Restricted Isometry Property on ${\cal H}-{\cal H}$

In terms of signal recovery, we will need a stronger RIP over the Minkowski difference of the signal set ${\cal H}$ defined by ${\cal H}-{\cal H}:=\{{\bf H}-{\bf H}^{\prime}:{\bf H},{\bf H}^{\prime}\in{\cal H}\}$. In this section, our goal is to develop a suitable notion of RIP over ${\cal H}-{\cal H}$. Our approach in this section follows from similar techniques in [21].

Let ${\bf B}$ be the $m\times n$ measurement matrix and let ${\mathbb{b}}={\rm vec}({\bf B})$ as before. Similar to the RIP over ${\cal H}$, for the RIP over ${\cal H}-{\cal H}$, it seems reasonable to impose the condition

[TABLE]

with some ${\delta}\in(0,1)$, to hold for all ${\bf H},{\bf H}^{\prime}\in{\cal H}$. As we will see in Section III-D, an RIP condition as in (12) is sufficient to guarantee a stable signal recovery for UOS. In this section, as in Section III-B, we attempt to construct a matrix ${\bf B}$ satisfying (12) via sampling the components of ${\bf B}$ randomly. Hence, we need to prove that, in a suitable regime of parameters $m,n,k$, and ${\delta}$, any realization of the matrix ${\bf B}$ and as a result ${\mathbb{b}}$ satisfies (12) for all ${\bf H},{\bf H}^{\prime}\in{\cal H}$ with a very high probability. Unfortunately, such a universal concentration result over ${\cal H}-{\cal H}$ does not immediately hold as can be seen from the following simple example111This should be contrasted with classical CS, where ${\cal H}$ is the class of $k$-sparse $n$-dim signals, thus, ${\cal H}-{\cal H}$ is a subset of $2k$-sparse signals, and deriving the RIP over ${\cal H}$ and ${\cal H}-{\cal H}$ requires quite similar techniques. In UOS, in contrast, although it is easy to derive an RIP over ${\cal H}$, obtaining a suitable notion of RIP over ${\cal H}-{\cal H}$ is more challenging..

Example 1

Let ${\bf B}$ be a random matrix with i.i.d. ${\sf N}(0,1)$ components. Consider the signals ${\bf H}={\bf y}^{\sf T}\otimes{\bf S}$ and ${\bf H}^{\prime}={{\bf y}^{\prime}}^{\sf T}\otimes{\bf S}^{\prime}$, where ${\bf y}={\bf y}^{\prime}$ and where all the rows of ${\bf S}$ and ${\bf S}^{\prime}$ are similar except the last row. For such a case, we have

[TABLE]

where $g=({\bf B}_{r_{m},:}-{\bf B}_{r^{\prime}_{m},:}){\bf y}$, where $r_{m},r^{\prime}_{m}\in[n]$ denote the indices of the last rows selected by ${\bf S}$ and ${\bf S}^{\prime}$ respectively. From (13), it results that for the given ${\bf H},{\bf H}^{\prime}$, $\|({\bf H}-{\bf H}^{\prime}){\mathbb{b}}\|^{2}=|g|^{2}$, where ${g\sim{\sf N}(0,2\|{\bf y}\|^{2})}$ from the independence of the rows of ${\bf B}$. However, $|g|^{2}$ does not concentrate very well around its mean ${\mathbb{E}}[|g|^{2}]$. This implies that the universal concentration result over ${\cal H}-{\cal H}$ in (12) can not hold for any ${\delta}\in(0,1)$. $\lozenge$

From Example 1, it is seen that we can not hope to construct, with a high probability, a matrix ${\bf B}$ satisfying RIP in (12) under the random sampling of components of ${\bf B}$. A direct inspection in the Example 1, however, reveals that the main obstacle to establishing the concentration (12) are those signals ${\bf H},{\bf H}^{\prime}$ where ${\bf y}\approx{\bf y}^{\prime}$ and ${\bf S}\approx{\bf S}^{\prime}$, thus, ${\bf H}\approx{\bf H}^{\prime}$. Since we use the RIP to guarantee a universal stable recovery for all the signals in ${\cal H}$, intuitively speaking, the troublesome cases ${\bf H}\approx{\bf H}^{\prime}$ for (12) are not problematic at all in terms of signal recovery. To take this into account, we develop a modified version of RIP over ${\cal H}-{\cal H}$ in (12) up to a fixed precision.

We need some notation first. We define the following metric over the signal set ${\cal H}$

[TABLE]

where we used $\|{\bf H}\|_{\sf F}^{2}=m\|{\bf y}\|^{2}$, $\|{\bf H}^{\prime}\|_{\sf F}^{2}=m\|{\bf y}^{\prime}\|^{2}$, and that

[TABLE]

where we also defined the similarity metric between selection matrices ${\bf S}$ and ${\bf S}^{\prime}$ by

[TABLE]

which gives the fraction of similar rows in ${\bf S}$ and ${\bf S}^{\prime}$. We define the Relaxed RIP over ${\cal H}-{\cal H}$ up to a precision ${\mu}$ as follows.

Definition 2** (R-RIP over ${\cal H}-{\cal H}$)**

Let ${\bf B}$ be a matrix and let ${\mathbb{b}}={\rm vec}({\bf B})$. The linear map ${\bf H}\mapsto{\bf H}{\mathbb{b}}$ induced by ${\mathbb{b}}$ satisfies the Relaxed RIP (R-RIP) over ${\cal H}-{\cal H}$ with a parameter ${\delta}$ and a precision ${\mu}$ if

[TABLE]

holds for all ${\bf H},{\bf H}^{\prime}$, with $d_{{\bf H},{\bf H}^{\prime}}\geq{\mu}\max\{\|{\bf H}\|_{\sf F},\|{\bf H}^{\prime}\|_{\sf F}\}$. $\lozenge$

The next proposition shows that, in a suitable regime of parameters $m,n,k$, and ${\delta},{\mu}$, a matrix ${\bf B}$ with ${\sf N}(0,1)$ elements satisfies R-RIP over ${\cal H}-{\cal H}$ with a high probability.

Proposition 2

Let ${\bf B}$ be a random matrix with i.i.d. ${\sf N}(0,1)$ elements and let ${\mathbb{b}}={\rm vec}({\bf B})$. There is a constant $c>0$ such that ${\bf B}$ satisfies R-RIP on ${\cal H}-{\cal H}$ with parameters ${\delta},{\mu}$ with a probability larger than $1-2(1+\frac{2}{{\delta}})^{2k}{n\choose m}^{2}e^{-cm{\delta}^{2}{\mu}^{2}}$. $\square$

Proof:

Proof in Appendix A-A.     \sqcap$$\sqcup

Remark 1

By taking the union bound over the concentration results derived in Proposition 1 and 2 and using the fact that the concentration in Proposition 1 is much sharper, we will assume in the rest of the paper that ${\bf B}$ satisfies both RIP over ${\cal H}$ with a parameter ${\epsilon}$ and R-RIP over ${\cal H}-{\cal H}$ with parameters ${\delta},{\mu}$, with a probability approximately given by $1-2(1+\frac{2}{{\delta}})^{2k}{n\choose m}^{2}e^{-cm{\delta}^{2}{\mu}^{2}}$. In particular, ${\epsilon}$ can always be taken much lower than ${\delta}$ (${\epsilon}\ll{\delta}$) without affecting this bound. $\lozenge$

For a suitably selected set of parameters $n,m,k,$ and ${\delta},{\mu}$, ${\bf B}$ satisfies R-RIP with a high probability if $2(1+\frac{2}{{\delta}})^{2k}{n\choose m}^{2}e^{-cm{\delta}^{2}{\mu}^{2}}$ is small. Asymptotically as $m,n,k\to\infty$, this is satisfied provided that

[TABLE]

In particular, if ${\frac{m}{n}\to\rho}$ and ${\frac{k}{n}\to{\kappa}}$, applying the Stirling’s approximation for ${n\choose m}$ [20], the condition (16) takes the form

[TABLE]

where ${\theta:=1-\rho}=\frac{n-m}{n}$ denotes the fractional sampling loss and where ${h(\theta)=-\theta\log(\theta)-(1-\theta)\log(1-\theta)}$ denotes the entropy function (computed with the natural logarithm). For a given R-RIP parameters ${\delta},{\mu}$, (17) is satisfied for a sufficiently small $\theta$ and ${\kappa}$ at the cost of incurring a large oversampling factor $\frac{n}{k}=\frac{1}{{\kappa}}$. When only $o(n)$ number of samples are missing, i.e., $m=n-o(n)$, thus, ${\theta\to 0}$, an oversampling of the order

[TABLE]

is sufficient to compensate for the missing labels.

III-D Guarantee for Signal Recovery

Let $\check{{\bf H}}=\check{{\bf y}}^{\sf T}\otimes\check{{\bf S}}\in{\cal H}$ and let $\check{{\bf x}}=\check{{\bf S}}{\bf B}\check{{\bf y}}+{\bf w}=\check{{\bf H}}{\mathbb{b}}+{\bf w}$ be the set of $m$ noisy linear measurements taken via ${\bf B}$, where ${\bf w}$ denotes the measurement noise. We define the measurement SNR by ${\mathsf{snr}}=\frac{\|\check{{\bf H}}{\mathbb{b}}\|^{2}}{\|{\bf w}\|^{2}}=\frac{\|\check{{\bf S}}{\bf B}\check{{\bf y}}\|^{2}}{\|{\bf w}\|^{2}}$ as before. We consider the following recovery problem

[TABLE]

Theorem 1

Let ${\bf B}$ be a random matrix that satisfies RIP over ${\cal H}$ with a parameter ${\epsilon}$ and R-RIP over ${\cal H}-{\cal H}$ with parameters ${\delta},{\mu}$. Let $\check{{\bf H}}$ and $\check{{\bf x}}$ be as before and let $\widehat{{\bf H}}$ be an estimate of $\check{{\bf H}}$ obtained from (19). Then, $\|\widehat{{\bf H}}-\check{{\bf H}}\|_{\sf F}\leq\chi\|\check{{\bf H}}\|_{\sf F}$, where \chi=\max\big{\{}{\mu}\frac{\sqrt{1+{\epsilon}}}{\sqrt{1-{\epsilon}}}(1+\frac{2}{\sqrt{{\mathsf{snr}}}}),\frac{2}{\sqrt{{\mathsf{snr}}(1-2{\delta})}}\big{\}}. $\square$

Proof:

Since $\widehat{{\bf H}}$ is the solution of (19) and $\check{{\bf H}}$ itself is feasible, we have that

[TABLE]

From triangle inequality, $\|(\widehat{{\bf H}}-\check{{\bf H}}){\mathbb{b}}\|\leq\|{\bf w}\|+\|{\bf w}\|=2\|{\bf w}\|$, thus, $\|\widehat{{\bf H}}{\mathbb{b}}\|\leq\|\check{{\bf H}}{\mathbb{b}}\|+2\|{\bf w}\|=\|\check{{\bf H}}{\mathbb{b}}\|(1+\frac{2}{\sqrt{{\mathsf{snr}}}})$. This yields

[TABLE]

where in $(i)$ and $(ii)$ we applied RIP over ${\cal H}$ to $\widehat{{\bf H}}$ and $\check{{\bf H}}$ respectively. From (22), it results that

[TABLE]

If $d_{\widehat{{\bf H}},\check{{\bf H}}}\leq{\mu}\max\{\|\widehat{{\bf H}}\|_{\sf F},\|\check{{\bf H}}\|_{\sf F}\}$, (23) yields

[TABLE]

Otherwise, if $d_{\widehat{{\bf H}},\check{{\bf H}}}\geq{\mu}\max\{\|\widehat{{\bf H}}\|_{\sf F},\|\check{{\bf H}}\|_{\sf F}\}$, R-RIP over ${\cal H}-{\cal H}$ holds for $\widehat{{\bf H}}$ and $\check{{\bf H}}$, where we obtain

[TABLE]

where in $(a)$ we applied R-RIP to $\check{{\bf H}}-\widehat{{\bf H}}$, and where in $(b)$ we used the RIP over ${\cal H}$ for $\check{{\bf H}}$. After simplification, we have

[TABLE]

where in $(i)$ we used ${\epsilon}\ll{\delta}$ (see Remark 1). Combining (24) and (25) completes the proof.     \sqcap$$\sqcup

Remark 2

Theorem 1 continues to hold for any other estimate $\widehat{{\bf H}}$ that merely satisfies the feasibility condition

[TABLE]

for any $\eta\geq 1$, where it yields $\|\check{{\bf H}}-\widehat{{\bf H}}\|_{\sf F}\leq\chi\|\check{{\bf H}}\|_{\sf F}$, with a \chi=\max\big{\{}{\mu}\frac{\sqrt{1+{\epsilon}}}{\sqrt{1-{\epsilon}}}(1+\frac{\eta+1}{\sqrt{{\mathsf{snr}}}}),\frac{\eta+1}{\sqrt{{\mathsf{snr}}(1-2{\delta})}}\big{\}}. $\lozenge$

Theorem 1 provides a recovery guarantee for any estimate $\widehat{{\bf H}}$ obtained from (19) or (26) when ${\bf B}$ satisfies R-RIP. An implicit requirement, however, is to have an efficient low-complexity algorithm to find such an estimate $\widehat{{\bf H}}$. In the next section, we design such a low-complexity algorithm via Alternating Minimization.

IV Recovery Algorithm

Let ${\check{{\bf y}}}$ be a signal and let ${\check{{\bf x}}}={\check{{\bf S}}}{\bf B}{\check{{\bf y}}}+{\bf w}$ be the vector of noisy unlabeled measurements where $\check{{\bf S}}\in{\cal S}$. We define the following cost function for the recovery of ${\check{{\bf y}}}$ as in (19):

[TABLE]

where ${{\cal Y}:={\mathbb{R}}^{k}}$. For a fixed ${\bf B}$ and i.i.d. Gaussian noise ${\bf w}$, the minimizer of $f({\bf S},{\bf y})$ denoted by $({\bf S}^{*},{\bf y}^{*})$ gives the maximum likelihood (ML) estimate of $({\check{{\bf S}}},{\check{{\bf y}}})$ and consequently the ML estimate ${\bf y}^{*}$ of the desired signal ${\check{{\bf y}}}$ [22]. Finding the ML estimate ${\bf y}^{*}$, however, requires a joint optimization over $({\bf S},{\bf y})$. This seems to require searching over all ${n\choose m}$ possible ${\bf S}\in{\cal S}$, which might be unfeasible for large $n$ and $m$. Here, instead of joint search over ${\bf S}$ and ${\bf y}$, we use an iterative alternating minimization with respect to ${\bf y}$ and ${\bf S}$ to reduce the complexity. Our proposed algorithm alternates between estimating a suitable selection matrix in ${\bf S}$ and a target signal ${\bf y}$ (Alternating Minimization) and resembles the Iterative Hard Thresholding (IHT) algorithm proposed for signal recovery in CS [23], which also alternates between finding a suitable support ${\bf L}$ and a suitable vector of nonzero coefficients ${\bf x}$ with the notation introduced here (see Fig. 2). This emphasizes further the underlying duality between UOS and CS. We call our algorithm $\mathtt{AltMin}$. We initialize $\mathtt{AltMin}$  with a random ${\bf S}^{1}\in{\cal S}$ and define $({\bf S}^{t},{\bf y}^{t})$ as the estimate of $({\check{{\bf S}}},\check{{\bf y}})$ obtained by $\mathtt{AltMin}$ at iteration $t=1,2,\dots$, via the following sequential projection operations

[TABLE]

where we used ${{\bf S}\mapsto{\bf y}}$ for ${\bf y}=\mathop{\rm arg\,min}_{{\bf y}\in{\cal Y}}f({\bf S},{\bf y})$ and ${\bf y}\mapsto{\bf S}$ for ${\bf S}=\mathop{\rm arg\,min}_{{\bf S}\in{\cal S}}f({\bf S},{\bf y})$ for the projection operators onto ${\cal Y}$ and ${\cal S}$ respectively.

IV-A Projection on ${\cal Y}$

For a fixed ${\bf S}^{t}$, finding ${\bf y}^{t}$ given by ${{\bf S}^{t}\mapsto{\bf y}^{t}}$ in (28) boils down to obtaining the least-squares solution of an over-complete set of linear equations (via the matrix ${\bf S}^{t}{\bf B}$). The optimal solution is given by ${\bf y}^{t}=({\bf S}^{t}{\bf B})^{\dagger}\check{{\bf x}}$, where † denotes the pseudo-inverse operator (for the tall matrix ${\bf S}^{t}{\bf B}$).

IV-B Projection on ${\cal S}$

Let ${\bf y}^{t}$ be the optimal solution obtained from (28) and set ${\bf z}^{t}={\bf B}{\bf y}^{t}$. Finding the optimal selection matrix ${\bf S}^{t+1}\in{\cal S}$ at (29) given by ${\bf y}^{t}\mapsto{\bf S}^{t+1}$ requires extracting a subvector of ${\bf z}^{t}$ of dimension $m$, while keeping the relative order of the components, that is closest to $m$-dim vector ${\check{{\bf x}}}$ in $l_{2}$ distance. We formulate this as a Dynamic Programming as follows. We define a 2D table of size $m\times n$ whose elements are labeled with $(r,c)\in[m]\times[n]$ and have the value ${\bf T}_{r,c}\in{\mathbb{R}}_{+}\cup\{\infty\}$ given by

[TABLE]

We initialize the diagonal elements of the table with {\bf T}_{i,i}=\|{\check{{\bf x}}}_{{\text{1$\,:\,$i}}}-{\bf z}^{t}_{{\text{1$\,:\,$i}}}\|^{2} since there is only one way to match the first $i$ elements of ${\bf z}^{t}$ with the $i$ elements in {\check{{\bf x}}}_{{\text{1$\,:\,$i}}}. We also initialize the elements in the first column of the table, i.e., ${\bf T}_{1,j}$ for $j\in[n]$, with ${\bf T}_{1,j}=\min_{j^{\prime}\in[j]}|{\check{{\bf x}}}_{1}-{\bf z}^{t}_{j^{\prime}}|^{2}$ since the single element ${\check{{\bf x}}}_{1}$ should be matched with the closest element in the subvector {\bf z}^{t}_{{\text{1$\,:\,$j}}} consisting of the first $j$ elements of ${\bf z}^{t}$.

To find the value of a generic element ${\bf T}_{r,c}$ in the table, we need to match {\check{{\bf x}}}_{{\text{1$\,:\,$r}}} with a suitable subsequence of {\bf z}^{t}_{{\text{1$\,:\,$c}}} of length $r$. In the optimal matching, the last component ${\check{{\bf x}}}_{r}$ is matched either with ${\bf z}^{t}_{c}$ or with ${\bf z}^{t}_{c^{\prime}}$ for some $c^{\prime}<c$. Thus, ${\bf T}_{r,c}$ can be computed from the already computed elements of the table as follows:

[TABLE]

With the already mentioned initialization and the recursion (30), we can complete the whole table by filling its $d$-th diagonal consisting of elements $\{{\bf T}_{1,d},{\bf T}_{2,d+1},{\bf T}_{3,d+2},\dots\}$ for $d=2,\dots,n$ one at a time. Overall, filling the whole table requires $O(mn)$ operations.

After filling the whole table, we can find the indices of those $m$ elements of ${\bf z}^{t}$ that are optimally matched to the elements of ${\check{{\bf x}}}$ as follows. We start from the element ${\bf T}_{m,n}$ located at the up-right corner of the table at location $(m,n)$. Note that by definition, $i\mapsto{\bf T}_{m,i}$, for $i\in[n]$, is a decreasing sequence of $i$ since by increasing $i$ the subvector {\bf z}^{t}_{{\text{1$\,:\,$i}}} becomes longer and provides more options to find a better subsequence thereof matching ${\check{{\bf x}}}$. The index of the last element in ${\bf z}^{t}$ that is matched to the last element ${\check{{\bf x}}}_{m}$ of ${\check{{\bf x}}}$ in the optimal matching is

[TABLE]

To find the next largest index $i_{m-1}$, we apply the same argument to the sub-table {\bf T}_{{\text{1$\,:\,$m-1}},{\text{1$\,:\,$i_{m}-1}}} and its up-right corner element located at $(m-1,i_{m}-1)$, where we obtain the following recursive formula for the remaining indices:

[TABLE]

where ${l\in[m-1]}$. Fig. 3 illustrates this for a vector ${\bf z}$ of dimension $7$ and a vector ${\bf x}$ of dimension $4$.

V Performance Analysis

In this section, we analyze the performance of $\mathtt{AltMin}$ under the assumption that the measurement matrix ${\bf B}$ satisfies R-RIP over ${\cal H}-{\cal H}$. Intuitively speaking, we show that $\mathtt{AltMin}$ can be seen as an iterative procedure for finding fixed points of an appropriate set-valued map. We use this to analyze the performance of $\mathtt{AltMin}$ rigorously.

V-A Set-valued Map

Let ${\cal X}$ and ${\cal Z}$ be two arbitrary sets and let ${\cal F}\subseteq{\cal X}\times{\cal Z}:=\{({\bf x},{\bf z}):{\bf x}\in{\cal X},{\bf z}\in{\cal Z}\}$. We define the set-valued or the multi-valued map corresponding to ${\cal F}$ by ${\cal F}:{\cal X}\to 2^{\cal Z}$, where $2^{\cal Z}$ denotes the power set of ${\cal Z}$, and where ${\cal F}$ assigns to each ${\bf x}\in{\cal X}$ a subset of ${\cal Z}$ given by

[TABLE]

Note that, for simplicity, we use the same symbol ${\cal F}$ both for a subset of ${\cal X}\times{\cal Z}$ and for the associated set-valued map. In the special case, where ${\cal F}({\bf x})$ is a singleton (has only a single element) for all ${\bf x}\in{\cal X}$, the set-valued map ${\cal F}$ corresponds to a single-valued function ${\cal F}:{\cal X}\to{\cal Z}$. We will use ${\cal F}:{\cal X}\rightrightarrows{\cal Z}$ for a set-valued map ${\cal F}$. We define the domain of ${\cal F}$ by ${\cal D}_{{\cal F}}:=\{{\bf x}\in{\cal X}:{\cal F}({\bf x})\neq\emptyset\}$ and the range of ${\cal F}$ by ${\cal R}_{{\cal F}}=\cup_{{\bf x}\in{\cal X}}{\cal F}({\bf x})$. In particular, ${\bf z}\in{\cal R}_{\cal F}$ if there is an ${\bf x}\in{\cal X}$ such that ${\bf z}\in{\cal F}({\bf x})$ or equivalently $({\bf x},{\bf z})\in{\cal F}$. We define the inverse of a set-valued map ${\cal F}$ as a mapping ${\cal F}^{-1}:{\cal Z}\rightrightarrows{\cal X}$, where ${\bf x}\in{\cal F}^{-1}({\bf z})$ if and only if ${\bf z}\in{\cal F}({\bf x})$. A simple example of a set-valued map is illustrated in Fig. 4.

We define the composition of two set-valued maps ${\cal F}:{\cal X}\rightrightarrows{\cal Z}$ and ${\cal G}:{\cal Z}\rightrightarrows{\cal W}$ by a set-valued map ${\cal J}={\cal G}\circ{\cal F}:{\cal X}\rightrightarrows{\cal W}$, where ${\bf w}\in{\cal J}({\bf x})$ if and only if there is a ${\bf z}\in{\cal Z}$ connecting ${\bf x}$ and ${\bf w}$, i.e., ${\bf z}\in{\cal F}({\bf x})$ and ${\bf w}\in{\cal G}({\bf z})$. This can also be written more compactly as ${\cal J}({\bf x})=\cup_{{\bf z}\in{\cal F}({\bf x})}{\cal G}({\bf z})$. We refer to [24, 25] for a more comprehensive introduction to the set-valued maps and additional related references.

V-B Decoding up to a Radius

As in Section IV, we denote the target signal by $\check{{\bf H}}=\check{{\bf y}}^{\sf T}\otimes\check{{\bf S}}$, the noisy samples by $\check{{\bf x}}=\check{{\bf H}}{\mathbb{b}}+{\bf w}=\check{{\bf S}}{\bf B}\check{{\bf y}}+{\bf w}$, and the sequence generated by $\mathtt{AltMin}$ by ${\bf S}^{1}\mapsto{\bf y}^{1}\mapsto{\bf S}^{2}\mapsto\cdots$. Let ${\bf H}^{i}={{\bf y}^{i}}^{\sf T}\otimes{\bf S}^{i}$ be the estimates produced by $\mathtt{AltMin}$ at the end of the $i$-th iteration. Since $\mathtt{AltMin}$ minimizes the cost function (27) at each iteration, we have that $\|{\bf H}^{i}{\mathbb{b}}-\check{{\bf x}}\|\leq\|\check{{\bf x}}\|$. Applying the triangle inequality and using the RIP over ${\cal H}$ for $\check{{\bf H}}$ and ${\bf H}^{i}$ (as in the proof of Theorem 1) yields

[TABLE]

where $(i)$ holds for a large ${\mathsf{snr}}$ and a small ${\epsilon}$. Let us consider the following condition for $({\bf H}^{i},\check{{\bf H}})$:

[TABLE]

where $d_{{\bf H}^{i},\check{{\bf H}}}=\|{\bf H}^{i}-\check{{\bf H}}\|_{{\sf F}}$ as in (14). Since ${\bf B}$ satisfies R-RIP over ${\cal H}-{\cal H}$, if (34) is violated then from R-RIP it results that

[TABLE]

where $(a)$ follows from (33). In words, the violation of (34) implies that the estimate ${\bf H}^{i}$ lies in a neighborhood of $\check{{\bf H}}$ of radius $2{\mu}\|\check{{\bf H}}\|_{\sf F}$, as illustrated in Fig. 5. For the analysis, we assume that (34) holds for all $({\bf H}^{i},\check{{\bf H}})$ and study $\mathtt{AltMin}$ under this condition. In particular, letting $\|{\bf H}^{i}-\check{{\bf H}}\|^{\text{Analysis}}_{\sf F}$ be the performance obtained under this additional condition, we have

[TABLE]

where \|{\bf H}^{i}-\check{{\bf H}}\|^{\text{{\mathtt{AltMin}}}}_{\sf F} denotes the true performance achieved by $\mathtt{AltMin}$ without the additional condition (34).

V-C Alternating Minimization as a Set-valued Map

We define the normalized cost function $g({\bf S},{\bf y}):=\frac{\sqrt{f({\bf S},{\bf y})}}{\sqrt{m}\|\check{{\bf y}}\|}$, where $f({\bf S},{\bf y})$ is as in (27). We introduce the normalized variables $\xi=\frac{\|{\bf y}\|^{2}}{\|\check{{\bf y}}\|^{2}}$ and $\vartheta=\frac{{\langle{\bf y},\check{{\bf y}}\rangle}}{\|\check{{\bf y}}\|^{2}}$, and ${\nu:=\nu_{{\bf S},\check{{\bf S}}}}$, where ${\nu_{{\bf S},\check{{\bf S}}}=\frac{\mathsf{tr}({\bf S}^{\sf T}\check{{\bf S}})}{m}}$ as in (15) denotes the fraction of common rows in ${\bf S}$ and $\check{{\bf S}}$. As ${\bf y}$ varies in ${{\cal Y}={\mathbb{R}}^{k}}$, the set of feasible $(\xi,\vartheta)$ is given by

[TABLE]

which is the area above the parabola $\xi={\vartheta}^{2}$ (see Fig. 6).

Applying the triangle and the reverse triangle inequality to $\sqrt{f({\bf S},{\bf y})}=\|{\bf S}{\bf B}{\bf y}-\check{{\bf S}}{\bf B}\check{{\bf y}}-{\bf w}\|$ and using (34) for ${{\bf H}={\bf y}^{\sf T}\otimes{\bf S}}$ and an arbitrary ${\check{{\bf H}}=\check{{\bf y}}^{\sf T}\otimes\check{{\bf S}}}$, we obtain the following upper and lower bounds for $g({\bf S},{\bf y})$:

[TABLE]

where $\zeta:=\frac{\|{\bf w}\|}{\|\check{{\bf y}}\|\sqrt{m}}\leq\frac{\|{\bf w}\|}{\|\check{{\bf S}}{\bf B}\check{{\bf y}}\|\sqrt{1-{\epsilon}}}=\frac{1}{\sqrt{{\mathsf{snr}}(1-{\epsilon})}}$, where ${\epsilon}$ is the parameter of RIP over ${\cal H}$, and where

[TABLE]

where in $(i)$ we used (14), and where $(\xi,{\vartheta})$ takes values in ${\cal E}$ as in (37). With this notation, we investigate the properties of the solution path ${\bf S}^{1}\mapsto{\bf y}^{1}\mapsto\cdots\mapsto{{\bf S}^{t}\mapsto{\bf y}^{t}}\mapsto{\bf S}^{t+1}\cdots$ generated by $\mathtt{AltMin}$.

V-C1 Projection on ${\cal Y}$

Consider an iteration $t$, where $\mathtt{AltMin}$  produces ${\bf y}^{t}$ as the projection ${{\bf S}^{t}\mapsto{\bf y}^{t}}$ of ${\bf S}^{t}$ onto ${\cal Y}$. From the upper bound in (38), we have that

[TABLE]

where ${\nu^{t}=\nu_{\check{{\bf S}},{\bf S}^{t}}}$. To derive (42), we used the fact that $\varphi(\xi,\vartheta,\nu^{t})$ in (41) is a decreasing function of $\xi$ for a fixed ${\vartheta}$ and $\nu^{t}\in[0,1]$, thus, the minimum of $\varphi(\xi,\vartheta,\nu^{t})$ over ${\cal E}$ is achieved at the boundary ${\{(\xi,{\vartheta}):\xi={\vartheta}^{2}\}}$ of ${\cal E}$ (see Fig. 6). After replacing $\xi={\vartheta}^{2}$ in $\varphi$, we need to minimize $\sqrt{1+{\vartheta}^{2}-2\nu^{t}{\vartheta}}$ with respect to ${\vartheta}\in{\mathbb{R}}$, where the minimum is achieved at ${\vartheta}=\nu^{t}$ and its value is given by $\sqrt{1-(\nu^{t})^{2}}$ as in (42). We could obtain this result also directly by minimizing (40) with respect to ${\bf y}$ (rather than $(\xi,{\vartheta})\in{\cal E}$), where we obtain the optimal solution ${\bf y}=\nu^{t}\check{{\bf y}}$ and (42) after replacement.

Let us denote by $\xi^{t}=\frac{\|{\bf y}^{t}\|^{2}}{\|\check{{\bf y}}\|^{2}}$ and ${\vartheta}^{t}=\frac{{\langle{\bf y}^{t},\check{{\bf y}}\rangle}}{\|\check{{\bf y}}\|^{2}}$ the variables corresponding to ${\bf y}^{t}$. Using the lower bound (39), we have $g({\bf S}^{t},{\bf y}^{t})\geq\sqrt{1-{\delta}}\,\varphi(\xi^{t},{\vartheta}^{t},\nu^{t})-\zeta$, which with (42) yields

[TABLE]

where $\varsigma:=\sqrt{\frac{1+{\delta}}{1-{\delta}}}$, where $\varrho:=\frac{2\zeta}{\sqrt{1-{\delta}}}$, and where we defined

[TABLE]

We write (V-C1) as ${(\xi^{t},{\vartheta}^{t})\in{\cal F}_{\text{\textxi,\texttheta}}(\nu^{t})}$, where ${\cal F}_{\text{\textxi,\texttheta}}(\nu)$ is the set-valued map defined by

[TABLE]

where ${\varpi}(\nu)$ is as in (44). Note that ${\cal F}_{\text{\textxi,\texttheta}}(\nu)$ assigns to any $\nu\in[0,1]$ a subset of feasible points ${\cal E}$ that can be potentially produced by $\mathtt{AltMin}$. Before we proceed we make the following crucial condition.

Condition 1

There is a $\bar{\nu}\in(0,1]$ such that for $\nu\in({\bar{\nu},1]}$

[TABLE]

where ${\varsigma}=\sqrt{\frac{1+{\delta}}{1-{\delta}}}$ as before. Such a $\bar{\nu}$ exists if

[TABLE]

where $(i)$ follows from $\zeta\leq\frac{1}{\sqrt{{\mathsf{snr}}(1-{\epsilon})}}$. We define $\nu_{0}$ as the smallest such $\bar{\nu}$, which from (46) is given by

[TABLE]

Condition 1 is satisfied for a sufficiently large ${\mathsf{snr}}$ and a sufficiently small ${\delta}\in(0,1)$, and is easy to meet in practice.

Proposition 3

Let $\nu_{0}$ be as in Condition 1 and let ${\nu\in(\nu_{0},1]}$. Then, for any $(\xi,{\vartheta})\in{\cal F}_{\text{\textxi,\texttheta}}(\nu)$, ${\vartheta}\geq 0$. $\square$

Proof:

From (44), it is seen that ${\varpi}(\nu)$ is a decreasing function of $\nu$, thus, from ${\nu\in(\nu_{0},1]}$ it immediately results that ${{\varpi}(\nu)<{\varpi}(\nu_{0})=1}$. From (45), this implies that $\sqrt{1+\xi-2{\vartheta}\nu}<1$. Since $\nu>0$ and $\xi\geq 0$, this is satisfied only if ${\vartheta}\geq 0$. This completes the proof.     \sqcap$$\sqcup

V-C2 Projection on ${\cal S}$

Let $(\xi^{t},{\vartheta}^{t})$ be the output of $\mathtt{AltMin}$  corresponding to ${\bf y}^{t}$. The next step is the projection ${\bf y}^{t}\mapsto{\bf S}^{t+1}$ onto ${\cal S}$, where we have

[TABLE]

where in $(i)$ we used the upper bound (38), where in $(ii)$ we used (41), where we defined $(a)_{+}=\max\{a,0\}$ as the positive part of $a\in{\mathbb{R}}$, and where in $(iii)$ we used the fact that when ${\vartheta}^{t}>0$ the minimum is achieved at $\nu=1$ whereas for ${\vartheta}^{t}\leq 0$ it is achieved at $\nu=0$. Using (52) and the lower bound (39) at $({\bf S}^{t+1},{\bf y}^{t})$, we obtain

[TABLE]

Although (53) does not give the exact value of $\nu^{t+1}$ generated by $\mathtt{AltMin}$, it implies that ${\nu^{t+1}\in{\cal F}_{\text{\textnu}}(\xi^{t},{\vartheta}^{t})}$, where ${\cal F}_{\text{\textnu}}(\xi,{\vartheta})$ is the set-valued map defined by

[TABLE]

V-C3 Full Iteration

Combining the two steps of projection in $\mathtt{AltMin}$  at iteration $t$, we obtain that starting from $\nu^{t}\in(\nu_{0},1]$, the algorithm returns $\nu^{t+1}\in{\cal F}(\nu^{t})$, where ${\cal F}$ is the composition of the set-valued maps ${\cal F}_{\text{\textxi,\texttheta}}$ and ${\cal F}_{\text{\textnu}}$ given by

[TABLE]

In words, ${\cal F}(\nu)$ is the set of all $\nu^{\prime}$ in $[0,1]$ that are potentially reachable from $\nu$ at a single iteration of $\mathtt{AltMin}$.

V-D Structure of ${\cal F}(\nu)$ and its lower envelope

We define the lower envelope of the set-valued map ${\cal F}$ as the single-valued function

[TABLE]

where ${\sf F}_{0}(\nu)\in[0,1]$ for all $\nu\in[0,1]$. Let $\nu_{0}$ be as in Condition 1. We will study the behavior of ${\sf F}_{0}$ in two regions ${\nu\in(\nu_{0},1]}$ and $\nu\in[0,\nu_{0}]$ separately.

V-D1 First case: ${\nu\in(\nu_{0},1]}$

For this range of $\nu$, from Proposition 3 it results that for any $(\xi,{\vartheta})\in{\cal F}_{\text{\textxi,\texttheta}}(\nu)$, we have ${\vartheta}\geq 0$, thus, $({\vartheta})_{+}={\vartheta}$. For each such a $(\xi,{\vartheta})\in{\cal F}_{\text{\textxi,\texttheta}}(\nu)$ with ${\vartheta}\geq 0$, (54) yields a lower bound on the possible values of $\nu$ as follows

[TABLE]

However, to find a lower envelope of ${\cal F}(\nu)$, from (V-C3) we need to consider all those $(\xi,{\vartheta})\in{\cal F}_{\text{\textxi,\texttheta}}(\nu)$. Thus, the lower envelope ${\sf F}_{0}(\nu)$ is given by ${\sf F}_{0}(\nu)=\max\{{\sf F}(\nu),0\}$, where

[TABLE]

The region ${\cal F}_{\text{\textxi,\texttheta}}(\nu)$ is illustrated in Fig. 6 for a generic value of $\nu\in(\nu_{0},1]$. It is seen that ${\cal F}_{\text{\textxi,\texttheta}}(\nu)$ for any fixed $\nu\in(\nu_{0},1]$ consist of all $(\xi,{\vartheta})$ lying above the parabola $\xi={\vartheta}^{2}$ and below the line $\xi=2\nu{\vartheta}+{\varpi}^{2}-1$ in the range $[{\vartheta}_{\min},{\vartheta}_{\max}]$, where ${\vartheta}_{\min}$ and ${\vartheta}_{\max}$ denote the ${\vartheta}$ coordinate of the two intersection points of the line and the parabola and correspond to the roots of the quadratic equation obtained by setting $\xi={\vartheta}^{2}$ in the equation of the line $\xi=2\nu{\vartheta}+{\varpi}^{2}-1$ given by:

[TABLE]

Solving (59) we have

[TABLE]

where ${\varpi}(\nu)$ is as in (44). Note that the roots are well-defined (real-valued) for all $\nu\in[0,1]$ and in particular for ${\nu\in(\nu_{0},1]}$ because the discriminant

[TABLE]

since ${\varsigma}=\sqrt{\frac{1+{\delta}}{1-{\delta}}}>1$ and $1-\nu^{2}\geq 0$ for $\nu\in[0,1]$. From the root relation for a quadratic function, we also have that

[TABLE]

for ${\nu\in(\nu_{0},1]}$ as ${\varpi}(\nu)\in[0,1)$ for ${\nu\in(\nu_{0},1]}$ from Condition 1. As ${{{\vartheta}_{\max}>0}}$, (62) also yields ${{\vartheta}_{\min}>0}$, thus, both roots are positive for ${\nu\in(\nu_{0},1]}$. Moreover, as ${{\vartheta}_{\min}{\vartheta}_{\max}\leq 1}$, the minimum root ${\vartheta}_{\min}$ is always less than $1$.

V-D2 Second case: $\nu\in[0,\nu_{0}]$

In this region, we see from Condition 1 and (44) that ${\varpi}(\nu)\geq 1$, thus, $(\xi,{\vartheta})=(0,0)\in{\cal F}_{\text{\textxi,\texttheta}}(\nu)$ from (45). We can see from (54) and the fact that ${\varsigma}\in(1,\infty)$ and ${\varrho}\in[0,1)$ from (47) that for $(\xi,{\vartheta})=(0,0)$, ${\cal F}_{\text{\textnu}}(\xi,{\vartheta})=[0,1]$. From the composition of set-valued maps in (V-C3), it results that ${\cal F}(\nu)=[0,1]$ for $\nu\in[0,\nu_{0}]$, and in particular, ${\sf F}_{0}(\nu)=0$.

V-E Explicit formula for ${\sf F}_{0}(\nu)$

Suppose ${\nu\in(\nu_{0},1]}$ and write the minimization (58) as

[TABLE]

where ${\sf G}:a\to a-({\varsigma}\sqrt{a}+{\varrho})^{2}$ is a function defined over the region ${\cal D}_{\sf G}:=\{1+\xi-2{\vartheta}:(\xi,{\vartheta})\in{\cal F}_{\text{\textxi,\texttheta}}(\nu)\}\subset{\mathbb{R}}$. Note that as $\xi\geq{\vartheta}^{2}$ over the region ${\cal F}_{\text{\textxi,\texttheta}}(\nu)$, we have that

[TABLE]

Hence, ${\cal D}_{\sf G}\subseteq{\mathbb{R}}_{+}$ and ${\sf G}$ is well-defined over the whole ${\cal D}_{\sf G}$. Moreover, since ${\cal D}_{\sf G}$ is the image of the closed connected (and bounded) set ${\cal F}_{\text{\textxi,\texttheta}}(\nu)$ (see Fig. 6) under the continuous affine map $(\xi,{\vartheta})\mapsto 1+\xi-2{\vartheta}$, the set ${\cal D}_{\sf G}$ must be closed and connected [26], thus, a closed interval of ${\mathbb{R}}_{+}$.

Lemma 2

${\sf G}$* is a negative decreasing function over ${\mathbb{R}}_{+}$ and in particular ${\cal D}_{G}$. $\square$*

Proof:

By expanding the expression for ${\sf G}$, we have that ${{\sf G}(a)=(1-{\varsigma}^{2})a-2{\varsigma}{\varrho}\sqrt{a}-{\varrho}^{2}}$. As ${\varsigma}>1$ and $a\in{\mathbb{R}}_{+}$ (especially when $a\in{\cal D}_{G}$), it results that ${\sf G}(a)\leq 0$ for all $a\in{\mathbb{R}}_{+}$ ($a\in{\cal D}_{G}$). Also, $\frac{d}{da}{\sf G}(a)={\sf G}^{\prime}(a)=(1-{\varsigma}^{2})-\frac{{\varsigma}{\varrho}}{\sqrt{a}}<0$ for $a\in{\mathbb{R}}_{+}$ ($a\in{\cal D}_{G}$). Thus, ${\sf G}$ is a decreasing function.     \sqcap$$\sqcup

Lemma 3

Let ${\nu\in(\nu_{0},1]}$. Then, the minimum in (58) over ${\cal F}_{\text{\textxi,\texttheta}}(\nu)$ is achieved at corner point $(\xi,{\vartheta})=({\vartheta}_{\min},{\vartheta}_{\min}^{2})$. $\square$

Proof:

We first fix a ${\vartheta}>0$ and consider ${\sf G}(1+\xi-2{\vartheta})$ as a function of $\xi$ in the range ${\cal D}_{\text{\textxi}}({\vartheta}):=\{\xi:(\xi,{\vartheta})\in{\cal F}_{\text{\textxi,\texttheta}}(\nu)\}$, which is a closed interval (see Fig. 6). As $1+\xi-2{\vartheta}$ for a fixed ${\vartheta}$ is an affine increasing function of $\xi$, from Lemma 2 it results that ${\sf G}(1+\xi-2{\vartheta})$ is a decreasing function of $\xi$ in ${\cal D}_{\text{\textxi}}({\vartheta})$, thus, its minimum over ${\cal D}_{\text{\textxi}}({\vartheta})$ is achieved at $\xi_{\max}({\vartheta})=\max{\cal D}_{\text{\textxi}}({\vartheta})$. As a result, $(\xi_{\max}({\vartheta}),{\vartheta})$ lies on the boundary line $1+\xi-2\nu{\vartheta}={\varpi}^{2}$, thus, $\xi_{\max}({\vartheta})={\varpi}^{2}-1+2\nu{\vartheta}$. Replacing $\xi_{\max}({\vartheta})$ for $\xi$, the minimum in (63) is given by

[TABLE]

Since ${\varpi}^{2}+2{\vartheta}(\nu-1)$ is an affine decreasing function of ${\vartheta}$ for ${\nu\in(\nu_{0},1]}$, we have that

[TABLE]

where we used the fact that $\nu-1\leq 0$ and that ${\sf G}^{\prime}(.)\leq 0$ from Lemma 2. This implies that

[TABLE]

where we used (66), the fact that ${\vartheta}>0$, and that ${\sf G}({\varpi}^{2}+2{\vartheta}(\nu-1))\leq 0$ from Lemma 2. As a result, $\frac{{\sf G}({\varpi}^{2}+2{\vartheta}(\nu-1))}{2{\vartheta}}$ is an increasing function of ${\vartheta}$ for ${\vartheta}\in[{\vartheta}_{\min},{\vartheta}_{\max}]$ and achieves its minimum at ${\vartheta}_{\min}$. Overall, this implies that for ${\nu\in(\nu_{0},1]}$ the minimum in (63) is achieved at the corner point $({\vartheta}_{\min}^{2},{\vartheta}_{\min})\in{\cal F}_{\text{\textxi,\texttheta}}(\nu)$. This completes the proof.     \sqcap$$\sqcup

Applying Lemma 3, we obtain the following explicit formula for ${\sf F}(\nu)$ for $\nu\in(\nu_{0},1]$:

[TABLE]

where ${\vartheta}_{\min}$ is as in (60) and where we used the fact that ${\sqrt{1+{\vartheta}_{\min}^{2}-2{\vartheta}_{\min}}=1-{\vartheta}_{\min}}$ as ${{\vartheta}_{\min}\in[0,1]}$ for ${\nu\in(\nu_{0},1]}$. Thus, we obtain

[TABLE]

V-F Properties of ${\sf F}_{0}(\nu)$

We need some preliminary results first.

Lemma 4

Let ${\vartheta}_{\min}(\nu)$ be as in (60). Then, ${\vartheta}_{\min}(\nu)$ is an increasing function of $\nu$ for ${\nu\in(\nu_{0},1]}$. Moreover, ${\vartheta}_{\min}(\nu_{0})=0$ and ${\vartheta}_{\min}(1)=1-{\varrho}$. $\square$

Proof:

From Condition 1, ${\varpi}(\nu_{0})=1$. Thus, ${\vartheta}_{\min}(\nu_{0})=\nu_{0}-\sqrt{\nu_{0}^{2}-1+{\varpi}(\nu_{0})^{2}}=0$. Also, replacing $\nu=1$ and using the fact that ${{\varpi}(\nu)={\varsigma}\sqrt{1-\nu^{2}}+{\varrho}={\varrho}}$ at $\nu=1$ yields ${{\vartheta}_{\min}(1)=1-{\varrho}}$. Replacing ${\varpi}(\nu)={\varsigma}\sqrt{1-\nu^{2}}+{\varrho}$ in ${\vartheta}_{\min}(\nu)$, we can write

[TABLE]

where ${{\sf G}:a\to a-({\varsigma}\sqrt{a}+{\varrho})^{2}}$ is as before. From Lemma 2, it results that $-{\sf G}(.)$ is a positive increasing function over ${\mathbb{R}}_{+}$ and in particular $[0,1]$. Thus, $-\sqrt{-{\sf G}(.)}$ is a decreasing and $-\sqrt{-{\sf G}(1-\nu^{2})}$ is an increasing function of $\nu$. Hence, ${\vartheta}_{\min}(\nu)=\nu-\sqrt{-{\sf G}(1-\nu^{2})}$ is an increasing function.     \sqcap$$\sqcup

We first consider (67) and write ${\sf F}(\nu)={\sf M}\circ{\vartheta}_{\min}(\nu)$, where ${{\sf M}:{\mathbb{R}}_{+}\to{\mathbb{R}}}$ is defined by

[TABLE]

Lemma 5

${\sf M}$* is an increasing function over $(0,1]$. $\square$*

Proof:

Taking the derivative of ${\sf M}(a)$, we have

[TABLE]

Note that ${\sf M}^{\prime}(a)$ has a positive root at $a_{0}=\sqrt{\frac{({\varsigma}+{\varrho})^{2}-1}{{\varsigma}^{2}-1}}\geq 1$ as ${\varsigma}>1$ and ${\varrho}\geq 0$. Moreover, ${\sf M}^{\prime}(0^{+})=+\infty$ and ${\sf M}^{\prime}(a)\geq 0$ for $a\in(0,a_{0}]$ and in particular for all $a\in(0,1]$. This implies that ${\sf M}(a)$ is an increasing function of $a$ for $a\in(0,1]$.     \sqcap$$\sqcup

Proposition 4

${\sf F}_{0}(\nu)$* is an increasing function over $(\nu_{0},1]$.*

Proof:

We first prove that ${\sf F}(\nu)$ is an increasing function of $\nu$. Note that ${\sf F}(\nu)={\sf M}\circ{\vartheta}_{\min}(\nu)$, where for ${\nu\in(\nu_{0},1]}$, ${\vartheta}_{\min}(\nu)\in(0,1]$, over which ${\sf M}$ is an increasing function from Lemma 5. Moreover, ${\vartheta}_{\min}(\nu)$ is also an increasing function of $\nu$ for ${\nu\in(\nu_{0},1]}$ from Lemma 4. This implies the composition function ${\sf F}(\nu)={\sf M}\circ{\vartheta}_{\min}(\nu)$ is an increasing function of $\nu$ over $\nu\in(\nu_{0},1]$. From (70) and the fact that $a\to(a)_{+}$ is an increasing function of $a$, we obtain that ${\sf F}_{0}(\nu)$ is an increasing function of $\nu$ over $[0,1]$. This completes the proof.     \sqcap$$\sqcup

Proposition 4 implies that ${\sf F}_{0}(\nu)\in[0,({\sf F}_{\max})_{+}]$, where

[TABLE]

where in $(a)$ we used the fact that ${\vartheta}_{\min}(1)=1-{\varrho}$ from Lemma 4. From (76), it is seen that ${\sf F}_{\max}\leq 1$ as ${\varsigma}>1$ and ${\varrho}\in[0,1)$. Moreover, ${\sf F}_{\max}>0$ when ${\varrho}$ is sufficiently small (large ${\mathsf{snr}}$) and ${\varsigma}$ is not far from $1$ (small R-RIP parameter ${\delta}$), and $\lim_{{({\varsigma},{\varrho})\to(1,0)}}{\sf F}_{\max}({\varsigma},{\varrho})=1$. Since ${\sf F}(\nu)$ is an increasing function of $\nu$, provided that ${\sf F}_{\max}({\varsigma},{\varrho})>0$ there would exist a ${\nu_{1}({\varsigma},{\varrho})\in(\nu_{0},1]}$ such that ${\sf F}(\nu)>0$, thus, ${\sf F}_{0}(\nu)>0$ for $\nu\in(\nu_{1},1]$. A direct calculation by setting ${\sf F}(\nu_{1})={\sf M}\circ{\vartheta}_{\min}(\nu_{1})=0$ yields

[TABLE]

which using ${\vartheta}_{\min}(\nu)=\nu-\sqrt{\nu^{2}-1+({\varsigma}\sqrt{1-\nu^{2}}+{\varrho})^{2}}$ can be solved to obtain the value of $\nu_{1}({\varsigma},{\varrho})$ explicitly.

V-G Fixed points of $\sf{\sf F}_{0}(\nu)$

The crucial step in our analysis is based on the fixed points of ${\sf F}_{0}(\nu)$ defined by $\{\nu\in[\nu_{1},1]:{\sf F}_{0}(\nu)=\nu\}$. Note that the possible fixed points of ${\sf F}_{0}(\nu)$ in the range $\nu\in[0,1]$ are also fixed points of ${\sf F}(\nu)$. Form (67), these fixed points, provided that they exist, should satisfy the following equation

[TABLE]

A straightforward calculation yields

[TABLE]

Using the identity ${{\sf Q}({\vartheta}_{\min})={\vartheta}_{\min}^{2}-2{\vartheta}_{\min}\nu+1-{\varpi}^{2}=0}$ as in (59) and doing some simplification results in

[TABLE]

which can be simplified to ${\vartheta}_{\min}(\nu)=1-\sqrt{1-\nu^{2}}$ and written, using (44) and (60), more explicitly in terms of $\nu$ as follows

[TABLE]

By introducing the auxiliary variable ${\nu=\sin(\alpha)}$ for $\alpha\in[0,\frac{\pi}{2}]$, we can write this equivalently as

[TABLE]

We will consider the noiseless (${\varrho}=0$) and the noisy (${\varrho}\in(0,1)$) cases separately.

V-G1 Noiseless Case (${\varrho}=0$)

In the noiseless case, using the fact that ${\cos(\alpha)+\sin(\alpha)\geq 1}$ for $\alpha\in[0,\frac{\pi}{2}]$, (81) yields

[TABLE]

By introducing $\tan(\beta)=1-\sqrt{{\varsigma}^{2}-1}$ for $\beta\in(-\frac{\pi}{2},\frac{\pi}{2})$, we can write (82) as $\sin(\alpha+\beta)=\cos(\beta)$. One of the solutions is $\alpha_{\max}=\frac{\pi}{2}$ and corresponds to the largest fixed point $\nu_{\max}=\sin(\frac{\pi}{2})=1$. The second solution is given by $\alpha+\beta=\frac{\pi}{2}-\beta$, or $\alpha_{\min}=\frac{\pi}{2}-2\beta$ and lies in the allowed range $[0,\frac{\pi}{2}]$ provided that $\beta\in[0,\frac{\pi}{4}]$ or equivalently $\tan(\beta)\in[0,1]$. This restricts the range of permitted ${\varsigma}$ to ${\varsigma}\in(1,\sqrt{2}]$. From ${\varsigma}=\sqrt{\frac{1+{\delta}}{1-{\delta}}}$, this yields the bound ${\delta}\in(0,\frac{1}{3})$ on R-RIP parameter ${\delta}$. Moreover,

[TABLE]

It is seen that $\nu_{\min}({\varsigma})\in(0,1)$ and $\nu_{\min}({\varsigma})\to 0$ as ${\varsigma}\to 1$.

V-G2 Noisy Case (${\varrho}\in(0,1)$)

A full analysis of the possible roots of (81) is more involved in the noisy case. We first write (81) after factorizing the first term as \big{(}({\varsigma}-1)\cos(\alpha)+{\varrho}\big{)}\big{(}({\varsigma}-1)\cos(\alpha)+{\varrho}\big{)}=\big{(}\cos(\alpha)+\sin(\alpha)-1\big{)}^{2}. Since ${\varsigma}\in(1,\infty)$ and $\cos(\alpha)\in[0,1]$ for $\alpha\in[0,\frac{\pi}{2}]$, we have (({\varsigma}\pm 1)\cos(\alpha)+{\varrho}\big{)}>0. Thus, we can write (81) equivalently as $\upsilon(\alpha)=0$ where

[TABLE]

Proposition 5

Let $\upsilon(\alpha)$ be as in (84). For any ${\varsigma}\in(1,\infty)$ and ${\varrho}\in(0,1)$, $\upsilon(\alpha)$ is a convex function of $\alpha$ for $\alpha\in(0,\frac{\pi}{2})$, with $\upsilon(0^{+})=\upsilon(\frac{\pi}{2}^{-})=+\infty$. Moreover, $\upsilon(\alpha)$ has at most two roots in $(0,\frac{\pi}{2})$. $\square$

Proof:

Proof in Appendix A-B.     \sqcap$$\sqcup

Fig. 7 illustrates $\upsilon(\alpha)$ for ${({\varsigma},{\varrho})=(1.03,0.06)}$, where it is seen that it has two roots. We will consider only those $({\varsigma},{\varrho})$ in

[TABLE]

and will denote the roots by $\alpha_{\min}$ and $\alpha_{\max}$ or by $\alpha_{\min}({\varsigma},{\varrho})$ and $\alpha_{\max}({\varsigma},{\varrho})$ to emphasize the implicit dependence on $({\varsigma},{\varrho})$. We also denote the corresponding fixed points of ${\sf F}_{0}$ by $\nu_{\min}=\sin(\alpha_{\min})$ and $\nu_{\max}=\sin(\alpha_{\max})$. In view of Proposition 5, a simple sufficient condition for $({\varsigma},{\varrho})\in{\cal D}_{\text{\textvarsigma,\textrho}}$ is obtained by setting $\upsilon(\frac{\pi}{4})<0$, which can be simplified to

[TABLE]

which is satisfied for ${\varsigma}\to 1$ (small R-RIP parameter ${\delta}$) and ${\varrho}\to 0$ (large ${\mathsf{snr}}$). We have the following useful result.

Proposition 6

Let $\upsilon(\alpha)$ be as in (84) and let $\alpha_{\min}({\varsigma},{\varrho})$ and $\alpha_{\max}({\varsigma},{\varrho})$ be the roots of $\upsilon(\alpha)$ for $({\varsigma},{\varrho})\in{\cal D}_{\text{\textvarsigma,\textrho}}$. Then, $(i)$ ${\cal D}_{\text{\textvarsigma,\textrho}}$ is an open subset of $(1,\infty)\times(0,1)$.

$(ii)$* $\alpha_{\min}$ and $\alpha_{\max}$ are differentiable functions in ${\cal D}_{\text{\textvarsigma,\textrho}}$.*

$(iii)$* $\alpha_{\min}$ ($\alpha_{\max}$) is an increasing (decreasing) function of $({\varsigma},{\varrho})$, i.e., $\alpha_{\min}({\varsigma}^{\prime},{\varrho}^{\prime})\geq\alpha_{\min}({\varsigma},{\varrho})$ for ${\varsigma}^{\prime}\geq{\varsigma}$ and ${\varrho}^{\prime}\geq{\varrho}$.*

$(iv)$* $\alpha_{\min}({\varsigma},{\varrho})\to 0$ and $\alpha_{\max}({\varsigma},{\varrho})\to\frac{\pi}{2}$ as ${({\varsigma},{\varrho})\to(1,0)}$.*

Proof:

Proof in Appendix A-C.     \sqcap$$\sqcup

V-H Evolution Equation for Alternating Minimization

Let $\nu^{t}=\nu_{\check{{\bf S}},{\bf S}^{t}}$ be the similarity factor of the solution ${\bf S}^{t}$ to the target $\check{{\bf S}}$ obtained at iteration $t$ of $\mathtt{AltMin}$. Although we cannot control the value of $\nu^{t}$, we can guarantee that $\nu^{t+1}\in{\cal F}(\nu^{t})$. In particular, this implies that $\nu^{t+1}\geq{\sf F}_{0}(\nu^{t})$, where ${\sf F}_{0}(\nu)$ denotes the lower envelope of ${\cal F}(\nu)$ as before. Since ${\sf F}_{0}$ is an increasing function, by repeated application of ${\sf F}_{0}$, it results that $\nu^{t+1}\geq{\sf F}_{0}^{(t)}(\nu_{\check{{\bf S}},{\bf S}^{1}})$, where

[TABLE]

denotes the $t$-th order composition of ${\sf F}_{0}$.

Proposition 7** (Evolution Equation)**

Let $({\varsigma},{\varrho})\in{\cal D}_{\text{\textvarsigma,\textrho}}$ and let $\nu_{\min}({\varsigma},{\varrho})$ and $\nu_{\max}({\varsigma},{\varrho})$ be the two fixed points of ${\sf F}_{0}$. Let ${\bf S}^{1}$ be an initialization with $\nu^{1}=\nu_{{\bf S}^{1},\check{{\bf S}}}>\nu_{\min}$. Then, $\mathop{\rm lim\ inf}_{t\to\infty}\nu^{t}\geq\nu_{\max}$. $\square$

Proof:

Applying induction and using $\nu^{1}>\nu_{\min}$ and ${\sf F}_{0}(\nu_{\min})=\nu_{\min}$, we can show that $\nu^{t}\geq\nu_{1}>\nu_{\min}$. Taking the limit and denoting by $\nu_{\infty}=\mathop{\rm lim\ inf}_{t\to\infty}\nu^{t}$, we obtain $\nu_{\infty}\geq\nu^{1}>\nu_{\min}$, thus, $\nu_{\infty}>\nu_{\min}$. Moreover, we have

[TABLE]

where in $(i)$ we used the fact that ${\sf F}_{0}$ is an increasing function. Since $\nu^{\infty}>\nu_{\min}$, the only region where $\nu^{\infty}\geq{\sf F}_{0}(\nu^{\infty})$ is satisfied is when $\nu^{\infty}\in[\nu_{\max},1]$. This completes the proof.     \sqcap$$\sqcup

Fig. 8 illustrates the set-valued map ${\cal F}$ and its lower envelope ${\sf F}_{0}$ for ${\varsigma}=1.03$ and ${\varrho}=0.06$. It also illustrates the evolution of $\nu^{t}$ produced by $\mathtt{AltMin}$ across consecutive iterations. Using Proposition 7, we are finally in a position to analyze the performance of $\mathtt{AltMin}$.

Theorem 6** (Noiseless case)**

Let ${\bf B}$ be a matrix with an R-RIP parameter ${\delta}\in(0,\frac{1}{3})$, let ${\varsigma}=\sqrt{\frac{1+{\delta}}{1-{\delta}}}$, and let $\nu_{\min}({\varsigma})$ be as in (83), where $\nu_{\min}({\varsigma})\to 0$ as ${{\delta}\to 0}$ (${{\varsigma}\to 1}$). Let $\check{{\bf x}}=\check{{\bf S}}{\bf B}\check{{\bf y}}$ be the unlabeled samples taken from the signal $\check{{\bf y}}$ and let ${{\bf S}^{1}\mapsto{\bf y}^{1}\mapsto\cdots}$ be the sequence generated by $\mathtt{AltMin}$ for the input $\check{{\bf x}}$ starting from an initialization ${\bf S}^{1}$ with ${\nu_{{\bf S}^{1},\check{{\bf S}}}>\nu_{\min}({\varsigma})}$. Then, $\lim_{t\to\infty}\|\check{{\bf y}}-{\bf y}^{t}\|=0$. $\square$

Proof:

For the noiseless case, we have $\nu_{\max}=1$. Thus, under the stated conditions, starting from $\nu^{1}=\nu_{{\bf S}^{1},\check{{\bf S}}}$, we have that $\mathop{\rm lim\ inf}_{t\to\infty}\nu^{t}\geq\nu_{\max}=1$ from Proposition 7, which implies that $\lim_{t\to\infty}\nu^{t}=1$. From (V-C1) (${\varrho}=0$), this yields

[TABLE]

Using $\varphi(\xi^{t},{\vartheta}^{t},\nu^{t})=\sqrt{1+{\xi^{t}}-2{\vartheta}^{t}\nu^{t}}$ from (41) and the fact that $\xi^{t}\geq({\vartheta}^{t})^{2}$, (89) implies that ${\vartheta}^{t}\to 1$ and $\xi^{t}\to 1$ as $\nu^{t}\to 1$. This yields $\frac{\|{\bf y}^{t}\|^{2}}{\|\check{{\bf y}}\|^{2}}\to 1$ and $\frac{{\langle{\bf y}^{t},\check{{\bf y}}\rangle}}{\|\check{{\bf y}}\|^{2}}\to 1$, and gives the desired result $\lim_{t\to\infty}\|{\bf y}^{t}-\check{{\bf y}}\|=0$.     \sqcap$$\sqcup

We have also the following results for the noisy case.

Theorem 7** (Decoding up to the Noise Radius)**

Let ${\bf B}$ be a matrix as in Theorem 6. Suppose that $({\varsigma},{\varrho})$ are in the region ${\cal D}_{\text{\textvarsigma,\textrho}}$ where ${\sf F}_{0}(\nu)$ has two fixed points. Also, let $\check{{\bf x}}=\check{{\bf S}}{\bf B}\check{{\bf y}}+{\bf w}$ be the noisy unlabeled samples taken from the signal $\check{{\bf y}}$. Assume that ${{\bf S}^{1}\mapsto{\bf y}^{1}\mapsto\cdots}$ is the sequence generated by $\mathtt{AltMin}$ for the input $\check{{\bf x}}$ starting from an initialization ${\bf S}^{1}$ with $\nu_{\check{{\bf S}},{\bf S}^{1}}\geq\nu^{\prime}>0$. Then, $\mathop{\rm lim\ sup}_{{({\varsigma},{\varrho})\to(1,0)}}\mathop{\rm lim\ sup}_{t\to\infty}\frac{\|{\bf S}^{t}{\bf B}{\bf y}^{t}-\check{{\bf x}}\|}{\|{\bf w}\|}\leq 3$. $\square$

Proof:

From proposition 6, we have that $\alpha_{\min}({\varsigma},{\varrho})\to 0$, thus, ${\nu_{\min}({\varsigma},{\varrho})\to 0}$, as ${({\varsigma},{\varrho})\to(1,0)}$. Hence, there is a sufficiently small neighborhood of $(1,0)$ such that $\nu_{\min}<\nu^{\prime}$. Fixing a $({\varsigma},{\varrho})$ in this neighborhood, (V-C1) yields

[TABLE]

where in $(i)$ we used $\mathop{\rm lim\ inf}_{t\to\infty}\nu^{t}\geq\nu_{\max}$ from Proposition 7, and where in $(ii)$ we replaced ${\nu_{\max}=\sin(\alpha_{\max})}$. Using $g({\bf S}^{t},{\bf y}^{t})=\frac{\sqrt{f({\bf S}^{t},{\bf y}^{t})}}{\sqrt{m}\|\check{{\bf y}}\|}$ and $f({\bf S}^{t},{\bf y}^{t})=\|{\bf S}^{t}{\bf B}{\bf y}^{t}-\check{{\bf x}}\|^{2}$, we have

[TABLE]

where in $(a)$ we used (38) and replaced $\zeta=\frac{\|{\bf w}\|}{\|\check{{\bf y}}\|\sqrt{m}}$, and where in $(b)$ we used (92). Using ${\varrho}=\frac{2\zeta}{\sqrt{1-{\delta}}}$, taking the limit as ${({\varsigma},{\varrho})\to(1,0)}$ (thus, ${\delta}\to 0$), and using the fact that $\lim_{{({\varsigma},{\varrho})\to(1,0)}}\alpha_{\max}({\varsigma},{\varrho})=\frac{\pi}{2}$ from Proposition 6, we obtain

[TABLE]

which is the desired result. This completes the proof.     \sqcap$$\sqcup

V-I Summary of the Analysis of $\mathtt{AltMin}$

Theorem 7 proves that $\mathtt{AltMin}$ is essentially able to decode the target signal up to thrice the noise radius when ${\mathsf{snr}}$ is sufficiently large (${\varrho}\to 0$) and the R-RIP parameter ${\delta}$ is sufficiently small (${\varsigma}\to 1$). In practice, we can afford only a finite ${\mathsf{snr}}$. Moreover, obtaining smaller ${\delta}$ requires taking much more measurements scaling like $O(\frac{1}{{\delta}^{2}})$ as in (18). However, in view of Remark 2, a reasonable recovery is still possible by decoding the signal up to a multiple of noise radius $\eta\|{\bf w}\|$, where $\eta$ can be selected sufficiently large such that the recovery is still possible for a reasonable ${\mathsf{snr}}$ and ${\delta}$. From Theorem 6, the situation is much better in the noiseless case, where ${\delta}\in(0,\frac{1}{3})$ along with a good initialization ${\bf S}^{1}$ will guarantee a suitable recovery of $\check{{\bf y}}$.

Recall that for analyzing $\mathtt{AltMin}$, we made the assumption that condition (34) holds for all the solutions ${\bf H}^{i}={{\bf y}^{i}}^{\sf T}\otimes{\bf S}^{i}$ produced by $\mathtt{AltMin}$. In particular, from our discussion in Section V-B (see also Fig. 5 and especially (36)), it results that under a good initialization ${\bf S}^{1}$ $\mathtt{AltMin}$ fulfills the conditions of Theorem 1. In brief, $\mathtt{AltMin}$ is able to recover the desired signal $\check{{\bf H}}$ up to a relative precision $2{\mu}$ (see also Fig. 5).

Remark 3

A key assumption in our analysis is that ${\bf B}$ satisfies R-RIP over ${\cal H}-{\cal H}$. Our numerical simulations in Section VI, however, illustrate that $\mathtt{AltMin}$ still performs quite well even in a regime of parameters where R-RIP does not hold. Therefore, as in CS [3], R-RIP seems to be sufficient but not necessary for signal recovery. In contrast, our simulations evidently confirm that a good initialization of ${\bf S}^{1}$ plays a crucial role on the performance of $\mathtt{AltMin}$, which partly validates the results we obtained in this section via a fixed-point analysis of $\mathtt{AltMin}$ (even though our results were derived under R-RIP). $\lozenge$

VI Simulation Results

We run numerical simulation to assess the performance of $\mathtt{AltMin}$. For each simulation, we generate an $n\times k$ Gaussian matrix ${\bf B}$ and a $k$-dim signal $\check{{\bf y}}$, and take noisy unlabeled samples from $\check{{\bf y}}$ given by ${\check{{\bf x}}}=\check{{\bf S}}{\bf B}\check{{\bf y}}+{\bf w}$, where $\check{{\bf S}}$ selects $m$ out of $n$ elements in ${\bf B}\check{{\bf y}}$ randomly and where ${\bf w}$ is the additive Gaussian noise. We denote the SNR by ${\mathsf{snr}}=\frac{\|\check{{\bf S}}{\bf B}\check{{\bf y}}\|^{2}}{\|{\bf w}\|^{2}}$ as before. For simulations, we assume an SNR of $20$ dB.

VI-A Probability of Success of the Algorithm

We run $\mathtt{AltMin}$ with the noisy input ${\check{{\bf x}}}$ and a random initialization ${\bf S}^{1}\in{\cal S}$. To see the effect of the initialization, we repeat the simulation with a genie-aided initialization of ${\bf S}^{1}$ with $\nu_{{\bf S}^{1},\check{{\bf S}}}=0.2$, where $20\%$ of the rows of ${\bf S}^{1}$ are set equal to the corresponding rows of $\check{{\bf S}}$ while the remaining rows are selected completely randomly among the remaining possible rows. In both cases, we define the output of $\mathtt{AltMin}$ by ${{\bf S}^{1}\mapsto{\bf y}^{1}\mapsto{\bf S}^{2}\mapsto{\bf y}^{2}\cdots}$ and denote the final output produced by the algorithm by ${\bf y}^{\infty}$. We call the recovery successful if the relative error satisfies $\frac{\|\check{{\bf y}}-{\bf y}^{\infty}\|^{2}}{\|\check{{\bf y}}\|^{2}}\approx O(\frac{1}{{\mathsf{snr}}})$. For simulations, we set $n=1000$ and define parameters ${\kappa}=\frac{k}{n}$ as the measurement ratio and $\rho=\frac{m}{n}$ as the sampling ratio as before. For each ${\kappa}$ and $\rho$, we run simulations for $1000$ independent realizations of ${\bf B}$ and $\check{{\bf S}}$ to obtain an estimate of the success probability of $\mathtt{AltMin}$. Fig. 9 illustrates the success probability as a function of ${\kappa},\rho\in(0,1)$ for the fully random initialization in (9a) and for the genie-aided initialization in (9b). The results clearly indicate that a good initialization is crucial for the recovery performance of $\mathtt{AltMin}$, as also mentioned in Remark 3, where the genie-aided case undergoes a much sharper phase-transition in ${\kappa}-\rho$ plane.

VI-B Application to System Identification

In this section, as a practical signal processing problem, we study the estimation of the impulse response of a linear time-invariant system, classically known as System Identification [27]. This is illustrated in Fig. 10, where a known pre-designed training sequence ${\bf b}=(b_{0},\dots,b_{\tau-1})^{\sf T}$ of length $\tau$ is applied to the input of a linear system with an impulse response of length at most $k$ given by ${\bf y}=(y_{0},\dots,y_{k-1})^{\sf T}$. We assume that an estimate of the delay spread of the channel $k$ is a priori known. The output of the linear system is given by ${\bf z}={\bf b}\star{\bf y}$, where $\star$ denotes the convolution operation, where the output ${\bf z}=(z_{0},\dots,z_{n-1})^{\sf T}$ is given by $z_{l}=\sum_{t=0}^{k-1}y_{t}b_{l-t}$ for $l=0,1,\dots,n-1$, where $n=k+\tau-1$ denotes the length of the output ${\bf z}$ and where $b_{r}=0$ for $r<0$. We consider a scenario in which the output ${\bf z}$ is observed only through a noisy deletion channel, which deletes some of the output samples ${\bf z}$ but preserves their underlying order. Denoting by ${\bf y}=(y_{0},\dots,y_{k-1})^{\sf T}$, we can write ${\bf z}={\bf B}{\bf y}$ with a measurement matrix ${\bf B}$ given by

[TABLE]

where it is seen that ${\bf B}$ is an $n\times k$ matrix that depends on the training sequence ${\bf b}$. We denote the final set of $m$ samples, for some $m\leq n$, available for system identification by ${\bf x}=(x_{0},\dots,x_{m-1})^{\sf T}$, where ${\bf x}=\check{{\bf S}}{\bf z}+{\bf w}=\check{{\bf S}}{\bf B}{\bf y}+{\bf w}$, where $\check{{\bf S}}$ is a selection matrix representing the location of those samples in ${\bf z}$ that are not deleted by the deletion channel, and where ${\bf w}$ is the additive noise. It is seen that the system identification in the scenario illustrated in Fig. 10 boils down to the UOS problem (1) with a matrix ${\bf B}$ given by (94).

For simulation, we assume that the training sequence ${\bf b}=(b_{0},\dots,b_{\tau-1})^{\sf T}$ has i.i.d. ${\sf N}(0,1)$ samples and is known to the system identification algorithm (i.e., ${\bf B}$ in (94) is known). Note that although the rows of ${\bf B}$ in (94) still consist of Gaussian variables, due to the special structure of the convolution operation, they are highly correlated. Nevertheless, we can still run $\mathtt{AltMin}$ with ${\bf B}$ as in (94). Fig. 11 illustrates the simulation results for $n=1000$ and for an SNR of $20$ dB. We assume that $\mathtt{AltMin}$ is initialized with an ${\bf S}^{1}$ with $\nu_{{\bf S}^{1},\check{{\bf S}}}=0.2$. It is seen that, as expected, for a given delay spread $k$, the performance improves by increasing the length of the training sequence $\tau$ (equivalently $n$) and the number of unlabeled samples $m$. We also observe that, in comparison with Fig. 9b where ${\bf B}$ has i.i.d. components across different rows, the correlation among the rows of ${\bf B}$ degrades the performance of $\mathtt{AltMin}$ only slightly.

VII Discussion and further Remarks

Fig. 9 illustrates the success probability of a single round of $\mathtt{AltMin}$. When the success probability is quite small, $\mathtt{AltMin}$ hits a local minimum with a high probability and fails to find a suitable estimate $\widehat{{\bf H}}$. To improve the performance, we can run $\mathtt{AltMin}$ several times each time with a different random initialization ${\bf S}^{1}$ and terminate when a good estimate is found. However, this requires a procedure to certify whether $\mathtt{AltMin}$ succeeds to find a good estimate. When the measurement matrix satisfies the R-RIP over ${\cal H}-{\cal H}$, we can develop such a procedure by simply checking whether $\mathtt{AltMin}$ has been able to decode the target signal up to the noise radius, e.g., $\|{\check{{\bf x}}-\widehat{{\bf H}}{\mathbb{b}}}\|\leq\eta\|{\bf w}\|$ for some $\eta=O(1)$ as in Remark 2.

Another direction to improve the performance of $\mathtt{AltMin}$ is to find a good initialization for ${\bf S}^{1}$. The hope is that such an ${\bf S}^{1}$ lies in the basin of attraction of the desired signal $(\check{{\bf S}},\check{{\bf y}})$ under $\mathtt{AltMin}$. In this paper, we use a random initialization for ${\bf S}^{1}$. Assuming that ${\bf B}$ satisfies the R-RIP over ${\cal H}-{\cal H}$ and ${\nu_{{\bf S}^{1},\check{{\bf S}}}\geq\nu^{\prime}>0}$ for a sufficiently large $\nu^{\prime}$, Theorem 1 together with Theorem 6 and 7 guarantee a suitable recovery of the target signal $(\check{{\bf S}},\check{{\bf y}})$. The condition ${\nu_{{\bf S}^{1},\check{{\bf S}}}\geq\nu^{\prime}>0}$ is, however, quite difficult to meet for some $\check{{\bf S}}$ under the naive random initialization of ${\bf S}^{1}$. For example, if $\check{{\bf S}}$ is the selection matrix that samples the first $m$ measurements in ${\bf B}\check{{\bf y}}$, then ${\mathbb{P}}[\nu_{{\bf S}^{1},\check{{\bf S}}}\geq\gamma]=\frac{{n-m\gamma\choose m-m\gamma}}{{n\choose m}}$ for a uniformly randomly sampled ${\bf S}^{1}\in{\cal S}$. In the practically interesting regime where ${n,m\to\infty}$ and ${\frac{m}{n}\to\rho\in(0,1)}$, this scales like $e^{-n\varepsilon(\rho,\gamma)}$ with an exponent $\varepsilon(\rho,\gamma)=h(\rho)-(1-\rho\gamma)h(\frac{\rho-\rho\gamma}{1-\rho\gamma})$, where $h(\rho)$ is the entropy function introduced before. A direct calculation reveals that ${\varepsilon(\rho,\gamma)>0}$ for all $\gamma\in(0,1)$, thus, it is almost surely impossible to achieve any $\gamma\in(0,1)$ asymptotically. This implies that a good initialization method is necessary even when the RIP holds. Running a nonconvex optimization problem with a good initialization has recently been of interest in other problems in Compressed Sensing such as phase retrieval [28], blind deconvolution [29], and blind calibration [30]. We leave developing a good initialization scheme for our algorithm as a future work.

The unlabeled sensing problem $x={\bf S}{\bf B}{\bf y}+{\bf w}$ studied in this paper can be extended to cases in which the signal ${\bf y}$ belongs to a structured class of signals ${\cal Y}\subseteq{\mathbb{R}}^{k}$ and ${\bf S}$ resides in a family of selection matrices ${\cal S}$ other than the ordered sampling matrices studied in this paper. As in the UOS, in terms of signal recovery, we need to check two main requirements. The first is to develop an R-RIP over ${\cal H}-{\cal H}$ under a suitable metric (e.g., $l_{2}$ distance as in this paper). Following Proposition 2 and assuming that ${\bf B}$ has i.i.d. Gaussian components, such an R-RIP with a constant ${{\delta}\in(0,1)}$ can be generally derived with a probability larger than ${1-2|{\cal S}|^{2}(1+\frac{2}{{\delta}})^{2}e^{-cm{\delta}^{2}{\mu}^{2}}}$, where $|{\cal S}|$ denotes the cardinality of ${\cal S}$. This provides a theoretical lower bound on the the number of measurements $m$ for a given ${\delta}$. The second requirement is to develop an algorithm to recover the signal ${{\bf H}}={{\bf y}}^{\sf T}\otimes{{\bf S}}\in{\cal H}$ from the noisy unlabeled measurements ${{\bf x}}={{\bf H}}{\mathbb{b}}+{\bf w}={{\bf S}}{\bf B}{{\bf y}}+{\bf w}$. From Remark 2, under the R-RIP, such an algorithm needs to recover an estimate $\widehat{{\bf H}}\in{\cal H}$ satisfying $\|\widehat{{\bf H}}{\mathbb{b}}-{{\bf x}}\|\leq\eta\|{\bf w}\|$ for some $\eta=O(1)$, where $\|.\|$ here denotes the metric with respect to which the RIP is derived. The performance guarantee we derived for UOS  in Theorem 1 immediately applies to such an estimate. For the UOS studied in this paper, we used the alternating minimization over ${\bf y}$ and ${\bf S}$, where the latter minimization was done with a feasible complexity by using the ordered structure of the matrices in ${\cal S}$ and applying the Dynamic Programming. Deriving such an algorithm for a general signal set ${\cal H}$ (for a general signal in ${\cal Y}$ and unlabeled sampling structure in ${\cal S}$) requires exploiting the algebraic as well as the geometric structure of ${\cal H}$.

VIII Conclusion

In this paper, we studied the Unlabeled Ordered Sampling  (UOS) problem, where the goal was to recover a signal from a set of linear measurements taken via a fully known measurement matrix when the labels of the measurements are missing but their order is preserved. We identified a duality between UOS and the traditional Compressed Sensing (CS), where the unknown support (location of nonzero elements) of a sparse signal in CS corresponds in a natural way to the unknown indices of the measurements kept in UOS. Motivated by this duality, we developed a Restricted Isometry Property (RIP) similar to that in CS. We also designed a low-complexity Alternating Minimization algorithm to recover the target signal from the set of its noisy unlabeled samples. We analyzed the performance of our proposed algorithm for different signal dimensions and number of measurements theoretically under the established RIP. We also provided numerical simulations to validate the theoretical results.

Appendix A Proofs

A-A Proof of Proposition 2

Since R-RIP is scale-invariant, we first define the following normalized sets

[TABLE]

where ${\cal D}_{1}\subset{\cal D}_{2}$. We also define the following probability events (on the random realization of ${\mathbb{b}}$):

[TABLE]

where for $a,b\in{\mathbb{R}}$ we denoted by $a\vee b=\max\{a,b\}$. Note that ${\cal E}_{1}\subset{\cal E}_{2}$ since the subevents corresponding to $({\bf H},{\bf H}^{\prime})\in{\cal D}_{1}$ with $d_{{\bf H},{\bf H}^{\prime}}^{2}\geq m{\mu}^{2}$ (thus, $d_{{\bf H},{\bf H}^{\prime}}^{2}\vee m{\mu}^{2}=d_{{\bf H},{\bf H}^{\prime}}^{2}$) in ${\cal E}_{1}$ are also included in ${\cal E}_{2}$, where ${\cal E}_{2}$ contains in addition those subevents corresponding to $d_{{\bf H},{\bf H}^{\prime}}^{2}\leq m{\mu}^{2}$. To prove R-RIP result, we need to find an upper bound on ${\mathbb{P}}[{\cal E}_{1}]$. Since ${\cal E}_{1}\subset{\cal E}_{2}$, thus, ${\mathbb{P}}[{\cal E}_{1}]\leq{\mathbb{P}}[{\cal E}_{2}]$, we will do this by deriving an upper bound for ${\mathbb{P}}[{\cal E}_{2}]$. So, we focus on the event ${\cal E}_{2}$ in the sequel. It is seen that ${\cal E}_{2}$ consists of the union of a continuum of events labeled with ${\bf H},{\bf H}^{\prime}\in{\cal D}_{2}$. As in the proof of Proposition 1, we will first derive a concentration bound for a fixed ${\bf H},{\bf H}^{\prime}\in{\cal D}_{2}$ and then extend it to the whole set ${\cal D}_{2}$ via a net argument and applying the union bound.

Consider a fixed $({\bf H},{\bf H}^{\prime})\in{\cal D}_{2}$ where ${\bf H}={\bf y}^{\sf T}\otimes{\bf S}$ and ${\bf H}^{\prime}={{\bf y}^{\prime}}^{\sf T}\otimes{\bf S}^{\prime}$ with $\|{\bf y}\|=1$ and $\|{\bf y}^{\prime}\|\leq 1$. Let us define ${\mathscr{E}}_{{\bf H},{\bf H}^{\prime}}=({\bf H}-{\bf H}^{\prime}){\mathbb{b}}$. Note that ${\mathscr{E}}_{{\bf H},{\bf H}^{\prime}}$ is an $m$-dim Gaussian vector with a zero mean and a covariance matrix

[TABLE]

where ${\bf I}_{m}$ denotes the identity matrix of order $m$ and where in $(i)$ we used the fact that

[TABLE]

where in $(ii)$ we used the property of the Kronecker product, where in $(iii)$ we used the fact that ${\bf y}^{\sf T}{\bf y}=\|{\bf y}\|^{2}$ is just a number and dropped the Kronecker product, and where in $(iv)$ we used the fact that the each row of ${\bf S}$ has only one $1$ at a specific column, thus, different rows are orthogonal to each other. A similar derivation gives the second and the third term in (A-A). It is seen that $\Sigma$ is not in general a diagonal matrix, thus, ${\mathscr{E}}_{{\bf H},{\bf H}^{\prime}}$ consists of correlated Gaussian variables and the conventional concentration result for the i.i.d. Gaussian variables does not immediately apply. We first prove that although $\Sigma$ is not a diagonal matrix, its singular values are bounded and in particular do not grow with the dimension $m$. In words, this implies that the components of ${\mathscr{E}}_{{\bf H},{\bf H}^{\prime}}$ are not that correlated. We denote by \hbox{\boldmath\Sigma}={\bf U}\hbox{\boldmath\Lambda}{\bf U}^{\sf T} the Singular Value Decomposition (SVD) of $\Sigma$ where \hbox{\boldmath\Lambda}={\text{diag}}(\lambda_{1},\dots,\lambda_{m}) is the diagonal matrix consisting of the singular values ${\boldsymbol{\lambda}}=(\lambda_{1},\dots,\lambda_{m})^{\sf T}$. We use the convention that the singular values are ordered with $\lambda_{m}\leq\dots\leq\lambda_{1}$. We have the following result.

Lemma 8

Let $\Sigma$ and ${\boldsymbol{\lambda}}=(\lambda_{1},\dots,\lambda_{m})^{\sf T}$ be as before. Then, all the singular values satisfy

[TABLE]

Proof:

Let us denote by ${{\sf R}:=r_{1}<\cdots<r_{m}}$ and ${\sf R}^{\prime}:=r^{\prime}_{1}<\cdots<r^{\prime}_{m}$ the ordered sequences consisting of indices of those rows of ${\bf B}$ selected by ${\bf S}$ and ${\bf S}^{\prime}$, where $r_{i},r^{\prime}_{i}\in[n]$. Also, let ${\cal C}:=\{i\in[m]:r_{i}=r^{\prime}_{i}\}$ be the index set of similar elements in ${\sf R}$ and ${\sf R}^{\prime}$. Since ${\bf S}$ and ${\bf S}^{\prime}$ have only one $1$ in each row at column set ${\sf R}$ and ${\sf R}^{\prime}$, we can simply check that ${\bf S}{{\bf S}^{\prime}}^{\sf T}$ has at most one $1$ at each row, where $({\bf S}{{\bf S}^{\prime}}^{\sf T})_{ij}=1$ if and only if $r_{i}=r^{\prime}_{j}$. In particular, the only nonzero diagonal elements of ${\bf S}{{\bf S}^{\prime}}^{\sf T}$ lie on the rows belonging to ${\cal C}$. This implies that the symmetric matrix \hbox{\boldmath\Gamma}:={\bf S}{{\bf S}^{\prime}}^{\sf T}+{\bf S}^{\prime}{{\bf S}}^{\sf T} has the diagonal element $2$ and zero off-diagonal terms at the rows belonging to ${\cal C}$. Moreover, $\Gamma$ has at most two $1$’s in the other rows (not belonging to ${\cal C}$), where those $1$’s do not lie on the diagonal of $\Gamma$. Therefore, from (A-A), we have the following two cases. On the rows belonging to ${\cal C}$, $\Sigma$ has only a diagonal element $\|{\bf y}\|^{2}+\|{\bf y}^{\prime}\|^{2}-2{\langle{\bf y},{\bf y}^{\prime}\rangle}$, which is also a singular values of $\Sigma$. On the rows not belonging to ${\cal C}$, $\Sigma$ has a diagonal element $\|{\bf y}\|^{2}+\|{\bf y}^{\prime}\|^{2}$ plus at most two off-diagonal terms given by $-{\langle{\bf y},{\bf y}^{\prime}\rangle}$. Hence, form the Gershgorin disk theorem [31], all the singular values of $\Sigma$ should lie in the range $\|{\bf y}\|^{2}+\|{\bf y}^{\prime}\|^{2}\pm 2|{\langle{\bf y},{\bf y}^{\prime}\rangle}|$. This completes the proof.     \sqcap$$\sqcup

Let \hbox{\boldmath\Sigma}={\bf U}\hbox{\boldmath\Lambda}{\bf U}^{\sf T} be the SVD of $\Sigma$ as before and let {\bf e}=\hbox{\boldmath\Lambda}^{-\frac{1}{2}}{\bf U}^{\sf T}{\mathscr{E}}_{{\bf H},{\bf H}^{\prime}}. We can check that ${\bf e}=(e_{1},\dots,e_{m})^{\sf T}$ consists of i.i.d. ${\sf N}(0,1)$ variables. Since ${\bf U}$ is an orthogonal matrix, i.e., ${\bf U}{\bf U}^{\sf T}={\bf U}^{\sf T}{\bf U}={\bf I}_{m}$, we have $\|{\mathscr{E}}_{{\bf H},{\bf H}^{\prime}}\|^{2}=\|{\bf U}^{\sf T}{\mathscr{E}}_{{\bf H},{\bf H}^{\prime}}\|^{2}=\sum_{i=1}^{m}\lambda_{i}e_{i}^{2}$. We also have ${\mathbb{E}}[\|{\mathscr{E}}_{{\bf H},{\bf H}^{\prime}}\|^{2}]=\|{\bf H}-{\bf H}^{\prime}\|_{\sf F}^{2}=d_{{\bf H},{\bf H}^{\prime}}^{2}=\sum_{i=1}^{m}\lambda_{i}$. Defining $T:=d_{{\bf H},{\bf H}^{\prime}}^{2}+{\delta}(d_{{\bf H},{\bf H}^{\prime}}^{2}\vee m{\mu}^{2})$ and setting $\theta$ to be a positive variable, we obtain the following concentration bound [20]

[TABLE]

where $E(\theta,{\boldsymbol{\lambda}}):=\theta T+\frac{1}{2}\sum_{i=1}^{m}\log(1-2\theta\lambda_{i})$. We also used the fact that for a ${\sf N}(0,1)$ variables $e_{i}$, ${\mathbb{E}}[e^{se_{i}^{2}}]=\frac{1}{\sqrt{1-2s}}$ for any $s\in(0,\frac{1}{2})$, thus, the feasible range of $\theta$ in (98) is given by $(0,\frac{1}{2\lambda_{1}})$, where $\lambda_{1}=\max\{\lambda_{i}:i\in[m]\}$ is the the maximum singular value. From Lemma 8, it results that

[TABLE]

where in $(i)$, we used the fact that $\|{\bf y}\|=1$ and $\|{\bf y}^{\prime}\|\leq 1$ for any $({\bf H},{\bf H}^{\prime})\in{\cal D}_{2}$. This implies that, for all $({\bf H},{\bf H}^{\prime})\in{\cal D}_{2}$, the set of permitted values of $\theta$ at least contains $(0,\frac{1}{8})$. Now let us consider a fixed $\theta\in(0,\frac{1}{8})$. To derive a concentration bound for (98), we need to find a strictly positive lower bound on $E(\theta,{\boldsymbol{\lambda}})$ that is independent of the configuration of the singular values ${\boldsymbol{\lambda}}$. We have the following lemma.

Lemma 9

Let $E(\theta,{\boldsymbol{\lambda}})$ be as before. Then, for any $\theta\in(0,\frac{1}{8})$, the vector ${\boldsymbol{\lambda}}$ minimizing $E(\theta,{\boldsymbol{\lambda}})$ (i.e., the worst case singular value configuration) is given by the vector ${\boldsymbol{\lambda}}^{*}:=(4,\dots,4,0,\dots,0)^{\sf T}$ that has $4$ at its $\frac{d_{{\bf H},{\bf H}^{\prime}}^{2}}{4}$ components and is [math] elsewhere. $\lozenge$

Proof:

Note that the only constraint we put on ${\boldsymbol{\lambda}}$ is that it should belong to the set

[TABLE]

where the upper bound $4$ results from (99). Since $T=d_{{\bf H},{\bf H}^{\prime}}^{2}+{\delta}(d_{{\bf H},{\bf H}^{\prime}}^{2}\vee m{\mu}^{2})$ is independent of ${\boldsymbol{\lambda}}$, from $E(\theta,{\boldsymbol{\lambda}})=\theta T+\frac{1}{2}\sum_{i=1}^{m}\log(1-2\theta\lambda_{i})$ and the concavity of the Logarithm, it results that $E(\theta,{\boldsymbol{\lambda}})$ is a concave function of ${\boldsymbol{\lambda}}$ over the constraint set (100). Therefore, it achieves its minimum at the boundary [26] of the constraint set (100), which corresponds to ${\boldsymbol{\lambda}}^{*}={(4,\dots,4,0,\dots,0)^{\sf T}}$ in the statement of the lemma. This completes the proof.     \sqcap$$\sqcup

From Lemma 9, it results that $E(\theta,{\boldsymbol{\lambda}})$ for $\theta\in(0,\frac{1}{8})$ and for all valid configurations of the singular values ${\boldsymbol{\lambda}}$ is lower bounded by the following function

[TABLE]

where ${\boldsymbol{\lambda}}^{*}=(4,\dots,4,0,\dots,0)^{\sf T}$ is as in Lemma 9, and where ${T=d_{{\bf H},{\bf H}^{\prime}}^{2}+{\delta}(d_{{\bf H},{\bf H}^{\prime}}^{2}\vee m{\mu}^{2})}$ is as before. The optimal $\theta^{*}\in(0,\frac{1}{8})$ minimizing $E(\theta)$ is given by $\theta^{*}=\frac{1}{8}(1-\frac{d_{{\bf H},{\bf H}^{\prime}}^{2}}{T})$, where after replacing in $E(\theta)$, yields the following exponent

[TABLE]

For $d_{{\bf H},{\bf H}^{\prime}}^{2}\in[m{\mu}^{2},\infty)$, $E_{\min}=\frac{1}{8}({\delta}-\log(1+{\delta}))d_{{\bf H},{\bf H}^{\prime}}^{2}$, which is an increasing function of $d_{{\bf H},{\bf H}^{\prime}}^{2}$ as ${\delta}-\log(1+{\delta})\geq 0$ for all ${\delta}\in{\mathbb{R}}_{+}$ and in particular ${\delta}\in(0,1)$, thus, the minimum of $E_{\min}$ is achieved at $d_{{\bf H},{\bf H}^{\prime}}^{2}=m{\mu}^{2}$. Similarly, for $d_{{\bf H},{\bf H}^{\prime}}^{2}\in[0,m{\mu}^{2}]$, $E_{\min}=\frac{m{\mu}^{2}{\delta}}{8}(1-\frac{d_{{\bf H},{\bf H}^{\prime}}^{2}}{m{\mu}^{2}{\delta}}\log(1+\frac{m{\mu}^{2}{\delta}}{d_{{\bf H},{\bf H}^{\prime}}^{2}}))$. We can check that $x\mapsto 1-x\log(1+\frac{1}{x})$ is a decreasing function of $x$ for $x\in{\mathbb{R}}_{+}$. Thus, setting $x$ equal to $\frac{d_{{\bf H},{\bf H}^{\prime}}^{2}}{m{\mu}^{2}{\delta}}$, it results that the minimum of $E_{\min}$ for $d_{{\bf H},{\bf H}^{\prime}}^{2}\in[0,m{\mu}^{2}]$ is also achieved at $d_{{\bf H},{\bf H}^{\prime}}^{2}=m{\mu}^{2}$. Overall, by setting $d_{{\bf H},{\bf H}^{\prime}}^{2}=m{\mu}^{2}$, we obtain that $E_{\min}$ is lower bounded by

[TABLE]

where in $(i)$ we used the inequality $\log(1+{\delta})\leq{\delta}-\frac{{\delta}^{2}}{2}+\frac{{\delta}^{3}}{3}$ for ${\delta}\in{\mathbb{R}}_{+}$, and where in $(ii)$ we used $\frac{{\delta}^{3}}{3}\leq\frac{{\delta}^{2}}{3}$ for ${\delta}\in(0,1)$. This establishes a strictly positive exponent $E_{0}(\mu,{\delta})=\frac{m{\mu}^{2}{\delta}^{2}}{48}$ for (98), which holds for all $({\bf H},{\bf H}^{\prime})\in{\cal D}_{2}$. By following similar steps, we can extend the concentration bound in (98) to the reverse inequality, where overall we obtain that for any fixed $({\bf H},{\bf H}^{\prime})\in{\cal D}_{2}$

[TABLE]

The final step is to generalize (103) to derive a concentration bound for all $({\bf H},{\bf H}^{\prime})\in{\cal D}_{2}$. As in the proof of Proposition 1, we do this by quantizing ${\cal D}_{2}$ into a $\frac{{\delta}}{2}$-net with a minimal size and applying the union bound. We can build such an $\frac{{\delta}}{2}$-net by finding a joint net for ${\bf y}\in{\cal B}_{2}^{k}$ and ${\bf y}^{\prime}\in{\cal B}_{2}^{k}$, each consisting of at most $N=(1+\frac{2}{{\delta}})^{k}$ points (as in the proof of Proposition 1). Taking the union bound over this joint net and also all possible selection matrices ${\bf S}$ and ${\bf S}^{\prime}$, we obtain that

[TABLE]

where $c>0$ is a constant independent of $m,n,k$ and ${\delta},{\mu}$. This completes the proof.

A-B Proof of Proposition 5

We first need the following two lemmas. We refer to [32] for an introduction to (strict) convexity and (strictly) convex functions needed in this section.

Lemma 10

Let ${\phi(\alpha)=\frac{\cos(\alpha)+\sin(\alpha)-1}{\cos(\alpha)+b}}$, where $b\in{\mathbb{R}}_{+}$ is a constant. Then, $\phi(\alpha)$ is strictly concave over $\alpha\in[0,\frac{\pi}{2}]$. $\square$

Proof:

Taking the derivative and simplifying, we obtain

[TABLE]

Also, taking the second derivative and simplifying yields

[TABLE]

Note that for $\alpha\in[0,\frac{\pi}{2}]$ the denominator $(\cos(\alpha)+b)^{3}$ of both terms is always positive. The numerator of the first term can be written as $b\cos(\alpha)(\sin(\alpha)-1-b)$ which is negative since $\cos(\alpha)\geq 0$, $b\in{\mathbb{R}}_{+}$, and $\sin(\alpha)-1-b\leq 0$. Now let us consider the numerator of the second term. If $\alpha=0$, $\sin(\alpha)=0$ and the numerator is given by $-1-b$, which is negative. So, we consider $\alpha\in(0,\frac{\pi}{2}]$, where $\sin(\alpha)\in(0,1]$. For this case, the numerator of the second term simplifies to

[TABLE]

where in $(i)$ we used the fact that $b\in{\mathbb{R}}_{+}$, thus, $1+b>0$, and the fact that $\sin(\alpha)>0$, and applied the inequality $l+\frac{1}{l}\geq 2$ for $l\in{\mathbb{R}}_{+}$ by replacing $l=\sin(\alpha)$. This yields $\phi^{\prime\prime}(\alpha)\leq 0$ for $a\in[0,\frac{\pi}{2}]$. Moreover, we can check that $\phi^{\prime\prime}(\alpha)<0$ except at some single points $a$, thus, $\phi(\alpha)$ is strictly concave.     \sqcap$$\sqcup

Lemma 11

Let $\psi(\alpha)=\frac{\cos(\alpha)+b}{\cos(\alpha)+\sin(\alpha)-1}$, where $b\in{\mathbb{R}}_{+}$ is a constant. Then, $\psi(\alpha)$ is convex over $\alpha\in(0,\frac{\pi}{2})$. $\square$

Proof:

First note that $\cos(\alpha)+\sin(\alpha)-1\geq 0$ and $\cos(\alpha)+b\geq 0$ for $\alpha\in(0,\frac{\pi}{2})$, thus, $\psi(\alpha)\geq 0$. We can write ${-\psi(\alpha)=\iota\circ\phi(\alpha)}$ where ${\iota(x)=-\frac{1}{x}}$ and is a concave increasing function over $(0,\infty)$, which contains $\phi(\alpha)$ for $\alpha\in(0,\frac{\pi}{2})$ (since $\phi(\alpha)>0$ in this domain). Let $\alpha,\alpha^{\prime}\in(0,\frac{\pi}{2})$ and $\lambda\in[0,1]$, and set $\bar{\lambda}=1-\lambda$. We have

[TABLE]

where in $(i)$ we used the concavity of $\phi$ proved in Lemma 10 and the fact that $\iota$ is an increasing function, and where in $(ii)$ we used the concavity of $\iota$. This implies that $-\psi(\alpha)$ is a concave function, thus, $\psi(\alpha)$ is convex over $\alpha\in(0,\frac{\pi}{2})$.     \sqcap$$\sqcup

We can now prove Proposition 5. First note that from (84)

[TABLE]

As ${{\varsigma}-1>0}$ and ${\varrho}>0$, the first term is a strictly convex function of $\alpha\in(0,\frac{\pi}{2})$ from Lemma 11 (by setting $b=\frac{{\varrho}}{{\varsigma}-1}>0$). The second term is also a concave function of $\alpha$ from Lemma 10 (by setting $b=\frac{{\varrho}}{{\varsigma}+1}>0$). Therefore, $\upsilon(\alpha)$ is strictly convex in $(0,\frac{\pi}{2})$.

As ${\upsilon(0^{+})=\upsilon(\frac{\pi}{2}^{-})=+\infty}$, the roots of $\upsilon$ should lie in $(0,\frac{\pi}{2})$. Suppose that $\upsilon(\alpha)$ has more than two distinct roots, i.e., there are $\alpha_{1}<\alpha_{2}<\alpha_{3}$ in $(0,\frac{\pi}{2})$ with $\upsilon(\alpha_{i})=0$, for $i=1,2,3$. Then, setting $\lambda=\frac{\alpha_{2}-\alpha_{1}}{\alpha_{3}-\alpha_{1}}\in(0,1)$, a simple calculation shows that $\alpha_{2}=\lambda\alpha_{3}+(1-\lambda)\alpha_{1}$. Thus, we have

[TABLE]

where in $(i)$ we used the the fact that $\lambda\in(0,1)$ and that $\upsilon(\alpha)$ is strictly convex for $\alpha\in(0,1)$. Since (116) is a contradiction, $\upsilon(\alpha)$ cannot have more than two roots.

A-C Proof of Proposition 6

First note that $\upsilon(\alpha,{\varsigma},{\varrho})$ is a continuous and differentiable function of $(\alpha,{\varsigma},{\varrho})\in(0,\frac{\pi}{2})\times(1,\infty)\times(0,1)$ of any order.

To prove $(i)$, let $({\varsigma}_{0},{\varrho}_{0})\in{\cal D}_{\text{\textvarsigma,\textrho}}$ be an arbitrary point and let $\alpha_{\min}({\varsigma}_{0},{\varrho}_{0})<\alpha_{\max}({\varsigma}_{0},{\varrho}_{0})$ be the corresponding two roots. From the strict convexity of $\upsilon$ proved in Proposition 5, we have $\upsilon(\alpha_{0},{\varsigma}_{0},{\varrho}_{0})<0$ for $\alpha_{0}=\frac{\alpha_{\min}({\varsigma}_{0},{\varrho}_{0})+\alpha_{\max}({\varsigma}_{0},{\varrho}_{0})}{2}$. Thus, from the continuity of $\upsilon(\alpha_{0},{\varsigma},{\varrho})$ it results that there is an open set ${\cal D}_{0}:=({\varsigma}_{0}-\varepsilon,{\varsigma}_{0}+\varepsilon)\times({\varrho}_{0}-\varepsilon,{\varrho}_{0}+\varepsilon)$, with a sufficiently small $\varepsilon>0$, around $({\varsigma}_{0},{\varrho}_{0})$ over which $\upsilon(\alpha_{0},{\varsigma},{\varrho})<0$. For all $({\varsigma},{\varrho})\in{\cal D}_{0}$, $\upsilon(\alpha,{\varsigma},{\varrho})$ must have two roots from Proposition 5. This implies that ${{\cal D}_{0}\subset{\cal D}_{\text{\textvarsigma,\textrho}}}$, thus, ${\cal D}_{\text{\textvarsigma,\textrho}}$ is an open set.

To prove $(ii)$, note that $\upsilon(\alpha,{\varsigma},{\varrho})=0$ defines $\alpha_{\min}$ and $\alpha_{\max}$ implicitly as functions of $({\varsigma},{\varrho})\in{\cal D}_{\text{\textvarsigma,\textrho}}$. Since ${\cal D}_{\text{\textvarsigma,\textrho}}$ is an open set from $(i)$, Implicit Function Theorem [25] implies that $\alpha_{{\max/\min}}$ would be differentiable at a point $({\varsigma}_{0},{\varrho}_{0})\in{\cal D}_{\text{\textvarsigma,\textrho}}$ provided that $\frac{\partial}{\partial\alpha}\upsilon(\alpha,{\varsigma}_{0},{\varrho}_{0})=\upsilon^{\prime}(\alpha,{\varsigma}_{0},{\varrho}_{0})\neq 0$ at ${\alpha=\alpha_{{\max/\min}}({\varsigma}_{0},{\varrho}_{0})}$. So, we need to prove that this condition is satisfied (see, e.g., Fig. 7 where $\upsilon^{\prime}(\alpha,{\varsigma}_{0},{\varrho}_{0})$ is strictly positive at $\alpha_{\max}$ and strictly negative at $\alpha_{\min}$). Suppose, for example, that $\upsilon^{\prime}(\alpha_{\min},{\varsigma}_{0},{\varrho}_{0})=0$. Then, from the strict convexity of $\upsilon(\alpha,{\varsigma}_{0},{\varrho}_{0})$ proved in Proposition 5, it results that $\upsilon^{\prime}(\alpha,{\varsigma}_{0},{\varrho}_{0})>0$ for $\alpha\in(\alpha_{\min},\alpha_{\max})$, which from $\upsilon(\alpha_{\min},{\varsigma}_{0},{\varrho}_{0})=0$ implies that $\upsilon(\alpha_{\max},{\varsigma}_{0},{\varrho}_{0})>0$, thus, contradicting the fact that $\alpha_{\max}$ is another root of $\upsilon(\alpha,{\varsigma}_{0},{\varrho}_{0})$. A similar argument shows that $\upsilon^{\prime}(\alpha_{\max},{\varsigma}_{0},{\varrho}_{0})\neq 0$. Therefore, both $\alpha_{\min}$ and $\alpha_{\max}$ are differentiable in ${\cal D}_{\text{\textvarsigma,\textrho}}$.

To prove $(iii)$, we first check that $\upsilon({\varsigma},{\varrho})$ is an increasing function of ${\varsigma}$ and ${\varrho}$. This follows simply from the fact that $\cos(\alpha)\geq 0$ for $\alpha\in[0,\frac{\pi}{2}]$, thus, the numerator/denominator of the first/second term in (84) is an increasing function of $({\varsigma},{\varrho})$. This implies that $\frac{\partial}{\partial{\varsigma}}\upsilon(\alpha,{\varsigma},{\varrho})>0$ and $\frac{\partial}{\partial{\varrho}}\upsilon(\alpha,{\varsigma},{\varrho})>0$ for all $\alpha\in(0,\frac{\pi}{2})$ including $\alpha_{\min/\max}$. Using the Implicit Function Theorem [25] and taking the derivative with respect to ${\varsigma}$ from $\upsilon(\alpha_{\min}({\varsigma},{\varrho}),{\varsigma},{\varrho})=0$ yields

[TABLE]

From $\upsilon^{\prime}(\alpha_{\min},{\varsigma},{\varrho})<0$ established in $(ii)$ (see also Fig. 7), this implies that $\frac{\partial}{\partial{\varsigma}}\alpha_{\min}({\varsigma},{\varrho})>0$. Similarly, we obtain that $\frac{\partial}{\partial{\varrho}}\alpha_{\min}({\varsigma},{\varrho})>0$. Thus, $\alpha_{\min}({\varsigma},{\varrho})$ is an increasing function of $({\varsigma},{\varrho})$. The result for $\alpha_{\max}({\varsigma},{\varrho})$ follows similarly with the only difference that $\upsilon^{\prime}(\alpha_{\max},{\varsigma},{\varrho})>0$ (see Fig. 7), which implies that $\alpha_{\max}({\varsigma},{\varrho})$ is a decreasing function of $({\varsigma},{\varrho})$.

To prove $(iv)$, it is easier to use (81). Note that

[TABLE]

where the last term converges to [math] as ${({\varsigma},{\varrho})\to(1,0)}$. Let us define the open region ${\cal D}=\{({\varsigma},{\varrho}):({\varsigma}+{\varrho})^{2}-1<\varepsilon^{2}\}$ with $\varepsilon\in(0,\sqrt{2}-1)$. For any $({\varsigma},{\varrho})\in{\cal D}$, it is seen from (117) that (81) has a solution satisfying $\cos(\alpha)+\sin(\alpha)<1+\varepsilon$. Note that such a solution exists since $1+\varepsilon\in(1,\sqrt{2})$, which is a subset of the range of $\cos(\alpha)+\sin(\alpha)-1$ for $\alpha\in[0,\frac{\pi}{2}]$ given by $[1,\sqrt{2}]$. Thus, ${\cal D}\subseteq{\cal D}_{\text{\textvarsigma,\textrho}}$. In particular, ${\cal D}$ contains all the paths along which $({\varsigma},{\varrho})$ can approach to $(1,0)$ from the inside of ${\cal D}_{\text{\textvarsigma,\textrho}}$. From the identity $\cos(\alpha)+\sin(\alpha)=\sqrt{2}\cos(\alpha-\frac{\pi}{4})$, it results that for any $({\varsigma},{\varrho})\in{\cal D}$

[TABLE]

Taking the limit as $\varepsilon\to 0$ yields that any limit point of $\alpha_{\max}({\varsigma},{\varrho})$ as ${({\varsigma},{\varrho})\to(1,0)}$ should be larger than $\lim_{\varepsilon\to 0}\frac{\pi}{4}+\cos^{-1}(\frac{1+\varepsilon}{\sqrt{2}})=\frac{\pi}{2}$. Since ${\alpha\in[0,\frac{\pi}{2}]}$, this implies that $\lim_{{({\varsigma},{\varrho})\to(1,0)}}\alpha_{\max}({\varsigma},{\varrho})=\frac{\pi}{2}$. Similarly, it results that $\lim_{{({\varsigma},{\varrho})\to(1,0)}}\alpha_{\min}({\varsigma},{\varrho})=0$. This completes the proof.