PaperlayerView exact revision on arXiv

Source-linked AI summary

Solving (most) of a set of quadratic equalities: Composite optimization for robust phase retrieval

John C. Duchi, Feng Ruan

arXiv:1705.02356v2math.STcs.ITmath.OC

TL;DR

Phase retrieval requires recovering a signal from quadratic equalities despite a nonsmooth, nonconvex objective and potentially corrupted measurements. The paper applies prox-linear composite optimization, obtaining high-probability recovery with a constant-factor measurement-to-dimension ratio, local quadratic convergence, and essentially no tuning.

  • Problem

    Phase retrieval involves solving quadratic equalities whose ℓ1 penalty is nonsmooth and nonconvex, while corrupted measurements make recovery more challenging.

  • Method

    The paper linearizes the smooth component of a convex-composite objective and repeatedly minimizes regularized convex prox-linear models.

  • Results

    m/n greater than a numerical constant yields extremely high-probability local quadratic convergence and exact recovery under appropriate random measurement models, including a constant fraction of adversarial corruptions.

  • Takeaways & Limitations

    The method provides a tuneless sequence of convex problems that applies to broad measurement matrices and achieves practical real-valued recovery near m/n = 2.

  • Takeaways & Limitations

    The objective remains nonsmooth and nonconvex, and guarantees require stability, an operator-norm bound, and suitable initialization.

Abstract

from arXiv · show

We develop procedures, based on minimization of the composition $f(x) = h(c(x))$ of a convex function $h$ and smooth function $c$, for solving random collections of quadratic equalities, applying our methodology to phase retrieval problems. We show that the prox-linear algorithm we develop can solve phase retrieval problems---even with adversarially faulty measurements---with high probability as soon as the number of measurements $m$ is a constant factor larger than the dimension $n$ of the signal to be recovered. The algorithm requires essentially no tuning---it consists of solving a sequence of convex problems---and it is implementable without any particular assumptions on the measurements taken. We provide substantial experiments investigating our methods, indicating the practical effectiveness of the procedures and showing that they succeed with high probability as soon as $m / n \ge 2$ when the signal is real-valued.

1 Introduction

The paper formulates phase retrieval as composite optimization and develops a prox-linear method that handles nonsmooth, nonconvex objectives through convex subproblems. Its theory and experiments establish recovery with near-minimal measurement scaling, robustness to corrupted measurements, and practical success without tuning.

  • Problem: Phase retrieval seeks x ∈ C^n whose quadratic measurements match nonnegative observations b_i = |⟨a_i, x⟩|^2, but the worst-case formulation is NP-hard.Optical sensors motivate the problem because they measure squared amplitudes rather than diffraction phases.
  • Problem: An ℓ1 exact-penalty objective replaces quadratic equalities with residual costs, tolerates gross measurement errors, and recovers signals only up to phase or sign.The resulting objective is nonsmooth and nonconvex, including cases that are not locally convex near the solution.
  • Method: Linearizing only the smooth inner map c in f(x) = h(c(x)) produces a convex local surrogate, yielding the prox-linear algorithm for the composite objective.The method repeatedly minimizes regularized convex models, with complex-valued problems treated using Wirtinger calculus.
  • Method: A stepsize α ≤ 1/(βL) guarantees descent and approximate stationarity when h is L-Lipschitz and ∇c is β-Lipschitz.Each iteration is a convex quadratic program in n variables, solvable with standard convex programming methods.
  • Guarantees: m/n greater than a numerical constant suffices for extremely high-probability local quadratic convergence and exact recovery under appropriate random measurement models.The required scaling is optimal up to numerical constants because m ≳ n is necessary for exact recovery.
  • Context: Unlike SDP approaches that manipulate an n × n matrix and scale roughly as n^3/poly(ϵ), the method maintains only a vector but solves convex quadratic programs per iteration.This trades the computational burden of full matrix optimization for more expensive iterations than typical gradient methods.
  • Analysis: The convergence and recovery analysis relies primarily on measurement stability, an operator-norm bound, and suitable initialization across noiseless and corrupted settings.Small-ball techniques establish stability for broad classes of random matrices, without requiring light-tailed measurements.

