Signal Reconstruction using Blind Super-resolution with Arbitrary Sampling

TL;DR

This paper introduces a novel semi-definite programming approach using Prolate Spheroidal Wave Functions for blind super-resolution of sparse signals, enabling spike localization without magnitude recovery, and demonstrates its effectiveness through numerical simulations.

Contribution

It presents a new SDP formulation leveraging PSWFs for blind super-resolution that localizes spikes without needing magnitude information, advancing prior methods.

Findings

Effective spike localization without magnitude recovery.

Outperforms recent related methods in simulations.

Applicable to arbitrary sampling schemes.

Abstract

In this paper the problem of blind super-resolution of sparse signals using arbitrary sampling scheme and atomic lift is discussed. After comprehensive description on blind superresolution problem, it is shown that using Prolate Spheroidal Wave Functions (PSWFs), it is possible to derive a new SemiDefinite Program (SDP) for the blind super-resolution problem. Unlike the previous results, the newly proposed SDP can localize spikes without magnitude recovery. Several numerical simulations were conducted to compare the performance of the proposed method with the recent related research.

Tables

Table 1. TABLE I: MSE of spike localization ( K = 4 , L = 3 , M = 50 formulae-sequence 𝐾 4 formulae-sequence 𝐿 3 𝑀 50 {K=4,L=3,M=50} ).

SNR	5	10	15	20	25	30
MSE	0.3262	0.3129	0.2560	0.1986	0.0941	0.0082

Equations

$$x(t)=\sum _{k=1}^K{a}_k\delta (t-{\bar{\tau }}_k)$$

$$y(t)=x(t)\ast g(t)=\sum _{k=1}^K{a}_{k}g(t-{\bar{\tau }}_k)$$

$$Y(f)=X(f)G(f)=\big{(}\sum_{k=1}^{K}a_{k}e^{-j2\pi f\bar{\tau}_{k}}\big{)}G(f)$$

$$Y(f_{m})=\big{(}\sum_{k=1}^{K}a_{k}e^{-j2\pi f_{m}\bar{\tau}_{k}}\big{)}G(f_{m}),\>m=1,...,M.$$

$$Y(\bar{f}_{m})=\big{(}\sum_{k=1}^{K}a_{k}e^{-j2\pi\bar{f}_{m}\tau_{k}}\big{)}G(\bar{f}_{m}),\>m=1,...,M.$$

$$\mathbf{y}=diag(\mathbf{g})\mathbf{x}$$

$$y_m$$

$$\mathbf{y}={\rm Y}(\mathbf{Z})$$

$$\mathbf{Z}=\mathbf{x}{\mathbf{h}}^T=\sum _{k=1}^K{a}_k\mathbf{c}({\tau }_k){\mathbf{h}}^T\in C^{M\times L}$$

$$\mathbf{c}({\tau }_k)=[e^{-j2\pi f_1{\tau }_k},e^{-j2\pi f_2{\tau }_k},\dots ,e^{-j2\pi f_M{\tau }_k}{]}^T$$

$$\mathcal{A}=\{\mathbf{A}(\tau ,\mathbf{u})=\mathbf{c}(\tau ){\mathbf{u}}^H\in C^{M\times L}|\tau \in [-\frac{1}{2},\frac{1}{2}],\Vert \mathbf{u}{\Vert }_2=1\}.$$

$$\Vert \mathbf{Z}{\Vert }_{\mathcal{A}}$$

$$=\underset{{\mathbf{u}}_k\in C^L:\Vert {\mathbf{u}}_k{\Vert }_2=1}{\underset{{\tau }_k\in [-\frac{1}{2},\frac{1}{2}]}{\inf }}\{\sum _k{c}_k|\mathbf{Z}=\sum _k{c}_k\mathbf{A}({\tau }_k,{\mathbf{u}}_k),c_k\ge 0\}$$

$$\underset{\mathbf{Z}}{\min }\Vert \mathbf{Z}{\Vert }_{\mathcal{A}}s.t.\mathbf{y}={\rm Y}(\mathbf{Z})$$

$$\Vert \mathbf{Z}{\Vert }_{\mathcal{A}}=$$

