Stable optimizationless recovery from phaseless linear measurements

Introduction

This recovery problem is difficult because the set of real or complex numbers with a given magnitude is nonconvex. In the real case, there are $2^{m}$ possible assignments of sign to the $m$ phaseless measurements. Hence, exhaustive searching is infeasible. In the complex case, the situation is even worse, as there are a continuum of phase assignments to consider. A method based of alternated projections avoids an exhaustive search but does not always converge toward a solution .

In , the authors convexify the problem by lifting it to the space of $n\times n$ matrices, where $\mathbf{x}\mathbf{x}^{*}$ is a proxy for the vector $\mathbf{x}$. A key motivation for this lifting is that the nonconvex measurements on vectors become linear measurements on matrices . The rank-1 constraint is then relaxed to a trace minimization over the cone of positive semi-definite matrices, as is now standard in matrix completion . This convex program is called PhaseLift in , where it is shown that $\mathbf{x}_{0}$ can be found robustly in the case of random $\mathbf{z}_{i}$, if $m=O(n\log n)$. The matrix minimizer is unique, which in turn determines $\mathbf{x}_{0}$ up to a global phase.

The contribution of the present paper is to show that trace minimization is unnecessary in this lifting framework for the phaseless recovery problem. The vector $\mathbf{x}_{0}$ can be recovered robustly by an optimizationless convex problem: one of finding a positive semi-definite matrix that is consistent with linear measurements. We prove there is only one such matrix, provided that there are $O(n\log n)$ measurements. In other words, the phase recovery problem can be solved by intersecting two convex sets, without minimizing an objective. We show empirically that two algorithms converge linearly (exponentially fast) toward the solution. We remark that these methods are simpler than methods for PhaseLift because they require less or no parameter tuning. A result subsequent to the posting of this paper has improved the number of required measurements to $O(n)$ by considering an alternative construction of the dual certificate that allows tighter probabilistic bounds .

The determinacy of the recovery problem over $n\times n$ matrices may be unexpected because there are $n^{2}$ unknowns and only $O(n\log n)$ measurements. What compensates for the apparent lack of data is the fact that the matrix we seek has rank one and is thus on the edge of the cone of positive semi-definite matrices. Most perturbed matrices that are consistent with the measurements cease to remain positive semi-definite. In other words, the positive semi-definite cone $\mathbf{X}\succeq 0$ is “spiky” around a rank-1 matrix $\mathbf{X}_{0}$. That is, with high probability, particular random hyperplanes that contain $\mathbf{X}_{0}$ and have large enough codimension will have no other intersection with the cone.

The present paper does not advocate for fully abandoning trace minimization in the context of phase retrieval. The structure of the sensing matrices appears to affect the number of measurements required for recovery. Consider measurements of the form $\mathbf{x}_{0}^{*}\Phi\mathbf{x}_{0}$, for some $\Phi$. Numerical simulations (not shown) suggest that $O(n^{2})$ measurements are needed if $\Phi$ is a matrix with Gaussian i.i.d. entries. On the other hand, it was shown in that minimization of the nuclear norm constrained by Tr$(\mathbf{X}\Phi)=\mathbf{x}_{0}^{*}\Phi\mathbf{x}_{0}$ recovers $\mathbf{x}_{0}\mathbf{x}_{0}^{*}$ with high probability as soon as $m=O(n\log n)$. Other numerical observations (not shown) suggest that it is the symmetric, positive semi-definite character of $\Phi$ that allows for optimizationless recovery.

The present paper owes much to , as our analysis is very similar to theirs. We wish to also reference the papers , where phase recovery is cast as synchronization problem and solved via a semi-definite relaxation of max-cut type over the complex torus (i.e., the magnitude information is first factored out.) The idea of lifting and semi-definite relaxation was introduced very successfully for the max-cut problem in . The paper also introduces a fast and efficient method based on eigenvectors of the graph connection Laplacian for solving the angular synchronization problem. The performance of this latter method was further studied in .

