Sparsity Order Estimation from a Single Compressed Observation Vector

TL;DR

This paper introduces a novel method for estimating the sparsity order from a single compressed observation vector by designing specific measurement matrices that enable direct inference without sparse recovery.

Contribution

It demonstrates that measurement matrices based on Khatri-Rao and Vandermonde structures can reliably infer sparsity order from a single snapshot, with new designs for low-coherence Vandermonde matrices.

Findings

Effective rank of measurement matrix indicates sparsity order in multiple snapshots.

Khatri-Rao and Vandermonde matrices enable single snapshot sparsity estimation.

New Vandermonde matrix design achieves low coherence, improving measurement quality.

Abstract

We investigate the problem of estimating the unknown degree of sparsity from compressive measurements without the need to carry out a sparse recovery step. While the sparsity order can be directly inferred from the effective rank of the observation matrix in the multiple snapshot case, this appears to be impossible in the more challenging single snapshot case. We show that specially designed measurement matrices allow to rearrange the measurement vector into a matrix such that its effective rank coincides with the effective sparsity order. In fact, we prove that matrices which are composed of a Khatri-Rao product of smaller matrices generate measurements that allow to infer the sparsity order. Moreover, if some samples are used more than once, one of the matrices needs to be Vandermonde. These structural constraints reduce the degrees of freedom in choosing the measurement matrix which may incur in a degradation in the achievable coherence. We thus also address suitable choices of the measurement matrices. In particular, we analyze Khatri-Rao and Vandermonde matrices in terms of their coherence and provide a new design for Vandermonde matrices that achieves a low coherence.

Figures

Equations

$$\mathbf{b}=\mathbf{A}\cdot \mathbf{S}\cdot \mathbf{x}+\mathbf{n},$$

$${\mathbf{b}}_i={\left[\begin{array}{ccc}{b}_{1+p\cdot (i-1)},& b_{2+p\cdot (i-1)},& \dots ,b_{\ell +p\cdot (i-1)}\end{array}\right]}^{\mathrm{T}}\in {\mathrm{C}}^{\ell }$$

$${\mathbf{b}}_i={\mathbf{A}}_i\cdot \mathbf{x}.$$

$${\mathbf{b}}_i={\mathbf{A}}_i\cdot \mathbf{x}={\mathbf{a}}_{i,q}\cdot x_q\text{for}i\in [k].$$

$${\mathbf{a}}_{i,q}={\phi }_{i,q}\cdot {\mathbf{\psi }}_q\text{for}i\in [k],$$

$$\mathbf{B}=\mathbf{\Psi }\mathrm{diag}(\mathbf{x}){\mathbf{\Phi }}^{\mathrm{T}},$$

$$\mathbf{B}={\mathbf{\Psi }}_S\mathrm{diag}({\mathbf{x}}_S){\mathbf{\Phi }}_S^{\mathrm{T}}.$$

$$\mathbf{B}=x_q[{\mathbf{a}}_{1,q},\dots ,{\mathbf{a}}_{k,q}].$$

$${\mathbf{a}}_{i,q}=z_{i,q}{\mathbf{a}}_{i-1,q}\text{for}i\in \{2,\dots ,k\},$$

$${\mathbf{a}}_{i+1,q}=z_{i+1,q}{\mathbf{a}}_{i,q}\text{for}i\in \{1,\dots ,k-1\}.$$

$$({\mathbf{a}}_{i,q}{)}_j=({\mathbf{a}}_{i-1,q}{)}_{p+j}\text{for}j\in [\ell -p]+1$$

$$({\mathbf{a}}_{i+1,q}{)}_j=({\mathbf{a}}_{i,q}{)}_{p+j}\text{for}j\in [\ell -p]$$

$$z_{i+1,q}{\left({\mathbf{a}}_{i,q}\right)}_j=z_{i,q}{\left({\mathbf{a}}_{i-1,q}\right)}_{p+j}\text{for}j\in [\ell -p]$$

$$z_{i,q}=z_{i+1,q}=z_q\text{for}i\in [k-1]$$

$${\mathbf{a}}_{i,q}=z_q{\mathbf{a}}_{i-1,q}=z_q^2{\mathbf{a}}_{i-2,q}=\cdots =z_q^{i-1}{\mathbf{a}}_{1,q}.$$

$$\mathbf{V}=\left[\begin{array}{ccc}1& \dots & 1\\ z_1& \dots & z_N\\ ⋮& \ddots & ⋮\\ z_1^{\lceil m/p\rceil -1}& \dots & z_N^{\lceil m/p\rceil -1}\end{array}\right],$$

$$\mathbf{B}=\hat{\mathbf{\Psi }}\cdot \mathrm{diag}(\mathbf{x})\cdot {\mathbf{V}}_k^{\mathrm{T}},$$

$$\det (\mathbf{V})=\prod _{1⩽i<j⩽n}(z_i-z_j).$$

$$\mathbf{m}=\left(\mathbf{V}\diamond \mathbf{\Psi }\right)\cdot \mathbf{x}=0.$$

$$\mathbf{M}=\mathbf{\Psi }\mathrm{diag}(\mathbf{x}){\mathbf{V}}^{\mathrm{T}}=\tilde{\mathbf{\Psi }}{\mathbf{V}}^{\mathrm{T}},$$

$$\hat{\mathbf{B}}=\mathbf{B}+\mathbf{N}.$$

$${\mathbf{J}}_i=\left[\begin{array}{ccc}{0}_{\ell \times p(i-1)},& {\mathbf{I}}_{\ell \times \ell },& 0_{\ell \times (N-\ell -p(i-1))}\end{array}\right]\in {\mathrm{R}}^{\ell \times m},$$

$${\mathbf{b}}_i={\mathbf{J}}_i\cdot \mathbf{b}\text{and}{\mathbf{n}}_i={\mathbf{J}}_i\cdot \mathbf{n},$$

$$\mathbb{E}\left(\mathbf{N}\cdot {\mathbf{N}}^{\mathrm{H}}\right)$$

$$=\sum _{i=1}^k{\mathbf{J}}_i\cdot \mathbb{E}\left(\mathbf{n}\cdot {\mathbf{n}}^{\mathrm{H}}\right)\cdot {\mathbf{J}}_i^{\mathrm{H}}$$

$$=\sum _{i=1}^k{\mathbf{J}}_i\cdot {\mathbf{J}}_i^{\mathrm{H}}=k\cdot {\mathbf{I}}_{\ell }.$$

$$\mu (\mathbf{A})=\underset{1⩽i<j⩽m}{\max }\frac{⟨{\mathbf{a}}_i,{\mathbf{a}}_j⟩}{\Vert {\mathbf{a}}_i\Vert \Vert {\mathbf{a}}_j\Vert },$$

$$\mathrm{rk}{}^{\ast }\mathbf{X}\ge \frac{1}{\mu (\mathbf{X})}.$$

$$K_{\max }=\frac{1}{2}\left(1+\frac{1}{\mu (\mathbf{A})}\right).$$

$$\mu (\mathbf{A})=\mu (\mathbf{\Phi }\diamond \mathbf{\Psi })\le \mu (\mathbf{\Phi }\otimes \mathbf{\Psi })=\mu (\mathbf{\Phi })\mu (\mathbf{\Psi }).$$

