Skip to Main Content

Algebraic Theory of Phase Retrieval

Tamir Bendory
Dan Edidin

Communicated by Notices Associate Editor Steven Sam

Article cover

The phase retrieval problem first arose in the early 20th century from work on X-ray crystallography. In the last century, X-ray crystallography has developed into the leading method for elucidating the atomic structure of molecules, leading to enormous scientific breakthroughs: at least 25 Nobel Prizes have been awarded for work directly or indirectly involving crystallography. The phase retrieval problem also occurs in numerous other scientific and engineering applications, such as diffraction imaging, ptychography, ultra-short pulse characterization, speech processing, radar, and astronomy.

In its most general form, the phase retrieval problem can be written as:

where is the measurement vector, is a “sensing matrix,” defines the space of signals of interest, and the absolute value is taken entry-wise. Typically, the phaseless measurements are invariant under symmetry groups which depend on and . Thus, only the orbit of  under this intrinsic symmetry group can be recovered.

Without model error or noise, namely, when the measurement model 1 is accurate, phase retrieval is a problem of solving real quadratic equations. As such, it is naturally amenable to algebraic techniques, including commutative algebra, algebraic geometry, and invariant theory. In particular, algebraic methods are powerful techniques for proving that phase retrieval is theoretically possible; i.e., that a particular phase retrieval problem has a unique solution up to the action of its intrinsic symmetry group. Although these methods are typically not algorithmic, they can be used to provide theoretical validations of existing algorithms, and to unveil the fundamental limitations of different methods. The purpose of this article is to discuss recent advances in this growing field of research, and to publicize open problems that we believe will be of interest to mathematicians in general, and algebraists in particular.

In the rest of the article, we succinctly introduce five specific phase retrieval setups that we find important and intriguing, present known results, and delineate open mathematical questions. We conclude by discussing the limitations of the algebraic point of view, and how other mathematical fields, such as information theory, statistics, combinatorics, and optimization, can provide indispensable insights into the phase retrieval problem.

Phase retrieval with general linear measurements

The revival of interest in the mathematics of the phase retrieval problem in the last 15 years emerged from the study of phase retrieval models with a general sensing matrix. Specifically, in this setup or is a “generic” matrix (usually assumed to be random), or a frame,⁠Footnote1 with , and is either all of or . While measurements in practice are not random, and thus this line of work is of theoretical rather than applicable interest, it attracted the attention of the mathematical community and led to fascinating results in mathematics, statistics, and optimization; see for example CESV15.

1

For the purpose of this article, a frame is a collection of vectors which span or .

For a general real matrix , we can only expect to recover up to a global sign from the phaseless measurements . Namely, the intrinsic symmetric group is . Likewise, for a general complex matrix, we can recover only up to multiplication by a scalar , where is the circle group.

For real matrices, there is a remarkably elegant characterization of when every vector can be recovered (up to a global sign) from the phaseless measurements . To state the result, we introduce Definition 1.

Definition 1.

A matrix with row vectors , …, satisfies the complement property if, for every subset , the vectors or the vectors span .

Theorem 1.

BCE06 Every signal can be recovered, up to a sign, from if and only if satisfies the complement property.

Note that a necessary condition for an matrix to satisfy the complement property is . This immediately implies that if , then for any matrix t exist a pair of vectors with such that .

Example 1.

The matrix

has full rank but does not satisfy the complement property because the third row is the sum of the first two rows. In particular, if then neither nor span . Therefore, we know that not all vectors can be recovered up to a sign from the phaseless measurements . Indeed, if and , then and so , although .

For complex matrices, the situation is more nuanced.

Theorem 2.

CEHV15 For a generic complex matrix with , every signal can be recovered, up to multiplication by a scalar , from the phaseless measurements .

Here, we view the space of complex matrices as a real algebraic variety of dimension . The generic assertion in Theorem 2 means that the set of matrices for which the conclusion of the theorem holds contains a non-empty open set in the Zariski topology on . However, unlike the case for real matrices, the locus of complex matrices for which is injective is not itself open in the Zariski topology on because its complement is a semi-algebraic rather than algebraic subset.

Problem 1.

