1

PhaseMax: Convex Phase Retrieval via Basis Pursuit

arXiv:1610.07531v3 [] 30 Jan 2018

Tom Goldstein and Christoph Studer

Abstract—We consider the recovery of a (real- or complexvalued) signal from magnitude-only measurements, known as
phase retrieval. We formulate phase retrieval as a convex
optimization problem, which we call PhaseMax. Unlike other
convex methods that use semidefinite relaxation and lift the
phase retrieval problem to a higher dimension, PhaseMax is
a “non-lifting” relaxation that operates in the original signal
dimension. We show that the dual problem to PhaseMax is Basis
Pursuit, which implies that phase retrieval can be performed
using algorithms initially designed for sparse signal recovery.
We develop sharp lower bounds on the success probability of
PhaseMax for a broad range of random measurement ensembles,
and we analyze the impact of measurement noise on the solution
accuracy. We use numerical results to demonstrate the accuracy
of our recovery guarantees, and we showcase the efficacy and
limits of PhaseMax in practice.

I. I NTRODUCTION
Phase retrieval is concerned with the recovery of an ndimensional signal x0 ∈ Hn , with H either R or C, from
m ≥ n squared-magnitude, noisy measurements [2]
b2i = |hai , x0 i|2 + ηi ,

i = 1, 2, . . . , m,

(1)