$${\mathbf{\nu }}_n:{\mathrm{C}}^m\to {\mathcal{V}}^{n\times m}\text{with}\mathbf{z}\mapsto \left(\begin{array}{ccc}{z}_1& \dots & z_m\\ z_1^2& \dots & z_m^2\\ ⋮& ⋮& ⋮\\ z_1^n& \dots & z_m^n\end{array}\right)$$

$$\sum _{k=1}^n{q}^k=\frac{q-q^{n+1}}{1-q}=S_n(q)$$

$$\underset{q\to 1}{\lim }{S}_n(q)=\underset{q\to 1}{\lim }\frac{q-q^{n+1}}{1-q}=\underset{q\to 1}{\lim }\frac{1-(n+1)q^n}{-1}=n.$$

Full text

\DTLsetseparator

,

Sparsity Order Estimation from a Single Compressed Observation Vector

Sebastian Semper∗, , Florian Römer, , Thomas Hotz, Giovanni DelGaldo,

This work was partially supported by the Deutsche Forschungsgemeinschaft (DFG) project CoSMoS (grant GA 2062/2-1) and the Carl-Zeiss Foundation under the postdoctoral scholarship project “EMBiCoS”.

The authors S. Semper and T. Hotz are with Technische Universität Ilmenau, Institute of Mathematics, Ilmenau, Germany. The authors F. Roemer and G. Del Galdo are with Technische Universität Ilmenau, Institute for Information Technology, Ilmenau, Germany. The author G. Del Galdo is also with Fraunhofer Institute for Integrated Circuits (IIS), Ilmenau, Germany. $*$ corresponding author

Abstract —** We investigate the problem of estimating the unknown degree of sparsity from compressive measurements without the need to carry out a sparse recovery step. While the sparsity order can be directly inferred from the effective rank of the observation matrix in the multiple snapshot case, this appears to be impossible in the more challenging single snapshot case. We show that specially designed measurement matrices allow to rearrange the measurement vector into a matrix such that its effective rank coincides with the effective sparsity order. In fact, we prove that matrices which are composed of a Khatri-Rao product of smaller matrices generate measurements that allow to infer the sparsity order. Moreover, if some samples are used more than once, one of the matrices needs to be Vandermonde. These structural constraints reduce the degrees of freedom in choosing the measurement matrix which may incur in a degradation in the achievable coherence. We thus also address suitable choices of the measurement matrices. In particular, we analyze Khatri-Rao and Vandermonde matrices in terms of their coherence and provide a new design for Vandermonde matrices that achieves a low coherence. **

I Introduction

The Compressive Sensing (CS) paradigm states that a signal $\bm{x}\in\bm{\mathrm{C}}^{N}$ can be acquired using $M<N$ linear measurements, provided that it possesses a sparse representation with $K\ll N$ coefficients in a suitable basis $\bm{A}\in\bm{\mathrm{C}}^{N\times N}$ [3]. This observation has sparked considerable interest in many engineering disciplines devoted to signal acquisition, including imaging (e.g., MRI [13], X-Ray [23]), Radar [8, 1] and spectrum sensing [16, 28], among many others.

The number of non-zeros $K$ is also referred to as the sparsity order or sparsity degree. Even though the number of measurements required $M$ quite clearly depends on the sparsity order, the latter is typically not known when acquiring a signal. Therefore, one usually determines $M$ according to some upper bound $K\leq K_{\rm max}$ which must be known in order to design the measurement. However, since $K$ may vary quite significantly between observation periods, such an approach can lead to systems that are too conservative, i.e., they take more measurements than necessary. Therefore, being able to estimate and monitor the sparsity order of a signal would be an attractive feature of a measurement device. It would allow to adjust the number of measurements to the current sparsity of the signal which might vary in time (e.g., with the number of transmissions in cognitive radio [25] or with the complexity of an image in CS-based image acquisition [6]). Moreover, knowing the sparsity order allows to improve reconstruction algorithms by tuning algorithm-specific parameters such as the regularization parameter $\lambda$ in LASSO-type techniques (which is related to $K$ [26]) or stopping criteria in greedy techniques such as the Orthogonal Matching Pursuit [27].

I-A Related Work

Due to the prominent role the sparsity order plays in sparse signal recovery, the lack of knowledge of the sparsity order has been recognized as a fundamental gap between theory and practice [30, 15]. Still, existing literature on this subject is quite scarce and very recent. Early papers on this subject have proposed to employ sequential measurements [15] and cross-validation type techniques [2, 30] where sequential reconstructions of the signal are considered. Similarly, [29] shows that the sparsity order can be estimated from the reconstruction, stating bounds on the number of measurements that are required for this step. However, the bounds are only found numerically, and the reconstruction process involves cumbersome optimization problems. As we show, this can be avoided by estimating the sparsity order directly based on the compressed observations.

A different approach is taken in [12], where the authors show that a specifically tailored measurement procedure which consists of a Cauchy and a Gaussian distributed measurement matrix allows to estimate a continuous measure of sparsity given by the ratio of squared one- and two-norm of the signal. However, this measure is not equal to the sparsity order. In fact, it is continuous and hence needs to be rounded to an integer number (which is not discussed in [12]). Moreover, the measurement process is very restrictive since the distribution of the measurement matrices is pre-specified.

The authors in [19] propose to use sparse sensing matrices since these allow to infer the degree of sparsity of a signal from the degree of sparsity of the measurement. The resulting estimator has a very low complexity. However, it is only approaching the true sparsity in the large system limit and hence not applicable to lower-dimensional problems. Moreover, the proposed measurement matrices incur a certain performance degradation at the reconstruction stage due to the somewhat higher coherence, which is imposed by the restriction to sparse matrices.

A link between sparsity order estimation and rank estimation was put forward in [22] for the multiple measurement vector (MMV) problem, which we have also studied in our own prior work, both for the stationary case [9, 11] as well as the case of time-varying support for block-stationary signals [10]. However, these approaches require a certain stationarity in the support patterns as well as a temporal variation in the coefficients of the sparse representation to create linearly independent observations. This limits their applicability in many practical problems.

I-B Contributions

This paper introduces a method for estimating the sparsity order of a signal from a single snapshot of the compressive measurement. In particular, we propose to consider rearrangements of the observation vector into a matrix and show under which conditions the rank of this matrix coincides with the sparsity order of the unknown signal. Thereby, the sparsity order can be estimated by applying any known rank estimation scheme [24]. Since there exist many efficient algorithms which estimate the model order in presence of perturbations such as noise, this approach allows us to handle noisy measurements as well as the case of approximate sparsity as well.

As we show, the proposed approach only requires the measurement matrix to possess a Khatri-Rao structure, which leaves considerable room for optimizing their choice. Moreover, in the case of overlapping blocks, one of the factors needs to be a Vandermonde matrix. Therefore, the second part of the manuscript discusses the design of measurement matrices in presence of these structural constraints. To facilitate the sparsity order estimation as well as the sparse reconstruction, the factors need to possess a low coherence. As we show, in the presence of a Khatri-Rao structure, it is best to optimize the factors of the Khatri-Rao product independently. We therefore investigate the coherence of Vandermonde matrices and propose a new design algorithm, which efficiently constructs Vandermonde matrices with low coherence. Additionally, we derive simple upper and lower bounds for the resulting coherence of this algorithm.

Compared to the earlier conference version [20], this paper contains more rigorous versions of the main Theorems on the sparsity order estimation, along with their proofs (which had to be partially omitted in the conference version due to lack of space). Moreover, the optimization of the measurement matrix with respect to the structural constraints was added.

I-C Notation