2 Composite Optimization, Algorithm and Convergence Analysis

The paper analyzes prox-linear optimization for composite objectives and applies it to robust phase retrieval. Its convergence framework uses stability, quadratic model accuracy, and a sufficiently good initializer.

  • Composite optimization: The prox-linear method linearizes the smooth component while preserving convexity in the resulting local model.For phase retrieval, the iteration becomes a convex quadratic problem.
  • Composite optimization: The complex-valued formulation uses Wirtinger calculus because the real-valued function c lacks an ordinary complex derivative.The Hermitian-conjugate Jacobian treats c as a mapping from R^2n to R.
  • Phase retrieval application: The three required conditions are established for noiseless and adversarially perturbed phase-retrieval measurements.The noiseless objective uses f(x) = 1/m |||Ax|^2 − b||_1, with global minimizers representing the signal’s phase-equivalent solutions.
  • Convergence conditions: Quadratic convergence relies on stability of f near its global minimizers and a quadratic approximation bound for the linearized model.These are formalized as Conditions C1 and C2.
  • Convergence conditions: The analysis allows inaccurate solutions of the intermediate convex quadratic programs through additive accuracy parameters ϵ_k.Theorem 1 explicitly incorporates these errors into relative error measures.

1. Stability (Condition C1),

The convergence proof derives per-iteration progress and then converts it into quadratic convergence, while accounting for inexact subproblem solutions. The theorem’s guarantees depend on quadratic model accuracy, initialization, and controlled errors.

  • Required conditions: The algorithm requires an initializer satisfying the constant relative error bound dist(x_0, X⋆) ≤ dist(0, X⋆) / λ.This is listed as the third condition supporting the theorem’s phase-retrieval guarantees.
  • Practical implications: About five iterations suffice for machine-precision solutions when intermediate problems are solved to near machine precision.Quadratic convergence can also persist when accuracies decrease sufficiently quickly, such as ϵ_k ≤ 2^-2^k.
  • Proof strategy: The convergence proof first establishes a per-iteration progress guarantee and then uses it to prove quadratic convergence.This two-step structure separates local progress from the final rate argument.
  • Proof strategy: The exact prox-linear subproblem is strongly convex after adding the quadratic regularizer around the current iterate.The proof compares its exact minimizer x⋆,k+1 with the computed iterate x_k+1.
  • Required conditions: The quadratic approximation condition bounds the discrepancy between the original objective and its linearized model.For phase retrieval, the bound follows from the quadratic expansion of |⟨a, y⟩|^2.
  • Proof strategy: When subproblem errors vanish, recursively applying the error inequality yields the theorem’s quadratic-convergence statement.With nonzero errors, a sequence lemma controls how the ϵ_k terms propagate through iterations.

3 Noiseless Phase Retrieval Problem