$$\displaystyle\>\>\>\>\>\>s.t.\left[\begin{array}[]{cc}\bm{W}&\bm{Z}^{\ast}\\ \bm{Z}&\bm{Q}\\ \end{array}\right]\succeq 0.$$

$$c=\pi {\mathbf{B}}_{max},$$

$$(\xi r)(f)={\int }_{-1}^1{e}^{jc\zeta f}r(\zeta )d\zeta ,\forall f\in [-1,1],$$

$$\xi {\phi }_j={\lambda }_j{\phi }_j$$

$${\mathbf{Q}}_{jl}={\mathbf{h}}_{jl}^T{\Phi }^{-1}[v_d^{\ast },\dots ,v_1^{\ast },v_0,v_1,\dots ,v_d{]}^T,$$

$${\mathbf{h}}_{jl}(k)$$

$$\mathbf{T}:$$

$$\Psi (\mathbf{T}):$$

$$-({\mathbf{J}}_1-{\mathbf{J}}_2)\mathbf{T}({\mathbf{J}}_1-{\mathbf{J}}_2{)}^{\ast },$$

$$\mu (\tau )=\sum _{k=1}^K{c}_{k}u(\tau -{\tau }_k),K\le d,|c_k|>0$$

$$\mathbf{Q}={\int }_{-\frac{1}{2}}^{\frac{1}{2}}\mathbf{P}(\tau )d\mu (\tau ),$$

$$\mathbf{T}=\sum _{k=1}^K{c}_k\omega (2\pi {\tau }_k{B}_{max}/d){\omega }^{\ast }(2\pi {\tau }_k{B}_{max}/d),$$

$$\underset{\mathbf{W}\in {\mathbb{R}}^{L\times L},v_0\in \mathbb{R},v_1,\dots ,v_d\in \mathbb{C}}{\mathrm{minimize}}$$

$$s.t.$$

$${\mathbf{Q}}_{jl}=$$

$$\mathbf{T}:=$$

$$\Psi (\mathbf{T}):=$$

$$-({\mathbf{J}}_1-{\mathbf{J}}_2)\mathbb{T}({\mathbf{J}}_1-{\mathbf{J}}_2{)}^{\ast }⪰0.$$

$$\displaystyle\left[\begin{array}[]{cc}\bm{W}&\bm{Z}^{\ast}\\ \bm{Z}&\bm{Q}\\ \end{array}\right]$$

$$Tr(\mathbf{W})\ge Tr({\mathbf{Z}}^{\ast }{\mathbf{Q}}^{-1}\mathbf{Z})=min\{Tr({\bar{\mathbf{Z}}}^{\ast }\bar{\mathbf{Z}}):\mathbf{Z}=\bar{\mathbf{Q}}\bar{\mathbf{Z}}\}$$

Full text

I Introduction

Many applications consider the problem of recovering a signal ${x(t)}$ convolved with a function ${g(t)}$ known as Point Spread Function (PSF)[1]. Without any assumption on ${g(t)}$ and ${x(t)}$ the problem is ill posed. Therefore, many methods have been proposed based on different assumptions [2, 3]. In practice, it is common to assume that the signal of interest is sparse and can be modeled with a sum of spikes. Though sparsity is an important assumption for signal recovery, it is not enough to solve the blind deconvolution problem successfully. Another important assumption is on PSF. Inspired by [2], it is possible to assume that the PSF lies in a known subspace but its orientation in the subspace remains unknown. This assumption provides us a mathematical model for PSF which is previously used in blind deconvolution problem [1].

Conventional sparsity based algorithms, which are known as blind compressed sensing methods, assume that frequencies lie on a predefined grid which is originated from the dictionary[4]. However, in many applications true frequencies can take any value and are not confined to the points on a grid. Therefore, basis mismatch problem arises[4]. Candés et.al introduced a new attitude towards sparse recovery which yielded sparse grid-less methods[5]. These methods achieved higher precision in frequency estimation and signal recovery than grid-based methods provided that a minimum frequency separation condition is satisfied. However, most of the recent sparse grid-less methods assume a grid based sampling scheme. Mahata and Hyder developed a sparse arbitrary sampling based frequency estimation method to gain higher precision in recovery of signals which are bandlimited and confined to a time window using Prolate Spheroidal Wave Functions (PSWFs)[6].