Capital bold face letters denote matrices, whereas lowercase boldface letters stand for vectors. We use $[N]$ for $N\in\bm{\mathrm{N}}$ as a shorthand for the set $\{1,\dots,N\}$. The set $\operatorname*{supp}(\bm{x})$ is the set of indices where the vector $\bm{x}$ is non-zero, if $\operatorname*{supp}(\bm{x})=s$ we call $\bm{x}$ $s$-sparse. The vector $\bm{x}_{S}$ subscribed with an index set $S$ is the restriction of $\bm{x}$ to the indices $S$ whereas the matrix $\bm{X}_{S}$ denotes the matrix $X$ restricted to the columns indexed by $S$. Moreover, ${\left\|\cdot\right\|}_{p}$ stands for the $\ell_{p}$ norm where $p$ is omitted if it is clear from the context. The operators $\bm{A}\otimes\bm{B}$ and $\bm{A}\diamond\bm{B}$ refer to the Kronecker product and the Khatri-Rao (column-wise Kronecker) product, respectively. Finally, $\operatorname*{rk}^{*}\bm{A}$ represents the Kruskal rank of a matrix $\bm{A}$, i.e. the largest natural number $r$ so that every set of $r$ columns of $\bm{A}$ is linearly independent.

II Sparsity Order Estimation

II-A Data Model Setup

We consider the following discrete sparse recovery problem

[TABLE]

where $\bm{x}\in\bm{\mathrm{C}}^{N}$ is $K$-sparse, $\bm{S}\in\bm{\mathrm{C}}^{N\times N}$ represents the sparsifying basis, and $\bm{A}\in\bm{\mathrm{C}}^{m\times N}$ with $m\ll N$ is the measurement matrix. Since it is assumed that $\bm{S}$ is a basis (and hence invertible) we can limit ourselves to the case $\bm{S}=\bm{I}_{N}$ without loss of generality if $\bm{A}$ can be chosen freely since we can always replace $\bm{A}$ by $\bm{\bar{A}}=\bm{A}\cdot\bm{S}^{-1}$ to account for an $\bm{S}\neq\bm{I}_{N}$.

Moreover, for clarity we consider the noise-free case $\bm{n}=\bm{0}$ first. The role of additive noise is discussed in Section II-D.

In order to estimate the sparsity order $K$ directly from $\bm{b}$, we propose to consider rearrangements of $\bm{b}$ into a matrix. Specifically, we divide $\bm{b}$ into $k$ blocks $\bm{b}_{i}$ of length $\ell$ where the $i$-th block is given by

[TABLE]

for $i=1,2,\ldots,k$, where $p$ specifies by how many samples consecutive blocks are advancing. It is clear from (2) that for $p=\ell$ the blocks $\bm{b}_{i}$ do not overlap whereas for $p<\ell$ the overlap grows with decreasing $p$, up to the case of maximum overlap for $p=1$. Moreover, since $m$ samples are available, the parameters $\ell,p,k,m$ should satisfy $\ell+p\cdot(k-1)=m$.

The blocks $\bm{b}_{i}$ can be used to form the columns of a matrix $\bm{B}=\begin{bmatrix}\bm{b}_{1},&\bm{b}_{2},&\ldots,&\bm{b}_{\ell}\end{bmatrix}\in\bm{\mathrm{C}}^{\ell\times k}$. Moreover, we define submatrices $\bm{A}_{i}$ by selecting the corresponding rows from $\bm{A}$ such that

[TABLE]

The main idea of the proposed approach is to show that for a suitably chosen $\bm{A}$, we have $\operatorname*{rk}\bm{B}=K$ for any $K$-sparse $\bm{s}$ and therefore, the sparsity order can be inferred from the rank of $\bm{B}$. Since the required structure differs for the case of overlapping and non-overlapping blocks, we treat these two cases in separate subsections.

II-B Sparsity order estimation for non-overlapping blocks

The following theorem summarizes the conditions on $\bm{A}$ to facilitate the sparsity order estimation for the case of non-overlapping blocks ($p=\ell$).

Theorem II.1** ([20]).**

For $k$ blocks of the measurement $\bm{b}$ of length $\ell=m/k$ and any $r\leq\operatorname*{min}(k,\ell)$ the following statements are equivalent.

For all $s\leq r$ and all $\bm{x}$ with ${\left|\operatorname*{supp}\{\bm{x}\}\right|}=s$ we have $\operatorname*{rk}\bm{B}=s$. 2. 2)

$\bm{A}=\bm{\Phi}\diamond\bm{\Psi}$* for some $\bm{\Phi}\in\bm{\mathrm{C}}^{k\times N}$, $\bm{\Psi}\in\bm{\mathrm{C}}^{\ell\times N}$ with $\operatorname*{rk}^{*}\bm{\Phi}\geq r$ and $\operatorname*{rk}^{*}\bm{\Psi}\geq r$.*

Proof.

Assume 1) holds and consider an arbitrary but fixed $1$-sparse $\bm{x}$ with $\operatorname*{supp}(\bm{x})=\{q\}$ for some $q\in[N]$. Now 1) implies that $\operatorname*{rk}\bm{B}=1$. On the other hand, denoting the $q$-th column of $\bm{A}_{i}$ by $\bm{a}_{i,q}$, we have

[TABLE]

Because all $\bm{b}_{i}$ must be linearly dependent, we get

[TABLE]

for some non-zero $\bm{\psi}_{q}\in\bm{\mathrm{C}}^{\ell}$ and some $\varphi_{i,q}\in\bm{\mathrm{C}}$. In terms of the columns of $\bm{A}$ this means $\bm{a}_{q}=\operatorname*{vec}(\bm{\psi}_{q}\bm{\varphi}_{q}^{\mathrm{T}})=\bm{\varphi}_{q}\otimes\bm{\psi}_{q}$ for $\bm{\varphi}_{q}=(\varphi_{1,q},\dots,\varphi_{k,q})^{\mathrm{T}}\in\bm{\mathrm{C}}^{k}$. Since $q$ was arbitrary we conclude that $\bm{A}=\bm{\Phi}\diamond\bm{\Psi}$, where $\bm{\Phi}=(\bm{\varphi}_{i})_{i\in[N]}\in\bm{\mathrm{C}}^{k\times N}$ and $\bm{\Psi}=(\bm{\psi}_{i})_{i\in[N]}\in\bm{\mathrm{C}}^{\ell\times N}$.

Now consider any $s$-sparse $\bm{x}$ for $s\leq r$. Then from the structure of $\bm{A}$ and $\bm{B}$, we get that

[TABLE]

since $p=\ell$. This can also be restricted to $S=\operatorname*{supp}(\bm{x})$ and then reads as

[TABLE]

Now seeking a contradiction, assume that $\operatorname*{rk}^{*}\bm{\Phi}<r$. This means that there is a set $T\subset[N]$ of size $r$ such that $\bm{\Phi}_{T}$ is rank deficient. Now set $\bm{x}=\sum_{i\in T}\bm{e}_{i}$, which yields $\bm{B}=\bm{\Psi}_{T}\bm{\Phi}_{T}^{\mathrm{T}}$. The fact that $\bm{B}$ has rank strictly smaller than $r$ implies the desired contradiction. An analogue argument holds if we assume $\operatorname*{rk}^{*}\bm{\Psi}<r$.