The paper analyzes noiseless phase retrieval by establishing stability, quadratic approximability, and a suitable initializer under random measurement assumptions. These conditions yield high-probability recovery with measurement counts proportional to the signal dimension, including local quadratic convergence.

  • Overview: The noiseless analysis separately establishes stability, quadratic approximability, and a good initializer for phase retrieval.The measurement vectors are treated as independent, identically distributed random copies, and each condition is developed with progressively stronger assumptions.
  • Real case: Real-valued stability follows from a small-ball condition requiring measurement vectors to have sufficient support in every pair of directions.The condition needs only a positive probability that both inner products are bounded away from zero, without light-tail assumptions.
  • Real case: The real stability condition holds with high probability when m/n is at least a numerical constant, under the stated random measurement assumptions.The paper compares this guarantee with earlier results requiring subgaussian and isotropic measurements, whereas its small-ball result imposes weaker restrictions.
  • Complex case: For complex measurements, stability requires both a small-ball assumption and sufficient directional uniformity, with the guarantee applying when m is proportional to n.The ratio between the growth constant and sub-Gaussian constant is bounded for rotationally invariant measurement distributions.
  • Quadratic approximation: The quadratic approximation condition also holds with high probability when m/n is at least a numerical constant, under a sub-Gaussian measurement assumption.The operator-norm control needed for this step is obtained from light-tailed measurement vectors.
  • Recovery guarantees: The resulting prox-linear method achieves local quadratic convergence and exact recovery with high probability from an initializer within a constant-factor neighborhood.The method requires convex quadratic subproblems per iteration; experiments report better recovery properties but roughly four-times slower overall runtime than the best gradient-based methods.

4 Phase retrieval with outliers

The paper extends phase retrieval to settings with corrupted measurements, showing that prox-linear recovery remains effective under random or adversarial outliers. With suitable random measurements and initialization, the method converges quadratically with high probability once m/n exceeds a numerical constant.

  • Corruption models: The corruption framework includes a fixed failure fraction with randomly selected outlier indices, while corrupted values may depend adversarially on measurements under Model M2.Model M1 requires noise-measurement independence; Model M2 allows dependence for outlying observations.
  • Stability under outliers: The objective remains stable under arbitrary corruptions when measurement vectors are light-tailed, without requiring Gaussian measurements.The analysis retains the noiseless quadratic-approximation conditions and requires a suitable initializer.
  • Stability under outliers: For Gaussian measurements, stability holds when pfail < 1/π ≈ .318 and m/n exceeds a numerical constant.The resulting bound controls objective growth relative to distance from the solution set.
  • Initialization: The corrupted-data initializer estimates signal direction spectrally and estimates its norm by minimizing a robust one-dimensional objective.Selecting measurements with small observed values supports the directional estimate, while the norm estimate remains effective for non-Gaussian data.
  • Recovery guarantee: Under the stated assumptions, prox-linear iterations converge quadratically to the signal up to sign, with exponentially high probability, once m/n exceeds a numerical constant.The paper reports non-trivial degradation as the failure fraction grows.

5 Optimization methods

The implementation solves each prox-linear subproblem as a convex program, using interior-point methods at small scale and POGS or conjugate gradients for larger structured systems. Inner accuracy is tightened across outer iterations to match the algorithm’s quadratic convergence.

  • Subproblem formulation: Each prox-linear subproblem is rewritten as a convex optimization problem involving a linear constraint Bx = y.The reformulation uses B = DA after linearizing the phase-retrieval objective.
  • Small-scale solvers: Small-scale instances are solved with standard interior-point software through Mosek and Convex.jl.The paper uses Mosek for the resulting quadratic-program formulation.
  • Large-scale solvers: Large-scale instances use POGS, an ADMM variant that minimizes f(x)+g(y) subject to Bx = y.The authors report POGS as empirically best among the specialized methods they tested.
  • Large-scale solvers: POGS updates are inexpensive except for the graph-projection matrix inversion, whose factorization can be cached for repeated iterations.With caching, the repeated multiplication costs n^2 + nm time.
  • Accuracy control: The inner stopping tolerance is decreased during outer prox-linear iterations, reflecting the method’s quadratic convergence guarantees.The experiments use ϵ = 10^-5 initially and ϵ = 10^-8 later with POGS.
  • Structured measurements: For structured Fourier, Hadamard, or sparse matrices, conjugate gradients exploit fast matrix-vector products to solve the positive-definite system.The experiments use unpreconditioned conjugate gradients for structured random-sign measurements.

6 Experiments