Is there any characterization of an complex matrix which guarantees that the map is injective?

Due to the lack of characterization, there are relatively few explicit examples of complex matrices that guarantee unique recovery of all vectors, up to a global phase. When , a family of examples was constructed by Bodmann and Hammen in BH15. If for an integer , a necessary condition for unique recovery is . However, for other values of , the optimal bound on is unknown. To the authors’ knowledge, if , the only examples of complex matrices which guarantee unique recovery of every signal (up to an symmetry) are when and ; the first examples were constructed by Vinzant Vin15.

Problem 2.

Determine the optimal bound on such that, if , then no signal can be uniquely recovered (up to an symmetry) from .

One step towards solving this problem was made by Heinosaari, Mazarella and Wolf HMW13, who used results from topology on minimal embeddings of projective spaces to prove that , where is the number of ones in the binary expansion of . When , this implies , leaving open the possibility that Vinzant’s construction does not yield the lowest possible bound.

X-ray crystallography

X-ray crystallography—a prevalent technology for determining the three-dimensional atomic structure of molecules—is by far the largest phase retrieval application. In X-ray crystallography, the signal is the electron density function of the crystal — a periodic arrangement of a repeating, compactly supported unit

where is the repeated motif and is a large, but finite, subset of a lattice ; the dimension is usually two or three. The crystal is illuminated with a beam of X-rays producing a diffraction pattern, which is equal to

where and are, respectively, the Fourier transforms of the signal and a Dirac ensemble defined on . Figure 1 presents an illustration of an X-ray crystallography experiment.⁠Footnote2

Figure 1.

An illustration of an X-ray crystallography experiment. The crystalline structure causes an X-ray beam to diffract into many specific directions. The crystallographic phase retrieval problem is to recover a signal (e.g., a molecular structure) from its diffraction pattern.

Graphic without alt text

As the size of the set grows (the size of the crystal), the support of the function is more concentrated in the dual lattice . Thus, the diffraction pattern is approximately equal to a discrete set of samples of on . This implies that the acquired data are the Fourier magnitudes of a -periodic signal on (or, equivalently, a signal on ), defined by its Fourier series

This signal is supported only at the sparsely spread positions of atoms. Thus, the crystallographic phase retrieval entails finding a -sparse signal satisfying , where is the discrete Fourier transform (DFT) matrix and . The problem can be equivalently formulated as recovering a -sparse signal from its periodic auto-correlation function:

While this article does not focus on algorithmic questions, we mention that a set of benchmark problems for evaluating crystallographic phase retrieval algorithms was designed in ELB18.

In most models, the signal is assumed to be real, so the auto-correlation function is a real quadratic function. In this case, , so we typically only consider the entries , …, of . Besides a global sign change, the periodic auto-correlation is also invariant under cyclic shifts , , and reflection . Hence, signal recovery is possible only up to the action of the group , where is the dihedral group. In particular, for any element and any signal with support , the signals have support and the same periodic auto-correlation as .

We say that two subsets are equivalent if they lie in the same orbit of the dihedral group . For a binary signal (all entries are zeros or ones), the auto-correlation is determined by the cyclic difference multi-set

where each difference is counted with multiplicity. For example, if , then and the binary auto-correlation vector is .

The phase retrieval problem for binary sets is equivalent to the combinatorial question of whether two subsets with the same cyclic difference multi-sets are dihedrally equivalent. This question does not have an affirmative answer. For example, the subsets of , and both have cyclic difference multi-sets but are not equivalent. However, simple numerical experiments indicate that this phenomenon is quite rare.

Problem 3.

For a fixed ratio , prove that the proportion of non-equivalent sets with the same difference set is asymptotic to 0 as .

At the other extreme, we can consider a model where the non-zero entries of the sparse vector are assumed to be arbitrary. In this case, the goal is to prove that a suitably generic sparse vector is determined up to the action of from its periodic auto-correlation function.

In BE20, we conjectured that if , then a generic vector with support in is uniquely determined by its periodic auto-correlation up to the action of the group . (Here denotes the number of distinct differences in the multi-set defined above.) Verifying this conjecture for a given value of can be done computationally, although not efficiently, in two steps.