For the opposite direction, assume 2) and again consider equation (4) for some $s$-sparse $\bm{x}$ with support $S$. Knowing that $\bm{\Phi}$ and $\bm{\Psi}$ have Kruskal rank at least $s$, we also have $\operatorname*{rk}\left[\bm{\Psi}_{S}\operatorname*{diag}(\bm{x}_{S})\right]=s$ and $\operatorname*{rk}\bm{\Phi}_{S}=s$. This implies that $\operatorname*{rk}\bm{B}=s$, because the two matrices are a rank factorization of $\bm{B}$. ∎

In other words, Theorem II.1 states that in the case of non-overlapping blocks, we can obtain the sparsity order $K$ from the rank of $\bm{B}\in\bm{\mathrm{C}}^{\ell\times k}$ if and only if the measurement matrix $\bm{A}$ possesses a Khatri-Rao structure. Moreover, the theorem also shows that the highest sparsity order we can estimate is given by $K_{\rm max}={\rm min}(\operatorname*{rk}^{*}(\bm{\Phi}),\operatorname*{rk}^{*}(\bm{\Psi}))\leq\min(k,\ell)$. Note that in the case of non-overlapping blocks we have $k\cdot\ell=m$ which means $K_{\rm max}\leq\min(k,\frac{m}{k})\leq\sqrt{m}$ and we should choose $k$ and $\ell$ close to $\sqrt{m}$ to maximize $K_{\rm max}$.

II-C Sparsity order estimation for overlapping blocks

As we have seen for non-overlapping blocks, if we let $\bm{A}$ be a Khatri-Rao product of two matrices of equal size, the maximum sparsity order that can be estimated is $K_{\rm max}\leqslant\sqrt{m}$. This bound is tight, iff $m$ is a square number $m=q^{2}$; then we can set $k=\ell=q$. If we let blocks overlap, the size of the matrix $\bm{B}$ grows, which allows to estimate larger sparsity orders $K_{\rm max}$. However, depending on the overlap, additional constraints on $\bm{A}$ have to be posed. This is formalized by the following theorem.

Theorem II.2** ([20]).**

For $k$ overlapping blocks of length $\ell$, block advance $p$ and any $r\leqslant\operatorname*{min}(k,\ell)$ the following statements are equivalent:

For all $s\leqslant r$ and all $s$-sparse $\bm{x}$ we have $\operatorname*{rk}\bm{B}=s$. 2. 2)

$\bm{A}$* consists of the first $m$ rows of $\bm{V}\diamond\bm{\Psi}$ with $\bm{\Psi}\in\bm{\mathrm{C}}^{p\times N}$ being arbitrary, $\bm{V}\in\bm{\mathrm{C}}^{\lceil m/p\rceil\times N}$ being a Vandermonde matrix such that the matrix $\hat{\bm{\Psi}}$ originating from restricting $\bm{V}_{\lceil\ell/p\rceil}\diamond\bm{\Psi}$ to its first $\ell$ rows has Kruskal rank $r$ and the matrix $\bm{V}$ restricted to its first $k$ rows also has Kruskal rank $r$.*

Proof.

Given 1), we again consider some $\bm{x}\in\bm{\mathrm{C}}^{N}$ with $\operatorname*{supp}(\bm{x})=\{q\}$ such that

[TABLE]

Since $\bm{B}$ has rank $1$, we know that its columns, which are also the columns of the blocks $\bm{A}_{i}$, fulfill

[TABLE]

for some $z_{q}$ independent of $i$. Summarizing, we get

[TABLE]

This means that we can only choose $z_{q}$ and because of the overlap the first $p$ elements of $\bm{a}_{1,q}$ freely. Now setting $\left(\bm{a}_{1,q}\right)_{[p]}=\bm{\psi}_{q}$ with $\bm{\psi}_{q}\in\bm{\mathrm{C}}^{p}$ then implies the Khatri-Rao structure $\hat{\bm{A}}=\bm{V}\diamond\bm{\Psi}$, where we get $\bm{A}$ by restricting $\hat{\bm{A}}$ to its first $m$ rows and $\bm{\Psi}=[\bm{\psi}_{1},\dots,\bm{\psi}_{N}]$. The Vandermonde structure of the first factor reads as

[TABLE]

if we repeat the reasoning above for any $q\in[N]$.

For the Kruskal ranks of the involved matrices we deduce from the Khatri-Rao structure and the overlap that

[TABLE]

where $\hat{\bm{\Psi}}$ consists of the first $\ell$ rows of $\bm{V}_{\lceil\ell/p\rceil}\diamond\bm{\Psi}$ and $\bm{V}_{\lceil\ell/p\rceil}$ and $\bm{V}_{k}$ are the restrictions of $\bm{V}$ to its first $\lceil\ell/p\rceil$ and $k$ rows respectively. Now consider the above equation for any $s$-sparse $\bm{x}$. As in the last proof, we argue that the factors involved in the Khatri-Rao product must fulfill the Kruskal rank conditions imposed on them in the statement of 2) for $\bm{B}$ to have rank $s$.

Now given 2), we proceed similarly to the proof of Theorem II.1, but this time, we get for an $s$-sparse $\bm{x}$ with $s\leqslant r$ that

[TABLE]

where $\bm{V}_{\lceil\ell/p\rceil}$ and $\bm{V}_{k}$ are the restrictions of $\bm{V}$ to its first $\lceil\ell/p\rceil$ and $k$ rows respectively. Because of the conditions on $\hat{\bm{\Psi}}$ and $\bm{V}$ we deduce as before that $\bm{B}$ has rank $s$. ∎

The assumptions in 2) in Theorem II.2 on $\bm{V}$ and $\bm{\Psi}$ seem rather strict and technical. In the following we will show, that they are fulfilled given some simple criteria. If we consider the Vandermonde matrix $\bm{V}$ it is enough to require that the $z_{1},\dots,z_{N}$ are pairwise distinct for $\bm{V}$ to have maximal Kruskal rank. For a square Vandermonde matrix $\bm{V}\in\bm{\mathrm{C}}^{n\times n}$ with generating elements $z_{1},\dots,z_{n}$ reads as

[TABLE]

This means that a Vandermonde matrix with more rows than columns and pairwise distinct generating elements always has full column rank. Moreover, if the matrix has more columns than rows, and pairwise distinct generating elements, every square submatrix has full rank and as such the whole matrix has full Kruskal rank.

Now turning to the arbitrary factor $\bm{\Psi}\in\bm{\mathrm{C}}^{p\times N}$, we can choose it in general position and in the following we show that this is sufficient for the requirements of Theorem II.2. To this end, fix some $r$-sparse $\bm{x}\in\bm{\mathrm{C}}^{N}$ with support set $R$, a Vandermonde matrix $\bm{V}\in\bm{\mathrm{C}}^{m\times N}$ which has full Kruskal rank and consider the Khatri-Rao product $\bm{V}\diamond\bm{\Psi}\in\bm{\mathrm{C}}^{m\cdot p\times N}$. Now assume that

[TABLE]

If we reshape the vector $\bm{m}$ according to the block construction in Theorem II.2 to

[TABLE]

which is a rank decomposition of $\bm{M}$ and this implies that for non-zero $\bm{x}$ the matrix $\bm{M}$ is non-zero.

Note that the special case $p=1$ in Theorem II.2 implies that the entire sensing matrix $\bm{A}$ is a Vandermonde matrix with rescaled columns. In the context of harmonic retrieval, the mapping from $\bm{b}$ to $\bm{B}$ is also known as spatial smoothing [21] and is applied as a preprocessing step for subspace-based estimators in order to decorrelate coherent signals.