Experiments compare prox-linear composite optimization with established phase-retrieval methods across random signals, initializations, measurement ratios, outliers, and image data. The method achieves strong recovery near m/n = 2, handles outliers, and recovers real images using structured sensing.

  • Experimental design: The experiments vary measurement matrices, signal distributions, initializations, measurement ratios, outlier rates, and subproblem accuracy.Signals are Gaussian or random signed vectors, while experiments use Gaussian and structured sensing matrices.
  • Random-signal recovery: At m/n = 2, prox-linear recovery is around 0.6 for n = 400 and exceeds 78%, 88%, and 91% for n = 1000, 2000, and 3000.The larger-dimensional figures report 95% confidence intervals of ±4, ±6, and ±6 percentage points, respectively, for the latter three rates.
  • Random-signal recovery: For larger problems, prox-linear outperforms truncated amplitude flow, with only one trial favoring TAF when n = 3000 and m/n = 2.Both methods succeed in all reported trials when m/n is at least 2.2.
  • Outlying measurements: With outlying measurements, the composite objective yields better recovery than median-truncated Wirtinger flow across the tested initializations and dimensions.Figures 3 and 4 report 100 experiments per plot, using accurate or approximate prox-linear subproblem solves.
  • Computational behavior: In the high-measurement regime, pfail = 0 requires no more than 7 prox-linear iterations, while m/n = 8 uses precisely 5 iterations for pfail ≤ .15.The authors identify this fixed small iteration count as quadratic convergence.
  • Real-image recovery: The method successfully recovers all eight real color images to relative accuracy at least 10^-4, averaging 15100 matrix-vector multiplications.The images use random Hadamard sensing matrices and problems of up to n = 222; one RNA reconstruction improves after an accurate refinement step.

A.1 Proof of Proposition 1

The proof controls the relevant function and threshold-set classes through VC dimension, showing the needed complexity is linear in n. VC concentration then supplies the probabilistic bound used in Proposition 1.

  • The function class F is controlled through the VC dimension of threshold sets involving |⟨x,u⟩| ∧ |⟨x,v⟩|.
  • The threshold-set collection G has VC dimension at most Cn for a numerical constant C.
  • The proof embeds G in set operations involving half-plane classes, whose VC dimension is at most n + 2.
  • Standard VC-concentration inequalities are applied to the threshold classes to obtain the required uniform probabilistic control.

A.2 Proof of Proposition 2

The proof establishes a uniform bound over low-rank Hermitian matrices by relating rank-two matrix norms to signal distances, controlling moments, and applying covering and concentration arguments.

  • A rank-two matrix X = xx^H − yy^H is linked to phase-invariant distances between x and y through its Frobenius norm.
  • The argument reduces the theorem to a high-probability inequality over rank-at-most-two Hermitian matrices normalized in Frobenius norm.
  • Sub-Gaussian moment bounds control the fourth moments needed for the expectation and variance estimates.
  • An epsilon-cover of low-rank Hermitian matrices provides a finite approximation set for the subsequent covering argument.
  • The resulting concentration and covering bounds yield the proposition’s uniform estimate under an appropriate truncation choice.

A.3 Proof of Proposition 3: deterministic part

The deterministic proof decomposes the spectral initialization error and shows that a small operator-norm perturbation preserves a useful directional estimate and signal-scale estimate.

  • A spectral eigengap separates the target direction from the remaining eigenvalues, enabling eigenvector perturbation control.
  • The initialization matrix is decomposed into error components whose operator norms are bounded separately.
  • The selected index sets and associated events organize the bounds on truncation, conditioning, and directional errors.
  • The directional estimate is close to the target direction whenever the aggregate perturbation has sufficiently small operator norm.
  • The initialization scale satisfies (1 − ϵ)∥x⋆∥2^2 ≤ br^2 ≤ (1 + ϵ)∥x⋆∥2^2.

A.4 Proof of Proposition 3: high probability events