The problem of blind deconvolution of sparse signals in gridless compressed sensing regime via lifting was previously investigated[1]. In this paper we analyze the same problem with different sampling scheme and therefore different mathematical perspective. We arbitrarily sample the received signal and extend the mathematical theory of [6] to fit in blind gridless compressed sensing. Using the fact that the spikes of a sparse time signal are confined to a time window and the received signal is bandlimited due to PSF, we use PSWFs to express the variables of an optimization problem. This presentation leads to a new Semi-Definite Programming (SDP) to solve blind grid-less compressed sensing problem. Numerical results show that the proposed method achieves higher rate of successful spike recovery than recent methods.

The rest of the paper is as follow: In section II we formulate the blind sparse gridless compressed sensing problem. In section III we derive a SDP for the latter problem using arbitrary sampling scheme, PSWFs and atomic norm minimization. In section IV numerical results are presented to evaluate the performance of the new method.

II Problem Formulation

Consider signal ${x(t)}$ consisting of ${K}$ Dirac functions at positions ${\bar{\tau}_{k}\in[-\frac{T_{max}}{2},\frac{T_{max}}{2}]}$ with amplitudes ${a_{k}\in\mathbb{C}}$, where ${T_{max}}$ is signal duration in time. Therefore, we have

[TABLE]

where ${\delta(.)}$ is Dirac function. According to the fact that received signal is the result of convolution ${x(t)\ast g(t)}$, we have

[TABLE]

where ${g(t)}$ is the PSF. By taking the Fourier transform of (2) we have

[TABLE]

where ${X(f)}$ and ${G(f)}$ are Fourier transform of ${x(t)}$ and ${g(t)}$ respectively and ${j=\sqrt{-1}}$. Note that ${G(f)}$ is bandlimited to ${[\frac{-B_{max}}{2},\frac{B_{max}}{2}]}$. By sampling ${Y(f)}$ in arbitrary points ${f_{m}\in[\frac{-B_{max}}{2},\frac{B_{max}}{2}],m=1,\dots,M}$ we can write

[TABLE]

Consider ${\bar{f}_{m}=f_{m}T_{max}}$ and ${\tau_{k}=\bar{\tau}_{k}/T_{max}}$. Therefore,

[TABLE]

We can rewrite (5) in a compact form as

[TABLE]

where ${\bm{y}\in\mathbb{C}^{M\times 1}}$ is the observation vector, ${diag(\bm{g})}$ presents an ${M\times M}$ digonal matrix with ${\bm{g}}$ as its main diagonal and ${\bm{x}\in\mathbb{C}^{M\times 1}}$ is the Fourier transform of signal ${x(t)}$ at arbitrary sampled points ${\bar{f}_{m}}$.

In this paper, we assume that ${\bm{g}=\bm{S}\bm{h}}$ lies in a known low dimensional subspace ${\bm{S}\in\mathbb{C}^{M\times L}}$ with unknown orientation ${\bm{h}\in\mathbb{C}^{L}}$ and ${L\ll M}$. This assumption is not far from the reality and holds in many applications[1]. Consider the matrix ${\bm{S}^{T}=[\bm{s}_{0},...,\bm{s}_{M}]}$ where ${\bm{s}_{m}\in C^{L}}$ is the ${m}$-th column of the matrix ${\bm{S}^{T}}$. By defining ${x_{m}}$ and ${y_{m}}$ as the ${m}$-th elements of ${\bm{x}}$ and ${\bm{y}}$, we can use the lifting trick and write(6) as [2]

[TABLE]

where ${e_{m}}$ is the ${m-th}$ standard basis vector of ${\mathbb{R}^{M}}$. Therefore, we can consider ${\Upsilon:C^{M\times L}\rightarrow C^{M}}$ such that

[TABLE]

where

[TABLE]

and

[TABLE]