Remark 1**.**

Regarding the choice of the parameters we observe that $K_{\rm max}\leq\min(k,\ell)$ where $k$ and $\ell$ satisfy $(k-1)\cdot p+\ell=m$. Therefore, for a given block advance $p$, the number of blocks $k$ is equal to $k=\frac{m-\ell}{p}+1$. This means that to maximize $K_{\rm max}$ we should choose $\ell$ as the closest integer to $\frac{m+p}{p+1}$, which leads to $k\approx\ell$. Obviously, larger amounts of overlap (corresponding to smaller values of $p$) lead to a higher maximum sparsity order $K_{\rm max}$ where the maximum overlap case $p=1$ corresponds to $K_{\rm max}=\lfloor\frac{m+1}{2}\rfloor$. This shows that there is a fundamental tradeoff between the sparsity order estimation (SOE) stage and the sparse recovery stage: while a larger amount of overlaps improves the SOE capability, it leads to a more rigidly structured measurement matrix with a higher coherence, which is detrimental to the sparse recovery step. The achievable coherence is analyzed in more depth in Section III.

II-D Sparsity Order Estimation in the presence of noise

In the presence of additive, say Gaussian, noise as in equation (1) the entries of the matrix $\bm{B}$ are disturbed with a Gaussian noise matrix $\bm{N}$ as well, i.e.

[TABLE]

It is easy to see that the rank of $\hat{\bm{B}}$ is full with probability $1$. It is still possible to determine the “effective” rank for a model like (11) if the statistics of the additive noise are known. Interestingly, although $\bm{N}$ is a reshaped version of the vector $\bm{n}$ in (1) with possibly the same noise sample at multiple positions (depending on the overlap), we can show that if the noise samples in $\bm{n}$ are i.i.d., the noise matrix $\bm{N}$ is “white” in the sense that $\mathbb{E}\left(\bm{N}\bm{N}^{\mathrm{H}}\right)=C\cdot\bm{I}_{\ell}$. To this end, let us assume the elements in $\bm{n}$ have zero mean and variance one and let us define the selection matrices

[TABLE]

which satisfies $\bm{J}_{i}\cdot\bm{J}_{i}^{\mathrm{H}}=\bm{I}_{\ell\times\ell}$. Then one can rewrite the measurement and noise vectors as

[TABLE]

where $i=[k]$. Now we calculate

[TABLE]

As we have shown, we need to determine the rank of the matrix $\bm{B}$ in presence of noise according to (11) (where in our case, the additive noise is white). This task is known as model order selection and a number of efficient algorithms are available. Examples include information-theoretic criteria such as MDL, AIC, BIC (see [24] for a survey), the Eigenvalue Threshold Test (ETT)[9] or the Exponential Fitting Test (EFT) presented in [18]. We make use of the latter for our numerical experiments in section IV, because it is derived specifically for models disturbed by additive white Gaussian measurement noise.

III Sensing matrix design

As we have seen in the last section, sparsity order estimation can be achieved via rank estimation of a matrix obtained by rearranging the measurement vector, provided that the sensing matrix obeys certain structural constraints. Firstly, it has to be a Khatri-Rao product of two smaller matrices and secondly, in the case of overlapping blocks, one of the blocks has to be a Vandermonde matrix.

In this section, we analyze the implications this particular sensing matrix structure has for the design of the measurement matrix, with a particular focus on the Vandermonde matrices.

III-A Khatri-Rao structured measurement matrix optimization

As we have seen in Theorem II.1 in the case of no overlap one is able to recover the sparsity order by using a Khatri-Rao structured sensing matrix, whose factors have to fulfill the condition of having a high Kruskal rank. This condition is difficult to optimize for, since the Kruskal rank is hard to compute. However, the following inequality links it to the coherence

[TABLE]

namely for any matrix $\bm{X}\in\bm{\mathrm{C}}^{p\times q}$, $p<q$ we have

[TABLE]

Therefore, for the Kruskal rank of the factors to be high, we aim for matrices with low coherence. These are desirable also from the viewpoint of the subsequent sparse recovery step, which in general works better the lower the coherence of the sensing matrix is. In fact, it has been shown that the under-determined system of equations $\bm{b}=\bm{A}\cdot\bm{s}$ has a unique solution for $K$-sparse vectors $\bm{s}$ if $K<K_{\rm max}$ with

[TABLE]

Note that Basis Pursuit (BP) as well as Orthogonal Matching Pursuit (OMP) are able to achieve this bound in the sense that they can recover any $K$-sparse $\bm{s}$ in the noise-free case as long as $K<K_{\rm max}$. For this reason we minimize the coherence of the measurement matrix $\bm{A}$.

Since $\bm{A}$ is a Khatri-Rao product, it might seem advantageous to take this structure into account when optimizing the coherence of the sensing matrix using the two factors. In the real valued setting, one can use packing arguments in projective matrix spaces of rank $1$ to optimize the coherence of Khatri-Rao products. As a result one obtains that the best coherence of a Khatri-Rao product $\bm{A}=\bm{\Phi}\diamond\bm{\Psi}\in\bm{\mathrm{R}}^{m\times N}$ is achieved if $\bm{\Phi}\in\bm{\mathrm{R}}^{m_{1}\times N}$ and $\bm{\Psi}\in\bm{\mathrm{R}}^{m_{2}\times N}$ contain repeated columns according to $\bm{\Phi}=\hat{\bm{\Phi}}\otimes\bm{1}_{1\times N_{2}}$ and $\bm{\Psi}=\bm{1}_{1\times N_{1}}\otimes\hat{\bm{\Psi}}$ where $\bm{\hat{\Phi}}\in^{m_{1}\times N_{1}}$, $\bm{\hat{\Psi}}\in^{m_{2}\times N_{2}}$ with $m=m_{1}\cdot m_{2}$ and $N=N_{1}\cdot N_{2}$. In this case we get $\bm{A}=\bm{\Phi}\diamond\bm{\Psi}=\hat{\bm{\Phi}}\otimes\hat{\bm{\Psi}}$ so that $\bm{A}$ is actually Kronecker structured. Due to the repeating columns, this implies that both $\bm{\Phi}$ and $\bm{\Psi}$ have Kruskal rank one, i.e. Khatri-Rao products of minimal coherence have factors that are not suitable for sparsity order estimation. Instead, one should optimize the coherence of the factors $\bm{\Phi}$ and $\bm{\Psi}$ independently.

This also has a positive effect on the coherence of $\bm{A}$ since we have the trivial upper bound

[TABLE]

Therefore, the remainder of this section discusses the minimization of coherence of the matrices $\bm{\Phi}$ and $\bm{\Psi}$.

In the case of no overlap these matrices can be arbitrary as there are no further structural requirements from the sparsity order estimation method. We can therefore apply any method for coherence minimization from the literature. In fact, there are already many iterative construction methods like in [5] and in coding theory this problem arises in constructing good spherical codes. An overview of theory and algorithms can be found in [31].

This means that the only regime where we can improve upon existing results is in the case of overlap during the construction of the matrix $\bm{B}$ as in Theorem II.2. Here, one factor obeys a Vandermonde structure. As long as the generating elements of the matrix $\bm{V}$ are pairwise different it is well known that this matrix has maximal Kruskal rank. But since we also need a low coherence for efficient recovery to happen, the following section contains a more thorough analysis of this problem.

III-B Coherence minimization for Vandermonde matrices

