Tightness of the maximum likelihood semidefinite relaxation for angular synchronization
Afonso S. Bandeira, Nicolas Boumal, Amit Singer
Introduction
Recovery problems in statistics and many other fields are commonly solved under the paradigm of maximum likelihood estimation, partly due to the rich theory it enjoys. Unfortunately, in many important applications, the parameter space is exponentially large and non-convex, often rendering the computation of the maximum likelihood estimator (MLE) intractable. It is then common to settle for heuristics, such as expectation-maximization algorithms to name but one example. However, it is also common for such iterative heuristics to get trapped in local optima. Furthermore, even when these methods do attain a global optimum, there is in general no way to verify this.
A now classic alternative to these heuristics is the use of convex relaxations. The idea is to maximize the likelihood in a larger, convex set that contains the parameter space of interest, as (well-behaved) convex optimization problems are generally well understood and can be solved in polynomial time. The downside is that the solution obtained might not be in the original feasible (acceptable) set. One is then forced to take an extra, potentially suboptimal, rounding step.
This line of thought is the basis for a wealth of modern approximation algorithms . One preeminent example is Goemans and Williamson’s treatment of Max-Cut , an NP-hard combinatorial problem that involves segmenting a graph in two clusters to maximize the number of edges connecting the two clusters. They first show that Max-Cut can be formulated as a semidefinite program (SDP)—a convex optimization problem where the variable is a positive semidefinite matrix—with the additional, non-convex constraint that the sought matrix be of rank one. Then, they propose to solve this SDP while relaxing (removing) the rank constraint. They show that the obtained solution, despite typically being of rank strictly larger than one, can be rounded to a (suboptimal) rank-one solution, and that it provides a guaranteed approximation of the optimal value of the hard problem. Results of the same nature abound in the recent theoretical computer science literature .
In essence, approximation algorithms insist on solving all instances of a given NP-hard problem in polynomial time, which, unless P = NP, must come at the price of accepting some degree of sub-optimality. This worst-case approach hinges on the fact that a problem is NP-hard as soon as every efficient algorithm for it can be hindered by at least one pathological instance.
Alternatively, in a non-adversarial setting where “the data is not an enemy,” one may find that such pathological cases are not prominent. As a result, the applied mathematics community has been more interested in identifying regimes for which the convex relaxations are tight, that is, admit a solution that is also admissible for the hard problem. When this is the case, no rounding is necessary and a truly optimal solution is found in reasonable time, together with a certificate of optimality. This is sometimes achieved by positing a probability distribution on the instances of the problem and asserting tightness with high probability, as for example in compressed sensing , matrix completion , clustering and inverse problems . In this approach, one surrenders the hope to solve all instances of the hard problem, in exchange for true optimality with high probability.
The main contribution of the present paper is a proof that, even though the angular synchronization problem is NP-hard , its MLE in the face of Gaussian noise can often be computed (and certified) in polynomial time. This remains true even for entry-wise noise levels growing to infinity as the size of the problem (the number of phases) grows to infinity. The MLE is obtained as the solution of a semidefinite relaxation described in Section 2. This arguably striking phenomenon has been observed empirically before (see also Figure 2), but not explained.
Computing the MLE for angular synchronization is equivalent to solving a non-bipartite Grothendieck problem. Semidefinite relaxations for Grothendieck problems have been thoroughly studied in theoretical computer science from the point of view of approximation ratios. The name is inspired by its close relation to an inequality of Grothendieck . We direct readers to the survey by Pisier for a discussion.
The proposed result is qualitatively different from most tightness results available in the literature. Typical results establish either exact recovery of a planted signal (mostly in discrete settings), or exact recovery in the absence of noise, joint with stable (but not necessarily optimal) recovery when noise is present . In contrast, this paper shows optimal recovery even though exact recovery is not possible. In particular, Demanet and Jugnon showed stable recovery for angular synchronization via semidefinite programming, under adversarial noise . We complement this by showing tightness in a non-adversarial setting, meaning the actual MLE is computed.
A similar semidefinite relaxation was studied in a digital communications context, where the parameters to estimate are th roots of unity . There, non-asymptotic results show the relaxation approximates the MLE within some factor, with high probability. The present paper is related to the limit , although the noise model considered is different, and we focus on exact MLE computation.
Our proof relies on verifying that a certain candidate dual certificate is valid with high probability. The main difficulty comes from the fact that the dual certificate depends on the MLE, which does not coincide with the planted signal, and is a nontrivial function of the noise. We use necessary optimality conditions of the hard problem to both obtain an explicit expression for the candidate dual certificate, and to partly characterize the point whose optimality we aim to establish. This seems to be required since the MLE is not known in closed form.
In the context of sparse recovery, a result with similar flavor is support recovery guarantee , where the support of the estimated signal is shown to be contained in the support of the original signal. Due to the noise, exact recovery is also impossible in this setting. Another example is a recovery guarantee in the context of latent variable selection in graphical models .
Besides the relevance of angular synchronization in and of its own, we are confident this new insight will help uncover similar results in other applications where it has been observed that semidefinite relaxations can be tight even when the ground truth cannot be recovered. Notably, this appears to be the case for the Procrustes and multi-reference alignment problems ).
The crux of our argument concerns the rank of the solutions of an SDP. We mention in passing that there are many other deterministic results in the literature pertaining to the rank of solutions of SDP’s. For example, it has been shown that, in general, an SDP with only equality constraints admits a solution of rank at most (on the order of) the square root of the number of constraints, see . Furthermore, Sagnol shows that under some conditions (that are not fulfilled in our case), certain SDP’s related to packing problems always admit a rank-one solution. Sojoudi and Lavaei study a class of SDP’s on graphs which is related to ours and for which, under certain strong conditions on the topology of the graphs, the SDP’s admit rank-one solutions—see also applications to power flow optimization .
2 Notation
The Angular Synchronization problem
Further letting the noise be i.i.d. (complex) Gaussian variables for , it follows that an MLE for is any vector of phases minimizing . Equivalently, an MLE is a solution of the following quadratically constrained quadratic program (sometimes called the complex constant-modulus QP [41, Table 2] in the optimization literature, and non-bipartite Grothendieck problem in theoretical computer science):
where denotes the conjugate transpose of . This problem can only be solved up to a global phase, since only relative information is available. Indeed, given any solution , all vectors of the form are equivalent solutions, for arbitrary phase .
Such relaxations lift the problem to higher dimensional spaces. Indeed, the search space of (QP) has dimension (or , discounting the global phase) whereas the search space of (SDP) has real dimension . In general, increasing the dimension of an optimization problem may not be advisable. But in this case, the relaxed problem is a semidefinite program. Such optimization problems can be solved to global optimality up to arbitrary precision in polynomial time .
It is known that the solution of (SDP) can be rounded to an approximate solution of (QP), with a guaranteed approximation ratio [54, § 4]. But even better, when (SDP) admits an optimal solution of rank one, then no rounding is necessary: the leading eigenvector of is a global optimum of (QP), meaning we have solved the original problem exactly. Elucidating when the semidefinite program admits a solution of rank one, i.e., when the relaxation is tight, is the focus of the present paper.
Problem (QP) is posed over the complex numbers. As a result, the individual variables in (QP) live on a continuous search space (the unit circle). One effect of this is that even small noise on the data precludes exact recovery of the signal (in general). This is the root of most of the complications that will arise in the developments hereafter. In order to first illustrate some of the ideas of the proof in a simpler context, this section proposes to take a detour through the real case. Besides this expository rationale, the real case is interesting in and of itself. It notably relates to correlation clustering and the stochastic block model . The analysis proposed here also appears in .
Let be the signal to estimate and let contain the measurements, with a Wigner matrix: its above-diagonal entries are i.i.d. (real) standard normal random variables, and its diagonal entries are zero. Each entry is a noisy measurement of the relative sign . For example, could represent political preference of an agent (left or right wing) and could be a measurement of agreement between two agents’ views . Consider this MLE problem:
Thus, for each . The corresponding relaxation reads:
Strong duality holds, which implies that a given feasible is optimal if and only if there exists a dual feasible matrix such that , or, in other words, such that .Indeed, is diagonal and , hence . Since both and are positive semidefinite, this is equivalent to requiring (a condition known as complementary slackness), and hence requiring .
For ease of exposition, we now assume (without loss of generality) that . Tentatively, let —for the complex case, we will see how to obtain this candidate without guessing. Then, by construction, is diagonal and .Using , we get . It remains to determine under what conditions is positive semidefinite.
Define the diagonal matrix and the (Laplacian-like) matrix . The candidate dual certificate is a sum of two Laplacian-like matrices:
This is indeed compatible with the empirical observation of Figure 1.
2 Back to synchronization over SO(2)
We now return to the complex case, which is the focus of this paper. As was mentioned earlier, in the presence of even the slightest noise, one can no longer reasonably expect the true signal to be an optimal solution of (QP) (this can be further quantified using Cramér-Rao bounds ). Nevertheless, we set out to show that (under some assumptions on the noise) solutions of (QP) are close to and they can be computed via (SDP).
The proof follows that of the real case in spirit, but requires more sophisticated arguments because the solution is no longer known explicitly. This is important because the candidate dual certificate itself depends on that solution. With this in mind, the proof of the upcoming main lemma (Lemma 3.2) follows this reasoning:
For small enough noise levels , any optimal solution of (QP) is close to the sought signal (Lemmas 4.1 and 4.2).
Solutions , a fortiori, satisfy necessary optimality conditions for (QP). First-order conditions take up the form , where depends smoothly on (see (4.9)). Second-order conditions will also be used.
Remarkably, this can be used as a dual certificate for solutions of (SDP). Indeed, is optimal if and only if is positive semidefinite (Lemma 4.4). The solution is unique if (Lemma 4.3). Thus, it only remains to study the eigenvalues of .
In the absence of noise, is a Laplacian for a complete graph with unit weights (up to a unitary transformation), so that its eigenvalues are with multiplicity , and with multiplicity . Then, is always the unique solution.
Adding small noise, because of the first point, the solution will move only by a small amount, and hence so will . Thus, the large eigenvalues should be controllable into remaining positive (Section 4.4).
The crucial fact follows: because of the way is constructed (using first-order optimality conditions), the zero eigenvalue is “pinned down” (as long as is a local optimum of (QP)). Indeed, both and change as a result of adding noise, but the property remains valid. Thus, there is no risk that the zero eigenvalue from the noiseless scenario would become negative when noise is added.
Following this road map, most of the work in the proof below consists in bounding how far away can be from , and in using that to control the large eigenvalues of . This constructive way of identifying the dual certificate (third point in the roadmap) already appears explicitly in work by Journée et al. , who considered a different family of real, semidefinite programs which also admit a smooth geometry when the rank is constrained. This points to smoothness of (QP) (and non-degeneracy of (SDP) ) as a principal ingredient in our analysis: the KKT conditions of (QP) are a subset of the KKT conditions of its relaxation. It is because the former are explicit (rather than existential as in Lemma 4.3) that they help in identifying .
Our main theorem follows. In a nutshell, it guarantees that: under (complex) Wigner noise , with high probability, solutions of (QP) are close to , and, assuming the noise level is smaller than (on the order of) , (SDP) admits a unique solution, it is of rank one and identifies the solution of (QP) (unique, up to a global phase shift).
then the semidefinite program (SDP), given by
has, as its unique solution, the rank-one matrix .
The numerical experiments (Figure 2) suggest it should be possible to allow to grow at a rate of (as in the real case), but we were not able to establish that (see Remark 4.6). Nevertheless, we do show that can grow unbounded with . To the best of our knowledge, this is the first result of this kind. We hope it might inspire similar results in other problems where the same phenomenon has been observed .
Main result
In this section we present our main technical result and show how it can be used to prove Theorem 2.1. We begin with a central definition in this paper. Intuitively, this definition characterizes non-adversarial noise matrices .A similar but different definition appeared in a previous version of this paper.
The next lemma is the main technical contribution of this paper. It is a deterministic, non-asymptotic statement.
Furthermore, if , then the semidefinite program (SDP) has, as its unique solution, the rank-one matrix .
We defer the proof of Lemma 3.2 to Section 4. The following proposition, whose proof we defer to Appendix A, shows how this lemma can be used to prove Theorem 2.1.
The latter result is not surprising. Indeed, the definition of -discordance requires two elements. Namely, (i) that be not too large as an operator, and (ii) that no row of be too aligned with . For a Wigner matrix independent of , those are indeed expected to hold. In fact, for Gaussian noise, the constants 3 can be replaced by for any , provided is large enough.
The definition of -discordance is not tightly adjusted to Wigner noise. As a result, it is expected that Lemma 3.2 will be applicable to show tightness of semidefinite relaxations for a larger span of noise models.
The proof
In this section, we prove Lemma 3.2. See Section 2.2 for an outline of the proof.
Without loss of generality, assume the global phase of is such that , i.e., and are optimally aligned. Expand the inequality using the data model to obtain
Since , divide both sides by to obtain
Combine this and (4.1) with (4.2) to obtain . ∎
Note that the constant 12 is pessimistic, because we used even though it is closer to . Assuming and , the argument above can be bootstrapped to reduce the constant. One obtains for all , with and . For small enough, converges arbitrarily close to 6. For example, if , then .
The next lemma establishes a bound on the largest individual error, , after proper global phase alignment of and . Interestingly, for , the bound shows that individual errors decay, uniformly, as increases.
and combine with the assumption to obtain
Invoke Lemma 4.1, namely, , to get
Since , it finally comes that
We discuss the problem of bounding in more details later on. For now, we simply use the suboptimal bound (4.11) (obtained independently of the present lemma), i.e., . Then, for all , using ,
(For , the constant could be replaced by one arbitrarily close to .) ∎
2 Optimality conditions for (SDP)
The global optimizers of the semidefinite program (SDP) can be characterized completely via the Karush-Kuhn-Tucker (KKT) conditions.
If, furthermore, , then has rank one and is the unique global optimizer of (SDP).
Certainly, if (SDP) admits a rank-one solution, it has to be of the form , with an optimal solution of the original problem (QP). Based on this consideration, our proof of Lemma 3.2 goes as follows. We let denote a global optimizer of (QP) and we consider as a candidate solution for (SDP). Using the optimality of and assumptions on the noise, we then construct and verify a dual certificate matrix as required per Lemma 4.3. In such proofs, one of the nontrivial parts is to guess an analytical form for given a candidate solution . We achieve this by inspecting the first-order optimality conditions of (QP) (which necessarily satisfies). The main difficulty is then to show the suitability of the candidate , as it depends nonlinearly on the global optimum , which itself is a complicated function of the noise . We show feasibility of via a program of inequalities, relying on -discordance of the noise (see Definition 3.1).
3 Construction of the dual certificate S𝑆S
Every global optimizer of the combinatorial problem (QP) must, a fortiori, satisfy first-order necessary optimality conditions. We derive those now.
over the smooth Riemannian manifold . Therefore, the first-order necessary optimality conditions for (QP) (i.e., the KKT conditions) can be stated simply as , where is the Riemannian gradient of at . This gradient is given by the orthogonal projection of the Euclidean (the classical) gradient of onto the tangent space of at [3, eq. (3.37)],
The projector and the Euclidean gradient are given respectively by:
Note that is Hermitian and . Referring to the KKT conditions in Lemma 4.3, it follows immediately that is feasible for (SDP) (conditions 1 and 2); that (condition 3); and that is a diagonal matrix (condition 4). It thus only remains to show that is also positive semidefinite and has rank . If such is the case, then is the unique global optimizer of (SDP). Note the special role of the first-order necessary optimality conditions: they guarantee complementary slackness, without requiring further work.
We will also use second-order necessary optimality conditions, namely, that the (Riemannian) Hessian of at an optimizer is positive semidefinite on the tangent space (4.7). The action of the Riemannian Hessian is obtained by projecting the directional derivatives of the Riemannian gradient [3, eq. (5.15)]. Explicitly, for any tangent vector , one can compute that (with denoting the directional derivative of at along )
Thus, if is a (local) optimizer, then, since is self-adjoint,
for all : a necessary but insufficient step towards making positive semidefinite. We will later use this condition along selected directions.
The following lemma further shows that is the right candidate dual certificate. More precisely, for a critical point of (QP), it is necessary and sufficient for to be positive semidefinite in order for to be optimal for (SDP). This confirms that, by studying , nothing is lost with respect to the original question. See also .
A feasible (of any rank) for (SDP) is optimal if and only if (4.9) is positive semidefinite. There exists no other certificate.
The if part follows from Lemma 4.3: set and observe that, since by construction, and since , it follows that . We now show the only if part. Assume is optimal. Then, by Lemma 4.3, there exists which satisfies and , where is diagonal. Thus, and . Consequently, . ∎
4 A sufficient condition for rank recovery
Knowing which certificate (4.9) to verify for a given , it remains to characterize the point . Of course, is the global optimizer of (QP), but this is not a convenient property to exploit. Instead, we let be a second-order critical pointA second-order critical point satisfies first- and second-order necessary optimality conditions, namely, the gradient is zero and the Hessian is positive semidefinite (if minimizing) [26, §3.2.1]. for (QP) which outperforms the planted signal (certainly, the MLE is one such point). In effect, we prove the following, which immediately yields Lemma 3.2.
To this end, we prove that, at such points , the certificate is positive semidefinite of rank . As a first observation, note that being critical () implies that for all . Thus, is real, and it follows that
Furthermore, since is second-order critical, the Hessian of at is positive semidefinite, implying for all (4.7). In particular, for all , is a tangent vector at (where and is the th canonical basis vector), and it follows that
Since by construction,As before, is independent of . Assuming involves no loss of generality. it follows in particular that is real and positive. Thus:
As a result, a sufficient condition for to be positive semidefinite with rank is to have
Combining with (4.10) yields the following sufficient condition:
This condition involves only and , and as such is the best result we prove. The bottleneck is the term in , which leads to a bound of the type . Further inspection of (4.12) shows is indeed a sufficient condition, for . This concludes the proof of Lemma 3.2.
One avenue to improve the current bound might be to expand around 0, to first order. Some computations (omitted here) show that such a development would lead to . Then,
Conclusions
Appendix A Wigner matrices are discordant
A matrix is -discordant if and only if is -discordant. Since has the same distribution as (owing to complex normal random variables having uniformly random phase), we may without loss of generality assume in the remainder of the proof.
Although tail bounds for the real version of this are well-known (see for example ) and they mostly hold verbatim in the complex case, for the sake of completeness we include a classical argument, based on Slepian’s comparison theorem and Gaussian concentration, for a tail bound in the complex valued case.
This means that we can use Slepian’s comparison theorem (see for example [39, Cor. 3.12]) to get
.
The random vector given by is jointly Gaussian where the marginal of each entry is a standard complex Gaussian. By a suboptimal union bound argument, the maximum absolute value among standard complex Gaussian random variables (not necessarily independent) is larger than with probability at most . Hence,
To support the discussion following Proposition 3.3, we further argue that
It is easy to see that is a real Gaussian random variable with zero mean and variance . This implies that:
Acknowledgments
A. S. Bandeira was supported by AFOSR Grant No. FA9550-12-1-0317. Most of this work was done while he was with the Program for Applied and Computational Mathematics at Princeton University. N. Boumal was supported by a Belgian F.R.S.-FNRS fellowship while working at the Université catholique de Louvain (Belgium), by a Research in Paris fellowship at Inria and ENS, the “Fonds Spéciaux de Recherche” (FSR UCLouvain), the Chaire Havas “Chaire Economie et gestion des nouvelles données” and the ERC Starting Grant SIPA. A. Singer was partially supported by Award Number R01GM090200 from the NIGMS, by Award Numbers FA9550-12-1-0317 and FA9550-13-1-0076 from AFOSR, by Award Number LTR DTD 06-05-2012 from the Simons Foundation, and by the Moore Foundation.