Problem (1) can be convexified by lifting it to a matrix recovery problem. Let $\mathcal{A}$ and its adjoint be the linear operators

where $\mathcal{H}^{n\times n}$ is the space of $n\times n$ Hermitian matrices. Observe that $\mathcal{A}(\mathbf{x}\mathbf{x}^{*})=A(\mathbf{x})$ for all vectors $\mathbf{x}$. Letting $\mathbf{X}_{0}=\mathbf{x}_{0}\mathbf{x}_{0}^{*}$, we note that $\mathcal{A}(\mathbf{X}_{0})=b$. We emphasize that $\mathcal{A}$ is linear in $\mathbf{X}$ whereas $A$ is nonlinear in $\mathbf{x}$.

The matrix recovery problem we consider is

Without the positivity constraint, there would be multiple solutions whenever $m<\frac{(n+1)n}{2}$. We include the constraint in order to allow for recovery in this classically underdetermined regime.

Our main result is that the matrix recovery problem (2) has a unique solution when there are $O(n\log n)$ measurements.

As a result, the phaseless recovery problem has a unique solution, up to a global phase, with $O(n\log n)$ measurements. In the real-valued case, the problem is determined up to a minus sign.

Theorem 1 suggests ways of recovering $\mathbf{x}_{0}$. If an $\mathbf{X}\in\{\mathbf{X}\succeq 0\}\cap\{\mathbf{X}\mid\mathcal{A}(\mathbf{X})=\mathbf{b}\}$ can be found, $\mathbf{x}_{0}$ is given by the leading eigenvector of $\mathbf{X}$. See Section 6 for more details on how to find $\mathbf{X}$.

2 Stability result

In practical applications, measurements are contaminated by noise. To show stability of optimizationless recovery, we consider the model

We note that all three terms in (3) scale quadratically in $\mathbf{x}$ or $\mathbf{x}_{0}$.

Problem (3) can be convexified by lifting it to the space of matrices. The noisy matrix recovery problem is

We show that all feasible $\mathbf{X}$ are within an $O(\varepsilon)$ ball of $\mathbf{X}_{0}$ provided there are $O(n\log n)$ measurements.

for some $C>0$. This probability is at least $1-e^{-\gamma\frac{m}{n}}$, for some $\gamma>0$.

As a result, the phaseless recovery problem is stable with $O(n\log n)$ measurements.

for some $\phi\in[0,2\pi)$, and for some $C>0$. This probability is at least $1-e^{-\gamma\frac{m}{n}}$, for some $\gamma>0$.

Theorem 3 ensures that numerical methods can be used to find $\mathbf{X}$. See Section 6 for ways of finding $\mathbf{X}\in\{\mathbf{X}\succeq 0\}\cap\{\mathcal{A}(\mathbf{X})\approx\mathbf{b}\}$. As the recovered matrix may have large rank, we approximate $\mathbf{x}_{0}$ with the leading eigenvector of $\mathbf{X}$.

3 Organization of this paper

4 Notation

We use boldface for variables representing vectors or matrices. We use normal typeface for scalar quantities. Let $z_{i,k}$ denote the $k$th entry of the vector $\mathbf{z}_{i}$. For two matrices, let $\langle\mathbf{X},\mathbf{Y}\rangle=\text{Tr}(\mathbf{Y}^{*}\mathbf{X})$ be the Hilbert-Schmidt inner product. Let $\sigma_{i}$ be the singular values of the matrix $\mathbf{X}$. We define the norms

In particular, we write the Frobenius norm of $\mathbf{X}$ as $\|\mathbf{X}\|_{2}$. We write the spectral norm of $\mathbf{X}$ as $\|\mathbf{X}\|$.

We let $\mathbf{I}$ be the $n\times n$ identity matrix. We denote the range of $\mathcal{A}^{*}$ by $\mathcal{R}(\mathcal{A}^{*})$.

Proof of Main Result