Vandermonde matrices are very rigidly structured, i.e. they can be generated from the knowledge of their first row.

Definition III.1**.**

Let $\mathcal{V}^{n\times m}(\bm{\mathrm{C}})$ be the space of $n\times m$ Vandermonde matrices then it holds that the non-linear mapping

[TABLE]

is bijective. So whenever we think about a Vandermonde matrix $\bm{V}$ it suffices to consider its first row.

In the following, we will call the $z_{i}$ the generating elements of $\bm{V}$.

This strong requirement on the structure of elements in $\mathcal{V}^{n\times m}$ will allow us to derive an explicit term for the inner product of two columns of a Vandermonde matrix that does not involve any summation and just depends on the amplitudes and phases of the generating elements of the respective columns. This is key, because of the definition of the coherence $\mu(\bm{A})$ of a matrix given in (13).

Further on we will be using the lemma below a couple of times, which calculates the value of the discontinued geometric series.

Lemma III.1**.**

For any $q\in\bm{\mathrm{C}}\setminus\{1\}$ it holds that

[TABLE]

In the case $q=1$ the above formula yields

[TABLE]

Proof.

This follows from a straightforward application of L’Hôpital’s rule. ∎

The following theorem contains all information needed to derive an algorithm for Vandermonde matrices with low coherence and bounds thereof. We treat the general case for two arbitrary Vandermonde columns but also the special case where the first row of the matrix is a subset of the complex unit circle. Furthermore it defines upper and lower envelope functions (see Figure 1) for the inner product of two Vandermonde columns depending only on the phase shift between two generating elements for fixed amplitudes.

Theorem III.1**.**

If we consider two vectors

[TABLE]

and

[TABLE]

where $c_{1}$, $c_{2}$ are positive real numbers and $\phi_{1}\leqslant\phi_{2}$ with $\phi_{1},\phi_{2}\in[0,2\pi]$ and define the function

[TABLE]

with $\phi=\phi_{2}-\phi_{1}$, then the following statements hold.

1. i)

[TABLE] 2. ii)

[TABLE] 3. iii)

We have symmetry along the unit circle in $\bm{\mathrm{C}}$. In other words

[TABLE] 4. iv)

For $c_{1}\neq 1$ and $c_{2}\neq 1$ we can bound $\lambda$ by

[TABLE]

where $\kappa$ and $\eta$ are defined as

[TABLE]

and

[TABLE]

Moreover, the bounds satisfy for $c_{1},c_{2}\in(0,1]$

[TABLE] 5. v)

For any $c>0$ and $c_{1}>0$ it holds that

[TABLE] 6. vi)

Roots of $\lambda(c_{1},c_{2},\cdot)$ exist if and only if $c_{1}c_{2}=1$ and these roots are $\phi=k\cdot(2\pi/n)$ for $k\in\{1,\dots,n-1\}$, i.e. there are exactly $n-1$ roots.

Proof.

1. i)

We calculate the inner product of $\bm{z}_{1}$ and $\bm{z}_{2}$ with the law of cosines and the formula in (16), which yields

[TABLE]

In this case we would have to consider L’Hôpital’s rule for $c_{1}c_{2}=1$ and $\phi=0$, which would yield that $c_{1}c_{2}e^{i\phi}=1$ and thus ${\left|{\left\langle\bm{z}_{1},\bm{z}_{2}\right\rangle}\right|}^{2}=n^{2}$. Because we want to derive an expression for the coherence of a matrix, we still need to divide by the norm of the vectors $\bm{z}_{1}$ and $\bm{z}_{2}$. Now using (16) again we get

[TABLE]

where we have to also consider L’Hôpital’s rule for $c_{1}=1$ or $c_{2}=1$. Then we would get ${\left\|\bm{z}_{1}\right\|}^{2}=n$ or ${\left\|\bm{z}_{2}\right\|}^{2}=n$ respectively. 2. ii)

If we make use of i) and take the appropriate limits, this yields

[TABLE] 3. iii)

We simply manipulate $\lambda(c_{1}^{-1},c_{2}^{-1},\phi)$ by factoring out the proper powers of $c_{1}$ and $c_{2}$ and get $\lambda(c_{1},c_{2},\phi)$. 4. iv)

The statement follows from the fact that $-1\leqslant\cos x\leqslant 1$ for all $x\in\bm{\mathrm{R}}$ and the locations where ${\left|\cos x\right|}=1$. 5. v)

If we assume that $c_{1}\cdot c_{2}=c<1$ then we can suppose without loss of generality that $c_{1}\leqslant c_{2}$. If we hold $c$ fixed, we can rewrite $\lambda$ as

[TABLE]

[TABLE]

[TABLE]

If we now take the limit for growing $c_{2}$ and fixed $\phi$, we get:

[TABLE] 6. vi)

Let us first consider the case $c_{1}\cdot c_{2}<1$. First we rewrite $\lambda$ in another form with the trigonometric identity $\sin^{2}x+\cos^{2}x=1$ by which we obtain

[TABLE]

From (21), we know that $\lambda$ is zero if and only if

[TABLE]

Let us assume that $\sin(n\phi)=0$, which implies that $\phi=k\cdot(\pi/n)$. Now we know that ${\left|\cos(k\pi)\right|}=1$ and the only case that $\cos(n\phi)-(c_{1}c_{2})^{n}=0$ is when $c_{1}c_{2}=1$ what we explicitly excluded. In other words for $c_{1}c_{2}$ different than $1$, we can never find two vectors $z_{1}$ and $z_{2}$ which are orthogonal to each other.

∎

Some results from above are visualized in Figure 1. The last Theorem allows us to derive a construction for orthogonal Vandermonde matrices, which can be considered a generalization of the fact that the Fourier matrix is orthogonal. The Fourier matrix is, depending on the definition, a scaled Vandermonde matrix, where the generating elements are placed on a regular grid on the unit circle, i.e. $\exp(2\pi ik/n)$, $k\in[n-1]$. The corollary below additionally allows a uniform phase shift of the generating elements.

Corollary III.1**.**

Let $z_{1},\dots,z_{n}$ be chosen with absolute value $1$ such that

[TABLE]

Then for the matrix $\bm{V}=\nu_{n}(z_{1},\dots,z_{n})$ it follows that

[TABLE]

Proof.

This is a direct consequence of vi) in Theorem III.1. ∎

But in the scenario typical for compressed sensing the involved matrices have more columns than rows and thus it is impossible for them to be orthogonal. Hence we are left with the objective to minimize the coherence as far as possible. The following algorithm, which is also one of the central novelties in this work, makes use of Theorem III.1 to achieve this.

The general geometric idea is twofold. First, the upper envelope $\kappa$ suggests that a large angular distance between generating elements yields a lower inner product, since it is monotonically decreasing on $(0,\pi)$. To maximize the mutual angular difference between each pair of generating elements, we need to place them on a regular angular grid on $[0,2\pi]$. Second, facts iii) and v) suggest that the inner product of two columns is minimized if their generating elements have amplitudes that are reziprocals of each other or, geometrically speaking, are reflections on the unit circle in $\bm{\mathrm{C}}$ of each other.

Combining the above two observations yields the algorithm as stated below. Loosely speaking it places the angularly nearest neighbors among the generating elements on different sides of the complex unit circle and the ratio between the amplitudes is the parameter $c$.

Algorithm III.1** (Vandermonde matrices).**

Input:

• Dimensions $n$ and $m$
• Positive constant $c$

