|
The aim of noisy phase retrieval is to estimate a signal x 0 ∈ C d from m noisy intensity measurements b j = | 〈 a j , x 0 〉 | 2 + η j , j = 1 , … , m , where a j ∈ C d are known measurement vectors and η = ( η 1 , … , η m ) ⊤ ∈ R m is a noise vector. A commonly used estimator for x 0 is to minimize the intensity-based loss function, i.e., x ˆ : = argmin x ∈ C d ∑ j = 1 m ( | 〈 a j , x 〉 | 2 − b j ) 2 . Although many algorithms for solving the intensity-based estimator have been developed, there are very few results about its estimation performance. In this paper, we focus on the performance of the intensity-based estimator and prove that the error bound satisfies min θ ∈ R ‖ x ˆ − e i θ x 0 ‖ 2 ≲ min { ‖ η ‖ 2 m 1 / 4 , ‖ η ‖ 2 ‖ x 0 ‖ 2 ⋅ m } under the assumption of m ≳ d and a j ∈ C d , j = 1 , … , m , being complex Gaussian random vectors. We also show that the error bound is rate optimal when m ≳ d log m . In the case where x 0 is an s -sparse signal, we present a similar result under the assumption of m ≳ s log ( e d / s ) . To the best of our knowledge, our results provide the first theoretical guarantees for both the intensity-based estimator and its sparse version. Introduction Assume that b j : = | 〈 a j , x 0 〉 | 2 + η j , j = 1 , … , m where a j ∈ C d are known sensing vectors and η : = ( η 1 , … , η m ) ⊤ ∈ R m is a noise vector. Throughout this paper, we make the assumption that the noise vector, denoted as η , is either fixed or random, and is independent of the measurement vectors a j , where j = 1 , … , m . To estimate x 0 ∈ C d from b : = ( b 1 , … , b m ) ⊤ ∈ R m is referred to as phase retrieval . Due to the physical limitations, optical sensors can record only the modulus of Fraunhofer diffraction pattern while losing the phase information, and hence phase retrieval has many applications in fields of physical sciences and engineering, which includes X-ray crystallography [19], [28], microscopy [27], astronomy [10], coherent diffractive imaging [33], [17] and optics [39] etc. Despite its simple mathematical form, it has been shown that to reconstruct a finite-dimensional discrete signal from its Fourier transform magnitudes is generally NP-complete [32]. Based on the least squares criterion, one can employ the following intensity-based empirical loss to estimate x 0 : min x ∈ C d ∑ j = 1 m ( | 〈 a j , x 〉 | 2 − b j ) 2 . For the case where x 0 is sparse, the following Lasso-type program can be utilized to estimate x 0 : min x ∈ C d ∑ j = 1 m ( | 〈 a j , x 〉 | 2 − b j ) 2 s.t. ‖ x ‖ 1 ≤ R , where R is a parameter which specifies a desired sparsity level of the solution. An advantage of the intensity-based estimator (1) is that the objective function is differentiable based on the Wirtinger derivatives. One can therefore try to find its minimum with high-order algorithms such as trust-region and Gauss-Newton methods. Though the objective functions are non-convex, the strategy of spectral initialization plus local gradient descent can be adopted to solve (1) and (2) efficiently under Gaussian random measurements. For instance, it has been proved that when m ≳ d and ‖ η ‖ ∞ ≲ ‖ x 0 ‖ 2 2 , with high probability the truncated spectral method given in [8] can return an initial guess z 0 which is close to the target signal in the real case, namely, min θ ∈ [ 0 , 2 π ) ‖ z 0 − e i θ x 0 ‖ 2 ≤ δ ‖ x 0 ‖ 2 for any fixed relative error tolerance δ . With this in place, the update rules such as Wirtinger Flow [5], Trust-Region [34] and Gauss-Newton [16] methods could find a global solution to (1) at least in the noiseless case. In the noiseless case, i.e., b j = | 〈 a j , x 0 〉 | 2 , j = 1 , … , m , the solution to (1) is exactly x 0 (up to a unimodular constant) if m ≥ 4 d − 4 and a j , j = 1 , … , m , are generic vectors in C d [7]. However, in the noisy case, the distance between the solution to (1) and the true signal x 0 is still unknown. The aim of this paper is to investigate the theoretical performance of (1) and (2). For the last two decades, numerous algorithms have been designed for phase retrieval, particularly in the noiseless case. These algorithms can be classified into two categories: convex and non-convex methods. The convex methods rely on the “matrix-lifting” technique which lifts the quadratic system to a linear rank-one positive semi-definite program. More specifically, a rank one matrix X = x x ⁎ is introduced to linearize the quadratic constrains and then a nuclear norm minimization is adopted as a convex surrogate of the rank constraint. Such methods include PhaseLift [6], [4], PhaseCut [38] etc. While convex methods have strong theoretical guarantees, they require solving a semi-definite program in the “lifted” space C d 2 instead of C d , where d is the signal dimension. As a result, the computational complexity and memory requirements can become very high. The non-convex methods operate directly on the original space, which achieves significantly improved computational performance. The oldest non-convex algorithms for phase retrieval are based on alternating projection including Gerchberg-Saxton [17] and Fienup [13], but lack of theoretical guarantees. The first non-convex algorithm with theoretical guarantees was given by Netrapalli et al. who showed that the AltMinPhase [29] algorithm converges linearly to the true solution up to a global phase with O ( d log 3 d ) resampling Gaussian random measurements. In [5], Candès, Li and Soltanolkotabi developed the Wirtinger Flow (WF) to solve (1) and proved WF algorithm can achieve the linear convergence with O ( d log d ) Gaussian random measurements. Lately, Chen and Candès improved the result to O ( d ) Gaussian random measurements by Truncated Wirtinger Flow (TWF) [8]. In [16], Gao and Xu proposed a Gauss-Newton algorithm to solve (1) and proved the Gauss-Newton method can achieve quadratic convergence for the real-valued signals with O ( d log d ) Gaussian random measurements. In [34], Sun, Qu, and Wright proved that, for O ( d log 3 d ) Gaussian random measurements, the objective function of (1) has a benign geometric landscape: (1) all local minimizers are global, and (2) the objective function has a negative curvature around each saddle point.1 They also developed the Trust-Region method to find a global solution. Another alternative approach for phase retrieval is to solve the following amplitude-based empirical loss: min x ∈ C d ∑ j = 1 m ( | 〈 a j , x 〉 | − ψ j ) 2 , where ψ j : = b j , j = 1 , … , m . For Gaussian random measurements, through an appropriate initialization, many algorithms can be used to solve (3) successfully such as Truncated Amplitude Flow (TAF) [40], Reshaped Wirtinger Flow (RWF) [44] and Perturbed Amplitude Flow (PAF) [15]. It has been proved that TAF, RWF and PAF algorithms converge linearly to the true solution up to a global phase under O ( d ) Gaussian random measurements. For sparse phase retrieval, a standard ℓ 1 relaxation technique leads to the corresponding sparse intensity-based estimator (2). It has been shown that when the noises are independent centered sub-exponential random variables with maximum sub-exponential norm σ and m ≳ α σ 2 s 2 log ( m d ) , in the real case, the initialization procedure given in [2] can return an initial guess z 0 which satisfies supp ( z 0 ) ⊂ supp ( x 0 ) and is very close to the target signal with high probability, where s is the sparsity level and α σ is a constant related to σ . Next, the projection gradient method [12] could be used to find a global minimizer to (2). Such two-step procedure has also been used in various signal processing and machine learning problems, such as Blind Deconvolution [26], matrix completion [35] and sparse recovery [31]. We refer the reader to papers [41], [18], [43] for accounts of recent developments in the algorithms of sparse phase retrieval. The theoretical results concerning the injectivity of sparse phase retrieval can be found in [42], [22]. Throughout this paper, we assume the measurement vectors a j ∈ C d , j = 1 , … , m are i.i.d. complex Gaussian random vectors. Here we say a ∈ C d is a complex Gaussian random vector if a ∼ 1 / 2 ⋅ N ( 0 , I d ) + i / 2 ⋅ N ( 0 , I d ) . We write z ∈ S C d − 1 if z ∈ C d and ‖ z ‖ 2 = 1 . We use the notations ‖ ⋅ ‖ 2 and ‖ ⋅ ‖ ⁎ to denote the operator norm and nuclear norm of a matrix, respectively. For any A , B ∈ R , we use A ≲ B to denote A ≤ C 0 B where C 0 ∈ R + is an absolute constant. The notion ≳ can be defined similarly. Moreover, A ≍ B means that there exist constants C 1 , C 2 > 0 such that C 1 A ≤ B ≤ C 2 A . In this paper, we use C , c and the subscript (superscript) form of them to denote universal constants whose values vary with the context. As stated earlier, the two-step strategy of spectral initialization plus local gradient descent can be used to solve the intensity-based estimator (1). To our knowledge, there is no result concerning the reconstruction error of (1) for noisy phase retrieval from the theoretical viewpoint. The goal of this paper is to study the estimation performance of the intensity-based estimator (1) and its sparse version (2). Our first result shows that the estimation error is quite small and bounded by the average noise per measurement, as stated below. We emphasize that this theorem does not assume any particular structure on the noise η . Theorem 1.1 Suppose that the measurement vectors a j ∼ 1 / 2 ⋅ N ( 0 , I d ) + i / 2 ⋅ N ( 0 , I d ) are i.i.d. complex Gaussian random vectors and the measurement number m ≳ d . Then the following holds with probability at least 1 − exp ( − c m ) : For all x 0 ∈ C d , any solution x ˆ ∈ C d to (1) with b j = | 〈 a j , x 0 〉 | 2 + η j , j = 1 , … , m , satisfies min θ ∈ [ 0 , 2 π ) ‖ x ˆ − e i θ x 0 ‖ 2 ≤ C min { ‖ η ‖ 2 m 1 / 4 , ‖ η ‖ 2 ‖ x 0 ‖ 2 ⋅ m } . Here, C and c are positive absolute constants. According to Theorem 1.1, the following holds with high probability: min θ ∈ [ 0 , 2 π ) ‖ x ˆ − e i θ x 0 ‖ 2 ≲ min { ‖ η ‖ 2 m 1 / 4 , ‖ η ‖ 2 ‖ x 0 ‖ 2 ⋅ m } = { ‖ η ‖ 2 m 1 / 4 if ‖ x 0 ‖ 2 < ‖ η ‖ 2 m 1 / 4 ‖ η ‖ 2 ‖ x 0 ‖ 2 ⋅ m if ‖ x 0 ‖ 2 ≥ ‖ η ‖ 2 m 1 / 4 . The following theorem presents a lower bound for the estimation error of any fixed x 0 , assuming m ≳ d log m and noise with the following structures: ‖ η ‖ 2 ≍ m ‖ x 0 ‖ 2 2 , ‖ η ‖ ∞ ≲ log m m ‖ η ‖ 2 , m ‖ η ‖ 2 ≲ | ∑ j = 1 m η j | ≲ m ‖ x 0 ‖ 2 2 . The result demonstrates that the estimation error presented in Theorem 1.1 is rate optimal for certain η , provided m ≳ d log m . Theorem 1.2 Suppose that the measurement vectors a j ∼ 1 / 2 ⋅ N ( 0 , I d ) + i / 2 ⋅ N ( 0 , I d ) are i.i.d. complex Gaussian random vectors. For any fixed vector x 0 ∈ C d and any fixed 0 < ε 0 ≤ 1 / ( 2 C ) , assume that η ∈ R m is a noise vector satisfying ε 0 m ‖ x 0 ‖ 2 2 ≤ ‖ η ‖ 2 ≤ m ‖ x 0 ‖ 2 2 2 C , and ‖ η ‖ ∞ ≲ log m m ‖ η ‖ 2 , m ‖ η ‖ 2 ≲ | ∑ j = 1 m η j | ≲ m ‖ x 0 ‖ 2 2 . If m ≥ C ( ε 0 ) d log m then the following holds with probability at least 1 − c 1 ( ε 0 ) m − 1 − c 2 exp ( − c 3 ( ε 0 ) d ) : any solution x ˆ ∈ C d to (1) with b j = | 〈 a j , x 0 〉 | 2 + η j , j = 1 , … , m , satisfies min θ ∈ [ 0 , 2 π ) ‖ x ˆ − e i θ x 0 ‖ 2 ≳ min { ‖ η ‖ 2 m 1 / 4 , ‖ η ‖ 2 ‖ x 0 ‖ 2 ⋅ m } . Here, C ( ε 0 ) , c 1 ( ε 0 ) , c 3 ( ε 0 ) are positive constants depending only on ε 0 , c 2 is a universal constant, and C is the constant in Theorem 1.1 . Remark 1.3 In Theorem 1.2, we require η satisfies the conditions ε 0 m ‖ x 0 ‖ 2 2 ≤ ‖ η ‖ 2 ≤ m ‖ x 0 ‖ 2 2 / 2 C , and ‖ η ‖ ∞ ≲ log m m ‖ η ‖ 2 , m ‖ η ‖ 2 ≲ | ∑ j = 1 m η j | ≲ m ‖ x 0 ‖ 2 2 . In fact, there exist many noises satisfying them. For instance, if η j are generated independently according to the biased Gaussian distribution, i.e., η j ∼ ‖ x 0 ‖ 2 2 ⋅ N ( μ , σ 2 ) for some constant 0 < μ ≤ 1 / ( 2 C ) , then the noise vector η satisfies those conditions with high probability. For the case where the noise vector η does not satisfy those conditions, it is possible to improve the lower bound given in Theorem 1.2. We next turn to the phase retrieval for sparse signals. This is motivated by the signal x 0 ∈ C d admitting a sparse representation under some linear transformation in many applications. Without loss of generality, we assume that x 0 ∈ C d is an s -sparse vector and wish to estimate x 0 from b = ( b 1 , … , b m ) ⊤ by solving min x ∈ C d ∑ j = 1 m ( | 〈 a j , x 〉 | 2 − b j ) 2 s.t. ‖ x ‖ 1 ≤ R . The estimation performance of (4) is stated as follows. Theorem 1.4 Suppose that the measurement vectors a j ∼ 1 / 2 ⋅ N ( 0 , I d ) + i / 2 ⋅ N ( 0 , I d ) are i.i.d. complex Gaussian random vectors and the measurement number m ≳ s log ( e d / s ) . Then the following holds with probability at least 1 − 5 exp ( − c m ) where c is a constant: for any s-sparse vector x 0 ∈ C d , min θ ∈ [ 0 , 2 π ) ‖ x ˆ − e i θ x 0 ‖ 2 ≲ min { ‖ η ‖ 2 m 1 / 4 , ‖ η ‖ 2 ‖ x 0 ‖ 2 ⋅ m } , where x ˆ ∈ C d is a solution to (4) with parameter R : = ‖ x 0 ‖ 1 . Remark 1.5 In [2], Cai, Li, and Ma establish the estimation error of the Thresholded Wirtinger Flow algorithm for real-valued signals under the assumption of centered sub-exponential noise. In short, they show that, with probability at least 1 − 47 / m − 10 / e s , the estimator x ˆ given by the Thresholded Wirtinger Flow algorithm obeys min { ‖ x ˆ + x 0 ‖ 2 , ‖ x ˆ − x 0 ‖ 2 } ≲ σ ‖ x 0 ‖ 2 s log d m provided m ≥ O ( s 2 log ( m d ) ) , where σ : = max 1 ≤ j ≤ m ‖ η j ‖ ψ 1 . Although the estimation error is slightly better than the upper bound given in Theorem 1.4, our result holds for any noise structure and complex-valued signals. Moreover, the probability of failure in Theorem 1.4 is exponentially small in the number of measurements. We first introduce the estimation performance of PhaseLift for noisy phase retrieval. In [4], Candès and Li suggest using the following empirical loss function to estimate x 0 : min X ∈ C d × d ∑ j = 1 m | a j ⁎ X a j − b j | s.t. X ⪰ 0 . They prove that the solution X ˆ to (5) obeys ‖ X ˆ − x 0 x 0 ⁎ ‖ F ≲ ‖ η ‖ 1 m with high probability provided m ≳ d and a j ∈ C d , j = 1 , … , m , are complex Gaussian random vectors. Although (5) is a convex optimization problem, it needs to be solved in a “lifted” space C d 2 . The computational cost required for this approach typically far exceeds the order of d 3 , which makes it unsuitable for large-dimensional data. However, for the intensity-based estimator, there is no need to lift the problem into higher dimensions. One can simply operate on the original space C d . As shown before, the amplitude-based empirical loss (3) is an alternative estimator for phase retrieval. In [20], Huang and Xu studied the estimation performance of the amplitude-based estimator (3) for real-valued signals. They prove that the solution x ˆ to (3) satisfies min { ‖ x ˆ + x 0 ‖ 2 , ‖ x ˆ − x 0 ‖ 2 } ≲ ‖ η ‖ 2 m with high probability provided m ≳ d and a j ∈ R d , j = 1 , … , m , are Gaussian random vectors. They also prove that the reconstruction error ‖ η ‖ 2 / m is sharp. Furthermore, in [20], Huang and Xu consider the following constrained nonlinear Lasso to estimate s -sparse signals x 0 : min x ∈ R d ∑ j = 1 m ( | 〈 a j , x 〉 | − ψ j ) 2 s.t. ‖ x ‖ 1 ≤ R . They show that any global solution x ˆ to (6) with R : = ‖ x 0 ‖ 1 obeys min { ‖ x ˆ + x 0 ‖ 2 , ‖ x ˆ − x 0 ‖ 2 } ≲ ‖ η ‖ 2 m with high probability provided m ≳ s log ( e d / s ) and a j ∈ R d , j = 1 , … , m , are Gaussian random vectors. The results from [20] hold only for real-valued signals and it seems highly nontrivial to extend them to the complex case. However, the results in this paper hold for both real-valued and complex-valued signals. Furthermore, the intensity-based estimator considered in this paper is differentiable comparing with amplitude-based estimator, which admits high-order algorithms. In [8], Chen and Candès consider the Poisson log-likelihood function ℓ j ( x , b j ) = b j log ( | a j ⊤ x | 2 ) − | a j ⊤ x | 2 and use min x ∈ R d − ∑ j = 1 m ℓ j ( x , b j ) to estimate real-valued signals x 0 . They establish stability estimates using Truncated Wirtinger Flow approach and show that if ‖ η ‖ ∞ ≲ ‖ x 0 ‖ 2 2 then the solution x ˆ to (7) satisfies min { ‖ x ˆ + x 0 ‖ 2 , ‖ x ˆ − x 0 ‖ 2 } ≲ ‖ η ‖ 2 ‖ x 0 ‖ 2 m with high probability provided m ≳ d and a j ∈ R d , j = 1 , … , m , are Gaussian random vectors. Furthermore, a lower bound on the minimax estimation error is also derived under the real-valued Poisson noise model, namely, with high probability inf x ˆ sup x 0 ∈ γ ( K ) E min { ‖ x ˆ + x 0 ‖ 2 , ‖ x ˆ − x 0 ‖ 2 } ≳ d m , where γ ( K ) : = { x 0 ∈ R d : ‖ x 0 ‖ 2 ∈ ( 1 ± 0.1 ) K } for any K ≥ log 1.5 m . The result (8) presents a lower bound on the expectation of the minimax estimation error under Poisson noise structure, whereas the result in Theorem 1.2 gives a tail bound with a reasonably large measurements and holds for a wide range of noises. Moreover, all the results in [8] are for real-valued signals, while ours hold for complex-value ones. In this subsection, we report some numerical experiments to verify that the global solutions to (1) and (2) can be obtained efficiently and the results given in Subsection 1.4 are rate optimal. In our experiments, the target signal x 0 and the measurement vectors a 1 , … , a m are independent standard complex Gaussian random vectors, whereas the noise vector η ∈ R m is a real Gaussian random vector with entries η j ∼ N ( 1 , 1 ) . Example 1.6 In this example, we verify the estimation error presented in Theorem 1.1 is rate optimal. We consider the case where d = 500 and vary m within the range [ 4 d , 50 d ] . The estimator (1) is solved by Algorithm 1, where we use the truncated spectral method proposed in [8] to obtain a good initial guess and then refine it by Wirtinger Flow [5]. Fig. 1 depicts the ratio ρ m against the number of measurements m , when averaged over 100 times independent trials. Here, the ratio ρ m is defined as ρ m : = min θ ∈ [ 0 , 2 π ) ‖ x ˆ − e i θ x 0 ‖ 2 ‖ η ‖ 2 / ( ‖ x 0 ‖ 2 ⋅ m ) . Numerical results show that ρ m tends to be a constant around 0.37, which verifies the estimation error ‖ η ‖ 2 ‖ x 0 ‖ 2 ⋅ m presented in Theorem 1.1 is rate optimal. Example 1.7 The purpose of this numerical experiment is to verify the estimation bound given in Theorem 1.4 is rate optimal when m = O ( s log ( e d / s ) ) . We choose d = 1000 and take the sparsity level s = 100 . The support of x 0 is uniformly distributed at random. The non-zero entries of x 0 are chosen randomly according to a standard normal distribution. We vary m between ⌈ 6 s log ( e d / s ) ⌉ and ⌈ 20 s log ( e d / s ) ⌉ . For each fixed m , we run 100 times trials and calculate the average ratio ρ m defined in (9). The constrained optimization problem (4) is solved by combining the initialization method introduced in [2] and the projection gradient descent onto the ℓ 1 - ball [12]. The result is plotted in Fig. 2. We can see that ρ m tends to be a constant around 0.72, which verifies the estimation error ‖ η ‖ 2 ‖ x 0 ‖ 2 ⋅ m presented in Theorem 1.4 is rate optimal for sparse signals. The paper is organized as follows. In Section 2, after introducing some definitions, we study the recovery of low-rank matrices from rank-one measurements, which plays a key role in the proofs of main results. We also believe that the results in Section 2 are of independent interest. Combining the Mendelson's small-ball method and the results in Section 2, we present the proofs of Theorem 1.1 and Theorem 1.2 in Section 3. The proof of Theorem 1.4 is given in Section 4. A brief discussion is presented in Section 5. Appendix collects the technical lemmas needed in the proof. Section snippets The recovery of low-rank matrices from rank-one measurements For convenience, we let A : H d × d ( C ) → R m be a linear map which is defined as A ( X ) : = ( a 1 ⁎ X a 1 , a 2 ⁎ X a 2 , … , a m ⁎ X a m ) ⊤ , where H d × d ( C ) : = { X ∈ C d × d : X ⁎ = X } . Its dual operator A ⁎ : R m → H d × d ( C ) is given by A ⁎ ( z ) = ∑ j = 1 m z j a j a j ⁎ . In this section, we focus on the following minimization problem: min X ∈ H d × d ( C ) ‖ A ( X ) − b ‖ 2 2 s.t. rank ( X ) ≤ r . Set H r d × d ( C ) : = { X ∈ C d × d : X ⁎ = X , rank ( X ) ≤ r } . A simple observation is that x ˆ is a solution to (1) if and only if X ˆ : = x ˆ x ˆ ⁎ is a solution to (12) with r = 1 . Hence, (12) can be regarded as a lifted version of Proofs of Theorem 1.1 and Theorem 1.2 Motivated by the observation b j = | 〈 a j , x 0 〉 | 2 + η j = a j ⁎ X 0 a j + η j , j = 1 , … , m , where X 0 : = x 0 x 0 ⁎ , we can employ Corollary 2.4 to prove Theorem 1.1. Proof of Theorem 1.1 Let X 0 : = x 0 x 0 ⁎ . Since x ˆ is the global solution to (1), X ˆ : = x ˆ x ˆ ⁎ is the global solution to (12) with r = 1 . From Corollary 2.4, we obtain that with probability at least 1 − exp ( − c m ) , it holds ‖ x ˆ x ˆ ⁎ − x 0 x 0 ⁎ ‖ F ≲ ‖ η ‖ 2 m provided m ≳ d . We claim that for any u , v ∈ C d , we have min θ ∈ [ 0 , 2 π ) ‖ u − e i θ v ‖ 2 ≤ 2 ‖ u u ⁎ − v v ⁎ ‖ F ‖ v ‖ 2 . Indeed, choosing θ : = Phase ( v ⁎ u ) and setting v ¯ : = e i θ v , then 〈 u , v ¯ 〉 ≥ 0 . Proof of Theorem 1.4 In this section, we present the proof of Theorem 1.4. Before proceeding, we introduce several notations and technical lemmas which are useful in the proof. For convenience, we set S d , s : = { x ∈ C d : ‖ x ‖ 2 ≤ 1 , ‖ x ‖ 0 ≤ s } , and K d , s : = { x ∈ C d : ‖ x ‖ 2 ≤ 1 , ‖ x ‖ 1 ≤ s } . The relationship of the two sets is characterized by the following lemma: Lemma 4.1 [30, Lemma 3.1] It holds that conv ( S d , s ) ⊂ K d , s ⊂ 2 conv ( S d , s ) , where conv ( S d , s ) denotes the convex hull of S d , s . The following lemma presents an upper bound for the spectral norm of A ⁎ ( ϵ ) Discussions This paper focuses on the evaluation of the performance of intensity-based estimators for both phase retrieval and its sparse variation. Upper and lower bounds are derived under complex Gaussian random measurements. There are some interesting problems that warrant further research. Firstly, in the presence of noises, many numerical experiments show that gradient descent algorithms can solve estimators (1) and (2). However, it would be of great practical value to provide theoretical guarantees Funding Meng Huang was supported by NSFC grant ( 12201022 ) and the Fundamental Research Funds for the Central Universities (Grant No. YWF-22-T-204 ). Zhiqiang Xu was supported by the National Science Fund for Distinguished Young Scholars ( 12025108 ) and the National Nature Science Foundation of China ( 12021001 , 12288201 ). We thank the anonymous reviewers for their careful reading of the paper and their constructive comments which helped us to substantially improve the paper. References (44) A. Conca et al. An algebraic characterization of injectivity in phase retrieval Appl. Comput. Harmon. Anal. (2015) R. Kueng et al. Low rank matrix recovery from rank one measurements Appl. Comput. Harmon. Anal. (2017) Y. Wang et al. Phase retrieval for sparse signals Appl. Comput. Harmon. Anal. (2014) A. Bourrier et al. Fundamental performance limits for ideal decoders in high-dimensional linear inverse problems IEEE Trans. Inf. Theory (2014) T.T. Cai et al. Optimal rates of convergence for noisy sparse phase retrieval via thresholded Wirtinger flow Ann. Stat. (2016) T.T. Cai et al. ROP: matrix recovery via rank-one projections Ann. Stat. (2015) E.J. Candès et al. Solving quadratic equations via PhaseLift when there are about as many equations as unknowns Found. Comput. Math. (2014) E.J. Candès et al. Phase retrieval via Wirtinger flow: theory and algorithms IEEE Trans. Inf. Theory (2015) E.J. Candès et al. Phaselift: exact and stable signal recovery from magnitude measurements via convex programming Commun. Pure Appl. Math. (2013) Y. Chen et al. Solving random quadratic systems of equations is nearly as easy as solving linear systems Commun. Pure Appl. Math. (2017) Y. Chen et al. Exact and stable covariance estimation from quadratic sampling via convex programming IEEE Trans. Inf. Theory (2015) J.C. Dainty et al. Phase retrieval and image reconstruction for astronomy V. De la Pena et al. Decoupling: from Dependence to Independence (2012) J. Duchi et al. Efficient projections onto the ℓ -1-ball for learning in high dimensions J.R. Fienup Phase retrieval algorithms: a comparison Appl. Opt. (1982) S. Foucart et al. A mathematical introduction to compressive sensing Bull. Am. Math. Soc. (2017) B. Gao et al. Perturbed amplitude flow for phase retrieval IEEE Trans. Signal Process. (2020) B. Gao et al. Phaseless recovery using the Gauss–Newton method IEEE Trans. Signal Process. (2017) R.W. Gerchberg et al. A practical algorithm for the determination of the phase from image and diffraction plane pictures Optik (1972) P. Hand et al. Compressed sensing from phaseless Gaussian measurements via linear programming in the natural parameter space R.W. Harrison Phase problem in crystallography JOSA A (1993) M. Huang et al. The estimation performance of nonlinear least squares for phase retrieval IEEE Trans. Inf. Theory (2020) View more references Cited by (0) Recommended articles (6) Research article An improvement on the upper bounds of the magnitudes of derivatives of rational triangular Bézier surfaces Computer Aided Geometric Design, Volume 31, Issue 5, 2014, pp. 265-276 Show abstract New bounds on the magnitudes of the first- and second-order partial derivatives of rational triangular Bézier surfaces are presented. Theoretical analysis shows that the proposed bounds are tighter than the existing ones. The superiority of the proposed new bounds is also illustrated by numerical tests. Research article Nonlinear matrix recovery using optimization on the Grassmann manifold Applied and Computational Harmonic Analysis, Volume 62, 2023, pp. 498-542 Show abstract We investigate the problem of recovering a partially observed high-rank matrix whose columns obey a nonlinear structure such as a union of subspaces, an algebraic variety or grouped in clusters. The recovery problem is formulated as the rank minimization of a nonlinear feature map applied to the original matrix, which is then further approximated by a constrained non-convex optimization problem involving the Grassmann manifold. We propose two sets of algorithms, one arising from Riemannian optimization and the other as an alternating minimization scheme, both of which include first- and second-order variants. Both sets of algorithms have theoretical guarantees. In particular, for the alternating minimization, we establish global convergence and worst-case complexity bounds. Additionally, using the Kurdyka-Lojasiewicz property, we show that the alternating minimization converges to a unique limit point. We provide extensive numerical results for the recovery of union of subspaces and clustering under entry sampling and dense Gaussian sampling. Our methods are competitive with existing approaches and, in particular, high accuracy is achieved in the recovery using Riemannian second-order methods. Research article The springback penalty for robust signal recovery Applied and Computational Harmonic Analysis, Volume 61, 2022, pp. 319-346 Show abstract We propose a new penalty, the springback penalty, for constructing models to recover an unknown signal from incomplete and inaccurate measurements. Mathematically, the springback penalty is a weakly convex function. It bears various theoretical and computational advantages of both the benchmark convex ℓ 1 penalty and many of its non-convex surrogates that have been well studied in the literature. We establish the exact and stable recovery theory for the recovery model using the springback penalty for both sparse and nearly sparse signals, respectively, and derive an easily implementable difference-of-convex algorithm. In particular, we show its theoretical superiority to some existing models with a sharper recovery bound for some scenarios where the level of measurement noise is large or the amount of measurements is limited. We also demonstrate its numerical robustness regardless of the varying coherence of the sensing matrix. The springback penalty is particularly favorable for the scenario where the incomplete and inaccurate measurements are collected by coherence-hidden or -static sensing hardware due to its theoretical guarantee of recovery with severe measurements, computational tractability, and numerical robustness for ill-conditioned sensing matrices. Research article Dissipative solitons characterization in singly resonant optical parametric oscillators: A variational formalism Optics Communications, 2023, Article 129820 Show abstract In this work, the emergence of single-peak temporal dissipative solitons in singly-resonant degenerate optical parametric oscillators is investigated analytically. Applying the Kantarovitch optimization method, through a Lagrangian variational formalism, an approximate analytical soliton solution is computed using a parameter-dependent ansatz. This permits to obtain analytical estimations for the DS energy, peak power, and boundaries of the DS existence region, which are of great value for experimentalist. To confirm the validity of this procedure, these analytical results are compared with a numerical study performed in the context of pure quadratic systems, showing a good agreement. Research article On sparse recovery algorithms in unions of orthonormal bases Journal of Approximation Theory, Volume 289, 2023, Article 105886 Show abstract We are concerned with the problem of recovering sparse vectors exactly by two commonly used methods, orthogonal matching pursuit (OMP) and basis pursuit (BP). Suppose that overcomplete dictionaries are unions of several orthonormal bases. Gribonval and Nielsen have provided a sufficient mutual coherence condition for BP to recover every sparse vector exactly. Then Tropp proved that it is a sharp sufficient condition for OMP by existence proofs. We establish here the sharpness of Gribonval and Nielsen’s sufficient condition for both OMP and BP by constructing several families of counterexamples. Our main tools are Hadamard matrices and mutually unbiased bases from quantum information theory. Research article Synthesis-based time-scale transforms for non-stationary signals Applied and Computational Harmonic Analysis, Volume 65, 2023, pp. 112-136 Show abstract This paper deals with the modeling of non-stationary signals, from the point of view of signal synthesis. A class of random, non-stationary signals, generated by synthesis from a random time-scale representation, is introduced and studied. Non-stationarity is implemented in the time-scale representation through a prior distribution which models the action of time warping on a stationary signal. A main originality of the approach is that models directly a time-scale representation from which signals can be synthesized, instead of post-processing a pre-computed time-scale transform. A maximum a posteriori estimator is proposed for the time warping parameters and the power spectrum of an underlying stationary signal, together with an iterative algorithm, called JEFAS-S, for the estimation, based upon the Expectation Maximization approach. Numerical results show the ability of JEFAS-S to estimate accurately time warping and power spectrum. This is in particular true when time warping involves fast variations, where a similar approach called JEFAS, proposed earlier, fails. In addition, as a by-product, the approach is able to yield extremely sharp time-scale representations, also in the case of fast varying non-stationarity, where standard approaches such as synchrosqueezing fail. View full text © 2023 Elsevier Inc. All rights reserved. About ScienceDirect Remote access Shopping cart Advertise Contact and support Terms and conditions Privacy policy We use cookies to help provide and enhance our service and tailor content and ads. By continuing you agree to the use of cookies . Copyright © 2023 Elsevier B.V. or its licensors or contributors. ScienceDirect® is a registered trademark of Elsevier B.V. ScienceDirect® is a registered trademark of Elsevier B.V.
|