As motivation for the introduction of an inexact dual certificate in the next section, observe that if $\mathcal{A}$ is injective on $T$, and if there exists a (exact) dual certificate $\mathbf{Y}\in\mathcal{R}(\mathcal{A}^{*})$ such that

then $\mathbf{X}_{0}$ is the only solution to $\mathcal{A}(\mathbf{X})=\mathbf{b}$. This is because

where the first equality is because $\mathbf{Y}\in\mathcal{R}(\mathcal{A}^{*})$ and $\mathcal{A}(\mathbf{X})=\mathcal{A}(\mathbf{X}_{0})$. The last implication follows from injectivity on $T$.

Conceptually, $\mathbf{Y}$ arises as a Lagrange multiplier, dual to the constraint $\mathbf{X}\succeq 0$ in the feasibility problem

Dual feasibility requires $\mathbf{Y}\succeq 0$. As visualized in Figure 1a, $\mathbf{Y}$ acts as a vector normal to a codimension-1 hyperplane that separates the lower-dimensional space of solutions $\{\mathcal{A}(\mathbf{X})=b\}$ from the positive matrices not in $T$. The condition $\mathbf{Y}_{T^{\perp}}\succ 0$ is further needed to ensure that this hyperplane only intersects the cone along $T$, ensuring uniqueness of the solution.

The nullspace condition $\mathbf{Y}_{T}=0$ is what makes the certificate exact. As $\mathbf{Y}\in\mathcal{R}(\mathcal{A}^{*})$, $\mathbf{Y}$ must be of the form $\sum_{i}\lambda_{i}\mathbf{z}_{i}\mathbf{z}_{i}^{*}$. The strict requirement that $\mathbf{Y}_{T}=0$ would force the $\lambda_{i}$ to be complicated (at best algebraic) functions of all the $\mathbf{z}_{j}$, $j=1,\ldots,m$. We follow in constructing instead an inexact dual certificate, such that $\mathbf{Y}_{T}$ is close to but not equal to , and for which the $\lambda_{i}$ are more tractable (quadratic) polynomials in the $\mathbf{z}_{i}$. A careful inspection of the injectivity properties of $\mathcal{A}$, in the form of the RIP-like condition in , is what allows the relaxation of the nullspace condition on $\mathbf{Y}$.

2 Central Lemma on Inexact Dual Certificates

for some $\delta\leq 1/9$. Suppose that there exists $\bar{\mathbf{Y}}\in\mathcal{R}(\mathcal{A}^{*})$ satisfying

Then, $\mathbf{X}_{0}$ is the unique solution to (2).

Because $\mathcal{A}(\mathbf{H})=0$ and $\bar{\mathbf{Y}}\in\mathcal{R}(\mathcal{A}^{*}),$

where (9) and (10) follow from (7) and (8), respectively. Because the constant in (10) is positive, we conclude $\mathbf{H}_{T}=0$. Then, (9) establishes $\mathbf{H}_{T^{\perp}}=0$. ∎

3 Proof of Theorem 1 and Corollary 2

We use Lemma 1 to prove Theorem 1 for real-valued signals.

We need to show that (5)–(7) hold with high probability if $m>cn\log n$ for some $c$. Lemmas 3.1 and 3.2 in show that (5) and (6) both hold with probability of at least $1-3e^{-\gamma_{1}m}$ provided $m>c_{1}n$ for some $c_{1}$. In section 3, we construct $\bar{\mathbf{Y}}\in\mathcal{R}({\mathcal{A}^{*}})$. As per Lemma 2, $\|\bar{\mathbf{Y}}_{T}\|_{1}\leq 1/2$ with probability at least $1-e^{-\gamma_{2}m/n}$ if $m>c_{2}n$. As per Lemma 3, $\|\bar{\mathbf{Y}}_{T^{\perp}}-2\mathbf{I}_{T^{\perp}}\|\leq 1$ with probability at least $1-2e^{-\gamma_{2}m/\log n}$ if $m>c_{3}n\log n$. Hence, $\bar{\mathbf{Y}}_{T^{\perp}}\succeq\mathbf{I}_{T^{\perp}}$ with at least the same probability. Hence, all of the conditions of Lemma 1 hold with probability at least $1-e^{-\gamma m/n}$ if $m>cn\log n$ for some $c$ and $\gamma$. ∎

