Sparse Signal Recovery from Quadratic Measurements via Convex Programming

Introduction

provided $k\leq O(m/\log(n/m))$. Another example is a recently proposed semidefinite programming framework for phase retrieval, called PhaseLift , by which a signal can be exactly recovered-up to a multiplicative constant- from quadratic measurements. The SDP is a combination of trace minimization and Shor’s SDP-relaxation for quadratic constraints. We review the results in below:

PhaseLift

If we assume $m\geq C_{0}n$ for some numerical constant $C_{0}$, then with high probability, $\bm{x}\bm{x}^{T}$ is the unique solution to the following convex optimization problem:

Notice that $\bm{x}\bm{x}^{T}$ is feasible since $\bm{x}\bm{x}^{T}\succeq\bm{0}$ and

There is an inherent ambiguity to the solution of (1.2), since multiplying by a phase factor ($\pm 1$ in the real case) does not change measurements. From now on, we only consider solutions modulo multiplication by phase.

In this paper , we consider model (1.2) in the case that $m<<n$. In this regime, (1.2) does not yield injective measurements. In fact, each equation in (1.2) is the union of two linear equations by assigning different signs, so generally we have $2^{m}$ solutions. However, if we assume that the unknown vector $\bm{x}$ is $k$-sparse, then under some mild conditions on the number of measurements, system (1.2) becomes well-posed:

where $\bm{v}_{T}$ means the restriction of $\bm{v}$ on the support $T$. The genericity of $\bm{b_{i}},~{}i=1,...,m_{2}$ implies the genericity of ${\bm{b_{i}}}_{T},~{}i=1,...,m$. Then since $m_{2}\geq 4(2k)-2=8k-2$ we have $\bm{y}_{T}=e^{i\psi}\bm{y^{\prime}}_{T}$ for some real number $\psi$ by Theorem 3.1 in . Therefore $\bm{y}=e^{i\psi}\bm{y^{\prime}}$.

Injectivity of the measurements of course doesn’t imply that efficient recovery is possible. Yet, inspired by the success of convex relaxations in compressed sensing and phase retrieval, it is natural to leverage the sparsity assumption to try to efficiently recover signals from fewer than $n$ intensity measurements. A convex formulation in this direction, which, to the best of our knowledge, was first proposed in to solve (1.2), is the following program:

The next theorem shows that when $\bm{z_{j}}$ are IID standard normal random vectors, the solution to (1.4) for an appropriate choice of $\lambda$, is exactly $\bm{x}\bm{x}^{T}$, provided that $k\leq O(\sqrt{m\over{\log n}})$.

Remark 1: By choosing $\lambda=\sqrt{m\over{4C_{0}\log n}}$, we have exact recovery with probability at least $1-(2\log n+3)(4e^{-\gamma\frac{m}{2\log(n)+3}}+{1\over{n^{3}}})-(5+2n^{2})e^{-\gamma m}$ if the number of measurements obeys $m\geq O(\|\bm{x}\|_{1}^{2}k\log n)$. Moreover, by choosing $\bm{x}$ to be a k-sparse vector with components $x_{i}=\pm{1\over{\sqrt{k}}}$, this reads $m\geq O(k^{2}\log n)$. Remark 2: In , the authors operate under an assumption that the sampling operator satisfies a generalization of the Restricted Isometric Property and mutual coherence, while in Theorem 1.2 of our paper we assume the $\bm{z_{j}}$’s are IID standard Normal vectors. In our setting the mutual coherence of the sampling operator defined in will be on the order of $O(1)$, since the diagonal entries of $\bm{z_{j}}\bm{z_{j}}^{T}$ are always $\chi^{2}$ random variables. Applying the result in we get $k=O(1)$ in our setting, which is a much smaller range of sparsity than considered in the result of the above theorem. The conclusion of Theorem 1.2 is far more restrictive than that of Theorem 1.1, so one may ponder whether 1.2 is optimal. The following result shows that indeed there is a substantial gap between solving (1.2) and (1.4).

Remark: Taking $\bm{x}$ to be a k-sparse vector with components $x_{i}=\pm{1\over{\sqrt{k}}}$, this reads $m\geq O(k^{2}/\log^{2}n)$. This theorem obtains sharp theoretical results on the performance of (1.4) in the Gaussian quadratic measurement setting, which may be surprising since it implies that there is a substantial gap between the sufficient number of measurements for injectivity and the necessary number of measurements for recovery via a class of natural convex relaxations.

2 Definitions and notations

The proof of Theorem 1.2

In this section we will prove Theorem 1.2. First we will cite and prove some supporting lemmas. Then we prove that it suffices to construct an approximate dual certificate matrix to the primal convex optimization problem. Finally we use a modification of the golfing scheme to construct such an approximate dual certificate with high probability. Both the idea of the approximate dual certificate and the golfing scheme are originally due to David Gross’ work in Matrix completion.