for ${\tau_{k}\in[-\frac{1}{2},\frac{1}{2}]}$. ${\bm{x}}$ and ${\bm{g}}$ can be recovered up to a scaling factor using left and right singular vectors of ${\bm{Z}}$. Hence, now the problem is to recover ${\bm{Z}^{*}}$ from ${\bm{y}}$. This problem is ill-posed because the number of unknowns is much more than the number of measurements. It is clear that ${\bm{Z}}$ columns are composed of ${K}$ complex sinusoids with the same frequencies. We use this joint spectral sparsity of ${\bm{Z}}$ columns to recover ${\bm{Z}^{*}}$ from ${\bm{y}}$ by minimizing the atomic norm for joint spectrally sparse signals[1].

III Atomic Norm Minimization and Arbitrary Sampling

To recover ${\bm{Z}^{*}}$ from ${\bm{y}}$, one should define the respective atomic set ${\mathcal{A}}$ as

[TABLE]

Atomic norm minimization is the problem of finding minimum number of atoms in ${\mathcal{A}}$ to define ${\bm{Z}^{*}}$ and can be expressed as [7]

[TABLE]

where ${Conv(\mathcal{A})}$ denotes the convex hull of the set ${\mathcal{A}}$. Therefore, we should solve the following problem

[TABLE]

where ${\bm{y}\in C^{M}}$ is the measurement vector and ${\Upsilon}$ is previously defined. There are many SDP characterizations of atomic norm minimization problem for different scenarios[7, 6, 8]. Most of the previous SDPs consider a uniform grid-based sampling scheme[8, 7]. As presented by [6], it is possible to derive a SDP characterization of atomic norm minimization problem based on arbitrary sampling scheme and PSWFs, which can provide us an orthogonal basis on ${\mathbb{L}^{2}}$ as the set of all square integrable functions on ${[-1/2,1/2]}$. We present the following lemma to cast (9) in a tractable manner.

Lemma 1

Consider ${\mathbb{H}}$ to be the set of all ${M\times M}$ Hermitian matrices and ${\bm{P}:[-\frac{1}{2},\frac{1}{2}]\rightarrow{\mathbb{H}}}$ be such that ${\bm{P}(\tau)=\bm{A}(\tau,\bm{u})\bm{A}^{\ast}(\tau,\bm{u})}$. Let ${\mathbb{K}}$ be closed conic hull of the set ${\{\bm{P}(\tau)\in\mathbb{H}:\tau\in[-\frac{1}{2},\frac{1}{2}]\}.}$ Then,

[TABLE]

Also, if ${\bm{Z}=\sum_{k}\hat{c}_{k}\bm{A}(\hat{\tau}_{k},\hat{\bm{u}}_{k})}$ is such that ${\|\bm{Z}\|_{\mathcal{A}}=\sum_{k}\hat{c}_{k}}$, then ${\bm{W}=\sum_{k}\hat{c}_{k}\hat{\bm{u}}_{k}\hat{\bm{u}}_{k}^{\ast}}$, ${\bm{Q}=\sum_{k}\hat{c}_{k}\bm{P}(\hat{\tau}_{k})}$ is a solution to (1). Conversely if, ${\bm{Q}=\sum_{k}\hat{c}_{k}\bm{P}(\hat{\tau}_{k})}$ is a solution to (1) for some ${\hat{c}_{k}>0}$, then the corresponding optimum value of ${\bm{W}}$ is ${\sum_{k}\hat{c}_{k}\hat{\bm{u}}_{k}\hat{\bm{u}}_{k}^{\ast}}$, ${\bm{Q}=\sum_{k}\hat{c}_{k}\bm{P}(\hat{\tau}_{k})}$, therefore ${Tr(\bm{W})=\sum_{k}\hat{c}_{k}=\|\bm{Z}\|_{\mathcal{A}}}$ and ${\bm{Z}=\sum_{k}\hat{c}_{k}\bm{c}(\hat{\tau}_{k})\hat{\bm{u}}_{k}^{\ast}.}$