The proof of Corollary 2 is immediate because, with high probability, Theorem 1 implies

Existence of Inexact Dual Certificate

To use Lemma 1 in the proof of Theorem 1, we need to show that there exists an inexact dual certificate satisfying (7) with high probability. Our inexact dual certificate vector is different from that in , but we use identical tools for its construction and analysis. We also adopt similar notation.

We note that $\mathcal{A}^{*}\mathcal{A}(\mathbf{X})=\sum_{i}\langle\mathbf{X},\mathbf{z}_{i}\mathbf{z}_{i}^{*}\rangle\mathbf{z}_{i}\mathbf{z}_{i}^{*}$, which can alternatively be written as

Alternatively, we can write the inexact dual certificate vector as

For ease of understanding, we first consider a candidate dual certificate given by

where $\mathbf{Y}_{i}$ is an independent sample of the random matrix

where $\mathbf{z}\sim\mathcal{N}(0,\mathbf{I})$. Because the vector Bernstein inequality requires bounded vectors, we truncate the dual certificate in the same manner as . That is, we consider $1_{E_{i}}\mathbf{Y}_{i}$, completing the derivation of (12).

2 Bounds on 𝐘¯¯𝐘\bar{\mathbf{Y}}

We now present two lemmas that establish that $\bar{\mathbf{Y}}$ is approximately $2(\mathbf{I}-\mathbf{e}_{1}\mathbf{e}_{1}^{*})$, and is thus an inexact dual certificate satisfying (7).

Let $\bar{\mathbf{Y}}$ be given by (12). There exists positive $\gamma$ and $c$ such that for sufficiently large $n$

Let $\bar{\mathbf{Y}}$ be given by (12). There exists positive $\gamma$ and $c$ such that for sufficiently large $n$

3 Proof of Lemma 2: 𝐘¯¯𝐘\bar{\mathbf{Y}} on T𝑇T

We prove Lemma 2 in a way that parallels the corresponding proof in . Observe that

where the first inequality follows because $\bar{\mathbf{Y}}_{T}$ has rank at most 2, and the second inequality follows because $\bar{\mathbf{Y}}_{T}$ can be nonzero only in its first row and column. We can write

where $\bar{\mathbf{y}}_{i}=\mathbf{y}_{i}1_{E_{i}}$, and $\mathbf{y}_{i}$ are independent samples of

First, we compute $\max\|\bar{\mathbf{y}}\|_{2}$. On the event $E$, $|z_{1}|\leq\sqrt{2\beta\log n}$ and $\|\mathbf{z}\|_{2}\leq\sqrt{3n}$. If $n$ is large enough that $2\beta\log n\geq 9$, then $|\xi|\leq 2\beta\log n$. Thus,

By symmetry, every entry of $\bar{\mathbf{y}}$ has zero mean except the first. Hence,

Applying the vector Bernstein inequality with $V=m(8n+16)$, we have that for all $t\leq(8n+16)/[\sqrt{24n}(\beta\log n)^{3/2}]$,

Using the triangle inequality and (20), we get

Lemma 2 follows by choosing $t,\beta$, and $m\geq cn$ where $n$ and $c$ are large enough that

We prove Lemma 3 in a way that parallels the corresponding proof in . We write

where $\mathbf{W}_{i}$ are independent samples of

We decompose $\mathbf{W}$ into the three terms

Letting $\bar{\mathbf{W}}^{(k)}_{i}=\mathbf{W}^{(k)}_{i}1_{E_{i}}$, it suffices to show that with high probability

We show that $m^{-1}\|\sum_{i}\mathbf{I}_{T^{\perp}}1_{E_{i}^{c}}\|=m^{-1}\sum_{i}1_{E_{i}^{c}}$ is small with probability at least $1-2e^{-\gamma m}$ for some constant $\gamma>0$. To do this, we use the scalar Bernstein inequality.