The first step is support recovery; that is, verifying that the auto-correlation determines the support of a generic signal up to dihedral equivalence. Note that for a generic signal with support the support of is the set of distinct differences in the cyclic difference set . With this observation in hand, support recovery can be verified from the following conjecture.

Conjecture 1.

Let be two subsets of of size , and let and be the subspaces of vectors with support in and , respectively. If are non-equivalent subsets with , then the incidence variety has dimension strictly less than , where is defined in 5.

The reason that Conjecture 1 implies support recovery follows from the fact that the image of under the projection onto the first factor is the set of such that there exists a vector with . If , then we know that there is a non-empty Zariski open set of vectors for which there is no vector with the same auto-correlation function. Since the number of possible subsets is finite, affirming Conjecture 1 implies that for a generic vector any vector with the same auto-correlation must have an equivalent support set.

The second step is signal recovery for a signal with known support. Let be the subgroup of that preserves a subset . Signal recovery follows from the Conjecture 2.

Conjecture 2.

If then the incidence variety has dimension and degree .

To see that Conjecture 2 implies signal recovery, note that always contains linear subspaces of dimension consisting of pairs . Thus, if then these linear subspaces must necessarily be irreducible components of of maximal dimension . If in addition , then these are the only irreducible components of dimension . It follows that for a generic vector , the only vectors with the same auto-correlation as are the vectors of the form for . This is illustrated further in Example 3.

Conjectures 1 and 2 are similar to questions about the dimensions and degrees of secant varieties that arise in the literature on tensor decompositions. Indeed, if denotes the cyclic translation operator on , the symmetric two-tensor encodes the same information as the auto-correlation vector . Thus, the phase retrieval problem can also be viewed as a tensor decomposition problem.

Although we do not yet know how to prove either Conjecture 1 or Conjecture 2, they can be easily verified for small values of using a computer algebra system BE20. Despite the fact that solving systems of quadratic equations is known to be NP hard, the run time for the examples below was less than one second using a simple Macaulay2 script.

Example 2.

The following example verifies a single case of Conjecture 1. Let and be subsets of . Let

By computing the auto-correlations and explicitly, it can be shown that if and only if the following five equations are satisfied:

Thus, the incidence

is the algebraic subset of defined by the set of equations 6. Therefore, the generators of the ideal of are the five polynomials in the left-hand side of 6 included in . The Hilbert polynomial of this ideal is , where denotes the Hilbert polynomial of a polynomial ring in variables. This means that is a 3-dimensional affine algebraic subset of and therefore its image under to has dimension at most 3. It follows that for a generic vector in the 4-dimensional vector space , there is no corresponding vector with the same auto-correlation as .

Example 3.

We give an example which verifies a case of Conjecture 2. Let and let be subspace of with support in . The set is preserved by the element of order two, which is a reflection composed with a shift by two. Thus, for this subset. If and , then if and only if the following equations are satisfied:

Let be the ideal in generated by the five polynomials in the left-hand side of 7. The equations 7 are clearly satisfied if or . Thus, the 4-dimensional linear subspaces and are in , where denotes the algebraic subset of defined by the ideal . In addition, are also solutions to equations 7, so there are two additional 4-dimensional linear subspaces and in .

We calculated the Hilbert polynomial of the ideal  to be . Since contains four linear subspaces , then , namely,

has dimension at most 3. Hence, for generic , if then for some .

Interestingly, using tools from harmonic analysis and information theory, it was recently proven that -sparse, symmetric signals are determined uniquely from their periodic auto-correlation for for large enough GR21.

Future directions. Thus far we have discussed binary and generic models for X-ray crystallography. In practice, however, the model should account for sparse signals whose non-zero entries are taken from a finite (small) alphabet; this alphabet models the relevant type of atoms, such as hydrogen, oxygen, carbon, nitrogen, and so on. Any analysis of this model would be probabilistic, but we expect that in the case of fixed finite alphabet, the probability that a sparse signal can be recovered from its auto-correlation is asymptotic to one as the signal length .

Fourier phase retrieval

In this section, we consider the problem of recovering a signal from its aperiodic auto-correlation:

where if and