The proof of Lemma 1 appears in appendix A. To present a tractable algorithm solving (1), a finite parametrization of ${\mathbb{K}}$ is needed. For ${\bm{Q}=\sum_{k}\hat{a}_{k}\bm{P}(\hat{\tau}_{k})}$ as a solution to (1), it is crucial that ${\bm{Q}\in\mathbb{K}}$. According to the fact that our signal is limited to a time window and also band-limited because of PSF, it is possible to use PSWFs as basis to express ${\bm{Q}\in\mathbb{K}}$. First, a brief description on PSWFs is given and then we use Lemma 2 from [6] to characterize ${\bm{Q}\in\mathbb{K}}$. Consider

[TABLE]

and let ${\mathbb{L}^{2}}$ be as previously defined. For any ${r\in\mathbb{L}^{2}}$, PSWFs are the eigenfunctions of the linear map ${\xi:\mathbb{L}^{2}\rightarrow\mathbb{L}^{2}}$ such that

[TABLE]

and therefore, the PSWF ${\varphi_{j}}$ should satisfy

[TABLE]

where ${\lambda_{j}}$ is the ${j-th}$ eigenvalue of ${\xi}$. Note that ${\varphi_{j}}$s are produced according to ${|\lambda_{j}|}$ where ${|.|}$ denotes the absolute value operator. Specifically, ${|\lambda_{j}|}$s reduce to zero when ${j}$ increases beyond ${2c/\pi}$. This means that we can limit the number of needed PSWFs according to our working precision. One can find ${d}$ by ${d=\{min\>\>\>\>\frac{j}{2}\>|\>\>\>\>|\lambda_{j}|<\epsilon\}}$ where ${\epsilon}$ is the working precision. More information on PSWFs can be found in [9].

Lemma 2** ([6])**

${\bm{Q}\in\mathbb{K}}$* if and only if there exists ${v_{0}\in\mathbb{R}}$ and ${v_{k}\in\mathbb{C},\>k=1,\dots,d}$ such that ${\bm{T}\succeq 0,\Psi(\bm{T})\succeq 0}$ satisfying*

[TABLE]

where

[TABLE]

Also, ${\theta_{0}=c/d,c=\pi B_{max}}$, ${\bm{I}_{d}}$ is a ${d\times d}$ identity matrix and ${\bm{0}_{d\times 1}}$ is a ${d\times 1}$ zero vector. In addition, if there exists a singular measure

[TABLE]

with ${u(\tau)}$ as the unit step function, satisfying

[TABLE]

then it is unique and ${\bm{T}}$ is singular with the Vandermonde decomposition

[TABLE]

where ${\omega(\theta):=[1,e^{-j\theta},\ldots,e^{-jd\theta}]^{\ast}}$.

Using Lemma 1 and Lemma 2, a new algorithm to solve (9) is obtained as

[TABLE]

According to the fact that ${T}$ has a toeplitz structure and is singular with rank ${K<d+1}$, one can recover ${\{\bar{\tau}_{k}\}_{1}^{K}}$ via Proney’s method[10]. Having ${\{\bar{\tau}_{k}\}}$s, ${\bm{x}}$ and ${\bm{g}}$ are estimated up to a scaling factor using left and right singular vectors of ${\bm{Z}}$. In the case of noisy measurements, one should add ${\frac{1}{\gamma}\|\hat{\bm{y}}-\bm{\Upsilon}(\bm{Z})\|_{2}^{2}}$ to the cost function of (15) and delete ${\bm{y}=\Upsilon(\bm{Z})}$. In this case, ${\hat{\bm{y}}=\bm{y}+\bm{n}}$ where ${\bm{n}}$ is generated using ${\mathcal{C}\mathcal{N}(0,\sigma^{2}\bm{I})}$, ${\gamma=\sigma\sqrt{Mln(B_{max})}}$ [10] and ${\bm{Z}}$ is added to the optimization variables. Note that our method does not utilize dual problem to localize spikes. However, the authors verify that spike recovery using dual polynomial derived from (15) is possible. According to the fact that sampling is not conducted on a grid, we are not dealing with orthogonal Fourier basis. Thus, dual certificate proofs proposed before [1] are not applicable anymore. This problem is due to the squared Fejér kernel’s high leakage when arbitrary sampling scheme is applied. Finding an appropriate kernel is left for future studies.