Let $\{X_{i}\}$ be a finite sequence of independent random variables. Suppose that there exists $V$ and $c$ such that for all $X_{i}$ and all $k\geq 3$,

Using the triangle inequality and taking $t$ and $\beta$ such that $\pi(\beta)+t\leq 1/8$ for sufficiently large $n$, we get

We show $m^{-1}\|\sum_{i}\bar{\mathbf{X}}^{(0)}\|$ is small with probability at least $1-2\exp(-\gamma/\log n)$. We write this norm as a supremum over all unit vector perpendicular to $\mathbf{e}_{1}$:

To control the supremum, we follow the same reasoning as in . We bound $\sum_{i}\langle\mathbf{u},\bar{\mathbf{W}}^{(0)}_{i}\mathbf{u}\rangle$ for fixed $\mathbf{u}$ and apply a covering argument over the sphere of $\mathbf{u}$’s. We write

where $\eta_{i}$ are independent samples of

Observing that $\langle\mathbf{z},\mathbf{u}\rangle$ is a chi-squared variable with one degree of freedom, we have

Applying the scalar Bernstein inequality with $V=16m$ and $c_{0}=4\beta\log n$, we get

Taking $t,\beta,m\geq c_{1}n$ with $n$ large enough so that $t+2\sqrt{\pi(\beta)}\leq 1/8$, we have

for some $\gamma^{\prime}>0$. To complete the bound on (29), we use Lemma 4 in :

where $\mathcal{N}_{1/4}$ is a 1/4-net of the unit sphere of vectors $\mathbf{u}\perp\mathbf{e}_{1}$. As $|\mathcal{N}_{1/4}|\leq 9^{n}$, a union bound gives

The bound for the $\|\sum_{i}\bar{\mathbf{W}}^{(1)}\|$ term is similar. We write

where $\eta_{i}$ are independent samples of

The rest of the bound is similar to that of $\|\sum_{i}\bar{\mathbf{X}}^{(0)}\|$ above.

Finally, we also bound $\|\sum_{i}\bar{\mathbf{W}}^{(2)}\|$ similarly. We write

where $\eta_{i}$ are independent samples of

we apply the scalar Bernstein inequality with $c_{0}=4$ and $V=32m$, giving

Stability

Suppose that $\mathcal{A}$ satisfies (5) – (6) and there exists $\mathbf{Y}=\mathcal{A}^{*}\lambda$ satisfying (7) and $\|\lambda\|_{1}\leq 5$. Then,

As before, we take $\mathbf{x}_{0}=\mathbf{e}_{1}$ and $\mathbf{X}_{0}=\mathbf{e}_{1}\mathbf{e}_{1}^{*}$ without loss of generality. Consider any $\mathbf{X}\succeq 0$ such that $\|\mathcal{A}(\mathbf{X})-\mathbf{b}\|_{2}\leq\varepsilon$, and let $\mathbf{H}=\mathbf{X}-\mathbf{X}_{0}$. Whereas $\mathcal{A}(\mathbf{H})=0$ in the noiseless case, it is now of order $\varepsilon$ because

Similarly, $|\langle\mathbf{H},\mathbf{Y}\rangle|$ is also of order $\varepsilon$ because

Analogous to the proof of Lemma 1, we use (7) to compute that

for some $C_{0},C_{1}>0$. Recalling that $\mathbf{H}_{T}$ has rank at most 2,

It remains to show $\|\lambda\|_{1}\leq 5$ for $\bar{\mathbf{Y}}=\mathcal{A}^{*}\lambda$. From (15), we identify $\lambda=m^{-1}(\mathbf{1}_{E}\circ\mathcal{A}\mathcal{S}^{-1}2(\mathbf{I}-\mathbf{e}_{1}\mathbf{e}_{1}^{*}))$. Computing,