Output:

• A Vandermonde matrix $\bm{V}\in\bm{\mathrm{C}}^{n\times m}$

Procedure:

Set $\hat{m}=2\lceil m/2\rceil$, $c_{1}=c$ and $c_{2}=1/cy$. 2. 2)

Set $\phi_{k}=4\pi\frac{k-1}{\hat{m}}$ and $z_{k}=c_{1}\cdot e^{i\cdot\phi_{k}}$ for $k\in[\hat{m}/2]$. 3. 3)

Set $\phi_{k}=4\pi\frac{k-1}{\hat{m}}+2\pi\frac{1}{\hat{m}}$ and $z_{k}=c_{2}\cdot e^{i\cdot\phi_{k}}$ for $k\in\{\hat{m}/2+1,\dots,\hat{m}\}$. 4. 4)

Return $\bm{V}=\nu(z_{1},\dots,z_{m})$.

The algorithm above still depends on the input parameter $c$, but the optimal value can be computed with a simple bisection algorithm and depends only on $n$ and $m$. This allows it to be stored in a persistent lookup table.

Additionally we would like to derive bounds on the coherence. This allows efficient estimation of the algorithm’s performance before its application. The bounds exploit the fact that the above algorithms output is deterministic and thus the upper and lower envelope $\kappa$ and $\eta$ can be applied. The following two Theorems deal with the upper and lower bounds respectively.

A lower bound on the coherence can be used together with equation (14) to derive lower bounds on the number of measurements required in order to achieve a recovery guarantee for a given sparsity order $K_{\rm max}$.

Theorem III.2**.**

Given $n$, $m>n$ and $c\in(0,1)$ and the output $\bm{V}$ of III.1 for inputs $n$, $m$ and $c$ it holds that

[TABLE]

Proof.

From Theorem III.1 we know that $\lambda(c_{1},c_{2},\phi)\geqslant\kappa(c_{1},c_{2},\phi)$ for all $c_{1}$, $c_{2}$ and $\phi$. Algorithm III.1 produces a matrix $\bm{V}$ so that the minimum angle between two elements in its first row is $2\pi/m$. Moreover those elements are located either on the circle with radius $c$ or the circle with radius $1/c$. If two elements share the same circle, their minimum angle between each other is $4\pi/m$. If not they enclose the angle $2\pi/m$. This yields the estimate

[TABLE]

which becomes what we asserted when we consider the symmetry of $\kappa$ as stated in fact iv) in Theorem III.1 and that $\lambda$ is the square of the inner product of two vectors with Vandermonde structure. ∎

Conversely we can use equation (14) and upper coherence bounds to derive upper bounds on the number of measurements required to guarantee successful recovery in scenarios with a certain level of sparsity.

Theorem III.3**.**

If we define the function $u\hskip 1.42262pt:\hskip 2.84526pt(0,1]\times\bm{\mathrm{N}}\rightarrow\bm{\mathrm{R}}$ by

[TABLE]

and call Algorithm III.1 with parameters $n$, $m$ and $c$, then for its output $\bm{V}$ it follows that:

[TABLE]

Proof.

Suppose we have an output $\bm{V}$ of Algorithm III.1 for given $c$, $n$ and $m$. We remind ourselves that the first row of $\bm{V}$ is placed on two circles with radii $c$ and $1/c$. Then we know that the angle between two elements in the first row of $\bm{V}$ which are not on the same circle is an integer multiple of $2\pi/m$. Elements on the same circle enclose an angle which is an integer multiple of $4\pi/m$.

Case 1: Suppose that $m<2n$. From the monotony of $\eta$ and its symmetry

[TABLE]

we obtain

[TABLE]

for $k\in\{1,\dots,\lfloor m/2\rfloor\}$ and

[TABLE]

for $k\in\{1,\dots,m\}$. Since the above estimates address all possible values $\lambda$ can take for $\bm{V}$ and the maximum over all values $\lambda$ attains the square of the coherence of $\bm{V}$, we have concluded the proof for $m<2n$.

Case 2: Suppose now that $4n\geqslant m\geqslant 2n$. From the monotony of $\eta$ we can again follow that

[TABLE]

Moreover, from $4n>m\geqslant 2n$ we get

[TABLE]

and estimate

[TABLE]

These estimates address all possible values for $\lambda$ and so the statement follows in this case.

Case 3: Suppose now that $m>4n$. From the monotony of $\lambda$ on $[0,2\pi/n]$ and the fact that $4\pi/m$ and $2\pi/m$ are left of the first touching point of $\lambda$ and $\eta$, namely $\pi/n$, we can estimate

[TABLE]

for $k_{1}\in\{1,\dots,\lfloor m/(4n)\rfloor\}$ and $k_{2}\in\{\lfloor m/(4n)\rfloor+1,\dots,\lfloor m/2\rfloor\}$. With a similar argument we estimate

[TABLE]

for $k_{1}\in\{1,\dots,\lfloor m/(2n)\rfloor\}$ and $k_{2}\in\{\lfloor m/(2n)\rfloor+1,\dots,m\}$. This again estimates all values $\lambda$ ever takes and we conclude the proof. ∎

IV Simulations

\DTLloaddb

vand_paramsmatlab/data/simu_vand_params.csv \DTLloaddbsoe_paramsmatlab/data/simu_soe_params.csv This chapter is dedicated to empirical investigations for showcasing the performance of the proposed methods. To this end, we consider measurement scenarios as in (1). In order to generate the ground truths $\bm{x}$, we draw the amplitudes on the support of $\bm{x}\in\bm{\mathrm{C}}^{N}$ for $N=512$ i.i.d. from the set $\{\pm 1\pm i,\pm 1\mp i\}$ according to a uniform distribution. We carry out all our simulations with $m$ noisy measurements where the components of the additive noise vector $\bm{n}\in\bm{\mathrm{C}}^{m}$ were drawn independently from a zero mean circularly symmetric Gaussian distribution with variance $\sigma^{2}$, where the total SNR is then defined as $\sigma^{-2}$. Depending on the situation at hand, the construction of $\bm{A}$ changes accordingly, which is chosen once for each scenario and kept fixed when sampling the ground truths for a certain scenario type and size. We carry out $2000$ trials for each scenario and level of SNR.

IV-A Sparsity Order Estimation

To depict the performance of the proposed method for SOE, we simulate the procedure for various combinations of parameters $m$ and $p$. Here, we set the sparsity order of the ground truth $\bm{x}\in\bm{\mathrm{C}}^{N}$ to $8$. After setting the number of measurements $m$ and the overlap $p$, we select $k$ and $\ell$ in order to maximize the size of the reshaped matrix $\bm{B}$ for optimal performance according to Remark 1. The Vandermonde factors in the columnwise Kronecker product are constructed by Algorithm III.1 and the unstructured factors are drawn once from a Gaussian ensemble, which means that all elements are drawn independently from a Gaussian distribution with variance $1$. Both factors in each of the occuring instances of a columnwise Kronecker product have normalized columns, resulting in the product having normalized columns as well.

To study the influence of the overlap we vary the parameter $p$ for fixed signal size $N$ and number of measurements $m$ while adapting $k$ and $\ell$ appropriately. For the estimation of the effective rank of $\bm{B}$, we apply the Empirical Eigenvalue-Threshold Test (ETT) [9] to our scenario with a target false rejection rate of $0.005$. The training stage of this model order selection method is the computationally most time consuming part, but this procedure only has to be done once in advance.