In this section we establish some useful properties of $\mathcal{A}$.

There is an event $E$ of probability at least $1-5e^{-\gamma_{0}m}$ such that on $E$, any positive symmetric matrix obeys

There is an event $E$ of probability at least $1-2n^{2}e^{-\gamma_{0}m}$ such that on $E$, any symmetric matrix obeys

Since $|z_{ja}{z_{jb}}|,~{}j=1...,m$ are IID sub-exponential variables with expectation $1$ or ${2\over\pi}$ and have finite $\psi_{1}$-norm. By Proposition 5.16 of , we have

with probability at least $1-2n^{2}e^{-\gamma_{0}m}$. On this event we have $m^{-1}\|\mathcal{A}(\bm{X})\|_{1}\leq{9\over 8}\|\bm{X}\|_{1}$.

2 Exact recovery by the existence of an approximate dual certificate.

In the classical theory of semidefinite programming, the existence of an exact dual certificate can be used to prove that a specific point is the solution to the primal problem. By using an idea in , in order to prove Theorem 1.2, it suffices to prove the existence of an approximate dual certificate.

Denote $\bm{X_{0}}=\lambda\bm{x}\bm{x}^{T}+\mathcal{P}_{T}(\operatorname{sgn}(\bm{x})\operatorname{sgn}(\bm{x})^{T})$. Suppose there exists $\bm{Y}=v_{1}\bm{z_{1}}\bm{z_{1}}^{T}+...+v_{m}\bm{z_{m}}\bm{z_{m}}^{T}$ for some real numbers $v_{1},...,v_{m}$ satisfying $\|\bm{Y}_{T\cap\Omega}-X_{0}\|_{F}\leq{{\|\bm{X_{0}}\|_{F}}\over{6n^{2}}}$, $\|\bm{Y}_{T^{\perp}\cap\Omega}\|\leq{{\|\bm{X_{0}}\|_{F}}\over 5}$ and $\|\bm{Y}_{\Omega^{\perp}}\|_{\infty}\leq{{C\sqrt{\log n}}\over{\sqrt{m}}}\|\bm{X_{0}}\|_{F}$, with some numerical constant $C$. Then assuming that $\mathcal{A}$ satisfies properties (2.1), (2.2) and (2.3), we have that $\bm{x}\bm{x}^{T}$ is the unique solution to the convex program (1.4), provided that $\lambda>\sqrt{k}\|\bm{x}\|_{1}+1$, $\lambda<{{n^{2}}\over{4}}$ and $m>64C^{2}\lambda^{2}\log n$.

Proof Let $\bm{\hat{X}}$ be the solution to the convex program (1.4) and let $\bm{H}=\bm{\hat{X}}-\bm{x}\bm{x}^{T}$. Then by the feasibility condition of the convex program (1.4) , we have

By equality (2.4), we have $\mathcal{A}(\bm{H}_{T\cap\Omega})=\mathcal{A}(\bm{H}_{T^{\perp}\cup\Omega^{\perp}})$. Then by (2.1), (2.2), (2.3) and (2.6), we have

Since $\operatorname{rank}(\bm{H}_{T\cap\Omega})\leq 2$, we have

Now let’s see what inequalities about $\bm{H}$ we can get from the objective function. Since both $\bm{\hat{X}}$ and $\bm{x}\bm{x}^{T}$ are feasible and $\bm{\hat{X}}$ is the minimizer, we have

It is easy to see that $\mathcal{P}_{T^{\perp}}(\operatorname{sgn}(\bm{x})\operatorname{sgn}(\bm{x})^{T})$ is positive semidefinite and combining with (2.6), we get

Notice that $\operatorname{Tr}(\bm{H}_{T^{\perp}})=\operatorname{Tr}(\bm{H}_{T^{\perp}\cap\Omega})+\operatorname{Tr}(\bm{H}_{B})$. By (2.6) and $\lambda\geq 0$, we have

By the construction of the approximate dual certificate $\bm{Y}$, we know $\bm{Y}=\mathcal{A}^{*}(\bm{v})$, which implies $\langle\bm{H},\bm{Y}\rangle=\langle\mathcal{A}(\bm{H}),\bm{v}\rangle=0$. Then we have

By the assumed properties of $\bm{Y}$, we have

Then together with the assumptions of $\lambda>\sqrt{k}\|\bm{x}\|_{1}+1$, $\lambda<{{n^{2}}\over{4}}$ and $m>64C^{2}\lambda^{2}\log n$, we have

by direct calculation. Therefore, by (2.9)

Equations (2.7) and (2.10) give $\bm{H}_{T\cap\Omega}=0$, and then by (2.10), we have $\bm{H}_{T^{\perp}\cap\Omega}=0$ and $\bm{H}_{\Omega^{\perp}}=0$. Hence $\bm{H}=0$, which implies $\bm{x}\bm{x}^{T}$ is the unique minimizer of the convex program (1.4).