2 Proof of Corollary 4

Now we prove Corollary 4, showing that stability of the lifted problem (4) implies stability of the unlifted problem (3). As before, we take $\mathbf{x}_{0}=\mathbf{e}_{1}$ without loss of generality. Hence $\|\mathbf{X}_{0}\|_{2}=1$. Lemma 4 establishes that $\|\mathbf{X}-\mathbf{X}_{0}\|\leq C_{0}\varepsilon$. Recall that $\mathbf{X}_{0}=\mathbf{x}_{0}\mathbf{x}_{0}^{*}$. Decompose $X=\sum_{j}\lambda_{j}\mathbf{v}_{j}\mathbf{v}_{j}^{t}$ with unit-normalized eigenvectors $\mathbf{v}_{j}$ sorted by decreasing eigenvalue. By Weyl’s perturbation theorem,

we use the triangle inequality to form the spectral bound

Complex Case

for some $\delta\leq 3/13$. Suppose that there exists $\bar{\mathbf{Y}}\in\mathcal{R}(\mathcal{A}^{*})$ satisfying

Then, $\mathbf{X}_{0}$ is the unique solution to (2).

The proof of this lemma is identical to the real-valued case. The conditions of the lemma are satisfied with high probability, as before.

The construction of the inexact dual certificate is slightly different because $\mathcal{S}(\mathbf{X})=\mathbf{X}+\text{Tr}(\mathbf{X})\mathbf{I}$ and $\mathcal{S}^{-1}(\mathbf{X})=\mathbf{X}-\frac{1}{n+1}\text{Tr}(\mathbf{X})\mathbf{I}$. As a result

The remaining modifications are identical to those in , and we refer interested readers there for details.

Numerical Simulations

In this section, we show that the optimizationless perspective allows for additional numerical algorithms that are unavailable for PhaseLift directly. These methods give rise to simpler algorithms with less or no parameter tuning. We demonstrate successful recovery under Douglas-Rachford and Nesterov algorithms, and we empirically show that the convergence of these algorithms is linear.

From the perspective of nonsmooth optimization, PhaseLift and the optimizationless feasibility problem can be viewed as a two-term minimization problem

See, for example, the introduction to . Numerical methods based on this splitting include Forward-Backward, ISTA, FISTA, and Douglas-Rachford . If $F$ is smooth, it enables a forward step based on a gradient descent. Nonsmooth terms admit backward steps involving proximal operators. We recall that the proximal operator for a function $G$ is given by

and we note that the proximal operator for a convex indicator function is the projector onto the indicated set.

PhaseLift can be put in this two-term form by softly enforcing the data fit. That gives the minimization problem

where $\iota_{\mathbf{X}\succeq 0}$ is the indicator function that is zero on the positive semidefinite cone and infinite otherwise, and where $\lambda$ is small and positive. If $\lambda=0$, (41) reduces to the optimizationless feasibility problem. The smoothness of $F$ enables methods that are forward on $F$ and backward on $G$. As a representative of this class of methods, we will consider a Nesterov iteration for our simulations below.

The optimizationless view suggests the splitting

where the data fit term is enforced in a hard manner by the indicator function $\iota_{\mathcal{A}(\mathbf{X})=\mathbf{b}}$. Because of the lack of smoothness, we can only use the proximal operators for $F$ and $G$. These operators are projectors on to the affine space $\mathcal{A}(\mathbf{X})=\mathbf{b}$ and $X\succeq 0$, which we denote by $\mathcal{P}_{\mathcal{A}(\mathbf{X})=\mathbf{b}}$ and $\mathcal{P}_{\text{psd}}$, respectively.

The simplest method for (42) is Projection onto Convex Sets (POCS), which is given by the backward-backward iteration $\mathbf{X}_{n+1}=\mathcal{P}_{\text{psd}}\mathcal{P}_{\mathcal{A}(\mathbf{X})=\mathbf{b}}\mathbf{X}_{n}.$ Douglas-Rachford iteration often gives superior performance than POCS, so we consider it as a representative of this class of backward-backward methods.