The high-probability proof verifies several concentration events governing the selected measurements, truncated sums, and conditional sub-Gaussian behavior.

  • The proof reduces the proposition to showing that events E1 through E5 each occur with high probability.
  • Event E1 satisfies P(E1) ≥ 1 − exp(−cmϵ^2/σ^4) when m/n ≥ σ^4/(cϵ^2).
  • Event E2 satisfies P(E2) ≥ 1 − exp(−2ϵ^2κ^2m).
  • Hoeffding’s inequality controls the cardinality difference between the selected index sets I−ϵ and I+ϵ.
  • Conditioning on I−ϵ preserves independence and yields a conditional sub-Gaussian law with parameter σ^2 log(e/p0(d⋆)).

B.1 Proof of Proposition 4

The proof establishes concentration bounds for quadratic measurement quantities, uniformly over x, with high probability when m/n is sufficiently large. These bounds support the proposition’s control of the relevant matrix and error terms.

  • Inequality (36) shows that the relevant matrix quantity is well concentrated.
  • The resulting bound holds for all x ∈ R^n with probability at least 1 − e^(−mt), after choosing t so that n/m + t ≤ 1.
  • A concentration lemma for independent sub-Gaussian measurement vectors holds with probability at least 1 − e^(−cmt^2) − e^(−cm) when m/n ≥ C.
  • The proof combines the quadratic difference identity with a preceding lemma to obtain a uniform bound over x.

B.2 Proof of Proposition 5

The proof first bounds the minimizing radius deterministically, then uses matrix concentration to show the random terms are small with high probability. Convexity confines every minimizer to a narrow interval around the true radius.

  • All minimizers of G lie in [1 ± δ]∥x⋆∥² whenever δ ≤ 1.
  • Matrix concentration controls both quadratic forms and the operator norm of A^T A with probability at least 1 − e^(−mt).
  • With probability at least 1 − 2e^(−mt), the radius error parameter δ satisfies the stated concentration bound.
  • Choosing t proportional to (1 − 2p_fail)^2/σ^4 yields the desired high-probability result.
  • Convexity implies that every minimizer of the comparison function lies in the derived region, which transfers to the radius minimizers.

B.3 Proof of Proposition 6

The proof analyzes the spectral initialization matrix by decomposing it into four terms and controlling each under independent and adversarial noise models. Eigenvector perturbation then yields recovery guarantees.

  • Term control: Z_2 is approximately a scaled multiple of I − d⋆d⋆^T, while Z_1 is small because it sums products of independent vectors.
  • Noise models: The outlier contribution Z_3 receives separate bounds under independent noise M1 and adversarial noise M2.
  • Conclusion: The proposition follows after combining the four lemmas with the eigenvector perturbation lemma.
  • Spectral perturbation: When z_2 > z_0, the comparison matrix has minimal eigenvector d⋆ and eigengap z_2 − z_0, enabling eigenvector perturbation analysis.
  • Spectral separation: Under the stated probability events, the proof obtains a positive separation z_2 ≥ z_0 + (1 − 2p_fail)δ_q/2.

C.2 Proof of Lemma B.1

The proof of Lemma B.1 uses Bernstein concentration, Gaussian matrix concentration, and an epsilon-cover of the sphere to obtain uniform control. Auxiliary Gaussian calculations establish the needed moment relation.

  • Uniform concentration: The proof uses a standard concentration-and-covering-number argument with minor technical refinements.
  • Concentration ingredients: Bernstein’s inequality controls products |⟨u, a_i⟩⟨v, a_i⟩|, while Lemma 3.1 controls the operator norm of A^T A.
  • Uniform concentration: An epsilon-cover of S^(n−1) extends pointwise bounds uniformly over pairs of sphere vectors.
  • Parameter choice: Choosing ε = 3√(n/m) gives the stated regime when m/n is sufficiently large.
  • Auxiliary moment bound: The auxiliary Gaussian argument relates conditional second moments below and above a threshold c to prove the required lemma.
Loading 1705.02356v2…