IV Numerical Results

To evaluate the performance of the proposed method numerical simulations were conducted. In these simulations we show that the proposed approach has a higher rate of successful spike recovery than the recent method of [1]. The simulations are divided to three sections. First, we show how the rate of successful spike recovery is affected by varying the number of spikes (${K}$). In sub-section IV-B, the dimension of the unknown PSF (${L}$) and in sub-section IV-C Signal to Noise Ratio (SNR) is changed to investigate the behaviour of the rate of successful spike recovery under different circumstances. In order to sample arbitrarily, ${f_{1}=0,f_{M}=B_{max}=64}$ and each ${f_{m},m=1,\ldots,M}$ is obtained uniformly from the interval ${(0,64)}$ while ${M<B_{max}}$. Spikes are separated by at least ${1/M}$. The subspace ${\bm{S}\in\mathbb{R}^{M\times L}}$ and the vector ${\bm{h}\in\mathbb{R}^{L}}$ are drawn from i.i.d standard Gaussian distribution. Considering ${\hat{\bm{Z}}}$ as the estimate of ${\bm{Z}}$, normalized mean square error (NMSE) ${\frac{\|\hat{\bm{Z}}-\bm{Z}^{\ast}\|_{F}^{2}}{\|\bm{Z}^{\ast}\|^{2}_{F}}}$ is used to compare the results. Also, mean square error (MSE) is used to compare spike localization accuracy. All of the results are compared to the state of the art [1]. The simulations were conducted in MATLAB using CVX toolbox[11]. In each part 100 Monte-Carlo simulations were conducted.

IV-A Number of the Spikes

In this part, the number of the spikes is varied from ${1\text{ to }12}$. The dimension of the PSF subspace is set to (${L=3}$). The number of the measurements is varied from 20 to 63. Figure 1 depicts the NMSE of matrix ${\bm{Z}}$ recovery for ${M=20,40,63}$ using the method of [1] and the proposed one. It is clear that the proposed method performs better than the analogous method. Note that the recovery of matrix ${\bm{Z}}$ is specifically important for estimating ${\bm{x}\text{ and }\bm{h}}$ after the spike localization. To demonstrate the performance of the proposed method in spike localization, phase transition diagrams for both methods are presented in figure 2. A successful recovery is considered as ${max\>\>\>|\hat{\tau}_{k}-\tau_{k}|<3\times 10^{-4},\>k=1,\ldots,K}$ where ${\hat{\tau}_{k}}$ is the estimation of ${\tau_{k}}$. As it can be seen in figure 2, the probability of successful spike localization is slightly higher when arbitrary sampling scheme is utilized.

IV-B Dimension of the Subspace of PSF

In this subsection, we investigate the effect of the dimension of ${\bm{h}\>(L)}$ on matrix ${\bm{Z}}$ recovery and spike localization. In the simulations, the subspace dimension is varied from ${1\text{ to }10}$ for ${M=30,50,63}$ and ${K=5}$. Figure 3 depicts the results. As it can be seen, the performance of the proposed method is slightly better than the method of [1]. To compare the rate of successful spike localization, figure 4 is depicted. This figure clearly suggests that the proposed method has higher success rate in spike localization than the method of [1].

IV-C Spike Localization under Noisy Environment

Considering more realistic scenarios, we investigated the performance of the proposed method in spike localization while the received signal is subjected to Gaussian noise ${\mathcal{N}(0,\sigma^{2}\bm{I})}$. The SNR is defined as ${10log(\|\Upsilon(\bm{Z})\|_{2}^{2}/(M\sigma^{2}))}$dB. Table 1 provides the MSE ${\frac{1}{N}\sum_{n=0}^{N-1}\|\hat{{\bm{\tau}}}_{n}-{\bm{\tau}}_{n}\|_{2}^{2}}$ where ${\hat{\bm{\tau}}_{n}}$ and ${{\bm{\tau}}_{n}}$ are the vectors of estimated and true spikes at ${n}$-th Monte Carlo simulation. As SNR increases, the MSE of spike localization of the proposed method decreases. Note that in noisy environments, the method of [1] can not localize spikes uniquely without magnitude recovery.