A strength of the optimizationless perspective is that it does not require as much parameter tuning as PhaseLift. For example, formulation (41) requires a numerical choice for $\lambda$. Nonzero $\lambda$ will generally change the minimizer. It is possible to consider a sequence of problems with varying $\lambda$, or perhaps to create a schedule of $\lambda$ within a problem, but these considerations are unnecessary because the optimizationless perspective says we can take $\lambda=0$. In particular, formulation (42) has the further strength of requiring no parameters at all.

We note that PhaseLift could alternatively give rise to the two-term splitting

where the data fit term is enforced in a hard manner. An iterative approach with this splitting would have an inner loop which approximates the proximal operator of $G$. This inner iteration is equivalent to solving the optimizationless problem.

2 Numerical Results

First, we present a Douglas-Rachford approach for finding $\mathbf{X}\in\{\mathbf{X}\succeq 0\}\cap\{\mathcal{A}(\mathbf{X})\approx\mathbf{b}\}$ by the splitting (42). It is given by the iteration

where $\mathcal{P}_{\text{psd}}$ is the projector onto the positive semi-definite cone of matrices, and $\mathcal{P}_{\mathcal{A}(\mathbf{X})=\mathbf{b}}$ is the projector onto the affine space of solutions to $\mathcal{A}(\mathbf{X})=\mathbf{b}$. In the classically underdetermined case, $m<\frac{(n+1)n}{2}$, we can write

In the case that $m\geq\frac{(n+1)n}{2}$, we interpret $\mathcal{P}_{\mathcal{A}(\mathbf{X})=\mathbf{b}}$ as the least squares solution to $\mathcal{A}(\mathbf{X})=b$.

Second, we present a Nesterov gradient-based method for solving the problem (41). Letting ${g(\mathbf{X})=\frac{1}{2}\|\mathcal{A}(\mathbf{X})-\mathbf{b}\|^{2}+\lambda\ \text{tr}(\mathbf{X})}$, we consider the following Nesterov iteration with constant step size $\alpha$:

Figure 2 shows the average recovery error for the optimizationless problem under the Douglas-Rachford method and the Nesterov method over a range of values of $n$ and $m$. For the Nesterov method, we consider the optimizationless case of $\lambda=0$, and we let the step size parameter $\alpha=2\cdot 10^{-4}$. Each pair of values was independently sampled 10 times, and both methods were run for 1000 iterations. The plot shows that the number of measurements needed for recovery is approximately linear in $n$, significantly lower than the amount for which there are an equal number of measurements as unknowns. The artifacts around the curve $m=\frac{n(n+1)}{2}$ appear because the problem is critically determined, and the only solution to the noisy $\mathcal{A}(\mathbf{X})=\mathbf{b}$ is not necessarily positive in that case.

Figure 3 shows recovery error versus iteration number under the Douglas-Rachford method, the Nesterov method for $\lambda=0$ and the Nesterov method for $\lambda=10^{-5}$. For the Nesterov methods, we let the step size parameter be $\alpha=10^{-4}$. For noisy data, convergence is initially linear until it tapers off around the noise level. For noiseless data, convergence for feasibility problem is linear under both the Douglas-Rachford and Nesterov methods. The Nesterov implementation of PhaseLift shows initial linear convergence until it tapers off. Because any nonzero $\lambda$ allows for some data misfit in exchange for a smaller trace, the computed minimum is not $\mathbf{X}_{0}$ and the procedure converges to some nearby matrix. The convergence rates of the Nesterov method could probably be improved by tuning the step-sizes in a more complicated way. Nonetheless, we observe that the Douglas-Rachford method exhibits a favorable convergence rate while requiring no parameter tuning.

We would like to remark that work subsequent to this paper shows that the number of measurements needed by the optimizationless feasibility problem is about the same as the number needed by PhaseLift . That is, the phase transition in Figure 2 occurs in about the same place for both problems.

References