3 Key lemma

The following lemma will be essential for the construction of a desirable dual certificate:

For any fixed $\bm{X}\in T\cap\Omega$, we have $\operatorname{rank}(\bm{X})\leq 2$. Consider an eigenvalue decomposition $\bm{X}=\lambda_{1}\bm{u_{1}}\bm{u_{1}}^{T}+\lambda_{2}\bm{u_{2}}\bm{u_{2}}^{T}$, where $\|\bm{u_{1}}\|=\|\bm{u_{2}}\|=1$, $\bm{u_{1}}^{T}\bm{u_{2}}=0$ and both $\bm{u_{1}}$ and $\bm{u_{2}}$ are supported on $G$. Define

provided $m\geq C_{1}k$. Here $\gamma$, $C_{0}$ and $C_{1}$ are numerical constants.

Before proving Lemma 2.4, we need to prove the following supporting lemma:

with probability at least $1-2e^{-\gamma m}$ provided $m\geq C_{0}n$.

Proof By rotational invariance, we can assume $\bm{u}=\bm{e_{1}}$. Define a matrix $\bm{D}=\operatorname{diag}({1\over{\sqrt{\beta_{4}}}},{1\over{\sqrt{\beta_{2}}}},...,{1\over{\sqrt{\beta_{2}}}})$. Define $\bm{w_{j}}=\bm{D}|z_{j1}1_{\{|z_{j1}|\leq 3\}}|\bm{z_{j}}$. It is immediate to check that the $\bm{w_{j}}$’s are IID copies of a zero-mean, isotropic and sub-Gaussian random vector $\bm{w}$. Standard results about random matrices with sub-gaussian rows—e.g. Theorem 5.39 in —give

with probability at least $1-2e^{-\gamma(\epsilon)m}$ provided that $m\geq C_{0}(\epsilon)n$, where $C_{0}$ is sufficiently large.

with probability at least $1-2e^{-\gamma m}$ provided $m\geq C_{1}n$. Similarly, since ${1\over m}\sum_{j=1}^{m}{\bm{z_{j}}}_{G}{\bm{z_{j}}}_{G}^{T}$ is Wishart when restricted on $\Omega$, standard results in random matrix theory— e.g. Corollary 5.35 in —assert that

with probability at least $1-2e^{-\gamma m}$ provided $m\geq C_{1}n$. Then Denote

We have with probability at least $1-4e^{\gamma m}$, $\|\bm{W_{a}}\|\leq{1\over{20}}$ provided $m\geq C_{1}k$. This actually gives us the conclusion by noticing that

For any fixed $a,b\in[n]$, $a>k$ or $b>k$, we know $Y_{ab}=\bm{e_{a}}^{T}\bm{Y}\bm{e_{b}}$ is the arithmetic mean of m IID centered sub-exponential random variables, whose $\psi_{1}-$ norm is bounded by $K(|\lambda_{1}|+|\lambda_{2}|)$ with a numerical constant $K$. Then by Proposition 5.16 in , we have

with probability at least $1-1/{n^{5}}$, which implies our claim.

4 Adaptation of the golfing scheme

In this section we will construct the dual certificate satisfying all the properties in Lemma 2.3 by using the golfing scheme.

Proof of Theorem 1.2: It suffices to construct $\bm{Y}$ satisfying all the properties in Lemma 2.3 with high probability. We divide the group of IID random vectors $\{\bm{z_{1}},...,\bm{z_{m}}\}$ into $l:=\lfloor 2\log(n)\rfloor+3$ groups

This implies that $m_{1}+...+m_{l}=m$. We use the same definition of $\bm{X_{0}}$ in Lemma (2.3). For i=1,..,l, as in Lemma 2.4, we define the eigenvalue decomposition

Moreover, we define $\bm{X_{i}}=\bm{X_{i-1}}-\mathcal{P}_{T\cap\Omega}({\bm{Y_{i}}})$, and $\bm{Y}=\sum_{i=1}^{l}\bm{Y_{i}}$. By definition we have $\bm{X_{i}}$’s are in $T\cap\Omega$, so $\bm{Y_{i}}$ is well-defined. By Lemma (2.4), with probability at least $1-l(4e^{-\gamma m_{i}}+1/{n^{3}})$, we have for $i=1,...,l$

provided $m_{1}\geq C_{1}k,...,m_{l}\geq C_{1}k$. Therefore, $\bm{Y}=v_{1}\bm{z_{1}}\bm{z_{1}}^{T}+...+v_{m}\bm{z_{m}}\bm{z_{m}}^{T}$ and