V Conclusion

We have proposed a new method for the problem of blind super-resolution using an arbitrary sampling scheme. Incorporating arbitrary sampling scheme in the problem generalized the work of [1]. More specifically, lemma 1 made it possible to use arbitrary samples in the problem. After characterizing the problem in a tractable manner, a new method for the blind super-resolution problem was derived. Simulation results revealed superior performance of the proposed method than that of [1]. According to the leakage of squared Fejér kernel while using arbitrary samples, the proof of the dual certificate was left for future work.

Appendix A

Proof of lemma 1: Consider the following linear matrix inequality (LMI)

[TABLE]

For any ${\bar{\bm{Q}}}$ satisfying ${\bm{Q}=\bar{\bm{Q}}\bar{\bm{Q}^{\ast}}}$ there exists a matrix ${\bar{\bm{Z}}}$ so that ${\bm{Z}=\bar{\bm{Q}}\bar{\bm{Z}}}$. Also, ${Tr(\bm{W})}$ constrained to (21) is such that[12]

[TABLE]

Considering ${\bm{W}_{\ast}}$ and ${\bm{Q}_{\ast}}$ as solutions to (15), we should have ${\bm{Q}_{\ast}\in\mathbb{K}}$. Thus,

[TABLE]

Note that this decomposition is not necessarily unique. Using the above description, there exist complex numbers ${\{\beta_{k}\}_{k=1}^{P}}$ such that

[TABLE]

According to (22) it is expected that ${|\beta_{k}|=\sqrt{\hat{c}_{k}}}$. Now we prove this claim. Consider ${|\beta_{k}|\neq\sqrt{\hat{c}_{k}}}$ and define

[TABLE]

Note that ${\|\bm{u}\|_{2}=1}$. Now it is easy to verify (21) for ${\hat{\bm{W}}}$ and ${\hat{\bm{Q}}}$. Thus ${\hat{\bm{W}}}$ and ${\hat{\bm{Q}}}$ are in the feasible set of (1). We have

[TABLE]

which contradicts the fact that ${\bm{W}_{\ast}}$ and ${\bm{Q}_{\ast}}$ are solutions to (1). Thus, ${|\beta_{k}|=\sqrt{\hat{c}_{k}}}$ and ${Tr(\bm{W}_{\ast})=\bm{e}^{\ast}\bm{Q}_{\ast}\bm{e}}$ so that ${\frac{1}{2}(Tr(\bm{W}_{\ast})+\bm{e}^{\ast}\bm{Q}_{\ast}\bm{e})=\sum_{k=1}^{P}\hat{c}_{k}}$. Now, (23) is an atomic decomposition of ${\bm{Z}}$ since ${\bm{A}(\tau_{k},\bm{u})=\bm{c}(\tau_{k})\bm{u}^{\ast}}$. According to (8), we conclude that

[TABLE]

To prove the other side of the inequality, we take a solution to (8). Considering ${\tilde{P}}$ positive numbers ${\{\tilde{c}(\tau_{k})\}_{k=1}^{\tilde{P}}}$ and spikes ${\{\tau\}_{k=1}^{\tilde{P}}}$ as the solution of (15), we have

[TABLE]

By taking ${\tilde{\bm{W}}=\sum_{k=1}^{\tilde{P}}\tilde{c}_{k}\bm{u}\bm{u}^{\ast}}$ and ${\tilde{\bm{Q}}=\sum_{k=1}^{\tilde{P}}\tilde{c}_{k}\bm{c}(\tilde{\tau}_{k})\bm{c}^{\ast}(\tilde{\tau}_{k})}$, one can verify that ${\tilde{\bm{Q}}\in\mathbb{K}}$ and (21) is positive semi-definite for ${\tilde{\bm{W}},\bm{Z}}$ and ${\tilde{\bm{Q}}}$. Thus,

[TABLE]

Using (24) and the above inequality, ${\|\bm{Z}\|_{\mathcal{A}}=\frac{1}{2}(Tr(\bm{W}_{\ast})+\bm{e}\bm{Q}_{\ast}\bm{e})}$.