Figure 2 displays the empirical mean of the estimated sparsity order across all trials for a varying levels of the SNR. As expected for the pure task of SOE the case $p=1$ displays the best performance, because $\operatorname*{min}\{k,\ell\}$ is maximal. For increasing $p\geqslant 2$, including the case $p=\ell$, the phase transitions happens at an increasing level of SNR and the transition itself is not as sharp as in the cases of overlap.

To study a more realistic processing pipeline where the true sparsity order is not known in advance but just a more or less sharp upper bound, we combine the process of SOE with the reconstruction using OMP by feeding the estimated sparsity order into the algorithm as a runtime parameter. This is compared to an “unguided” OMP reconstruction, where we assume to only know upper bounds $K_{\mathrm{Max}}=40$ and $K_{\mathrm{Max}}=20$. Then OMP is run for $K$ steps.

As one can see in Figure 3, where the success rate of SOE and the resulting $l_{2}$ error after reconstruction are displayed in a combined plot, as soon as SOE starts to work, the reconstruction error decreases significantly below the one of the “unguided” reconstruction. So in this case the additional computational effort for SOE results in a significant increase in reconstruction precision.

As a last study, we investigate the reconstruction performance of the matrices that obey the structural constraints imposed by our means of SOE. To this end, we simulate a scenario with a priori known sparsity order to remove the influence of the SOE procedure on the reconstruction process. Here we expect differing performance since for decreasing $p$, the highly structured Vandermonde block in the sensing matrix becomes more dominant thus increasing the overall coherence of $\bm{A}$. The results are depicted in Figure 4. Here one can see that the case of no overlap $p=\ell$ comes close to the Gaussian measurement matrices despite the imposed columnwise Kronecker structure. Moreover one notices that the Vandermonde factor in the columnwise Kronecker product for the case of $p<\ell$ is the reason for a worse reconstruction performance.

IV-B Vandermonde Matrices

As we have seen before, in the case of block overlap our method of SOE requires that one of the factors has a Vandermonde structure and motivated by that, we derived Algorithm III.1 in order to get good reconstruction performance even with these rigid structural constraints. To conclude our numerical investigations we present simulation results that display the gain achieved through Algorithm III.1 compared to several other methods for constructing Vandermonde structured sensing matrices and also compared to unstructured matrices drawn from a Gaussian ensemble as a baseline, since they are well known for their good performance during reconstruction [4]. As a means of reconstruction, we present Orthogonal Matching Pursuit (OMP) [17], because it is a well known and well studied algorithm. Moreover, it is known that the coherence of the involved sensing matrix is involved in performance and stability results of OMP.

For comparison, we consider Gaussian sensing matrices constructed by drawing each entry identically and independently from a complex zero mean circular symmetric distribution. Then each column is normalized to unit length. We compared these Gaussian matrices with three types of Vandermonde matrices: (a) independently drawing their generating elements $z_{1},\dots,z_{N}$ from $Y\sim\exp(i2\pi U)$, where $U\sim\operatorname*{Unif}[0,1]$; (b) deterministically constructing the generating elements by placing them on a regular grid on the complex unit circle, so $z_{k}=\exp(i2\pi(k-1)/N)$; (c) according to Algorithm III.1. As in the case of the Gaussian matrices, we rescale each column of the Vandermonde matrices to unit length.

We also restrict the support sets of the underlying signals $\bm{x}$ such that the minimum distance between nonzero elements is larger than a certain quantity. To motivate this structural assumption, note that Vandermonde matrices have highly correlated columns, if their generating elements’ projections on the complex unit circle are closely located. This has already been established in Theorem III.1, where we show that

[TABLE]

Clearly, if one increases the number of generating elements for the mentioned methods (b) or (c) and fixes one arbitrary generating element, the distance with respect to the adjacent generating elements’ arguments decreases and thus the coherence of the resulting $\bm{V}$ gets closer to $1$. In terms of recovering a sparse vector $\bm{x}$ in this context means that if two distinct $i,j\in\operatorname*{supp}(\bm{x})$ are close to each other in the sense that ${\left|i-j\right|}/m\leqslant\varepsilon(n)$, then stably recovering this vector becomes impossible. This drawback originates from the structural impositions on $\bm{V}$ and as such cannot be circumvened.

In the context where the argument of the generating elements correspond to a parameter $\vartheta\in[0,2\pi)$ and the columns of $\bm{V}$ are atoms described by this parameter, e.g. [14], we can think of the number of columns $m$ of $\bm{V}$ as a granularity on $[0,2\pi)$. So increasing $m$ yields a finer grid on the parameter set, or in terms of the function $\lambda$ a finer discrete sampling. However, this does not decrease $\varepsilon(n)$ from above and as such does not reduce the distance between two parameters that are present in $\bm{x}$ and one still can resolve during reconstruction.

To include these observations in our simulations, we restrict the set of possible supports of the underlying signals $\bm{x}$. This is done by considering a Vandermonde matrix $\bm{V}\in\bm{\mathrm{C}}^{n\times m}$ with its $m$ generating elements regularly placed on the unit circle. Then we just sample support patterns such that the resulting ground truth $\bm{x}$ only contains columns such that their generating elements’ distances with respect to their arguments exceeds $2\pi/n$. This corresponds to making an assumption about the structure of the occuring signals and is known under the notion structured sparsity [7].

The results in Figure 5 display the effect of Algorithm III.1 on the rate of correct support detection during reconstruction. This is a valid performance measure for reconstruction schemes that produce exactly sparse results, like the OMP algorithm does, since it chooses the amplitudes on the detected support from the solution to the least squares problem on the span of the selected columns. The proposed method yields matrices that improve on the performance of the Vandermonde matrices with elements placed on a uniform grid on the complex unit circle. It is also clear why not all supports are correctly detected, since some of them have non separable elements due to the small distance of the corresponding generating elements of the Vandermonde matrices, which in the presence of noise results in reconstruction errors.

As a second performance metric we also measure the error in the $l_{2}$-norm between reconstructions and ground truths. These results can be found in Figure 6. Here, the numerical findings correspond to those above and the proposed method displays a lower reconstruction error than the other two means of constructing Vandermonde matrices. As expected, no instance of these highly structured matrices can compete with the performance displayed by matrices from the Gaussian ensemble, which do not have to fulfill any structural requirements.

V Conclusions

In this paper we investigate the problem of estimating the unknown sparsity order from compressive measurements without the need to carry out a sparse recovery stage. We focus on the more challenging case of a single compressed observation vector. As we show, the sparsity order estimation problem can be transformed into a rank estimation problem of a rearranged matrix version of the observation vector if and only if the compressed sensing measurement matrix has a Khatri-Rao structure. Moreover, if the columns of the matrix contain partially overlapping signal blocks, one of the matrices needs to possess a Vandermonde structure. A larger amount of overlap allows a higher sparsity order to be estimated but also leads to more stringent structural constraints on the measurement matrix.

The Khatri-Rao and Vandermonde constraints on the sompressed sensing measurement motivate us to analyze the suitable choice of the measurement matrices that allows the proposed sparsity order estimation and yet achieves a low coherence. In particular, we analyze the achievable coherence of Vandermonde matrices and propose a design that yields Vandermonde matrices with a low coherence.

Our numerical results demonstrate the trade-off between the maximal sparsity order that can be estimated and the coherence of the corresponding measurement matrices. Moreover, they clearly show the benefit of the proposed low-coherence Vandermonde matrix design.