When $m\geq(2\log n+3)C_{1}k$, we can always make such a division of $\{\bm{z_{1}},...,\bm{z_{m}}\}$, so the proof is complete.

The proof of Theorem 1.3

Then we can further assume $(\bm{v_{1}},...,\bm{v_{m_{1}}})$ only depend on $(\bm{a_{1}},...,\bm{a_{m_{1}}})$ and are independent of $(\bm{b_{1}},...,\bm{b_{m_{2}}})$. Then we have

Since $\bm{b_{j}}$ are IID $\mathcal{N}(\bm{0},\bm{I})$ random vectors, and are independent from the orthonormal vectors $\bm{v_{1}},...,\bm{v_{m_{1}}}$, we have

By the Chernoff upper bound for the $\chi^{2}$ distribution, we have

with probability $1-m_{2}e^{-0.09(N-m_{1})}$. On the other hand, we have

We start by defining the event $E=E(\bm{z_{1}},...,\bm{z_{m}})$. First, we define an event

By the assumption that $\|\bm{x}\|_{2}=1$ and $\bm{z_{j_{G}}}\sim\mathcal{N}(\bm{0},\bm{I}_{k\times k})$, we have

Hereafter all our discussions will be on the event $E$. We now come back to derive the necessary condition for $\bm{x}\bm{x}^{T}$ to be an optimal point of (1.4). By section 5.9.2 of , the condition is

which, using the definition of the subgradient, is equivalent to

One can verify that $\bm{S}\preceq\bm{0}$ and $\left<\bm{S},\bm{x}\bm{x}^{T}\right>=0$ is equivalent to $\bm{S}\preceq\bm{0}$ and $\mathcal{P}_{T}(\bm{S})=0$. Thus the necessary condition for $\bm{x}\bm{x}^{T}$ to be a minimizer of this program is the existence of a dual certificate $\bm{Y}$ with the following properties:

Project both sides of (3.1) on $\Gamma$, we have

It is also obvious that $\|\bm{L}_{\Gamma}\|_{\infty}\leq\|\bm{L}_{\Omega^{\perp}}\|_{\infty}\leq 1$, which implies

On the other hand, project both sides of (3.1) on $T$, we have

Case 1: λ<−k2𝜆𝑘2\lambda<-{k\over 2}.

By the assumption $k\leq m\leq{n\over{40\log n}}$, we can assume the eigenvalue decomposition

where $\{\bm{u_{1}},...,\bm{u_{n-k}}\}$ is an orthogonal basis of $\operatorname{span}(\bm{e_{k+1}},...,\bm{e_{n}})$. Then by (3.4), we have

Since $\{\bm{u_{1}},...,\bm{u_{n-k}}\}$ is an orthogonal basis of $\operatorname{span}(\bm{e_{k+1}},...,\bm{e_{n}})$, we have

By (3.8) and the assumption $4\leq k\leq m\leq{n\over{40\log n}}$, we have

By (3.9) and (3.10), we have $\leq(n-k)\sqrt{m}\geq(n-k-m){k\over 2}-(n-k)$ which implies

Case 2: λ≥−k2𝜆𝑘2\lambda\geq-{k\over 2}.

Let $I_{+}=\{k\in\{1,2\ldots,m\};c_{k}\geq 0\}$ and $I_{-}=\{k\in\{1,2,\ldots,m\};c_{k}<0\}$. By (3.7) and the definition of $E\subset E_{0}$, we have

By the definition of $E$ and Lemma 3.1, we have

Notice that $\|\bm{L_{\Gamma}}\|_{F}\leq(n-k)\|\bm{L_{\Gamma}}\|_{\infty}\leq n-k$. By (3.11) and (3.12),

By the assumption that $k\leq m\leq{n\over{40\log n}}$ and $\lambda\geq-{k\over 2}$, we have

Therefore, by putting Case 1 and Case 2 together, we have

Discussion

We provide theoretical guarantees on the recovery of a sparse signal from quadratic Gaussian measurements via convex programming and show that our results are sharp for a class of recently proposed convex relaxations. For this model, unlike classical compressed sensing, compressive phase retrieval imposes a stricter limitation on the number of measurements needed for recovery via naive convex relaxation than is needed for well-posedness. This leads to a natural open question: can we narrow the gap by using other convex programs besides (1.4)?

Theorem 1.3 shows the limitations of (1.4) in the sense of exact recovery, since we only need to recover the support of the unknown vector to recover $\bm{x}$ by using the PhaseLift algorithm to solve the resulting overdetermined system of quadratic equations. Mathematically, recovering the support is at least as easy as exact recovery. Can we do better than (1.4) by formulating the right support recovery problem? We leave these considerations for future research.

Acknowledgements

We are thankful for fruitful discussion with Emmanuel Candès and also Mahdi Soltanolkotabi, who generously provided us with the results of his numerical experiments on sparse recovery.

References