where ai ∈ Hn , i = 1, 2, . . . , m, are the (known) measurement
vectors and ηi ∈ R, i = 1, 2, . . . , m, models measurement
noise. Let x̂ ∈ Hn be an approximation vector1 to the true
signal x0 . We recover the signal x0 by solving the following
convex problem we call PhaseMax:
(
maximize
hx, x̂i<
x∈Hn
(PM)
subject to |hai , xi| ≤ bi , i = 1, 2, . . . , m.
Here, hx, x̂i< denotes the real-part of the inner product between
the vectors x and x̂. The main idea behind PhaseMax is to
find the vector x that is most aligned with the approximation
vector x̂ and satisfies a convex relaxation of the measurement
constraints in (1).
Our main goal is to develop sharp lower bounds on the
probability with which PhaseMax succeeds in recovering the
true signal x0 , up to an arbitrary phase ambiguity that does
not affect the measurement constraints in (1). By assuming
noiseless measurements, one of our main results is as follows.
T. Goldstein is with the Department of Computer Science, University of
Maryland, College Park, MD (e-mail: tomg@cs.umd.edu).
C. Studer is with the School of Electrical and Computer Engineering, Cornell
University, Ithaca, NY (e-mail: studer@cornell.edu).
This paper was presented in part at the 32th International Conference on
Machine Learning (ICML) [1].
The work of T. Goldstein was supported in part by the US National Science
Foundation (NSF) under grant CCF-1535902, the US Office of Naval Research
under grant N00014-17-1-2078, and by the Sloan Foundation. The work of C.
Studer was supported in part by Xilinx, Inc. and by the US NSF under grants
ECCS-1408006, CCF-1535897, CNS-1717559, and CAREER CCF-1652065.
1 Approximation vectors can be obtained via a variety of algorithms, or can
even be chosen at random. See Section VI.

Theorem 1. Consider the case of recovering a complex-valued
signal x ∈ Cn from m noiseless measurements of the form
(1) with measurement vectors ai , i = 1, 2, . . . , m, sampled
independently and uniformly from the unit sphere. Let


hx0 , x̂i<
angle(x0 , x̂) = arccos
kx0 k2 kx̂k2

be the angle between the true vector x0 and the approximation x̂, and define the constant
α=1−

2
π

angle(x0 , x̂)

that measures the approximation accuracy. Then, the probability that PhaseMax recovers the true signal x0 , denoted by
pC (m, n), is bounded from below as follows:
!
2
(αm − 4n)
pC (m, n) ≥ 1 − exp −
(2)
2m
whenever αm > 4n.
In words, if m > 4n/α and α > 0, then PhaseMax will
succeed with non-zero probability. Furthermore, for a fixed
signal dimension n and an arbitrary approximation vector x̂
that satisfies angle(x0 , x̂) < π2 , i.e., one that is not orthogonal
to the vector x0 , we can make the success probability of
PhaseMax arbitrarily close to one by increasing the number
of measurements m. As we shall see, our recovery guarantees
are sharp and accurately predict the performance of PhaseMax
in practice.
We emphasize that the convex formulation (PM) has been
studied by several other authors, including Bahmani and
Romberg [3], whose work appeared shortly before our own.
Related work will be discussed in detail in Section I-C.

A. Convex Phase Retrieval via Basis Pursuit
It is quite intriguing that the following Basis Pursuit problem
[4], [5]
(
minimize
kzk1
z∈Hm
(BP)
subject to x̂ = AB−1 z,
with B = diag(b1 , b2 , . . . , bm ) and A = [a1 , a2 , . . . , am ]
is the dual problem to (PM); see, e.g., [6, Lem. 1]. As a
consequence, if PhaseMax succeeds, then the phases of the
solution vector z ∈ Hm to (BP) are exactly the phases that
were lost in the measurement process in (1), i.e., we have
yi = phase(zi )bi = hai , x0 i,

i = 1, 2, . . . , m,

with phase(z) = z/|z| for z 6= 0 and phase(0) = 1. This
observation not only reveals a fundamental connection between
phase retrieval and sparse signal recovery, but also implies that

2

Basis Pursuit solvers can be used to recover the signal from
the phase-less measurements in (1).
B. A Brief History of Phase Retrieval

where β = m/n is constant and m → ∞. The authors of
[41] derive an exact asymptotic bound of the performance of
PhaseMax, and a related non-rigorous analysis was given in
[42]. It was also shown that better signal recovery guarantees
can be obtained by iteratively applying PhaseMax. The resulting
method, called PhaseLamp, was analyzed in [41].

Phase retrieval is a well-studied problem with a long
history [7], [8] and enjoys widespread use in applications
such as X-ray crystallography [9]–[11], microscopy [12], [13],
imaging [14], and many more [15]–[18]. Early algorithms, D. Contributions and Paper Outline
such as the Gerchberg-Saxton [7] or Fienup [8] algorithms,
In contrast to algorithms relying on semidefinite relaxation
rely on alternating projection to recover complex-valued signals or non-convex problem formulations, we propose PhaseMax, a
from magnitude-only measurements. The papers [2], [19], [20] novel, convex method for phase retrieval that directly operates
sparked new interest in the phase retrieval problem by showing in the original signal dimension. In Section II, we establish
that it can be relaxed to a semidefinite program. Prominent a deterministic condition that guarantees uniqueness of the
convex relaxations include PhaseLift [2] and PhaseCut [21]. solution to the (PM) problem. We borrow methods from
These methods come with recovery guarantees, but require geometric probability to derive sharp lower bounds on the
the problem to be lifted to a higher dimensional space, which success probability for real- and complex-valued systems in
prevents their use for large-scale problems. Several authors Section III. Section V generalizes our results to a broader
recently addressed this problem by presenting efficient methods range of random measurement ensembles and to systems with
that can solve large-scale semidefinite programs without the measurement noise. We show in Section VI that randomly
complexity of representing the full-scale (or lifted) matrix of chosen approximation vectors are sufficient to ensure faithful
unknowns. These approaches include gauge duality methods recovery, given a sufficiently large number of measurements.
[22], [23] and sketching methods [24].
We numerically demonstrate the sharpness of our recovery
More recently, a number of non-convex algorithms have guarantees and showcase the practical limits of PhaseMax in
been proposed (see e.g., [25]–[31]) that directly operate in Section VII. We conclude in Section VIII.
the original signal dimension and exhibit excellent empirical
performance. The algorithms in [25], [27]–[29], [32] come with E. Notation
recovery guarantees that mainly rely on accurate initializers,
Lowercase and uppercase boldface letters stand for column
such as the (truncated) spectral initializer [25], [28], the Null
vectors and matrices, respectively. For a complex-valued
initializer [33], the orthogonality-promoting method [29], or
matrix A, we denote its transpose and Hermitian transpose by
more recent initializers that guarantee optimal performance [34],
AT and A∗ , respectively; the real and imaginary parts are A<
[35] (see Section VI for additional details). These initializers
and A= . The ith column of the matrix A is denoted by ai and
enable non-convex phase retrieval algorithms to succeed, given
the kth entry of the ith vector ai is [ai ]k ; for a vector a without
a sufficiently large number of measurements; see [36] for more
index, we simply denote the kth entry by ak . We define the
details on the geometry of such non-convex problems.
inner product between two complex-valued vectors a and b
as ha, bi = a∗ b. We use j to denote the imaginary unit.
C. Related Work on Non-Lifting Convex Methods
The `2 -norm and `1 -norm of the vector a are kak2 and kak1 ,
Because of the extreme computational burden of traditional respectively.
convex relaxations (which require “lifting” to a higher dimension), a number of authors have recently been interested in
II. U NIQUENESS C ONDITION
non-lifting convex relaxations for phase retrieval. Shortly before
The measurement constraints in (1) do not uniquely define
the appearance of this work, the formulation (PM) was studied a vector. If x is a vector that satisfies (1), then any vector
by Bahmani and Romberg [3]. Using methods from machine x0 = ejφ x for φ ∈ [0, 2π) also satisfies the constraints. In
learning theory, the authors derived bounds on the recovery of contrast, if x is a solution to (PM), then ejφ x with φ 6= 0
signals with and without noise. The results in [3] are stronger will not be another solution. In fact, consider any vector x
than those presented here in that they are uniform with respect in the feasible set of (PM) with hx, x̂i 6= 0. By choosing
=
to the approximation vector x̂, but weaker in the sense that ω = phase(hx, x̂i), we have
they require significantly more measurements.
hωx, x̂i< = |hx, x̂i| > hx, x̂i< ,
Several authors have studied PhaseMax after the initial
appearance of this work. An alternative proof of accurate signal
which implies that given such a vector x, one can always
recovery was derived using measure concentration bounds in
increase the objective function of (PM) simply by aligning
[37]. Compressed sensing methods for sparse signals [38], and
x with the approximation x̂ (i.e., modifying its phase so that
corruption-robust methods for noisy signals [39] have also been
hx, x̂i is real valued). The following definition makes this
presented. The closely related non-lifting relaxation BranchHull
observation rigorous.
has also been proposed for recovering signals from entry-wise
Definition 1. A vector x is said to be aligned with another
product measurements [40].
Finally, we note that recent works have proved tight vector x̂, if the inner product hx, x̂i is real-valued and nonasymptotic bounds for PhaseMax in the asymptotic limit, i.e., negative.

3

From all the vectors that satisfy the measurement constraints
in (1), there is only one that is a candidate solution to the convex
problem (PM), which is also the solution that is aligned with x̂.
For this reason, we adopt the following important convention
throughout the rest of this paper.
The true vector x0 denotes a solution to (1) that is
aligned with the approximation vector x̂.

x0 = x0 .
Theorem 2 has an intuitive geometrical interpretation. If
x0 is an optimal point and δ is an ascent direction, then one
cannot move in the direction of δ starting at x0 without leaving
the feasible set. This condition is met if there is an ai such
that x0 and δ both lie on the same side of the plane through
the origin orthogonal to the measurement vector ai .
III. P RELIMINARIES : CLASSICAL SPHERE COVERING

PROBLEMS AND GEOMETRIC PROBABILITY
Remark 1. There is an interesting relation between the convex
formulation of PhaseMax and the semidefinite relaxation
To derive sharp conditions on the success probability of
method PhaseLift [2], [19], [20]. Recall that the set of solutions PhaseMax, we require a set of tools from geometric probability.
to any convex problem is always convex. However, the solution Many classical problems in geometric probability involve
set of the measurement constraints (1) is invariant under phase calculating the likelihood of a sphere being covered by random
rotations, and thus non-convex (provided it is non-zero). It is “caps,” or semi-spheres, which we define below.
therefore impossible to design a convex problem that yields
n−1
n
this set of solutions. PhaseMax and PhaseLift differ in how Definition 2. Consider the set SH = {x ∈ H | kxk2 = 1},
n
n
they remove the phase ambiguity from the problem to enable the unit sphere embedded in H . Given a vector a ∈ H , the
a convex formulation. Rather than trying to identify the true cap centered at a with central angle θ is defined as
n−1
vector x0 , PhaseLift reformulates the problem in terms of the
CH (a, θ) = {δ ∈ SH
| ha, δi< > cos(θ)}.
(3)
quantity x0 (x0 )H , which is unaffected by phase rotations in
x0 . Hence, PhaseLift removes the rotation symmetry from the This cap contains all vectors that form an angle with a of less
solution set, yielding a problem with a convex set of solutions. than θ radians. When θ = π/2, we have a semisphere centered
PhaseMax does something much simpler: it pins down the at a, which is simply denoted by
phase of the solution to an arbitrary quantity, thus removing
n−1
CH (a) = CH (a, π/2) = {δ ∈ SH
| ha, δi< > 0}. (4)
the phase ambiguity and restoring convexity to the solution set.
This arbitrary phase choice is made when selecting the phase
We say that a collection of caps covers the entire sphere
of the approximation x̂.
if the sphere is contained in the union of the caps. Before
we can say anything useful about when a collection of caps
We are now ready to state a deterministic condition under
covers the sphere, we will need the following classical result,
which PhaseMax succeeds in recovering the true vector x0 . The
which is often attributed to Schläfli [43]. Proofs that use simple
result applies to the noiseless case, i.e., ηi = 0, i = 1, 2, . . . , m.
induction methods can be found in [44]–[46]. For completeness,
In this case, all inequality constraints in (PM) are active at x0 .
we briefly include a short proof in Appendix A.
The noisy case will be discussed in Section V-B.
Lemma 1. Consider a sphere S n−1 ⊂ Rn . Suppose we slice
Theorem 2. The true vector x0 is the unique maximizer of the sphere with k planes through Rthe origin. Suppose also that
(PM) if, for any unit vector δ ∈ Hn that is aligned with the every subset of n planes have linearly independent normal
approximation x̂,
vectors. These planes divide the sphere into
∃i, [hai , x0 i∗ hai , δi]< > 0.
n−1
X k − 1
r(n, k) = 2
(5)
i
Proof. Suppose the conditions of this theorem hold, and
i=0
consider some candidate solution x0 in the feasible set for
regions.
(PM) with hx0 , x̂i ≥ hx0 , x̂i. Without loss of generality, we
0
0
0
assume x to be aligned with x̂. Then, the vector ∆ = x − x
Classical results in geometric probability study the likelihood
is also aligned with x̂, and satisfies
of a sphere being covered by random caps with centers chosen
independently and uniformly from the sphere’s surface. For our
h∆, x̂i = hx0 , x̂i − hx0 , x̂i ≥ 0.
purposes, we need to study the more specific case in which caps
0
are only chosen from a subset of the sphere. While calculating
Since x is a feasible solution for (PM), we have
this probability is hard in general, it is quite simple when the
|hai , x0 + ∆i|2 =
set obeys the following symmetry condition.
|hai , x0 i|2 + 2[hai , x0 i∗ hai , ∆i]< + |hai , ∆i|2 ≤ b2i , ∀i.
Definition 3. We say that the set A is symmetric if, for all
But |aTi x0 |2 = b2i , and so
x ∈ A, we also have −x ∈ A.

[hai , x0 i∗ hai , ∆i]< ≤ − 12 |hai , ∆i|2 ≤ 0, ∀i.

Now, if k∆k2 > 0, then the unit-length vector δ = ∆/k∆k2
satisfies [hai , x0 i∗ hai , δi]< ≤ 0 for all i, which contradicts
the hypothesis of the theorem. It follows that k∆k2 = 0 and

We are now ready to prove a general result that states when
the sphere is covered by random caps.
Lemma 2. Consider some symmetric set A ⊂ SRn−1 of
positive (n − 1 dimensional) measure. Choose some set of

4

A
mA measurements {ai }m
i=1 uniformly from A. Then, the caps
{CR (ai )} cover the sphere SRn−1 with probability
n−1
X mA − 1
1
pcover (mA , n) = 1 − mA −1
.
2
k

k=0

This is the probability of turning up n or more heads when
flipping mA − 1 fair coins.
Proof. Consider the following two-step process for constructing
the set {ai }. First, we sample mA vectors {a0i } independently
and uniformly from A. Note that, because A has positive
measure, it holds with probability 1 that any subset of n
vectors will be linearly independent.
Second, we define ai = ci a0i , where {ci } are i.i.d. Bernoulli
variables that take value +1 or −1 with probability 1/2. We
can think of this second step as randomly “flipping” a subset
of uniform random vectors. Since A is symmetric and {a0i } is
sampled independently and uniformly, the random vectors {ai }
also have an independent and uniform distribution over A. This
construction may seem superfluous since both {ai } and {a0i }
have the same distribution, but we will see below that this
becomes useful.
Given a particular set of coin flips {ci }, we can write the
set of points that are not covered by the caps {CR (ai )} as
\
\
CR (−ai ) =
CR (−ci a0i ).
(6)
i

i

mA

Fig. 1. (left) A diagram showing the construction of x, y, and x̃ in the proof
of Lemma 3. (right) The reflections defined in (7-9) map vectors ai lying in
T
the half-sphere aT
i x < 0 onto the half-sphere a x > 0.

Lemma 3. Consider two vectors x, y ⊂ SRn−1 , and the caps
CR (x) and CR (y). Let α = 1 − π2 angle(x, y) be a measure
of the similarity between the vectors x and y. Draw some
collection {ai ∈ CR (x)}m
i=1 of m vectors uniformly from CR (x)
so that m > 2n/α. Then
[
CR (y) ⊂
CR (ai )
i

with probability at least



Note that there are 2
such intersections that can be formed,
(αm − 2n)2
pcover (m, n; x, y) ≥ 1 − exp −
.
one for each choice of the sequence {ci }. The caps {CR (ai )}
2m
cover the sphere whenever
 the intersection (6) is empty.
Consider the set of planes {x | hai , xi = 0} . From Lemma 1, Proof. Due to rotational symmetry, we assume y =
we know that mA planes with a common intersection point [1, 0, . . . , 0]T without loss of generality. Consider the point
divide the sphere into r(n, mA ) non-empty regions. Each of x̃ = [x1 , −x2 , . . . , −xn ]T . This is the reflection of x over
these regions corresponds to the intersection (6) for one possible y (see Figure 1). Suppose we have some collection {ai }
choice of {ci }. Therefore, of the 2mA possible intersections, independently and uniformly distributed on the entire sphere.
only r(n, mA ) of them are non-empty. Since the sequence {ci } Consider the collection of vectors

is random, each intersection is equally likely to be chosen, and
if hai , xi ≥ 0
(7)

 ai ,
so the probability of covering the sphere is
0
ai =
ai − 2hai , yiy
if hai , xi < 0, hai , x̃i < 0 (8)


r(n, mA )
−a
if hai , xi < 0, hai , x̃i ≥ 0. (9)
i
pcover (mA , n) = 1 −
.
2mA
The mapping ai → a0i maps the lower half sphere {a | ha, xi <
0} onto the upper half sphere {a | ha, xi > 0} using a combination of reflections and translations (depicted in Figure 1).
Remark 2. Several papers have studied the probability of
Indeed, for all i we have ha0i , xi ≥ 0. This is clearly true in
covering the sphere using points independently and uniformly
case (7) and (9). In case (8), we use the definition of y and x̃
chosen over the entire sphere. The only aspect that is unusual
to write
about Lemma 2 is the observation that this probability remains
the same if we restrict our choices to the set A, provided
ha0i , xi = hai , xi − 2hai , yihy, xi
A is symmetric. We note that this result was observed by
= hai , xi − 2[ai ]1 x1 = −hai , x̃i ≥ 0.
Gilbert [45] in the case n = 3, and we generalize it to any
Because the mapping ai → a0i is onto and (piecewise)
n > 1 using a similar argument.
isometric, {a0i } will be uniformly distributed over the half
We now present a somewhat more complicated covering sphere {a | ha, x0 i > 0} whenever {ai } are independently and
theorem. The next result considers the case where the measure- uniformly distributed over the entire sphere.
ment vectors are drawn only from a semisphere. We consider
Consider the “hourglass” shaped, symmetric set
the question of whether these vectors cover enough area to
contain not only their home semisphere, but another nearby A = {a | ha, xi > 0, ha, x̃i > 0}∪{a | ha, xi < 0, ha, x̃i < 0}.
semisphere as well.
S
0
We now make the following claim: CR (y) ⊂
i CR (ai )

5

which is only valid for αm > 2n.

whenever
SRn−1

⊂

[
ai ∈A

CR (ai ).

(10)

In words, if the caps defined by the subset of {ai } in A cover
the entire sphere, then the caps {CR (a0i )} (which have centers
in CR (x)) not only cover CR (x), but also cover its neighbor cap
CR (y). To justify this claim, suppose that (10) holds. Choose
some δ ∈ CR (y). This point is covered by some cap CR (ai )
with ai ∈ A. If hai , xi > 0 and hai , x̃i > 0, then ai = a0i and
δ is covered by CR (a0i ). Otherwise, we have hai , xi < 0 and
hai , x̃i < 0, then
hδ, a0i i = hδ, ai − 2hai , yiyi

= hδ, ai i − 2hai , yihδ, yi ≥ hδ, ai i ≥ 0.

Note we have used the fact that hδ, yi is real and non-negative
because δ ∈ CR (y). We have also used hai , yi = [ai ]1 =
1
2 (hai , xi + hai , x̃i) < 0, which follows from the definition of
x̃ and the definition of A. In either case, hδ, a0i i > 0, and we
have δ ∈ CR (a0i ), which proves our claim.
S
We can now see that the probability that CR (y) ⊂ i CR (a0i )
is at least as high as the probability that (10) holds. Let
pcover (m, n; x, y | mA ) denote the probability of covering C(y)
conditioned on the number mA of points lying in A. From
Lemma 2, we know that
pcover (m, n; x, y | mA ) ≥ pcover (mA , n) > pcover (mA +1, n+1).
As noted in Lemma 2, the expression on the right is the chance
of turning up more than n heads when flipping mA fair coins.
The probability pcover (m, n; x, y) is then given by
pcover (m, n; x, y) = EmA [pcover (m, n; x, y | mA )]
≥ EmA [pcover (mA , n)]

Remark 3. In the proof of Lemma 3, we obtained a bound on
EmA [pcover (mA + 1, n + 1)] using an intuitive argument about
coin flipping probabilities. This expectation could have been
obtained more rigorously (but with considerably more pain)
using the method of probability generating functions.
Lemma 3 contains most of the machinery needed for the
proofs that follow. In the sequel, we prove a number of exact
reconstruction theorems for (PM). Most of the results rely on
short arguments followed by the invocation of Lemma 3.
We finally state a result that bounds from below the
probability of covering the sphere with caps of small central
angle. The following Lemma is a direct corollary of the results
of Burgisser, Cucker, and Lotz in [47]. A derivation that uses
their results is given in Appendix B.
Lemma 4. Let n ≥ 9, and m > 2n. Then the probability of
covering the sphere SRn−1 with independent uniformly sampled
caps of central angle φ ≤ π/2 is lower bounded by
pcover (m, n, φ) ≥
√


(em)n n − 1
sinn−1 (φ)(m − n)
√
1−
cos(φ)
exp −
(2n)n−1
8n


(m − 2n + 1)2
− exp −
.
2m − 2
IV. R ECOVERY G UARANTEES
Using the uniqueness condition provided by Theorem 2 and
the tools derived in Section III, we now develop sharp lower
bounds on the success probability of PhaseMax for noiseless
real- and complex-valued systems. The noisy case will be
discussed in Section V-B.

> EmA [pcover (mA + 1, n + 1)].
The expression on the right hand side is the probability of
getting more than n heads when one fair coin is flipped for
every measurement ai that lies in A.
Let’s evaluate how often this coin-flipping event occurs. The
region A is defined by two planes that intersect at an angle of
β = angle(x, x̃) = 2 angle(x, y). The probability of a random
point ai lying in A is given by α = 2π−2β
= 1 − πβ , which is
2π
the fraction of the unit sphere that lies either above or below
both planes. The probability of a measurement ai contributing
to the heads count is half the probability of it lying in A,
or 12 α. The probability of turning up more than n heads is
therefore given by
k 
m−k  
n 
X
1
1
m
1−
α
1− α
.
2
2
k
k=0

Using Hoeffding’s inequality, we obtain the following lower
bound
k 
m−k  
n 
X
1
1
m
pcover (m, n; x, y) ≥ 1 −
α
1− α
2
2
k
k=0


2
−(αm − 2n)
≥ 1 − exp
,
2m

A. The Real Case
We now study problem (PM) in the case that the unknown
signal and measurement vectors are real valued. Consider some
collection of measurement vectors {ai } drawn independently
and uniformly from SRn−1 . For simplicity, we also consider the
collection {ãi } = {phase(hai , x0 i)ai } of aligned vectors that
satisfy hãi , x0 i ≥ 0 for all i. Using this notation, Theorem 2
can be rephrased as a simple geometric condition.
Corollary 1. Consider the set {ãi } = {phase(hai , x0 i)ai } of
aligned measurement vectors. Define the half sphere of aligned
ascent directions
DR = CR (x̂) = {δ ∈ SRn−1 | hδ, x̂i ∈ R ≥ 0}.
The true vector x0 will be the unique maximizer of (PM) if
[
DR ⊂
CR (ãi ).
i

Proof. Choose some ascent direction δ ∈ DR . If the assumptions of this Corollary hold, then there is some i with
δ ∈ CR (ãi ), and so hãi , δi ≥ 0. Since this is true for any
δ ∈ DR , the conditions of Theorem 2 are satisfied and exact
reconstruction holds.

6

Using this observation, we can develop the following lower
bound on the success probability of PhaseMax for real-valued
systems.

Form the vector δ 0 = δ + jhδ, x0 i= x0 , which is the projection
of δ onto B. We can verify that δ 0 ∈ B by writing
hδ 0 , x0 i = hδ, x0 i + hjhδ, x0 i= x0 , x0 i

Theorem 3. Consider the case of recovering a real-valued
= hδ, x0 i − jhδ, x0 i= hx0 , x0 i
signal x0 ∈ Rn from m noiseless measurements of the form
= hδ, x0 i − jhδ, x0 i= = hδ, x0 i< ,
(1) with measurement vectors ai , i = 1, 2, . . . , m, sampled
independently and uniformly from the unit sphere SRn−1 . Then, which is real valued. Furthermore, δ 0 ∈ CC (x̂) because
the probability that PhaseMax recovers the true signal x0 ,
hδ 0 , x̂i = hδ, x̂i+hjhδ, x0 i= x0 , x̂i = hδ, x̂i−jhδ, x0 i= hx0 , x̂i.
denoted by pR (m, n), is bounded from below as follows:


The first term on the right is real-valued and non-negative
−(αm − 2n)2
pR (m, n) ≥ 1 − exp
,
(because
δ ∈ DC ), and the second term is complex valued
2m
(because x0 is assumed to be aligned with x̂). It follows that
where α = 1 − π2 angle(x0 , x̂) and m > 2n/α.
hδ 0 , x̂i< ≥ 0 and δ 0 ∈ CC (x̂). Since we already showed that
δ 0 ∈ B, we have δ 0 ∈ CC (x̂) ∩ B. Suppose now that (12) holds.
Proof. Consider the set of m independent and uniformly
The claim will be proved if we can show that δ ∈ D is covered
sampled measurements {ai ∈ SRn−1 }m
i=1 . The aligned vectors
by one of the CC (ãi ). Since δ 0 ∈ CC (x̂) ∩ B, there is some i
0
{ãi = phase(hai , x i)ai } are uniformly distributed over the
with δ 0 ∈ CC (ãi ). But then
half sphere CR (x0 ). Exact reconstruction happens when the
condition in Corollary 1 holds. To obtain a lower bound on the
0 ≤ hδ 0 , ãi i< = hδ, ãi i< + hjhδ, x0 i= x0 , ãi i< = hδ, ãi i< .
probability of this occurrence, we can simply invoke Lemma 3
(13)
with x = x0 and y = x̂.
We see that δ ∈ CC (ãi ), and the claim is proved.
We now know that exact reconstruction happens whenever
condition (12) holds. We can put a bound on the frequency of
B. The Complex Case
this using Lemma 3. Note that the sphere SCn−1 is isomorphic
2n−1
, and the set B is isomorphic to the sphere SR2n−2 .
We now prove Theorem 1 given in Section I, which charac- to SR
{ãi } are uniformly distributed over a half
terizes the success probability of PhaseMax for phase retrieval The aligned vectors
0
sphere
in
C
(x
)
∩
B,
which is isomorphic to the upper half
C
in complex-valued systems. For clarity, we restate our result
2n−2
sphere
in
S
.
The
probability
of these vectors covering the
R
in shorter form.
cap CC (x̂) ∩ B is thus given by pcover (m, 2n − 1; x0 , x̂) from
0
Theorem 1. Consider the case of recovering a complex-valued Lemma 3. We instead use the bound for pcover (m, 2n; x , x̂),
0
n
signal x ∈ C from m noiseless measurements of the form which is slightly weaker.
(1), with {ai }m
i=1 sampled independently and uniformly from
Remark 4. Theorems 1 and 3 guarantee exact recovery for
the unit sphere SCn−1 . Then, the probability that PhaseMax a sufficiently large number of measurements m provided that
recovers the true signal x0 is bounded from below as follows: angle(x0 , x̂) < π . In the case angle(x0 , x̂) > π , our theorems
2
2
!
2
guarantee convergence to −x0 (which is also a valid solution)
(αm − 4n)
,
pC (m, n) ≥ 1 − exp −
for sufficiently large m. Our theorems only fail for large m
2m
if arccos(x0 , x̂) = π/2, which happens with probability zero
when the approximation vector x̂ is generated at random. See
where α = 1 − π2 angle(x0 , x̂) and m > 4n/α.
Section VI for more details.
Proof. Consider the set {ãi } = {phase(hai , x0 i)ai } of aligned
measurement vectors. Define the half sphere of aligned ascent
directions
DC = {δ ∈ SCn−1 | hδ, x̂i< ∈ R+
0 }.
By Theorem 2, the true signal x0 will be the unique maximizer
of (PM) if
[
DC ⊂
CC (ãi ).
(11)
i

Let us bound the probability of this event. Consider the set B =
{δ | hδ, x0 i= = 0}. We now claim that (11) holds whenever
[
CC (x̂) ∩ B ⊂
CC (ãi ).
(12)
i

To prove this claim, consider some δ ∈ DC . To keep notation
light, we will assume without loss of generality that kx0 k2 = 1.

V. G ENERALIZATIONS
Our theory thus far addressed the idealistic case in which
the measurement vectors are independently and uniformly
sampled from a unit sphere and for noiseless measurements. We
now extend our results to more general random measurement
ensembles and to noisy measurements.
A. Generalized Measurement Ensembles
The theorems of Section IV require the measurement vectors
{ai } to be drawn independently and uniformly from the surface
of the unit sphere. This condition can easily be generalized to
other sampling ensembles. In particular, our results still hold
for all rotationally symmetric distributions. A distribution D is
rotationally symmetric if the distribution of a/kak2 is uniform
over the sphere when a ∼ D. For such a distribution, one
can make the change of variables a ← a/kak2 , and then

7

apply Theorems 1 and 3 to the resulting problem. Note that
this change of variables does not change the feasible set for
(PM), and thus does not change the solution. Consequently,
the same recovery guarantees apply to the original problem
without explicitly implementing this change of variables. We
thus have the following simple corollary.
Corollary 2. The results of Theorem 1 and Theorem 3 still hold
if the samples {ai } are drawn from a multivariate Gaussian
distribution with independent and identically distributed (i.i.d.)
entries.
Proof. A multivariate Gaussian distribution with i.i.d. entries
is rotationally symmetric, and thus the change of variables
a ← a/kak2 yields an equivalent problem with measurements
sampled uniformly from the unit sphere.
What happens when the distribution is not spherically
symmetric? In this case, we can still guarantee recovery, but
we require a larger number of measurements. The following
result is, analogous to Theorem 1, for the noiseless complex
case.
D
Theorem 4. Suppose that mD measurement vectors {ai }m
i=1
are drawn from the unit sphere with (possibly non-uniform)
probability density function D : SCn−1 → R. Let `D ≤
inf x∈S n−1 D(x) be a lower bound on D over the unit sphere
C
2π n
and let α = 1 − π2 angle(x0 , x̂) as above. We use sn = Γ(n)
to
n−1
denote the “surface area” of the complex sphere SC , and set
mU = bmD sn `D c. Then, exact reconstruction is guaranteed
with probability at least


(αmU − 4n)2
1 − exp −
2mU

whenever αmU > 4n and `D > 0. In other words, exact
recovery with mD non-uniform measurements happens at least
as often as with mU uniform measurements.
Proof. We compare two measurement models, a uniform
measurement model in which mU measurements are drawn
uniformly from a unit sphere, and a non-uniform measurement
model in which mD measurements are drawn from the
distribution D. Note that the sphere SCn−1 has surface area
2π n
sn = Γ(n)
, and the uniform density function U on this
sphere has constant value s−1
n . Consider some collection of
mU
measurements {aU
i }i=1 drawn from the uniform model. This
ensemble of measurements lies in Rn×mU , and the probability
density of sampling this measurement ensemble (in the space
Rn×mU ) is
U
mU !s−m
.
n

(14)

mD
Now consider a random ensemble {aD
i }i=1 drawn with density
U
U
D. Given {ai }, the event that {ai } ⊂ {aD
i } has density (in
the space Rn×mU )
m
U
Y

mD !
D(ai ).
(mD − mU )! i=1

(15)

The ratio of the non-uniform density (15) to the uniform density

(14) is

 mU



mU
mD Y
mD
mD sn `D
mU
sn D(ai ) ≥
(sn `D )
≥
,
mU
mU i=1
mU
(16)

k
D
where we have used the bound m
mU > (mD /mU ) to obtain
the estimate on the right hand side. The probability of exact
reconstruction using the non-uniform model will always be
at least as large as the probability under the uniform model,
provided the ratio (16) is one or higher. This holds whenever
mU ≤ mD sn `D . It follows that the probability of exact
recovery using the non-uniform measurements is at least the
probability of exact recovery from a uniform model with
mU = bmD sn `D c measurements. This probability is what
is given by Theorem 1.

B. Noisy Measurements
We now analyze the sensitivity of PhaseLift to the measurement noise {ηi }. For brevity, we focus only on the case of
complex-valued signals. To analyze the impact of noise, we
re-write the problem (PM) in the following equivalent form:
(
maximize
hx, x̂i<
n
x∈H

subject to |hai , xi|2 ≤ b̂2i + ηi , i = 1, 2, . . . , m.

(17)

Here, b̂2i = |hai , x0 i|2 is the (unknown) true magnitude
measurement and b2i = b̂2i + ηi . We are interested in bounding
the impact that these measurement errors have on the solution
to (PM). Note that the severity of a noise perturbation of size ηi
depends on the (arbitrary) magnitude of the measurement
vector ai . For this reason, we assume the vectors {ai } have
unit norm throughout this section.
We will begin by proving results only for the case of nonnegative noise. We will then generalize our analysis to the
case of arbitrary bounded noise. The following result gives a
geometric characterization of the reconstruction error.
Theorem 5. Suppose the vectors {ai ∈ Cn } in (17) are
normalized to have unit length, and the noise vector η is
non-negative. Let r be the maximum relative noise, defined by
 
ηi
r = max
,
(18)
i=1,2,...,m
b̂i
and let DC = {δ ∈ SCn−1 | hδ, x̂i< ≥ 0, hδ, x̂i= = 0} be the
set of aligned descent directions. Choose some error bound
ε > r/2, and define the angle θ = arccos(r/2ε). If the caps
{CC (ãi , θ)} cover DC , then the solution x? of (PM), and
equivalently of the problem in (17), satisfies the bound
kx? − x0 k2 ≤ ε.
Proof. We first reformulate the problem (17) as

hx0 + ∆, x̂i<

 maximize
n
∆∈C

subject to |hãi , (x0 + ∆)i|2 ≤ b̂2i + ηi ,


i = 1, 2, . . . , m,

(19)

8

where ∆ = x? − x0 is the recovery error vector and {ãi } =
{phase(hai , x0 i)ai } are aligned measurement vectors. In this
form, the recovery error vector ∆ appears explicitly. Because
we assume the errors {ηi } to be non-negative, the true signal
x0 is feasible for (17). It follows that the optimal objective
of the perturbed problem (17) must be at least as large as the
optimal value achieved by x0 , i.e., h∆, x̂i< ≥ 0. Furthermore,
the solution x0 + ∆ must be aligned with x̂, as is the true
signal x0 , and so h∆, x̂i ∈ R. For the reasons just described,
we know that the unit vector δ = ∆/k∆k2 ∈ DC .
Our goal is to put a bound on the magnitude of the recovery
error ∆. We start by reformulating the constraints in (19) to
get
|hãi , (x0 + ∆)i|2 =

|hãi , x0 i|2 + 2[hãi , x0 i∗ hãi , ∆i]< + |hãi , ∆i|2 ≤ b̂2i + ηi .

Proof. Define the following two sets:
D = {δ ∈ SCn−1 | hδ, x̂i ∈ R+
0}

D0 = {δ ∈ SCn−1 | hδ, x0 i ∈ R+
0 }.
We now claim that the conditions of Theorem 5 hold whenever
[
D0 ⊂
CC (ãi , φ)
(22)
i

0

where {ãi = phase(hai , x i)ai } is the set of aligned measurement vectors. To prove this claim, choose some δ ∈ D
and assume that (22) holds. Since the half-sphere D0 can be
obtained by rotating D by a principal angle of angle(x̂, x0 ),
there is some point δ 0 ∈ D0 with angle(δ, δ 0 ) ≤ angle(x̂, x0 ).
By property (22), there is some cap CC (ãi , φ) that contains δ 0 .
By the triangle inequality for spherical geometry it follows
that:

Subtracting |hãi , x0 i|2 = b̂2i from both sides yields

angle(δ, ãi ) ≤ angle(δ, δ 0 ) + angle(δ 0 , ãi )

2[hãi , x0 i∗ hãi , ∆i]< + |hãi , ∆i|2 ≤ ηi .

≤ angle(x0 , x̂) + φ ≤ θ.

Therefore, δ ∈ CC (ãi , θ), and the claim is proved.
It only remains to put a bound on the probability that (22)
ηi ≥ 2[hãi , x0 i∗ hãi , ∆i]< + |hãi , ∆i|2
occurs. Note that the aligned vectors {ãi } are uniformly
= 2hãi , x0 ihãi , ∆i< + |hãi , ∆i|2
distributed in D0 , which is isomorphic to a half-sphere in
2n−2
. The probability of covering the half sphere SR2n−2
≥ 2hãi , x0 ihãi , ∆i< + |hãi , ∆i< |2 .
(20) SR
with uniformly distributed caps drawn from that half sphere is
This inequality can hold only if hãi , ∆i< is sufficiently small. at least as great as the probability of covering the whole sphere
In particular we know
SR2n−2 with caps drawn uniformly from the entire sphere. This
p
probability is given by Lemma 4 as pcover (m, 2n − 1, φ), and
hãi , ∆i< ≤ −hãi , x0 i + (hãi , x0 i)2 + ηi
is lower bounded by pcover (m, 2n, φ).
ηi
r
ηi
=
≤ .
(21)
≤
0
2hãi , x i
2
2b̂i
We now consider the case of noise that takes on both
positive and negative values. In this case, we bound the error
Now suppose that k∆k2 > ε. From (21) we have
by converting the problem into an equivalent problem with
r
hãi , δi< ≤ hãi , ∆/k∆k2 i< < ,
non-negative noise, and then apply Theorem 6.
2ε
Therefore δ 6∈ CC (ãi , θ) where θ = arccos(r/2ε). This is Theorem 7. Suppose the vectors {ai } in (17) are normalized
2
2
a contradiction because the caps {CC (ãi , θ)} cover DC . It to have unit length, and that ηi > −b̂i for all i. Define the
following measures of the noise
follows that k∆k2 ≤ ε.
)
( )
(
b2i
b̂2i + ηi
Using this result, we can bound the reconstruction error in
2
= min
s = min
the noisy case. For brevity, we present results only for the
i=1,2,...,m
i=1,2,...,m
b̂2i
b̂2i
complex-valued case.
and


Theorem 6. Suppose the vectors {ai } in (17) are independenly
1
ηi
n−1
2
2
and uniformly distributed in SC , and the noise vector η is
r=
max
b̂ − s b̂i +
.
s i=1,2,...,m i
b̂i
non-negative. Let r be the maximum relative error defined in
(18). Choose some error bound ε > r/2, and define the angle Choose some error bound ε > r/2, and define the angle
φ = arccos(r/2ε) − angle(x0 , x̂). Then, the solution x? to φ = arccos(r/2ε) − angle(x0 , x̂). Then, we have the bound
(PM) satisfies
kx? − x0 k2 ≤ ε + (1 − s)kx0 k2 .
kx? − x0 k ≤ ε.
with probability at least
with probability at least
pcover (m, 2n, φ) ≥
pcover (m, 2n, φ) ≥
√


√


(em)2n 2n − 1
sin2n (φ)(m − n)
sin2n−1 (φ)(m − n)
(em)2n 2n − 1
√
1
−
exp
−
cos(φ)
√
1−
exp −
cos(φ)
(4n)2n−1
16n
(4n)2n−1
16n




(m − 4n + 1)2
(m − 4n + 1)2
− exp −
− exp −
2m − 2
2m − 2
Since ηi is non-negative, we have

when n ≥ 5 and m > 4n.

2 This

assumption is required so that b2i = b̂2i + ηi > 0.

9

when n ≥ 5 and m > 4n.

A. Random Initialization

Proof. Consider the “shrunk” version of problem (17)
(
maximize
hx, x̂i<
n

Consider the use of approximation vectors x̂ drawn randomly
from the unit sphere SRn−1 . Do we expect such approximation
vectors to be accurate enough to recover the unknown signal?
To find out, we analyze the inner product between two realvalued random vectors on the unit sphere. Note that we only
care about the magnitude of this inner product. If the inner
product is negative, then PhaseMax simply recovers −x0 rather
than x0 . Our analysis will make use of the following result.

x∈H

subject to |hai , xi|2 ≤ s2 b̂2i + ζi ,

i = 1, 2, . . . , m.
(23)

for some real-valued “shrink factor” s > 0. Clearly, if x0 is
aligned with x̂ and satisfies |hai , x0 i| = bi for all i, then sx0
is aligned with x̂ and satisfies |hai , sx0 i| = sbi . We can now
transform the noisy problem (17) into an equivalent problem
with non-negative noise by choosing
(
)
b̂2i + ηi
2
s = min
and ζi = b̂2i − s2 b̂2i + ηi ≥ 0.
i=1,2,...,m
b̂2i
We then have (sbi )2 + ζi = b2i + ηi2 , and so problem (23) is
equivalent to problem (17). However, the noise ζi in problem
(23) is non-negative, and thus we can apply Theorem 6. This
theorem requires the constant r for the shrunken problem,
which is now
 
ζi
rshrunk = max
i=1,2,...,m sb̂
i
n
o
1
b̂2i − s2 b̂i + ηi /b̂i .
=
max
s i=1,2,...,m
The solution to the shrunk problem (23) satisfies kx? −
sx0 k2 ≤ , with probability pcover (m, 2n, φ), where φ =
arccos(rshrunk /2) − angle(x0 , x̂). If this condition is fulfilled,
then we have
kx? − x0 k2 ≤ kx? − sx0 + sx0 − x0 k2

≤ kx? − sx0 k2 + ksx0 − x0 k2

≤  + (1 − s)kx0 k2 ,
which concludes the proof.

Remark 5. Theorem 7 requires the noise to be sufficiently small
so that the measurements are non-negative, and the PhaseMax
formulation is feasible. A natural extension that avoids this
caveat is to enforce constraints with a hinge penalty rather
than a hard constraint. In addition, it is possible to achieve
better noise robustness in this case. This direction was studied
in [39].

Lemma 5. Consider the angle β = angle(x, y) between
n−1
two random vectors x, y ∈ SH
sampled independently and
uniformly from the unit sphere. Then, the expected magnitude
of the cosine distance between the two random vectors satisfies
s
r
2
2
≤ E[| cos(β)|] ≤
, for H = R
(24)
πn
π(n − 12 )
s
r
1
4
≤ E[| cos(β)|] ≤
, for H = C. (25)
πn
π(4n − 1)
Proof. We first consider the real case. The quantity cos(β) =
hx, yi/(kxk2 kyk2 ) is simply the sample correlation between
two random vectors, whose distribution function is given
by [49]
n−3

f (z) =

(1 − z 2 ) 2

n−1 ,
2n−2 B n−1
2 , 2

(26)

where B is the beta function and z ∈ [−1, +1]. Hence, the
expectation of the magnitude of the inner product is given by
R1
n−3
z(1 − z 2 ) 2 dz
0
E[| cos(β)|] = 2 n−2 n−1 n−1  .
2
B 2 , 2
The integral in the numerator was studied in [50, Eq. 31] and
evaluates to 12 B(1, n−1
2 ). Plugging this expression into
√ (26),
and using the identity Γ(a)/Γ(a/2) = 2a−1 Γ( a+1
)/
π, we
2
get
Γ( n2 )
E[| cos(β)|] = √
.
π Γ( n+1
2 )
Finally, by using bounds on ratios of Gamma functions [51], we
obtain the bounds in (24) for real-valued vectors. The bounds
in (25) for complex-valued vectors are obtained by noting that
SCn−1 is isomorphic to SR2n−1 and by simply replacing n ← 2n
in the bounds for the real-valued case.

VI. H OW TO C OMPUTE A PPROXIMATION V ECTORS ?
There exist a variety of algorithms that compute approximation vectors3 , such as the (truncated) spectral initializer [25],
[28] or corresponding optimized variants [34], [35], the Null
initializer [33], the orthogonality-promoting method [29], or
least-squares methods [48]. We now show that even randomly
generated approximation vectors guarantee the success of
PhaseMax with high probability given a sufficiently large
number of measurements. We then show that more sophisticated
methods guarantee success with high probability if the number
of measurements depends linearly on n.
3 Approximation

vectors are also known as initialization or anchor vectors.

For such randomly-generated approximation vectors, we
now consider the approximation accuracy α that appears in
Theorems 1 and 3. Note that E[|β|] ≤ π2 − E[| cos(β)|], and,
thus
r
2
2
8
E[α] = 1 − E[|β|] ≥ E[| cos(β)|] ≥
π
π
π3 n
for the real case. Plugging this bound on the expected value
for α into Theorem 3, we see that, for an average
prandomlygenerated approximation vector (one with α > 8/(π 3 n)),
the probability of exact reconstruction goes to 1 rapidly as n
goes to infinity, provided that the number of measurements

10

p
satisfies m > cn3/2 for any c > π 3 /2. For complex-valued
√
signals and measurement matrices, this becomes c > 2 π 3 .
Our results indicate that the use of random approximation
vectors requires O(n3/2 ) measurements rather than the O(n)
required by other phase retrieval algorithms with other initialization methods, e.g., the ones proposed in [19], [29] (see also
Section VII). Hence, it may be more practical for PhaseMax
to use approximation vectors obtained from more sophisticated
initialization algorithms.
B. Spectral Initializers

of required measurements) as compared to PhaseLift, TWF,
and TAW, when used together with the truncated spectral
initializer [19]. While the constants c0 , c1 , and c2 in the
recovery guarantees for all of the other methods are generally
very large, our recovery guarantees contain no unspecified
constants, explicitly depend on the approximation factor α,
and are surprisingly sharp. We next demonstrate the accuracy
of our results via numerical simulations.
B. Tightness of Recovery Guarantees

We now investigate the tightness of the recovery guarantees
Spectral initializers were first used in phase retrieval by [25],
of PhaseMax using both experiments and theory. First, we
and enable the computation of an approximation vector x̂
compare the empirical success probability of PhaseMax in
that exhibits strong theoretical properties [28]. This class of
a noiseless and complex-valued scenario with measurement
initializers was analyzed in detail in [34], and [35] developed
vectors taken independently and uniformly from the unit sphere.
an optimal spectral initializer for a class of problems including
All experiments were obtained by using the implementations
phase retrieval.
provided in the software library PhasePack [52]. This software
Fix some δ > 1, and consider a measurement process with
library provides efficient implementations of phase retrieval
m = δn measurement vectors sampled from a Gaussian distrimethods using fast adaptive gradient solvers [53]. We declare
bution. Then, the spectral initializers in [34], [35] are known
signal recovery a success whenever the relative reconstruction
0
to deliver an initializer x̂ with limn→∞ angle(x̂, x ) >  for
error satisfies
some positive . Equivalently, there is some accuracy parameter
2
0
kx0 − xk22
α? = limn→∞ 1 − π angle(x̂, x ) > 0. By combining this
RRE =
< 10−5 .
(27)
kx0 k22
result with Theorem 1, we see that spectral initializers enable
PhaseMax to succeed for large n with high probability provided We compare empirical rates of success to the theoretical lower
that m > n/α∗ measurements are used for signal recovery4 ; a bound in Theorem 1. Figure 2 shows results for n = 100
number of measurements that is linear in n. Note that similar and n = 500 measurements, where we artificially generate
non-asymptotic guarantees can be achieved for finite n using an approximation x̂ for different angles β = angle(x̂, x0 )
other results on spectral initializers (e.g., Prop. 8 in [28]), measured in degrees. Clearly, our theoretical lower bound
although with less sharp bounds.
accurately predicts the real-world performance of PhaseMax.
In Section VII-B, we investigate the numerical value For large n and large β, the gap between theory and practice
of α∗ and make quantitative statements about how many becomes extremely tight. We furthermore observe a sharp
measurements are required when using a spectral initializer. phase transition between failure and success, with the transition
As we will see there, PhaseMax requires roughly m = 5.5n getting progressively sharper for larger dimensions n.
noiseless Gaussian measurements in theory5 , and somewhat
Next, we investigate the tightness of the recovery guarantees
fewer measurements in practice.
when an initial vector is chosen using the spectral initialization
method [34], [35]. Using analytical formulas, the authors of [35]
VII. D ISCUSSION
calculate the accuracy of their proposed spectral initializer for
This section briefly compares our theoretical results to that of real-valued signal estimation from Gaussian measurements with
existing algorithms. We furthermore demonstrate the sharpness large n. They find that, when 2n measurements are used for
of our recovery guarantees and show the practical limits of spectral estimation, the initializer and signal have a squared
cosine similarity of over 0.6. This corresponds to an accuracy
PhaseMax.
parameter of α > 0.55, and PhaseMax can recover the signal
exactly with m = n/α ≈ 3.5n measurements. This is within
A. Comparison with Existing Recovery Guarantees
a factor of 2 of the information-theoretic lower bound for
Table I compares our noiseless recovery guarantees in a real-valued signal recovery, which is 2n − 2 measurements.
complex system to that of PhaseLift [19], truncated Wirtinger
Note that our analysis of PhaseMax assumes that the
flow (TWF) [19], and truncated amplitude flow (TAW) [29]. We measurements used for recovery are statistically independent
also compare to the recovery guarantee provided for PhaseMax from those used by the initializer. For this reason, our theory
using classical machine learning methods in [3]. 6 We see does not allow us to perform signal recovery by “recycling”
that PhaseMax requires the same sample complexity (number the measurements used for initialization. While we do not
empirically observe any change in behavior of the method when
4 Our results assume the initial vector is independent of the measurement
vectors. The total number of measurements needed is thus m > n/α∗ + δn the initialization measurements are used for recovery, the above
for initialization and recovery combined.
reconstruction bounds formally require 5.5n measurements (2n
5 The constants that appear in the bounds [34], [35] were approximated using
for initialization, plus 3.5n for recovery). The bounds in [3]
stochastic numerical methods, and were not obtained exactly using analysis.
are uniform with respect to the initializer, and thus enable
6 Since AltMinPhase [25] requires an online measurement model that differs
significantly from the other algorithms considered here, we omit a comparison. measurement recycling (although this does not result in tighter

11

TABLE I
C OMPARISON OF T HEORETICAL R ECOVERY G UARANTEES FOR N OISELESS P HASE R ETRIEVAL

Algorithm

Sample complexity

Lower bound on pC (m, n)

PhaseMax

m > 4n/α

1 − e−(αm−4n)

PhaseLift [19]

m ≥ c0 n

TWF [19]

m ≥ c0 n

TAF [29]

m ≥ c0 n

Bahmani and Romberg [3]

m>

32
sin4 (α)

2

1 − c1 e−c2 m

1 − c1 e−c2 m

log

8e
sin4 (α)



2
1 − (m + 5)e−n/2 − c1 e−c2 m − 1/n
 

n

1 − 8e

− sin4 (α) M − sin32
4 (α) log

8e
sin4 (α)

N /16

1
Success probability pC (m, n)

Success probability pC (m, n)

1

0.8

0.6

0.4
β = 45◦
β = 36◦
β = 25◦

0.2

0

/(2m)

2

4

6
8
Oversampling ratio m/n

10

12

(a) n = 100

0.8

0.6

0.4
β = 45◦
β = 36◦
β = 25◦

0.2

0

2

4

6
8
Oversampling ratio m/n

10

12

(b) n = 500

Fig. 2. Comparison between the empirical success probability (solid lines) and our theoretical lower bound (dashed lines) for varying angles β between the
true signal and the approximation vector. Our theoretical results accurately characterize the empirical success probability of PhaseMax. Furthermore, PhaseMax
exhibits a sharp phase transition for larger dimensions.

bounds for the overall number of measurements because the
required constants are larger).
Finally, we would like to mention that asymptotically exact
performance bounds for PhaseMax have recently been derived
in [41], and tighter bounds for non-lifting phase retrieval have
been shown for the method PhaseLamp, which repeatedly uses
PhaseMax within an iterative process.
C. Experimental Results
We compare PhaseMax to other algorithms using synthetic
and empirical datasets. The implementations used here are
publicly available as part of the PhasePack software library [52].
This software library contains scripts for comparing different
phase retrieval methods, including scripts for reproducing the
experiments shown here.
1) Comparisons Using Synthetic Data: We briefly compare
PhaseMax to a select set of phase retrieval algorithms in terms
of relative reconstruction error with Gaussian measurements.
We emphasize that this comparison is by no means intended to
be exhaustive and serves the sole purpose of demonstrating the
efficacy and limits of PhaseMax (see, e.g., [17], [21] for more
extensive phase retrieval algorithm comparisons). We compare

the Gerchberg-Saxton algorithm [7], the Fienup algorithm [8],
truncated Wirtinger flow [28], and PhaseMax—all of these
methods use the truncated spectral initializer [28]. We also
run simulations using the semidefinite relaxation (SDR)-based
method PhaseLift [2] implemented via FASTA [53]; this is,
together with PhaseCut [21], the only convex alternative to
PhaseMax, but lifts the problem to a higher dimension.
Figure 3 reveals that PhaseMax requires larger oversampling
ratios m/n to enable faithful signal recovery compared to
non-convex phase-retrieval algorithms that operate in the
original signal dimension. This is because the truncated spectral
initializer requires oversampling ratios of about six or higher to
yield sufficiently accurate approximation vectors x̂ that enable
PhaseMax to succeed.
2) Comparisons Using Empirical Data: The PhasePack
library contains scripts for comparing phase retrieval algorithms
using synthetic measurements as well as publicly available real
datasets. The datasets contain measurements obtained using
the experimental setup described in [54], in which a binary
mask is imaged through a diffusive medium.
Figure 4 shows reconstructions from measurements obtained
using two different test images. Both images were acquired
using phaseless measurements from the same measurement

12

GS
Fienup
TWF
TAF
PhaseLift
PhaseMax

−20

−40

−60

−80

2

4

6
8
Oversampling ratio m/n

10

GS
Fienup
TWF
TAF
PhaseMax

0
Relative reconstruction error [dB]

Relative reconstruction error [dB]

0

12

−20

−40

−60

−80

2

4

(a) n = 100

6
8
Oversampling ratio m/n

10

12

(b) n = 500

Fig. 3. Comparison of the relative reconstruction error. We use the truncated spectral initializer for Gerchberg-Saxton (GS), Fienup, truncated Wirtinger
flow (TWF), truncated amplitude flow (TAF), and PhaseMax. PhaseMax does not achieve exact recovery for the lowest number of measurements among the
considered methods, but is convex, operates in the original dimension, and comes with sharp performance guarantees. PhaseLift only terminates in reasonable
computation time for n = 100.

Original

Gerchberg-Saxton
0.13019

Wirtinger Flow
0.13930

Trunc. Wirt. Flow
0.18017

PhaseMax
0.23459

PhaseLift
0.35452

0.12468

0.13386

0.17836

0.24947

0.35282

Fig. 4. Reconstruction of two 64 × 64 masks from empirical phaseless measurements obtained through a diffusive medium [54]. The numbers on top of each
of the recovered images denote the relative measurement error k|Ax| − bk2 /kbk2 achieved by each method.

operator A. Reconstructions are shown for the GerchbergSaxton [7], Wirtinger flow [27], and truncated Wirtinger
flow [28] methods. We also show results of the two convex
methods PhaseLift (which lifts the problem dimension) [2]
and PhaseMax (the proposed non-lifting relaxation). For each
algorithm, Figure 4 shows the recovered images and reports
the relative measurement error k|Ax| − bk2 /kbk2 . All phase
retrieval algorithms were initialized using a variant of the
spectral method proposed in [35].
We find that the Gerchberg-Saxton algorithm outperforms
all other considered methods by a small margin (in terms
of measurement error) followed by the Wirtinger flow. On
this dataset, measurements seem to be quite valuable, and
the truncated Wirtinger flow method (with default truncation
parameters as described in [28]) performs slightly worse

than the original Wirtinger flow, which uses the full dataset.
PhaseMax produces results that are visually comparable to that
of Wirtinger flow, but with slightly higher relative measurement
error.
D. Advantages of PhaseMax
While PhaseMax does not achieve exact reconstruction with
the lowest number of measurements, it is convex, operates
in the original signal dimension, can be implemented via
efficient solvers for Basis Pursuit, and comes with sharp
performance guarantees that do not sweep constants under
the rug (cf. Figure 2). The convexity of PhaseMax enables
a natural extension to sparse phase retrieval [55], [56] or
other signal priors (e.g., total variation or bounded infinity
norm) that can be formulated with convex functions. Such

13

non-differentiable priors cannot be efficiently minimized using
simple gradient descent methods (which form the basis of
Wirtinger or amplitude flow, and many other methods), but
can potentially be solved using standard convex solvers when
combined with the PhaseMax formulation.
VIII. C ONCLUSIONS
We have proposed a novel, convex phase retrieval algorithm,
which we call PhaseMax. We have provided accurate bounds
on the success probability that depend on the signal dimension,
the number of measurements, and the angle between the approximation vector and the true vector. Our analysis covers a broad
range of random measurement ensembles and characterizes the
impact of general measurement noise on the solution accuracy.
We have demonstrated the sharpness of our recovery guarantees
and studied the practical limits of PhaseMax via simulations.
There are many avenues for future work. Developing a
computationally efficient algorithm for solving PhaseMax (or
the related PhaseLamp procedure) accurately and for large
problems is a challenging problem. An accurate analysis of
more general noise models is left for future work.
A PPENDIX A
P ROOF OF L EMMA 1
The proof is by induction. As a base case, we note that
r(n, 1) = 2 for n ≥ 1 and r(2, k) = 2k for k ≥ 1. Now
suppose we have a sphere SRn−1 in n dimensions sliced by
k − 1 planes into r(n, k − 1) “original” regions. Consider the
effect of adding a kth plane, Pk . Every original region that is
intersected by Pk will split in two, and increase the number
of total regions by 1. The increase in the number of regions is
then the number of original regions that are intersected by Pk .
Equivalently, this is the number of regions formed inside Pk
by the original k − 1 planes. Any subset of k − 1 planes will
have normal vectors that remain linearly independent when
projected into Pk . By the induction hypothesis, the number
of regions formed inside Pk is then given by r(n − 1, k − 1).
Adding this to the number of original regions yields
r(n, k) = r(n, k − 1) + r(n − 1, k − 1).
We leave it to the reader to verify that (5) satisfies this
recurrence relation and base case.
A PPENDIX B
P ROOF OF L EMMA 4
In this section, we prove Lemma 4. This Lemma is a direct
corollary of the following result of Burgisser, Cucker, and Lotz
[47]. For a complete proof of this result, see Theorem 1.1
of [47], and the upper bound on the constant “C” given in
Proposition 5.5.
Theorem 8. Let m > n ≥ 2. Then the probability of covering
the sphere SRn−1 with independent and uniform random caps
of central angle φ ≤ π/2 is bounded by
pcover (m, n, φ) ≥
  Z 
2
m
1−
C
(1 − t2 )(n −2n−1)/2 (1 − λ(t))m−n dt
n
0

−

n−1
X

1
2m−1

k=0

where λ(t)

=

Vol(SRn−1 )

n/2

=

2π
Γ(n/2) ,


m−1
k

Vn−1
Vn

C=

R arccos(t)
0
√
n n−1
2n−1 ,

sinn−2 (φ) dφ, Vn

=

and  = cos(φ).

While Theorem 8 provides a bound on pcover (m, n, φ), the
formulation of this bound does not provide any intuition of
the scaling of pcover (m, n, φ) or its dependence on m and n.
For this reason, we derive Lemma 4, which is a weaker but
more intuitive result. We restate Lemma 4 here for clarity.
Lemma 4. Let n ≥ 9, and m > 2n. Then, the probability of
covering the sphere SRn−1 with caps of central angle φ ≤ π/2
is lower bounded by
pcover (m, n, φ) ≥
√


sin(n−1) (φ)(m − n)
(em)n n − 1
√
exp
−
cos(φ)
1−
(2n)n−1
8n


(m − 2n + 1)2
− exp −
.
2m − 2
Proof. Let us simplify the result of Theorem 8. If we assume
m > 2n, then Hoeffding’s inequality yields



n−1 
1 X m−1
(m − 2n + 1)2
≤ exp −
.
2m−1
k
2m − 2
k=0

Next, we derive a lower bound as follows:
Z arccos(t)
Γ(n/2)
√
sinn−2 (φ) dφ
λ(t) =
Γ((n − 1)/2) π 0
Z arccos(t)
p
≥ (n/2 − 1)/π
sinn−2 (φ) cos(φ) dφ
0

p
= (n/2 − 1)/π

1
sinn−1 arccos(t)
n−1

1
≥ √ (1 − t2 )(n−1)/2 .
8n
p
1
We have used the fact that (n/2 − 1)/π n−1
>
Γ(n/2)
4, and also the “Wallis ratio” bound Γ((n−1)/2)

[57], [58]. Finally, we plug in the inequality m
n

√1
8n
p

≥

for n ≥

n/2 − 1

(em)n
nn .

≤
now have
  Z 
2
m
C
(1 − t2 )(n −2n−1)/2 (1 − λ(t))m−n dt
n
0
√
Z
2
(em)n n − 1 
≤
(1 − t2 )(n −2n−1)/2
n−1
(2n)
0
m−n

1
dt.
× 1 − √ (1 − t2 )(n−1)/2
8n

We

Now we simplify the integral. Using the identity (1 − x)a <
e
, which holds for x ≤ 1, we can convert each term in the
integrand into an exponential. We do this first with x = t2 and
then with x = √18n (1 − t2 )(n−1)/2 to obtain
−ax

(1 − t2 )(n

2

−2n−1)/2



1
1 − √ (1 − t2 )(n−1)/2
8n

m−n

14

 2 2

t (n − 2n − 1) (1 − t2 )(n−1)/2 (m − n)
√
≤ exp −
.
−
2
8n
(28)
We then apply the Cauchy-Schwarz inequality to get
 2 2

Z 
t (n − 2n − 1) (1 − t2 )(n−1)/2 (m − n)
√
exp −
dt
−
2
8n
0
1/2
Z 

2 2
exp −t (n − 2n − 1) dt
≤
0


 1/2
(1 − t2 )(n−1)/2 (m − n)
√
exp −
dt
2n
0


1/2
(1 − 2 )(n−1)/2 (m − n)
1/2
√
≤ []
 exp −
2n


2 (n−1)/2
(1 −  )
(m − n)
√
=  exp −
.
8n
Z



×

Replacing the integral with this bound and using the definition
 = cos(φ) yields the result.

R EFERENCES
[1] T. Goldstein and C. Studer, “Convex phase retrieval without lifting via
PhaseMax,” in Proceedings of the 34th International Conference on
Machine Learning, ser. Proceedings of Machine Learning Research,
D. Precup and Y. W. Teh, Eds., vol. 70. International Convention
Centre, Sydney, Australia: PMLR, 06–11 Aug 2017, pp. 1273–1281.
[Online]. Available: http://proceedings.mlr.press/v70/goldstein17a.html
[2] E. J. Candès, T. Strohmer, and V. Voroninski, “PhaseLift: Exact and stable
signal recovery from magnitude measurements via convex programming,”
Commun. Pure Appl. Math., vol. 66, no. 8, pp. 1241–1274, 2013.
[3] S. Bahmani and J. Romberg, “Phase retrieval meets statistical learning
theory: A flexible convex relaxation,” in Intl. Conf. on Artificial
Intelligence and Statistics (AISTATS), May 2017, pp. 252–260.
[4] S. S. Chen, D. L. Donoho, and M. A. Saunders, “Atomic decomposition
by basis pursuit,” SIAM Rev., vol. 43, no. 1, pp. 129–159, Jul. 2001.
[5] S. Chen and D. Donoho, “Basis pursuit,” in Proc. Asilomar Conf. Signals,
Syst., Comput., vol. 1, Oct. 1994, pp. 41–44.
[6] C. Studer, W. Yin, and R. G. Baraniuk, “Signal representations with
minimum `∞ -norm,” in Proc. Allerton Conf. Commun., Contr., Comput.,
Oct. 2012, pp. 1270–1277.
[7] R. W. Gerchberg and W. O. Saxton, “A practical algorithm for the
determination of phase from image and diffraction plane pictures,” Optik,
vol. 35, pp. 237–246, Aug. 1972.
[8] J. R. Fienup, “Phase retrieval algorithms: a comparison,” Appl. Opt.,
vol. 21, no. 15, pp. 2758–2769, Aug. 1982.
[9] R. W. Harrison, “Phase problem in crystallography,” J. Opt. Soc. Am. A,
vol. 10, no. 5, pp. 1046–1055, May 1993.
[10] J. Miao, T. Ishikawa, Q. Shen, and T. Earnest, “Extending X-ray
crystallography to allow the imaging of noncrystalline materials, cells,
and single protein complexes,” Ann. Rev. Phys. Chem., vol. 59, pp.
387–410, Nov. 2008.
[11] F. Pfeiffer, T. Weitkamp, O. Bunk, and C. David, “Phase retrieval and
differential phase-contrast imaging with low-brilliance X-ray sources,”
Nat. Phys., vol. 2, no. 4, pp. 258–261, Apr. 2006.
[12] S. S. Kou, L. Waller, G. Barbastathis, and C. J. Sheppard, “Transportof-intensity approach to differential interference contrast (TI-DIC)
microscopy for quantitative phase imaging,” Opt. Lett., vol. 35, no. 3,
pp. 447–449, Feb. 2010.
[13] H. Faulkner and J. Rodenburg, “Movable aperture lensless transmission
microscopy: a novel phase retrieval algorithm,” Phys. Rev. Lett., vol. 93,
no. 2, Jul. 2004.
[14] J. Holloway, M. S. Asif, M. K. Sharma, N. Matsuda, R. Horstmeyer,
O. Cossairt, and A. Veeraraghavan, “Toward long-distance subdiffraction
imaging using coherent camera arrays,” IEEE Trans. Comput. Imag.,
vol. 2, no. 3, pp. 251–265, Sept. 2016.
[15] F. Fogel, I. Waldspurger, and A. d’Aspremont, “Phase retrieval for
imaging problems,” Math. Prog. Comp., vol. 8, no. 3, pp. 311–335, Sept.
2016.

[16] E. J. Candès, E.es, X. Li, and M. Soltanolkotabi, “Phase retrieval from
coded diffraction patterns,” Appl. Comput. Harm.. Anal., vol. 39, no. 2,
pp. 277–299, Sept. 2015.
[17] K. Jaganathan, Y. C. Eldar, and B. Hassibi, “Phase retrieval: An overview
of recent developments,” arXiv:1510.07713, Oct. 2015.
[18] L. Tian and L. Waller, “3D intensity and phase imaging from light field
measurements in an LED array microscope,” Optica, vol. 2, no. 2, pp.
104–111, Feb. 2015.
[19] E. J. Candès and X. Li, “Solving quadratic equations via phaselift when
there are about as many equations as unknowns,” Found. Comput. Math.,
vol. 14, no. 5, pp. 1017–1026, Oct. 2014.
[20] E. J. Candès, Y. C. Eldar, T. Strohmer, and V. Voroninski, “Phase retrieval
via matrix completion,” SIAM Rev., vol. 57, no. 2, pp. 225–251, Nov
2015.
[21] I. Waldspurger, A. d’Aspremont, and S. Mallat, “Phase recovery, maxcut
and complex semidefinite programming,” Math. Prog., vol. 149, no. 1-2,
pp. 47–81, Feb. 2015.
[22] M. P. Friedlander, I. Macedo, and T. K. Pong, “Gauge optimization and
duality,” SIAM Journal on Optimization, vol. 24, no. 4, pp. 1999–2022,
2014.
[23] A. Y. Aravkin, J. V. Burke, D. Drusvyatskiy, M. P. Friedlander, and
K. MacPhee, “Foundations of gauge and perspective duality,” arXiv
preprint arXiv:1702.08649, 2017.
[24] A. Yurtsever, M. Udell, J. A. Tropp, and V. Cevher, “Sketchy decisions:
Convex low-rank matrix optimization with optimal storage,” in Intl. Conf.
on Artificial Intelligence and Statistics (AISTATS), May 2017.
[25] P. Netrapalli, P. Jain, and S. Sanghavi, “Phase retrieval using alternating
minimization,” in Adv. Neural Inf. Process. Syst., 2013, pp. 2796–2804.
[26] P. Schniter and S. Rangan, “Compressive phase retrieval via generalized
approximate message passing,” IEEE Trans. Sig. Process., vol. 63, no. 4,
pp. 1043–1055, Feb. 2015.
[27] E. J. Candès, X. Li, and M. Soltanolkotabi, “Phase retrieval via Wirtinger
flow: Theory and algorithms,” IEEE Trans. Inf. Theory, vol. 61, no. 4,
pp. 1985–2007, Feb. 2015.
[28] Y. Chen and E. Candès, “Solving random quadratic systems of equations
is nearly as easy as solving linear systems,” in Adv. Neural Inf. Process.
Syst., 2015, pp. 739–747.
[29] G. Wang, G. B. Giannakis, and Y. C. Eldar, “Solving systems of random
quadratic equations via truncated amplitude flow,” arXiv: 1605.08285,
Jul. 2016.
[30] K. Wei, “Solving systems of phaseless equations via kaczmarz methods:
A proof of concept study,” Inverse Problems, vol. 31, no. 12, p. 125008,
2015.
[31] W.-J. Zeng and H. So, “Coordinate descent algorithms for phase retrieval,”
arXiv preprint arXiv:1706.03474, 2017.
[32] Z. Yuan and H. Wang, “Phase retrieval via reweighted Wirtinger flow,”
Applied optics, vol. 56, no. 9, p. 2418, 2017.
[33] P. Chen, A. Fannjiang, and G.-R. Liu, “Phase retrieval with one or
two diffraction patterns by alternating projections of the null vector,”
arXiv:1510.07379, Apr. 2015.
[34] Y. M. Lu and G. Li, “Phase transitions of spectral initialization for highdimensional nonconvex estimation,” arXiv preprint arXiv:1702.06435,
2017.
[35] M. Mondelli and A. Montanari, “Fundamental limits of weak recovery
with applications to phase retrieval,” arXiv preprint: 1708.05932, 2017.
[36] J. Sun, Q. Qu, and J. Wright, “A geometric analysis of phase retrieval,”
arXiv:1602.06664, Mar. 2016.
[37] P. Hand and V. Voroninski, “An elementary proof of convex phase
retrieval in the natural parameter space via the linear program phasemax,”
arXiv preprint: 1611.03935, 2016.
[38] ——, “Compressed sensing from phaseless Gaussian measurements via
linear programming in the natural parameter space,” arXiv preprint:
1611.05985, 2016.
[39] ——, “Corruption robust phase retrieval via linear programming,” arXiv
preprint: 1612.03547, 2016.
[40] A. Aghasi, A. Ahmed, and P. Hand, “BranchHull: convex bilinear
inversion from the entrywise product of signals with known signs,” arXiv
preprint: 1702.04342, 2017.
[41] O. Dhifallah, C. Thrampoulidis, and Y. M. Lu, “Phase retrieval via linear
programming: Fundamental limits and algorithmic improvements,” arXiv
preprint: 1710.05234, 2017.
[42] O. Dhifallah and Y. M. Lu, “Fundamental limits of phasemax for phase
retrieval: A replica analysis,” arXiv preprint arXiv:1708.03355, 2017.
[43] L. Schläfli, Gesammelte Mathematische Abhandlungen I. Springer Basel,
1953.
[44] J. G. Wendel, “A problem in geometric probability,” Math. Scand., vol. 11,
pp. 109–111, 1962.

15

[45] E. Gilbert, “The probability of covering a sphere with n circular caps,”
Biometrika, vol. 52, no. 3/4, pp. 323–330, Dec. 1965.
[46] Z. Füredi, “Random polytopes in the d-dimensional cube,” Disc. Comput.
Geom., vol. 1, no. 4, pp. 315–319, Dec. 1986.
[47] P. Bürgisser, F. Cucker, and M. Lotz, “Coverage processes on spheres
and condition numbers for linear programming,” Ann. Probab., vol. 38,
no. 2, pp. 570–604, 2010.
[48] T. Bendory and Y. C. Eldar, “Non-convex phase retrieval from STFT
measurements,” IEEE Trans. Inf. Theory, Aug. 2017.
[49] J. F. Kenney and E. Keeping, Mathematics of Statistics, Part 2. D. Van
Nostrand, 1951.
[50] L. Jacques, “A quantized Johnson–Lindenstrauss lemma: The finding of
Buffon’s needle,” IEEE Trans. Inf. Theory, vol. 61, no. 9, pp. 5012–5027,
Sept. 2015.
[51] F. Qi and Q.-M. Luo, “Bounds for the ratio of two gamma functions—
from Wendel’s and related inequalities to logarithmically completely
monotonic functions,” Banach J. Math. Anal, vol. 6, no. 2, pp. 132–158,
May. 2012.
[52] R. Chandra, Z. Zhong, J. Hontz, V. McCulloch, C. Studer, and
T. Goldstein, “PhasePack: A phase retrieval library,” arXiv preprint,
2017.
[53] T. Goldstein, C. Studer, and R. Baraniuk, “A field guide to forwardbackward splitting with a FASTA implementation,” arXiv:1411.3406,
Feb. 2014.
[54] C. A. Metzler, M. K. Sharma, S. Nagesh, R. G. Baraniuk, O. Cossairt,
and A. Veeraraghavan, “Coherent inverse scattering via transmission
matrices: Efficient phase retrieval algorithms and a public dataset,” in
Computational Photography (ICCP), 2017 IEEE International Conference
on. IEEE, 2017, pp. 1–16.
[55] K. Jaganathan, S. Oymak, and B. Hassibi, “Sparse phase retrieval: Convex
algorithms and limitations,” in Proc. IEEE Int. Symp. Inf. Theory (ISIT),
Jul. 2013, pp. 1022–1026.
[56] Y. Shechtman, A. Beck, and Y. C. Eldar, “GESPAR: efficient phase
retrieval of sparse signals,” IEEE Trans. Sig. Process., vol. 62, no. 4, pp.
928–938, Jan. 2014.
[57] C. Mortici, “New approximation formulas for evaluating the ratio of
gamma functions,” Math. Comput. Model., vol. 52, no. 1, pp. 425–433,
Jul. 2010.
[58] W. Gautschi, “Some elementary inequalities relating to the gamma and
incomplete gamma function,” J. Math. Phys., vol. 38, no. 1, pp. 77–81,
Apr. 1959.

Tom Goldstein is an Assistant Professor at the University of Maryland
in the Department of Computer Science. Before joining the faculty at
UMD, Tom completed his PhD at UCLA, and held research positions at
Stanford University and Rice University. Tom’s research focuses on efficient,
low complexity optimization routines. His work ranges from large-scale
computing on distributed architectures to inexpensive power-aware algorithms
for small-scale embedded systems. Applications of his work include scalable
machine learning, computer vision, and signal processing methods for wireless
communications. Tom has been the recipient of several awards, including
SIAM’s DiPrima Prize, and a Sloan Fellowship.

Christoph Studer (S’06–M’10–SM’14) received his Ph.D. degree in Electrical
Engineering from ETH Zurich in 2009. In 2005, he was a Visiting Researcher
with the Smart Antennas Research Group at Stanford University. From 2006 to
2009, he was a Research Assistant in both the Integrated Systems Laboratory
and the Communication Technology Laboratory (CTL) at ETH Zurich. From
2009 to 2012, Dr. Studer was a Postdoctoral Researcher at CTL, ETH Zurich,
and the Digital Signal Processing Group at Rice University. In 2013, he has
held the position of Research Scientist at Rice University. Since 2014, Dr.
Studer is an Assistant Professor at Cornell University and an Adjunct Assistant
Professor at Rice University.
Dr. Studer’s research interests include signal and information processing
as well as the design of digital very large-scale integration (VLSI) circuits.
His current research areas include applications in wireless communications,
nonlinear signal processing, optimization, and machine learning.
Dr. Studer received an ETH Medal for his M.S. thesis in 2006 and for his
Ph.D. thesis in 2009. He received a two-year Swiss National Science Foundation
fellowship for Advanced Researchers in 2011 and a US National Science
Foundation CAREER Award in 2017. In 2016, Dr. Studer won a Michael Tien
’72 Excellence in Teaching Award from the College of Engineering, Cornell
University. He shared the Swisscom/ICTnet Innovations Award in both 2010
and 2013. Dr. Studer was the winner of the Student Paper Contest of the 2007
Asilomar Conf. on Signals, Systems, and Computers, received a Best Student
Paper Award of the 2008 IEEE Int. Symp. on Circuits and Systems (ISCAS),
and shared the best Live Demonstration Award at the IEEE ISCAS in 2013.

