In super-resolution, we are interested in recovering the high-end spectral information of signals from observations of its low-end spectral components [5]. In one setting of super-resolution problems, one aims to recover a superposition of complex exponential functions from time-domain samples. In fact, many problems arising in science and engineering involve high-dimensional signals that can be modelled or approximated by a superposition of a few complex exponential functions. In particular, if we choose the exponential functions to be complex sinusoids, this superposition of complex exponentials models signals in acceleration of medical imaging [25], analog-to-digital conversion [43] and array signal processing [38]. Accelerated nuclear magnetic resonance (NMR) spectroscopy, which is a prerequisite for studying short-lived molecular systems and monitoring chemical reactions in real time, is another application where signals can be modelled or approximated by a superposition of complex exponential functions [19,22,34]. How to recover the superposition of complex exponential functions or parameters of these complex exponential functions is of critical importance in these applications.
In this paper, we consider the recovery of superpositions of complex exponentials from linear undersampled measurements. More specifically, let x ∈ C 2N-1 be an unknown vector (which we would like to recover) whose j-th element satisfies
where c k ∈ C are complex coefficients and z k ∈ C, k = 1, . . . , R, are some unknown complex numbers with 2R ≤ 2N -1 being a positive integer. In other words, x is a superposition of R complex exponential functions. When z k = e 2πıf k , with ı = √ -1 and k = 1, . . . , R, x is a superposition of complex sinusoids. When z k = e -τ k e 2πıf k , with k = 1, . . . , R, and τ k > 0, x can model signals acquired by NMR spectroscopy in monitoring real-time chemical reactions and studying short-lived molecular systems [3,4]. The results of this paper apply to general τ k > 0, but, to simplify presentations, we choose to restrict certain theorems to τ k = 0. More specifically, τ k > 0 is allowed in Theorem 3.1 (we do not require knowledge of τ k to establish recovery guarantees), while τ k 's are assumed to be equal to 0 (thus known) in Theorems 4.1 and 4.3. In our simulations, we also set τ k to be equal to 0, where we consider the superposition of complex sinusoids.
Since 2R ≤ 2N -1 and often 2R 2N -1, the degree of freedom to determine x is much less than the ambient dimension 2N -1. Therefore, one can still recover x by its undersampling [8,16]. In particular, in undersampling, we assume that x is unknown, and we consider recovering x from its linear measurements b = A(x), (1.2) where A is a linear subsampling map from C 2N-1 to C M , M < 2N -1. After x is recovered, we can use the single-snapshot MUSIC [24] or Prony's method to recover the parameters z k . The problem of recovering x from its linear measurements (1.2) can be solved using compressed sensing [8], by discretizing the dictionary of basis vectors into grid points corresponding to discrete values of z k . When the parameters f k 's in signals from spectral compressed sensing (or the parameters (f k , τ k )'s from signals in accelerated NMR spectroscopy) fall on the grid, compressed sensing is a powerful tool to recover those signals even when the number of samples is far below its ambient dimension (M 2N -1) [8,16]. Nevertheless, the parameters in our problem setting often take continuous values, leading to a continuous dictionary, and may not exactly fall on a grid. The basis mismatch problem between the continuously valued parameters and the grid-valued parameters degrades the performance of conventional compressed sensing [14].
In two seminal papers [5,41], the authors proposed using total variation (TV) minimization or atomic norm minimization to recover x or to recover z k , when z k = e ı2π f k with f k taking continuous values from [0, 1). In these papers, the authors show that TV minimization or atomic norm minimization can recover the continuously valued frequency f k 's in the noise free case. However, as shown in [5,40,41], for TV minimization or atomic norm minimization to recover spectrally sparse data or the associated frequencies correctly, it is required that adjacent frequencies be sufficiently separated from each other. For example, for complex exponentials with z k 's taking values on the complex unit circle, it is required that adjacent frequencies f k ∈ [0, 1]'s be at least 2 (2N-1)×2π apart [40]. This separation condition is necessary for the TV minimization, even if we observe all the (2N -1) data samples, and even if the observations are noiseless. Please notice that the tolerable noise in correctly identifying frequencies to a certain accuracy depends on the frequency separation according to information-theoretic study [28].
This raises a natural question, 'Can we super-resolve the superposition of complex exponentials with continuously valued parameters z k , without requiring frequency separations, from compressed measurements, in noiseless or low-noise measurements?' In this paper, we answer this question in the affirmative. More specifically, we show that a Hankel matrix recovery approach using nuclear norm minimization can super-resolve the superposition of complex exponentials with continuously valued parameters z k from compressed measurements, without requiring frequency separations. This separation-free super-resolution result holds even when we only compressively observe x over a subset M ⊆ {0, ..., 2N -2} with cardinality |M| = M.
In this paper, we give the worst-case and average-case performance guarantees of Hankel matrix recovery in recovering the superposition of complex exponentials. In establishing the worst-case performance guarantees, we establish conditions under which Hankel matrix recovery can recover the underlying complex exponentials, no matter what values the coefficients c k 's of the complex exponentials take. For the average-case performance guarantee, we assume that the phases of the coefficients c k are uniformly distributed over [0, 2π). For both the worst and average cases, we establish that Hankel matrix recovery can super-resolve complex exponentials with continuously valued parameters z k , no matter how close they are to each other. We further introduce a new concept of orthonormal atomic norm minimization (OANM), and discover that the success of Hankel matrix recovery in separation-free superresolution comes from the fact that the nuclear norm of the Hankel matrix is an orthonormal atomic norm. In particular, we show that, in traditional atomic norm minimization, for successful signal recovery, the underlying parameters must be well separated, if the atoms continuously depend on continuously valued parameters; however, it is possible the OANM can succeed even if the original atoms are arbitrarily close. As a byproduct of this research, we discover one interesting matrix-theoretic inequality of nuclear norm, and give its proof using the theory of compressed sensing.
We remark that some of our results can require high percentage of all the samples for successful recovery using Hankel matrix recovery. For example, in Theorem 3.1, we need |M| > 2N -1 -N 2R out of 2N -1 samples to successfully recover a signal consisting of R atoms. Similarly for Theorems 4.1 and 4.3, we assume the number of samples |M| = 2N -2. We are also only able to show the tightness of Theorem 3.1, for a very special sampling scheme. Even for these special sampling schemes, to the best of our knowledge, our results are the first to show that Hankel matrix recovery can successfully recover superposition of complex exponentials regardless of their frequency or parameter separations. We also want to point out that these theoretical bounds can be conservative, compared with numerical results.
For example, in the experiments demonstrating how frequency separation can influence signal recovery performance, even though we only use |M| ≈ 51% × (2N -1) samples, the Hankel matrix recovery approach can already recover faithfully the superposition of a large number of complex sinusoids which have small frequency separations.
Low-rank Hankel matrix recovery approaches have been used for recovering parsimonious models in system identifications, control and signal processing. In [26], Markovsky considered low-rank approximations for Hankel structured matrices with applications in signal processing, system identifications and control. Fazel et al. [17,18] introduced low-rank Hankel matrix recovery via nuclear norm minimization, motivated by applications including realizations and identification of linear time-invariant systems, inferring shapes (points on the complex plane) from moments estimation (which is related to super-resolution with z k from the complex plane), and moment matrix rank minimization for polynomial optimization. Fazel et al. [18] further designed optimization algorithms to solve the nuclear norm minimization problem for low-rank Hankel matrix recovery. Chen and Chi [11] proposed to use multi-fold Hankel matrix completion for spectral compressed sensing, studied the performance guarantees of spectral compressed sensing via structured multifold Hankel matrix completion and derived performance guarantees of structured matrix completion. However, the results in [11] require that the Dirichlet kernel associated with underlying frequencies satisfies certain incoherence conditions, and these conditions require the underlying frequencies to be well separated from each other. In [15,44], the authors derived performance guarantees for Hankel matrix completion in system identifications. However, the performance guarantees in [15,44] require a very specific sampling pattern of fully sampling the upper-triangular part of the Hankel matrix. Moreover, the performance guarantees in [15,44] require that the parameters z k be very small (or smaller than 1) in magnitude. In our earlier work [3], we established performance guarantees of Hankel matrix recovery for spectral compressed sensing under Gaussian measurements of x. By comparison, this paper considers direct observations of x over a set M ⊆ {0, 1, 2, ..., 2N -2}, which is a more relevant sampling model in many applications. In the single-snapshot MUSIC algorithm [24], the Prony's method [33] or the matrix pencil approach [21], one would need the full (2N -1) consecutive samples to perform frequency identifications, while the Hankel matrix recovery approach can work with compressed measurements. When prior information of the locations of the frequencies is available, one can use weighted atomic norm minimization to relax the separation conditions in successful signal recovery [27].
In [39], the authors considered super-resolution without separation using atomic norm minimization, but under the restriction that the coefficients are positive and for a particular set of atoms (point spread functions) ψ(s, t) = e -(s-t) 2 , where s represents the spatial location in an image and t is the location of a point source of light (t's are the parameters to be estimated). In contrast, the analysis in our paper deals with complex-numbered coefficients, and does not require the coefficients to be positive. In addition, the techniques used for establishing the recovery guarantees in [39] are quite different from ours. In [39], the main idea for establishing the recovery guarantee is constructing a dual certificate, by applying the machinery of Tchebycheff systems to Gaussian point spread functions. By comparison, in this paper, we exploit the null space condition for nuclear-norm-minimization-based signal recovery. This machinery of investigating the null space condition does not require constructing a dual certificate, thus making it possible to show that the Hankel matrix recovery can reconstruct the signals without imposing separation conditions on the underlying atoms.
Additionally, in previous research [30][31][32], Prony's method and related modified algorithms were studied to recover parameters of signals, which are expressed as sum of exponentials. Specifically, in [31], the authors considered the recovery of signal parameters with equispaced sampled and real-valued data through the approximate Prony's method. In [32], the authors considered the signal parameter estimation problem for the sum of non-increasing exponentials with equispaced sampled data, and introduced the connections among various parameter estimate methods including the Prony's method, the matrix pencil method and the ESPRIT method. In [30], the authors considered the Prony's method with equispaced and non-equispaced sampled data for the signal parameter estimation. In particular, the considered nonequispaced sampled data case in [30] is similar to our problem considered in this paper, which is about signal parameter estimation with missing data. A major difference between [30] and ours is the different methods used in recovering missing entries: in [30], the authors use the method based on interpolation, while we use Hankel matrix completion in this paper. In [28], the author established a phase transition on the cut-off frequency at which noisy super-resolution is possible, and showed that there exists pairs of frequency component hypothesis that cannot be distinguished when the noise is big. However, the results in [28] focus on recovering the frequencies, and do not apply to recovering missing data samples, the recovery of which is shown in this paper to be robust against noises (please see Theorem 3.2). Furthermore, the results in [28] focused on finding worst-case pairs of hard-to-distinguish frequency component hypotheses under worst-case noises, but they do not exclude the possibility of recovering missing data samples or frequencies from randomly distributed noises or magnitude noises with small enough magnitudes.
The rest of the paper is organized as follows. In Section 2, we present the problem model, and introduce the Hankel matrix recovery approach. In Section 3, we investigate the worst-case performance guarantees of recovering spectrally sparse signals regardless of frequency separation, using the Hankel matrix recovery approach. In Section 4, we study the Hankel matrix recovery's average-case performance guarantees of recovering spectrally sparse signals regardless of frequency separation. In Section 5, we show that in traditional atomic norm minimization, for successful signal recovery, the underlying parameters must be well separated, if the atoms depend on the continuously valued parameters. In Section 6, we introduce the concept of OANM, and show that it is possible that atomic norm minimization is successful even though the original atoms are arbitrarily close. In Section 7, as a byproduct of this research, we provide a new matrix-theoretic inequality of nuclear norm from the theory of compressed sensing. Numerical results are given in Section 8 to validate our theoretical predictions. We conclude our paper in Section 9. The proofs of the technical lemmas are in the appendix.
We denote the set of complex numbers and real number as C and R, respectively. We use calligraphic uppercase letters to represent index sets, and use | • | to represent the cardinality of a set. When we use an index set as the subscript of a vector, we refer to the part of the vector over the index set. For example, x Ω is the part of vector x over the index set Ω. We use C n 1 ×n 2 r to represent the set of matrices from C n 1 ×n 2 with rank r. We denote the trace of a matrix X by Tr(X), and denote the real and imaginary parts of a matrix X by Re(X) and Im(X), respectively. The superscripts T and * are used to represent transpose, and conjugate transpose of matrices or vectors. The Frobenius norm, nuclear norm and spectral norm of a matrix are denoted by • F , • * and • 2 (or • ), respectively. The notation • represents the spectral norm if its argument is a matrix, and represents the Euclidean norm if its argument is a vector.
The probability of an event S is denoted by P(S). A vector or matrix with all its entries zero will be denoted by 0. The inner product for matrices is defined as M, N = Tr(M T N) when M and N are real valued, or M, N = Re(Tr(M * N)) when M, N are complex valued.
We consider the signal model in (1.1) with
where f k ∈ [0, 1) are normalized frequencies, τ k > 0 are the damping terms and the observation set is
To estimate the continuous parameter f k in [5,41], the authors proposed using TV minimization or atomic norm minimization to recover x or to recover the parameter z k , when z k = e ı2π f k with f k taking continuous values from [0, 1). However, as shown in [5,40,41], in order for atomic norm minimization to recover spectrally sparse data or the associated frequencies correctly, it is necessary that adjacent frequencies be separated far enough from each other. As shown in [40], for complex exponentials with z k 's taking values on the complex unit circle, it is required that adjacent frequencies f k ∈ [0, 1) be at least 2 2π(2N-1) apart. This separation condition is necessary, even if we observe the full (2N -1) data samples, and even if the observations are noiseless. A formal statement of minimal separation of frequencies can be found in Definition 2.1. Definition 2.1. ( [5]) For a frequency subset F ⊂ [0, 1) with a group of points, the minimum separation is defined as smallest distance between two arbitrary different elements in F, i.e. dist(F) = inf
where d(f i , f l ) is the wrap around distance between two frequencies.
Following the idea of the matrix pencil method in [21] and Enhanced Matrix Completion (EMaC) in [11], we construct a Hankel matrix based on signal x. More specifically, define the Hankel matrix
3)
The expression (1.1) and (2.3) lead to a rank-R decomposition:
Instead of reconstructing x directly, we reconstruct the rank-R Hankel matrix H, subject to the observation constraints. Low-rank matrix recovery has been widely studied in recovering a matrix from incomplete observations [2,6,7,9,35]. It is well known that minimizing the nuclear norm can lead to a solution of low-rank matrices. We therefore use the nuclear norm minimization to recover the low-rank matrix H.
which can be easily solved with existing convex program solvers such as interior point algorithms.
After successfully recovering all the time samples, we can use the single-snapshot MUSIC algorithm (as discussed in [24]) to identify the underlying frequencies f k . For the recovered Hankel matrix H(x), let its singular value decomposition (SVD) be
and we define the vector φ N (f ) and imaging function J(f ) as
The single-snapshot MUSIC algorithm is given in Algorithm 1. In [24], the authors showed that the MUSIC algorithm can exactly recover all the frequencies by finding the local maximal of J(f ). Namely, for undamped signal (1.1) with the set of frequencies F and
While it is true that when R < N 2 , one can use MUSIC to find the frequencies from the first N -2 samples and then recover the signals. However, we focus on the nuclear norm for the following reasons: Algorithm 1 Single-Snapshot MUSIC algorithm 24 1: require: solution x, parameter R and N 2: form Hankel matrix Z ∈ C N×N
(1) MUSIC does not apply to the case where R > N 2 if we use the first N -2 samples. However, the nuclear norm minimization can be more powerful than MUSIC in that it can recover more frequencies. For example, our results in Theorem 4.1 show that nuclear norm minimization can recover close to N frequencies, while MUSIC can only recover up to N 2 if it uses only the first N -2 samples. (2) In addition, this paper focuses on using nuclear norm minimization, rather than MUSIC. So while the MUSIC can recover the same number ( N 2 ) of frequencies in the worst case by taking the first N -2 samples, our results can still be useful in understanding the power of nuclear norm minimization.
In this section, we provide the worst-case performance guarantees of Hankel matrix recovery for recovering the superposition of complex exponentials. Namely, we provide the conditions under which the Hankel matrix recovery can uniformly recover the superposition of every possible R complex exponentials. Our results show that the Hankel matrix recovery can achieve separation-free superresolution, even if we consider the criterion of worst-case performance guarantees. Later, we further show that our derived worst-performance guarantees are tight, namely we can find examples where nuclear norm minimization fails to recover the superposition of complex exponentials for a larger R. We start with introducing necessary notations and a few lemmas needed for proofs of our main results.
Let w i be the number of elements in the i-th anti-diagonal of any N × N matrix like H, namely,
Here, we call the (2N -1) anti-diagonals of H(x) from the left top to the right bottom as the first antidiagonal up to the (2N -1)-th anti-diagonal. For example, consider a Hankel matrix
then correspondingly, w 1 = 1, w 2 = 2, w 3 = 3, w 4 = 4, w 5 = 3, w 6 = 2 and w 7 = 1. We also define w min as
where M is the observation set. For a Hankel matrix
we define Ind(j) as the set of indices of rows which intersect with the j-th anti-diagonal of the Hankel matrix, and define {j : i ∈ Ind(j)} as the set of indices of all non-zero anti-diagonals which intersect with the i-th row. For example, consider the matrix in (3.2), we have With the setup above, we give the following Theorem 3.1 concerning the worst-case performance guarantee of Hankel matrix recovery, the proof of which is based on the strong null space conditions from [29,36]. We will specifically use the following Lemma 3.1 from [29].
Lemma 3.1 is called the strong null space condition, and it characterizes the conditions under which the nuclear norm minimization (3.5) can recover any matrix with rank at most R. We use this lemma to show below that the recovery of H(x) from the observation set M is guaranteed by a condition that does not mandate sufficient frequency separation. Lemma 3.1. ( [29]) Consider an N × N matrix X 0 of rank at most R, and a linear mapping A(X 0 ) = b. Then the nuclear norm minimization (3.5)
can uniquely and correctly recover every matrix X 0 with rank no more than R if, for all non-zero Z ∈ N(A),
where Z * R is the sum of the largest R singular values of Z, and N(A) is the null space of A.
When there is noise in the observation data, we have the following results guaranteeing the robustness of the nuclear norm minimization in recovering data. Lemma 3.2. Consider an N × N matrix X 0 of rank at most R, and a linear mapping B such that B(X 0 + E) = b, where E is a perturbation matrix with E * ≤ . Then the nuclear norm minimization (3.7)
will output a X such that X 0 -X * ≤ 3C+1 C-1 , if for all non-zero Z ∈ N(B),
where C > 1 is a constant, Z * R is the sum of the largest R singular values of Z, and N(B) is the null space of B.
Proof. We let X = X 0 + E + W, where W is a non-zero matrix from the null space of B. Then
where W is a matrix from the null space of B, σ i (•) represents the i-th largest singular value of a matrix and inequality (3.14) is due to X 0 's rank being at most R and Lemma 3.3 as shown below. Thus, we have 2
where σ i (X) and σ i (Y) are, respectively, the i-th largest singular values of X and Y, respectively.
Because of the two lemmas above, we arrive at the following lemma particularly about the condition for Hankel matrix recovery in super-resolution. Lemma 3.4. Consider the signal model (1.1) with z k defined in (2.1), and the observation set M ⊆ {0, 1, 2, ..., 2N -2}. Then the nuclear norm minimization (2.4) will uniquely recover H(x), regardless of the (frequency) separation between the R continuously valued (frequency) parameters, if for every non-zero element in the null space N(A) of the linear mapping A corresponding to the sampling set M satisfies the condition in Lemma 3.1. Namely, every non-zero matrix in the set
This lemma shows that the guarantee of recovering missing data via Hankel matrix recovery essentially does not require frequency separations of frequencies, but instead depends on the sampling pattern only.
We are now ready to introduce our main results. Theorem 3.1. Consider the signal model (1.1) with z k defined in (2.1), and the observation set M ⊆ {0, 1, 2, ..., 2N -2}. We further define w i of an N × N Hankel matrix H(x) as in (3.1), and define w min as in (3.3). Then the nuclear norm minimization (2.4) will uniquely recover H(x), regardless of the (frequency) separation between the R continuously valued (frequency) parameters, if
On the one hand, the performance guarantees given in Theorem 3.1 can be conservative: for averagecase performance guarantees, even when the number of complex exponentials R is bigger than predicted by Theorem 3.1, the Hankel matrix recovery can still recover the missing data, even though the sinusoids can be very close to each other. On the other hand, the bounds on recoverable sparsity level R given in Theorem 3.1 is tight for worst-case performance guarantees, as shown in the next section. We will show the tightness of recoverable sparsity guaranteed by Theorem 3.1, which is built on the developments in Section 4.1.
Proof of Theorem 3.1: We can change (2.4) to the following optimization problem:
We can see (3.18) as a nuclear norm minimization problem, where the null space of the linear operator (applied to (B,x)) in the constraints of (3.18) is given by (H(z), z) such that A(z) = 0. Using Lemma 3.1, we can see that (3.18) or (2.4) can correctly recover x as a superposition of R complex exponentials if 2 H(z) * R < H(z) * holds true for every non-zero z from the null space of A.
From the sampling set M ⊆ {0, 1, 2, ..., 2N -2}, the null space of the sampling operator A is composed of (2N -1) × 1 vectors z's such that z i+1 = 0 if i ∈ M. For such a vector z, let us denote the element across the i-th anti-diagonal of H(z) as a i , and thus Q = H(z) is a Hankel matrix with its (i + 1) anti-diagonal element equal to 0 if i ∈ M. Let σ 1 , • • • , σ N be the N singular values of H(z) arranged in a descending order. To verify the null space condition for nuclear norm minimization, we would like to find the largest R such that
for every non-zero z in the null space of A.
Towards this goal, we first obtain an upper bound for its largest singular value (for a matrix Q, we use Q i,: to denote its i-th row vector):
where Ind(j) is the set of indices of rows which intersect with the j-th anti-diagonal, {j : i ∈ Ind(j), a j = 0} is the set of all non-zero anti-diagonals intersecting with the i-th row and (3.20) is obtained by looking at the N rows. Note that in the calculations above, we only need to pay attention to the non-zero entries in the i-th row of Q. From the Cauchy-Schwarz inequality, we have
Because u 2 = 1 and, for each i, |u i | 2 appears in i∈Ind(j 1 ) j 2 ∈{j:i∈Ind(j),a j =0} |u j 2 -i+1 | 2 for no more than (2N -1 -M) times, we have
2N-1
Furthermore, summing up the energy of the matrix Q, we have
, then for any integer k ≤ N, we have
then for every non-zero vector z in the null space of A, and the corresponding Hankel matrix
It follows that, for any superposition of R < w min 2(2N-1-M) complex exponentials, we can correctly recover x over the whole set {0, 1, ..., 2N -2} using the incomplete sampling set M, regardless of the separations between different frequencies (or between the continuously valued parameters z k 's for damped complex exponentials).
We further have the following results of the robustness of the Hankel matrix recovery in super-resolution. Theorem 3.2. Consider a ground-truth signal model (1.1) with z k defined in (2.1), and the observation set M ⊆ {0, 1, 2, ..., 2N -2}. Suppose that the ground-truth signal is given by
and we denote the ground-truth vector x 0 = [x 0 0 , x 0 1 , ..., x 0 2N-2 ] T . We further represent the noise-corrupted signal by
where e = [e 0 , e 1 , ..., e 2N-2 ] T is the noise vector. We define a partial-observation error vector ẽ ∈ C 2N-1 such that ẽM = e M and ẽ{0,1,2,...,2N-2}\M = 0. We define = N ẽ 1 , and denote the measurement results by
Then the optimization formulation
will output a x such that
and
where C > 1, regardless of the (frequency) separation between the R continuously valued (frequency) parameters, if
Proof. The proof is similar to the proof of Theorem 3.1. Instead of using Lemma 3.1, we use Lemma 3.2.
First of all, we observe that the nuclear norm of the error matrix E = H(ẽ) is upper bounded by N × ẽ 1 . This is because the Hankel matrix, which has non-zero elements only in the i-th (1 ≤ i ≤ 2N -1) anti-diagonal and zeros elsewhere, has its nuclear norm bounded by the number of non-zero elements (which is upper bounded by N) times the magnitude of that non-zero element. By the triangle inequality, we have the nuclear norm of the error matrix E = H(ẽ) upper bounded by N × ẽ 1 .
We let X 0 = H(x 0 ) and E = H(ẽ). Then the optimization formulation
where B(•) is a linear mapping such that its null space is {H(x) :
Then we can apply Lemma 3.2 to obtain the error bounds as in this theorem. In fact, following similar proofs as in the proof of Theorem 3.1, we have, if
then for every non-zero vector z in the null space of A, and the corresponding Hankel matrix Q = H(z), we have
Applying Lemma 3.2, we immediately have
where • denotes the Frobenius norm, we have
Let us now consider the recovery of the underlying frequencies. In fact, the robustness of recovering missing data samples translates into the robustness of recovering the frequencies. When the observations are noiseless, for the first time, we have shown that the nuclear norm minimization method can recover the data correctly regardless of frequency separation, and moreover, using the recovered data samples and the Prony's method (or the single snapshot MUSIC algorithm), we can even exactly recover the frequency accurately regardless of frequency separation. When the observations are noisy, the Hankel matrix recovery method can guarantee recovering the missing data samples with robustness: regardless of frequency separations, the norms of the recovery error are nicely bounded by the norms of the noise vector, as proved by Theorem 3.2 in this paper. According to theorems 3 and 4 in [24], we know that when the error in the recovered data samples goes to zero, using the single-snapshot MUSIC algorithm, the error in the noise-space correlation function used in the single-snapshot algorithm also goes to zero; and, moreover, there exist local minimizers of the noise-space correlation function converging to the true underlying frequencies. Of course, the ratio of the magnitude of frequency recovery error to the magnitude of noise depends on frequency separations; however, at least this method is robust against noises in the sense that the magnitude of frequency recovery error goes to zero as the noise's magnitude goes to zero. This is very different from atomic norm minimization, which will identify wrong frequencies even with noiseless and complete data if the separation of two frequencies is below a certain threshold [40]. Our results in this paper demonstrate that Hankel matrix recovery method exhibits very different behaviours than the atomic norm minimization, even though both reduce to a semidefinite programming problem for minimizing the nuclear norm of a matrix.
In this section, we study the performance guarantees for Hankel matrix recovery, when the phases of the coefficients of the R sinusoids are i.i.d. and uniformly distributed over [0, 2π). For average-case performance guarantees, we show that we can recover the superposition of a larger number of complex exponentials than Theorem 3.1 offers. We introduce Lemmas 4.1 and 4.2 first, which are used in the derivations of the average-case recovery guarantee stated in Theorems 4.1 and 4.3.
and
then
, ∀t ≥ 0. (4.3) Lemma 4.1, called the Matrix Bernstein Inequality, characterizes the tail behaviour of the spectral norm of a random matrix. Lemma 4.2. Let X 0 be any M ×N matrix of rank R in C M×N , and we observe it through a linear mapping A(X 0 ) = b. We also assume that X 0 has an SVD X 0 = UΣV * , where U ∈ C M×R , V ∈ C N×R , and Σ ∈ C R×R is a diagonal matrix. Then the nuclear norm minimization (4.4)
correctly and uniquely recovers X 0 if, for every non-zero element
where Ū and V are such that [U Ū] and [V V] are unitary.
We use Lemma 4.2 as the condition for successful signal recovery through nuclear norm minimization. This lemma is an extension of lemma 13 in [29]. The key difference is that Lemma 4.2 deals with complex-numbered matrices. Moreover, Lemma 4.2 gets rid of the 'iff' claim for the null space condition in lemma 13 of [29], because we find that the condition in lemma 13 of [29] is a sufficient condition for the success of nuclear norm minimization, but not a necessary condition for the success of nuclear norm minimization [46]. We give the proof of Lemma 4.2 in Appendix.
Theorem 4.1. Let us consider the signal model of the superposition of R complex exponentials (1.1) with z k defined in (2.1) where τ k = 0. We assume that the R frequencies f 1 , f 2 ,..., and f R are such that the atoms (1, e ı2π f i , ..., e ı2πf i (N-1) ) T , 1 ≤ i ≤ R, are orthogonal to each other. We let the observation set be M = {0, 1, 2, ..., 2N -2} \ {N -1}. We assume the phases of coefficients c 1 , • • • , c R in signal model (1.1) are independent and uniformly distributed over [0, 2π). Then the nuclear norm minimization (2.4) will successfully and uniquely recover H(x) and x, with probability approaching 1 as
where c > 0 is a constant. We can think of (4.7) as a nuclear norm minimization, where the null space of the linear operator (applied to B and x) in the constraints of (3.18) is given by (H(z), z) such that A(z) = 0.
J. YI ET AL.
Then we can see that (4.7) or (2.4) can correctly recover x as a superposition of R complex exponentials if
hold true for any z which is a non-zero vector from the null space of A, where H(x) = UΣV * is the SVD of H(x) with U ∈ C N×R and V ∈ C N×R , and Ū and V are such that [U, Ū] and [V, V] are unitary. Without loss of generality, let
where s k 's are distinct integers between 0 and N -1. Then
and
and
where θ k 's (1 ≤ k ≤ R) are i.i.d. random variables uniformly distributed over [0, 2π). When the observation set M = {0, 1, 2, ..., 2N -2} \ {N -1}, any Q = H(z) with z from the null space of A takes the following form:
Notice that random variable e ıθ i e -ı2π s i (N-1) N are mutually independent random variables uniformly distributed over the complex unit circle. We will use the concentration lemma (see for example, [42]) to provide a concentration of measure result for the summation of e ıθ i e -ı2π s i (N-1) N .
Applying Lemma 4.1 to the 1×2 matrix composed of the real and imaginary parts of e ıθ i e -ı2π s i (N-1)
R+t/3 , ∀t > 0.
(4.16)
We further notice that, Ū's ( V' s) columns are normalized orthogonal frequency atoms (or their complex conjugates) with frequencies l i N , with integer l i 's different from the integers s k 's of those R complex exponentials. Thus,
since we can always find two unitary matrices Q 1 and Q 2 such that Q 1 QQ 2 becomes an identity matrix without changing the nuclear norm, i.e. Building on the discussions above, we now discuss some results on the tightness of the bounds in Theorem 3.1.
We first show that the recoverable sparsity R provided by Theorem 3.1 is tight for M = {0, 1, 2, ..., 2N -2} \ {N -1}. For such a sampling set, w min = N, and Theorem 3.1 provides a bound R < w min 2 = N 2 .
In fact, we can show that if R ≥ N 2 , we can construct signal examples where the Hankel matrix recovery approach cannot recover the original signal x. Consider the signal in Theorem 4.1. We choose the coefficients c i 's such that e ıθ i e -ı2π s i (N-1)
and we also have
with z in the null space of the sampling operator.
Thus, the Hankel matrix recovery cannot recover the ground truth x. This shows that the prediction by Theorem 3.1 is tight.
According to Theorem 3.1, in the worst case, we need a very large number of samples to guarantee successful recovery of signal comprising a small number of atoms. We have established the tightness of Theorem 3.1 for a special sampling scheme, i.e. M = {0, 1, • • • , 2N -1} \ {N -1}, and even under the restricted sampling patterns considered in this paper, these are the first results in the literature that show successful recovery is guaranteed regardless of frequency separation using nuclear norm minimization in the setting of noiseless measurements.
In fact, we have the following general results suggesting that the tightness of the results in Theorem 3.1 in terms of recovering a low-rank matrix. Then there is a non-zero vector u ∈ C 2N-1 ∈ such that A(u) = 0 and the corresponding Hankel matrix H(u) has a rank of w min ; moreover, there is a matrix Q ∈ C N×N of rank w min /2 such that
In other words, let B be linear mapping from C N×N to C |M| such that its null space is the same as the set {H(u) | A(u) = 0}, then there exists a low-rank matrix Q of rank w min /2 such that it is not the unique solution to the following nuclear norm minimization problem:
Proof. We consider the index i, 0 ≤ i ≤ (2N -2) such that w i+1 = w min . Then we take the low-rank matrix Q which has w min /2 1's in its (i + 1)-th anti-diagonal, and all the other elements in Q are 0. We take u ∈ C N such that u i = -1 at index i but equal to 0 at all the other indices. Then H(u) + Q will be a matrix having w minw min /2 '-1's at its (i + 1)-th diagonal, and have a no bigger nuclear norm than Q, proving this theorem.
We also want to emphasize that while results in our theorems are conservative, and the practical performance can still be quite encouraging. For example, in Theorem 4.3, the sampling scheme can hardly be regarded as undersampling since we need 2N-2 2N-1 = 126 127 ≈ 99% of the samples for successful recovery. But in the experiments, we actually used about 51% of samples. The experimental results actually imply the possibility of improving the bounds in the recovery guarantees.
In this section, we further show that the Hankel matrix recovery can successfully recover the superposition of complex exponentials with arbitrarily close frequency atoms. In particular, we give average performance guarantees on recoverable sparsity R when frequency atoms are arbitrarily close, and show that the Hankel matrix recovery can deal with much larger recoverable sparsity R for average-case signals than predicted by Theorem 3.1. For the tractability of our theoretical analysis, our results assume that certain frequency atoms are orthogonal to each other; however, in our results, there are always frequency atoms which are arbitrarily close to each other, enforcing small frequency separation in this research problem. We note that in the literature the separation requirement is often defined as the smallest circular separation between frequency atoms. Following this definition, the frequency separations used in this section are still consistent with the definition of small frequency separation, even though some of the frequency atoms are orthogonal to each other. Numerical results suggest that the assumptions that some of the frequency atoms are orthogonal to each other are not necessarily needed; however, removing such assumptions would require some future work. We consider a signal x composed of R complex exponentials. Among them, R-1 orthogonal complex exponentials (without loss of generality, we assume that these R -1 frequencies take values l i N , where l i 's are integers). The other complex exponential c cl e ı2π f cl j has a frequency arbitrarily close to one of the R -1 frequencies, i.e.
where f cl is arbitrarily close to one of the first (R -1) frequencies. For this set-up with arbitrarily close atoms, we have the following average-case performance guarantee. Remarks: 1. We can extend this result to cases where neither of the two arbitrarily close frequencies has on-the-grid frequencies, and the Hankel matrix recovery method can recover a similar sparsity; 2. Under a similar number of complex exponentials, with high probability, the Hankel matrix recovery can correctly recover the signal x, uniformly over every possible phases of the two complex exponentials with arbitrarily close frequencies; 3. We can extend our results to other sampling sets, but for clarity of presentations, we choose M = {0, 1, 2, ..., 2N -2} \ {N -1}.
The proof of Theorem 4.3 depends on Lemma 4.3, a perturbation result on the polar decomposition of a matrix, which we will present first before introducing the formal proof of Theorem 4.3. Let us consider two matrices X and X from C m×n r such that X = X + E, where E ∈ C m×n . Suppose that the matrix X has SVD
where
The polar decomposition of X is given by
and the matrix P is the unitary polar factor. Similarly for X, we have
To prove this theorem, we consider a perturbed signal x. The original signal x and the perturbed x satisfy
where x j is the same as in (4.21), and x is defined as
with f rm being a frequency such that its corresponding frequency atom is mutually orthogonal with the atoms corresponding to f 1 ,..., and f R-1 . We further define
Let us define X = H(x), X = H(x), and the error matrix E such that
where
Both X and X are rank-R matrices. For the error matrix E, we have
Following the derivations in Theorem 4.1, X has the following SVD
where Ũ1 ∈ C N×R , Σ1 ∈ C R×R and Ṽ1 ∈ C N×R are defined as in (4.9), (4.10) and (4.11). The polar decomposition of X using its SVD is given by
J. YI ET AL.
and the matrix P is called the unitary polar factor. Similarly, for X, we have
and its polar decomposition for X,
Suppose that we try to recover x through the following nuclear norm minimization:
As in the proof of Theorem 4.1, we analyse the null space condition for successful recovery using Hankel matrix recovery. Towards this, we first bound the unitary polar factor of X through the perturbation theory for polar decomposition.
For our problem, we have m = n = N, r = R. Let σ R be the R-th singular value of X. From (4.11) and (4.31), an explicit form for σ R is σ R = Nd min .
(4.40)
Then Lemma 4.3 implies that
From the null space condition for nuclear norm minimization (4.39), we can correctly recover x if Let us define ΔP = P -P = Ũ1 Ṽ * 1 -U 1 V * 1 , and we have
where we used the Cauchy-Schwartz inequality. Notice that E = X -X, and then it follows from (4.41) and (4.34) that
Thus if
namely,
where d rel is defined as
, then solving (4.39) will correctly recover x. The proof of Theorem 4.1 leads to the concentration inequality:
To have successful signal recovery with high probability, we let
where c > 0 is a constant. So we have
where
)
and c 1 = 2 3 c. Since
solving the quadratic equation leads to
Plugging in c 1 = 2 3 c, we have
In the previous sections, we have shown that the Hankel matrix recovery can recover the superposition of complex exponentials, even though the frequencies of the complex exponentials can be arbitrarily close. In this section, we show that, broadly, the atomic norm minimization must obey a non-trivial resolution limit. This is very different from the behaviour of Hankel matrix recovery. Our results in this section also greatly generalize necessity of resolution limit results in [40], to general continuously parametered dictionary, beyond the dictionary of frequency atoms. Moreover, this new lower bound can work for nonuniform sampling, where some samples at certain time indices can be missing, while the lower bound from [40] applies to uniform sampling. Our new lower bound further applies to samples obtained at noninteger time indices, while the lower bound from [40] does not apply to sampling at non-integer time indices. In terms of technical tools used, our analysis of using matrix analysis is very different from the derivations in [40], which used the Markov-Bernstein type inequalities for finite-degree polynomials. Theorem 5.1. Let us consider a dictionary with its atoms parameterized by a continuous-valued parameter τ ∈ C. We also assume that each atom a(τ ) belongs to C N , where N is a positive integer. We assume that the set of all the atoms spans a Q-dimensional subspace in C N .
Suppose the signal is the superposition of several atoms:
with non-zero c k ∈ C for each k, and R is a positive integer representing the number of active atoms. Consider any active atoms a(τ k 1 ) and a(τ k 2 ). With the other (R-2) active atoms and their coefficients fixed (this includes the case R = 2, namely there are only two active atoms in total), if the atomic norm minimization can always identify the two active atoms, and correctly recover the coefficients for a(τ k 1 ) and a(τ k 2 ), then the two atoms a(τ k 1 ) and a(τ k 2 ) must be well separated such that
where S is a positive integer, M S is the set of matrices with S columns and with each of these S columns corresponding to an atom, and σ min (•) is the smallest singular value of a matrix.
Proof of Theorem 5.1: Define the sign of a coefficient c k as
Then according to [41], a necessary condition for the atomic norm to identify the two active atoms, and correctly recover their coefficients is that there exists a dual vector q ∈ C N such that
We pick c k 1 and c k 2 such that |sign(c
Thus
Now we take a group of S atoms, denoted by a select,1 ,..., and a select,S , and use them to form the S columns of a matrix A. Then
where the last inequality comes form (5.4). It follows that √ S among all the choices for A and S, namely,
(5.9)
Table 1 Comparison between [40] and ours in lower bounds on frequency separations [40] 0.1592 0.0796 0.0531 0.0398 0.0318 0.0265 0.0227 0.0199 0.0177 Our bound 0.1667 0.0628 0.0351 0.0232 0.0167 0.0128 0.0102 0.0084 0.0071 Combining (5.6) and (5.8), we have
(5.10) proving this theorem. We would like to remark that, besides applications to frequency atoms, our bounds from Theorem 5.1 on atom separations are quite general, and can be applied to any dictionary of atoms with continuous parameters. The derivations in [40] apply only to dictionary of frequency atoms under complete (without missing samples) uniform sampling at integer time indices. When applied to frequency atoms, our bound from Theorem 5.1 works for non-uniform sampling, where some samples are missing at certain time indices, and works for sampling at non-integer time indices (for example, sampling at time indices j = √ 1, √ 2, √ 3,..., √ 2N -2; sampling at non-integer times indices may happen because of sampling time jittering), where the bound in [40] does not apply. When further restricted to uniformly sampled-atinteger-time-index frequency atoms, our calculations show that Theorem 5.1 provides similar or even tighter lower bounds on frequency separations than in [40]. In addition, our derivations are novel, and do not resort to the Markov-Bernstein-type inequalities for polynomials as used in [40]. It is possible to further tighten our method of deriving lower bound using matrices, and we will leave it for future work.
In the following calculations, we consider frequency atoms given in the following form (atoms will change for non-uniform sampling or sampling at non-integer time indices),
(5.11)
For two frequencies f a and f b , we call a(f a )a(f b ) the atom separation (Theorem 5.1 can directly produce a lower bound on atom separation, while [40] works with frequency separation), and call |f af b | wrap (wrap-around distance) the frequency separation. We note that, in [40], the author established a lower bound of 2 (2N-2) * 2π for frequency separation (see Theorem 3 in [40]).
In this case, the samples are observed at time indices j = 0, 1, 2, ..., (2N -2). For this sampling pattern, under small N, our bound can be tighter than the lower bound given in [40], and, under large N, the lower bound in [40] is tighter.
Table 1 shows comparison result between [40] and ours on lower bound on frequency separation for N = 2, ..., 10.
Lower bound from [40] when 2N -2 = 2: In this case (three time-domain samples in total), the lower bound on frequency separation in [40] is 0.159. We construct two atoms with frequencies f a = 0 and f b = 0.159, and then calculate the atom separation. This gives an atom separation of 1.94.
New lower bound when 2N -2 = 2: To calculate our own lower bound, we consider a set of orthogonal atoms (any set of linearly independent atoms will give us a valid lower bound). The frequencies of these atoms are taken to be f 1 = 0, f 2 = 1 3 and f 3 = 2 3 . Our Theorem 5.1 gives a lower bound of 2.00 for a(f a )a(f b ) . To obtain this lower bound based on Theorem 5.1, we take matrix A = [a(f 1 ) a(f 2 ) a(f 3 )], and S = 3 since the three atoms span a three-dimensional linear space. Using these parameters Theorem 5.1 gives a(f a )a(f b ) ≥ 2.00. This corresponds to a frequency separation |f af b | wrap = 0.167. In this case, the bound given by Theorem 5.1 is indeed tighter than the frequency separation bound (0.159) from [40].
Lower bound from [40] when 2N -2 = 62: When 2N -2 = 62 (totally 63 time-domain samples respectively), the frequency separation requirement in [40] will be 0.0051. We construct two atoms with frequencies f a = 0 and f b = 0.0051, and then calculate the distance. This gives an atom distance of 8.31.
New lower bound when 2N -2 = 62: To calculate our own lower bound, we consider a set of orthogonal atoms (any set of linearly independent atoms will give us a valid lower bound). The frequencies of these atoms are taken to be
Our Theorem 5.1 gives a bound of 2.00 for a(f a )a(f b ) . To obtain this lower bound based on Theorem 5.1, we take matrix A = [a(f 1 ) a(f 2 ) • • • a(f 63 )], and S = 63 since the 63 atoms span a 63-dimensional linear space. Theorem 5.1 gives a(f a )a(f b ) ≥ 2.00. This corresponds to a frequency separation |f af b | wrap = 1.1 × 10 -3 , which is less tight than, but comparable with the frequency separation bound from [40].
Lower bound from [40] when 2N -2 = 126: When 2N -2 = 126 (totally 2N -1 = 127 timedomain samples, respectively), the frequency separation required [40] is 0.0025. We simply construct two atoms with frequency f a = 0 and f b = 0.0025, and then calculate the distance between these two atoms. This gives atom distance of 11.78.
New lower bound when 2N -2 = 126: To calculate our own lower bound, we consider a set of orthogonal atoms (any set of linearly independent atoms will give us a valid lower bound). The frequencies of these atoms are taken to be f 1 = 0, f 2 = 1 127 , f 3 = 2 127 , • • • , f 127 = 126 127 . Our Theorem 5.1 gives a bound of 2.00 for a(f a )a(f b ) . We take matrix A = [a(f 1 ) a(f 2 ) • • • a(f 127 )], take S = 127 since the 127 atoms form a 127-dimensional space, and Theorem 5.1 gives a(f a )a(f b ) ≥ 2.00. This corresponds to a frequency separation |f af b | wrap = 4 × 10 -4 , which is less tight than, but still comparable with the frequency separation bound from [40].
Theorem 5.1 can still provide valid lower bounds on frequency or atom separations for non-uniform sampling. The lower bound from [40] does not apply to non-uniform sampling directly, but a simple modification can yield a lower bound 2 N max ×2π on frequency separation, where N max is the largest observed time index.
As one example, we consider N = 256 corresponding to 2N -1 = 511 total time indices, and we assume that samples are taken at only five time indices: j = 0, 1, 2, 3 and 124. In this case, the atom at a certain frequency f ∈ [0, 1) is given by a(f ) = [1 e i2π f e i4π f e i6π f e i2π f (124) ] T ∈ C 5 .
(5.12)
The modified frequency separation As we can see, for this non-uniform sampling example, our new lower bound is tighter than the lower bounds (both unmodified or modified for non-uniform sampling ) from [40] for N = 256.
We further consider the same set-up for N = 256 + 4 × k, where 1 ≤ k ≤ 50 is an integer, and the corresponding observation index is j = 0, 1, 2, 3 and 124 + 4 × k. In this case, the atom at a certain frequency f ∈ [0, 1) is given by a(f ) = [1 e i2π f e i4π f e i6π f e i2π f (124+4×k) ] T ∈ C 5 .
(5.13) Then Fig.
shows the lower bound on the frequency separation obtained from [40] and Theorem 5.1, and also shows the lower bound on the atom separation obtained from [40] and Theorem 5.1.
For general sampling at non-integer indices, the lower bound from [40] does not apply, since the proof of the lower bound therein depends on the Markov-Bernstein-type inequalities for (integer-exponent) polynomials. In that case, a trivial lower bound 0 exists for frequency separation.
On the other hand, our lower bound from Theorem 5.1 can still be easily calculated for sampling at non-integer indices. Let us consider N = 4, and we observe 2N -1 = 7 samples. The corresponding time index j for these samples are 0, 1, 2 1.5 ,..., and 6 1.5 . To use Theorem 5.1, instead of optimizing over seven frequencies, we simply randomly selected seven frequencies between 0 and 1 to be frequency atoms for matrix A, and obtained a lower bound 0.1947 for atom separation (notice any matrix A consisting of atoms can give a lower bound in Theorem 5.1). This corresponds to a non-trivial frequency separation of 0.0015. To the best of knowledge, this is the first approach in the literature being able to give a lower bound on frequency separation for such sampling patterns.
Our results in the previous sections naturally raise the following question: why can Hankel matrix recovery work without requiring separations between the underlying atoms while it is necessary for the atomic norm minimization to require separations between the underlying atoms? In this section, we introduce the concept of orthonormal atomic norm and its minimization, which explains the success of Hankel matrix recovery in recovering the superposition of complex exponentials regardless of the separations between frequency atoms.
Let us consider a vector w ∈ C N , where N is a positive integer. We denote the set of atoms by AMSET, and assume that each atom a(τ ) (parameterized by τ ) belongs to C N . Then the atomic norm is given by [5,10,41]
We say the atomic norm w ATOMIC is an orthonormal atomic norm if, for every w,
where w = s k=1 c k a(τ k ), a(τ k ) 2 = 1 for every k, and a(τ k )'s are mutually orthogonal to each other. In the Hankel matrix recovery, the atom set AMSET is composed of all the rank-1 matrices in the form uv * , where u and v are unit-norm vectors in C N [10,12]. Let us assume x ∈ C 2N-1 . We can see that the nuclear norm of a Hankel matrix is an orthonormal atomic norm of H(x):
where andv k is the k-th column of V. This is because the matrices u k v * k 's are orthogonal to each other and each of these rank-1 matrices has unit energy.
Let us now further assume that x ∈ C 2N-1 is the superposition of R complex exponentials with R ≤ N, as defined in (1.1). Then H(x) is a rank-R matrix, and can be written as
Even though the original R frequency atoms a(τ k )'s for x can be arbitrarily close, we can always write H(x) as a superposition of R orthonormal atoms u k v * k 's from the SVD of H(x). Because the original R frequency atoms a(τ k )'s for x can be arbitrarily close, they can violate the necessary separation condition set forth in (5.2). However, for H(x), its composing atoms can be R orthonormal atoms u k v * k 's from the SVD of H(x). These atoms u k v * k 's are of unit energy, and are orthogonal to each other. Thus, these atoms are well separated and have the opportunity of not violating the necessary separation condition set forth in (5.2). This explains why the Hankel matrix recovery approach can break free from the separation condition which is required for traditional atomic norm minimizations.
In this section, we present a new matrix-theoretic inequality of nuclear norms, and give a proof of it from the theory of compressed sensing (using nuclear norm minimization). To the best of our knowledge, we have not seen the statement of this inequality of nuclear norms, or its proof elsewhere in the literature.
Theorem 7.1. Let A ∈ C m×n , and t = min(m, n). Let σ 1 , σ 2 ,..., and σ t be the singular values of A arranged in descending order, namely
Let k be any integer such that
Then for any orthogonal projector P onto k-dimensional subspaces in C m , and any orthogonal projector
In particular,
where A 1:k,1:k is the submatrix of A with row indices between 1 and k and column indices between 1 and k, and A (k+1):m,(k+1):n is the submatrix of A with row indices between k + 1 and m and column indices between k + 1 and n.
Before giving the proof of Theorem 7.1, we will state Lemmas 7.1, 7.2 and 7.3 which will be used in the proof of Theorem 7.1. Lemma 7.1. ( [29]) Let G and H be two matrices of the same dimension in
Lemma 7.1 actually provides a lower bound for the nuclear norm of the difference between two matrices, and the lower bound is just the sum of the differences of singular values. Lemma 7.2. ([1; 20]) For arbitrary matrices X, Y and Z = X -Y ∈ C m×n . Let σ 1 , σ 2 ,..., and σ t (t = min{m, n}) be the singular values of a matrix A ∈ C m×n arranged in descending order, namely X,Y) be the distance between the i-th singular value of X and Y, namely,
J. YI ET AL.
Let s
where Z k is defined as k i=1 σ i (Z). Different from Lemma 7.1, Lemma 7.2 considers the partial sum of singular values of the difference between two matrices, and similarly, this partial sum can be lower bounded by the partial sum of the difference of singular values. Lemma 7.3. Suppose a rank-R matrix X ∈ C M×N admits an SVD X = UΣV * , where Σ ∈ R R×R is a diagonal matrix, and
Define S ⊆ C M×N as
and define
Then we have
Remarks: 1. Lemma 7.3 gives the subdifferential of the nuclear norm of a general complex-numbered matrix (including non-square matrix and non-symmetric matrix); 2. To find the subdifferential of X * , we need to derive
Re(X) Im(X) , for which we have Lemma 7.3. Note that in our earlier work [3], we have already shown one direction of (7.6), namely
Now we show the two sets are indeed equal. The proof of Lemma 7.3 is given in Appendix.
Proof of Theorem 7.1: We first consider the case where all the elements of A are real numbers. Without loss of generality, we consider
and
Then
and We will prove this theorem by contradiction. We first show that if σ 1 + • • • + σ k < σ k+1 + • • • + σ t holds for A, then for any matrix X ∈ C m×n with rank no more than k, for any positive number l > 0, we have
The proof of (7.12) follows similar arguments as in [29]. Since X has a rank at most k, then σ i (X) = 0, ∀i > k. Then we have
where, for the first inequality, we used the Lemma 7.2.
We next show that if A 1:k,1:k * > A (k+1):m,(k+1):n * , one can construct a matrix X with rank at most k such that X + lA * ≤ X * for a certain l > 0. We divide this construction into two cases: when A 1:k,1:k has rank equal to k, and when A 1:k,1:k has rank smaller than k.
When A 1:k,1:k has rank k, we denote its SVD as
Then the SVD of A 1:k,1:k 0 0 0 is given by
We now construct
then the subdifferential of • * at X is given by
where
), I 2 = Tr(M * A). (7.16) Let us partition the matrix A into four blocks:
where -k) . Then we have
Since M * X = 0, XM * = 0, M 2 ≤ 1 and X is a rank-k left top corner matrix, we must have
where where the last inequality is from the fact that the nuclear norm is the dual norm of the spectral norm, and M 22 2 ≤ 1.
J. YI ET AL.
Thus, we have
because we assume that A 1:k,1:k * > A (k+1):m,(k+1):n * . From Theorem 23.4 in [37], we have
where f (X; A) is defined as the one-sided directional derivative in [37] f
This means that, when A 11 has rank equal to k, we can always find l > 0 such that X+lA * -X * < 0, or we can reduce the nuclear norm along direction specified by A. Thus, there exists a positive number l > 0, such that
Let us suppose instead that A 11 has rank b < k. We can write the SVD of A 11 as
Then we construct
By going through similar arguments as above (except for taking care of extra terms involving U 3 and V 3 ), one can obtain that Z, A < 0 for every Z ∈ ∂ X * . In summary, no matter whether A 11 has rank equal to k or smaller to k, there always exists a positive number l > 0, such that X + lA * ≤ X * . However, this contradicts (7.12), and we conclude A 1:k,1:k * ≤ A (k+1):m,(k+1):n * , when A has real-numbered elements.
We further consider the case when A is a complex-numbered matrix. We first derive the subdifferential of X * for any complex-numbered matrix m × n X. For any α ∈ R m×n and any β ∈ R m×n , we define F : R 2m×n → R as
For a complex-numbered matrix A, we will similarly show that, if A 1:k,1:k * > A (k+1):m,(k+1):n * , one can construct a matrix X with rank at most k such that X + lA * ≤ X * for a certain l > 0. We divide this construction into two cases: when A 1:k,1:k has rank equal to k, and when A 1:k,1:k has rank smaller than k.
When A 1:k,1:k has rank k, we denote its SVD as
Then the SVD of A 1:k,1:k 0 0 0 is given by
We now construct
We denote
then by Lemma 7.3, the subdifferential of • * at X is given by
Since we define M, N = Re(Tr(M * N)) for complex valued matrices M, N, then for any Z ∈ ∂ X * , we have
where
Similar to the real-numbered case, let us partition the matrix A into four blocks:
where
So
Since M * X = 0, XM * = 0, M 2 ≤ 1 and X is a rank-k left top corner matrix, we must have
where where the last inequality is because the nuclear norm is the dual norm of the spectral norm, and M 22 2 ≤ 1. Thus, we have
because we assume that A 1:k,1:k * > A (k+1):m,(k+1):n * . Following the same arguments in the real numbered case, there exists a positive number l > 0, such that
Let us suppose instead that A 11 has rank b < k. We can write the SVD of A 11 as
Then we construct
By going through similar arguments as above (except for taking care of extra terms involving U 3 and V 3 ), one can obtain that Z, A < 0 for every Z ∈ ∂ X * . In summary, no matter whether complex-numbered A 11 has rank equal to k or smaller to k, there always exists a positive number l > 0, such that X + lA * ≤ X * . However, this contradicts (7.12), and we conclude A 1:k,1:k * ≤ A (k+1):m,(k+1):n * .
In this section, we perform numerical experiments to demonstrate the empirical performance of Hankel matrix recovery, and show its robustness to the separations between atoms. We use superpositions of complex sinusoids as test signals. But we remark that Hankel matrix recovery can also work for superpositions of complex exponentials. We consider the non-uniform sampling of entries studied in [11,41], where we uniformly randomly observe M entries (without replacement) of x from {0, 1, . . . , 2N -2}. By non-uniform sampling, we refer to the following procedure: we randomly select |M|=M distinct entries of x whose locations are uniformly distributed over {0, 1, • • • , 2N -2}. We note that the 'nonuniform sampling' is in contrast to 'uniform sampling' where every entry of x in {0, 1, • • • , 2N -2} is observed. We also consider two signal (frequency) reconstruction algorithms: the Hankel nuclear norm minimization and the atomic norm minimization.
We fix N = 64, i.e. the dimension of the ground truth signal x is 127. We conduct experiments under different M and R for different approaches. For each approach with a fixed M and R, we test 100 trials, where each trial is performed as follows. We first generate the true signal x = [x 0 , x 1 , . . . , x 126 ] T with x j = R k=1 c k e ı2πf k j for j = 0, 1, . . . , 126, where f k are frequencies drawn from the interval [0, 1) uniformly at random, and c k are complex coefficients satisfying the model c k = (1 + 10 0.5m k )e i2πθ k with m k and θ k uniformly randomly drawn from the interval [0, 1]. Let the reconstructed signal be represented by x. If x-x 2
x 2 ≤ 10 -3 , then we regard it as a successful reconstruction. We also provide the simulation results under the Gaussian measurements of x as in [3].
We plot in Fig. 2 the rate of successful reconstruction with respect to different M and R for different approaches. The black and white regions indicate a 0% and 100% of successful reconstruction, respectively, and the grey regions represent successful recovery rate between 0% and 100%. From this figure, we see that the atomic norm minimization still suffers from non-negligible probability of failure even if the number of measurements approach the full 127 samples. The reason is that, since the underlying frequencies are randomly chosen, there is a sizable probability that some frequencies are close to each other. When the frequencies are close to each other violating the atom separation condition, the atomic norm minimization can still fail even if we observe the full 127 samples. By comparison, the Hankel matrix recovery approach experiences a sharper phase transition, and is robust to the frequency separations. We also see that under both the Gaussian projection and the non-uniform sampling models, the atomic norm minimization and the Hankel matrix recovery approach have similar performance.
We further demonstrate the robustness of the Hankel matrix recovery approach to the separations between frequency atoms, as we vary the separations between frequencies. In our first set of experiments, we take N = 64, |M| = M = 65 (≈ 51% sampling rate) and R = 8, and consider noiseless measurements. Again we generate the magnitude of the coefficients as 1 + 10 0.5p , where p is uniformly randomly generated from [0, 1), and the realized magnitudes are 3.1800, 2.5894, 2.1941, 2.9080, 3.9831, 4.0175, 4.1259, 3.6182 in this experiment. The corresponding phases of the coefficients are randomly generated as 2π s, where s is uniformly randomly generated from [0, 1). In this experiment, the realized phases are 4.1097, 5.4612, 5.4272, 4.7873, 1.0384, 0.4994, 3.1975 and 0.5846. The first R -1 = 7 frequencies of exponentials are generated uniformly randomly over [0, 1), and then the last frequency is added in the proximity of the third frequency. In our six experiments, the eighth frequency is chosen such that the frequency separation between the eighth frequency and the third frequency is, respectively, 0.03, 0.01, 0.003, 0.001, 0.0003 and 0.0001. Specifically, in our six experiments, the locations of the eight frequencies are, respectively, {0.3923, 0.9988, 0.3437, 0.9086, 0.6977, 0.0298, 0.4813, 0.3737}, {0.3923, 0.9988, 0.3437, 0.9086, 0.6977, 0.0298, 0.4813, 0.3537}, {0.3923, 0.9988, 0.3437, 0.9086, 0.6977, 0.0298, 0.4813, 0.3467}, {0.3923, 0.9988, 0.3437, 0.9086, 0.6977, 0.0298, 0.4813, 0.3447}, {0.3923, 0.9988, 0.3437, 0.9086, 0.6977, 0.0298, 0.4813, 0.3440} and {0.3923, 0.9988, 0.3437, 0.9086, 0.6977, 0.0298, 0.4813, 0.3438}. Hankel matrix recovery approach gives relative error x-x 2
x 2 = 8.2556 × 10 -9 , 1.6709 × 10 -9 , 2.9204 × 10 -9 , 8.9444 × 10 -9 , 8.6000 × 10 -9 and 2.5026 × 10 -9 , respectively. With the recovered data x, we use the MUSIC algorithm to identify the frequencies. The recovered frequencies for these six cases are, respectively, {0.0298, 0.3437, 0.3737, 0.3923, 0.4813, 0.6977, 0.9086, 0.9988}, {0.0298, 0.3437, 0.3537, 0.3923, 0.4813, 0.6977, 0.9086, 0.9988}, {0.0298, 0.3437, 0.3467, 0.3923, 0.4813, 0.6977, 0.9086, 0.9988}, {0.0298, 0.3437, 0.3447, 0.3923, 0.4813, 0.6977, 0.9086, 0.9988}, {0.0298, 0.3437, 0.3440, 0.3923, 0.4813, 0.6977, 0.9086, 0.9988} and {0.0298, 0.3437, 0.3438, 0.3923, 0.4813, 0.6977, 0.9086, 0.9988}. We illustrate these six cases in Fig. 3, where the peaks of the imaging function J(f ) are the locations of the recovered frequencies and the red cross marks indicate where the true frequencies are. We can see the Hankel matrix recovery successfully recovers the missing data and correctly locates the frequencies.
Below, we demonstrate the performance of these two different approaches, i.e. Hankel matrix recovery approach and atomic norm minimization, alongside each other. In both the noiseless-observation case and noisy-observation case (other formulations for noisy cases yield similar results), the atomic norm minimization is modelled as
where • Atomic is the atomic norm.
Let us first consider noiseless observations. When the frequency separations are 10 -2 , 3 × 10 -2 , 10 -3 , 3 × 10 -3 , 10 -4 and 3 × 10 -4 , respectively, the performance of frequency identification is shown in Fig 3 and4. For atomic norm minimization, we have used the commonly used dual polynomial method to identify frequency atoms at the frequency locations where the magnitude of the dual polynomial is equal to 1. This means that we only need to identify the locations where the dual polynomial achieves a magnitude of 1 to recover the frequencies.
From the results in Fig. 4, we can see that when the frequency separation is big, e.g. 0.03 or 0.01, the atomic norm can identify the frequencies accurately, and the dual polynomial has exactly the same number of peaks as that of the underlying frequencies. Once the frequency separation is below a certain threshold, e.g. when the frequency separation is 0.003, the atomic norm minimization can fail to identify not only the frequencies that are close to each other, but also the frequencies which are far apart. In addition, there are many spurious peaks of the dual polynomial which do not identify the locations of the true frequencies. However, from the results in Fig. 3, the Hankel matrix recovery approach performs much better than the atomic norm minimization. The Hankel matrix recovery approach can still identify frequency atoms even if they are very close to each other, e.g. when the frequency separation is 0.0001 which is much smaller than 0.003.
We also provide quantitative results of signal recovery errors in Table 2, which clearly shows the advantages of Hankel matrix completion when frequency separation is small. We would also like to mention that when the separation is 3 × 10 -3 , the atomic norm minimization achieves a reconstruction error xx 2 of 0.2763 and a relative reconstruction error x-x 2 x 2 of 2.6958 × 10 -3 , while the Hankel matrix recovery can achieve a reconstruction error of 2.993 × 10 -7 and a relative reconstruction error of 2.9204 × 10 -9 .
We further demonstrate the performance of Hankel matrix recovery under noisy measurements. Again, we consider N = 64, |M| = 65 and R = 8. The magnitudes of the coefficients of the complex sinusoids is obtained by 1 + 10 0.5p , where p is uniformly randomly generated from [0, 1). The realized magnitudes are 3.1800, 2.5894, 2.1941, 2.9080, 3.9831, 4.0175, 4.1259 and 3.6182, respectively, in our experiment. The phase of the coefficients are obtained by 2π s, where s is uniformly randomly generated from [0, 1). In this example, the realized phases are 4.1097, 5.4612, 5.4272, 4.7873, 1.0384, 0.4994, 3.1975 and 0.5846, respectively. The first R-1 = 7 frequencies of exponentials are generated uniformly randomly from [0, 1), and then the last frequency is added with frequency separation 5 × 10 -3 from the third frequency. The locations of the realized frequencies are 0.3923, 0.9988, 0.3437, 0.3487, 0.9086, 0.6977, 0.0298 and 0.4813.
We generate the noise vector v ∈ C 2N-1 as s 1 + ıs 2 where each element of s 1 ∈ R 2N-1 and s 2 ∈ R 2N-1 is independently generated from the zero-mean unit-variance standard Gaussian distribution. And we further normalize v such that v 2 = 0.1. In this noisy case, we solve the problem to get the recovered signal x. For ANM, the relative reconstruction error is 5.7905 × 10 -3 , and the reconstruction error is 6.0706 × 10 -1 . For HMC, the relative reconstruction error is 1.0666 × 10 -3 , and the reconstruction error is 1.1182 × 10 -1 .
We illustrate the locations of the recovered frequencies in Fig. 5. We can see that the Hankel matrix recovery can also provide robust data and frequency recovery under noisy measurements. Similar to the noiseless situation, the separation 5 × 10 -3 can be too small for the atomic norm to identify the frequencies successfully, and there are many spurious peaks with magnitude 1. However, the Hankel matrix recovery approach can give relatively accurate identification of frequencies. For noise level v 2 = 0.5 and 1, the recovery performance of ANM and HMC is demonstrated in Fig. 5.
We now provide empirical evidence for the tightness of the proposed theorems. Although it is challenging to analytically prove the tightness of the bounds in these theorems at this moment, we provide empirical results to show the tightness or the sharpness of our results. We point out that even showing empirically the tightness of Theorem 3.1 (which establish worst-case performance guarantees) is difficult, since there are an infinite number of elements in the null space of the linear mapping for measurements. Hence, in this subsection, we only provide empirical evidences for the tightness of Theorems 4.1 and 4.3.
In this new set of experiments, we consider N = 64, and the number of atoms R is taken to be ranging from 2 to 64. The magnitude of the coefficient is 1 + 10 0.5p , where p is uniformly generated from [0, 1), and its phase is 2πs where s is generated uniformly from [0, 1). For each choice of R, we do 100 trials and calculate the corresponding successful recovery rate. A trial is regarded as a success if the decoding algorithm can achieve x-x 2 x 2 ≤ 10 -6 . We consider the case where only the middle sample of the signal vector of length (2N -1) is missing, and this is the sampling pattern which appears in Theorems 4.1 and 4.3.
In the first set of experiments, all the frequencies are chosen with frequency separations as non-zero integer multiples of 1 127 , and the corresponding frequency atoms are orthogonal. The decoder observes all the samples except for the middle one, i.e. |M| = 126. The results are presented in Part (a) of Fig. 6. As we can see, the HMC approach achieves far better performance than the ANM when the number of atoms increases and the separations between frequencies become smaller. When the number of frequency atoms is 30, the ANM achieves a successful recovery rate of almost 0, while the HMC approach achieves a successful rate of almost 1. Besides, the curve of HMC is quite sharp in the sense that the success rate drops quickly once the number of atoms exceeds a certain threshold. In Theorem 4.1, we predict that the recoverable number R of atoms approaches N, when N is large enough. The empirical results indeed show that the recoverable number R is around 60, which is very close to N = 64. This set of experimental results provide empirical evidence for the sharpness of Theorem 4.1.
In the second set of experiments, all but one frequencies are chosen such that the corresponding atoms are orthogonal, while the last frequency is chosen such that it is away from the first frequency with a distance of 10 -4 . All the other configurations are the same as in the experiments for orthogonal atoms. The numerical results for the second set of experiments are shown in Part (b) of Fig. 6. On the one hand, compared with numerical results for orthogonal atoms, we notice that the frequency separations greatly affect the performance of ANM: for example, when the number of frequency atoms is 10, the success rate of ANM under frequency separation 0.0001 is less than 0.5, while the success rate of ANM under separation 1 127 is about 1. On the other hand, we notice that frequency separations do not affect much the successful recovery rate of the HMC approach, as predicted by our theoretical analysis. In Theorem 3, we analytically showed that the number R of recoverable atoms approaches N asymptotically under the HMC approach, when only the middle sample of the signal vector is missing. The empirical results indeed show that the recoverable number R of frequencies is around 60, which is very close to N = 64, even when the frequency separation is small. This set of experimental results provide empirical evidence for the sharpness of Theorem 4.3.
In this paper, we have shown, theoretically and numerically, that the Hankel matrix recovery can be robust to frequency separations in super-resolving the superposition of complex exponentials. By comparison, the TV minimization and atomic norm minimization require the underlying frequencies to be wellseparated to super-resolve the superposition of complex exponentials even when the measurements are noiseless. In particular, we show that Hankel matrix recovery approach can super-resolve the R frequencies, regardless of how close the frequencies are to each other, from compressed non-uniform measurements. We presented a new concept of OANM, and showed that this concept helps us understand the success of Hankel matrix recovery in separation-free super-resolution. We further show that, in traditional atomic norm minimization, the underlying parameters must be well separated so that the signal can be successfully recovered if the atoms are changing continuously with respect to the continuously valued parameters; however, for OANM, it is possible the atomic norm minimization are successful even when the original atoms are arbitrarily close. As a byproduct of this research, we also provide one matrixtheoretic inequality of nuclear norm, and give its proof from the theory of compressed sensing. In future works, it would be interesting to extend the results in this paper to super-resolving the superposition of complex exponentials with higher dimensional frequency parameters [11,13,45]. From convex analysis and Ω = E Re(X) Im(X) , the subdifferential of F is given by
where E * is the adjoint of the linear operator E.
On the one hand, the adjoint E * is given by, for any Δ = We are now ready to show (7.6).
Firstly, we show that any element in H ≡ α β α + ıβ ∈ S must also be in ∂F Re(X) Im(X) , namely (9.8). In fact, for any W = Δ 1 + ıΔ 2 satisfying U * W = 0, WV = 0 and W 2 ≤ 1, we choose
. This choice of Δ satisfies the constraints on Δ in (9.8). Furthermore, UV * + W = (which comes from using the variational characterization of spectral norm). This concludes the proof of (9.10).
Combining (9.9) and (9.10), we arrive at (7.6).
Downloaded from https://academic.oup.com/imaiai/article/12/3/iaad033/7243208 by University of Iowa user on 29 August 2023
Z * R < Z * , (3.6) Downloaded from https://academic.oup.com/imaiai/article/12/3/iaad033/7243208 by University of Iowa user on 29 August 2023
Downloaded from https://academic.oup.com/imaiai/article/12/3/iaad033/7243208 by University of Iowa user on 29 August 2023 AN OANM APPROACH FOR SEPARATION-FREE SUPER-RESOLUTION
J. YI ET AL.
N max ×2π required by[40] is 2.6 * 10 -3 , where N max = 124. We simply construct
J. YI ET AL. Fig.3. Frequency identification performance with noiseless measurements: Hankel matrix recovery method is used.
Fig.4. Frequency identification performance with noiseless measurements: atomic norm minimization is used. According to[41], the frequencies identified by the atomic norm minimization are the locations where the dual polynomial achieves magnitude 1.
J. YI ET AL. Fig. 5. Frequency identification performance with noisy measurements. Left column: Hankel matrix recovery. Right column: atomic norm minimization.