Linear algebra

Why column-subset selection and curvature

The column-subset selection problem (CSSP) asks for $k$ columns of a real matrix $A$ whose span approximates all of $A$ as well as possible in squared Frobenius norm. It is a basic primitive of low-rank approximation, feature selection and numerical linear algebra, studied both for its exact complexity and for approximation algorithms (see Boutsidis, Drineas and Magdon-Ismail, doi:10.1137/12086755X, cited as [2] in the source paper).

Sviridenko, Vondrák and Ward study submodular maximization and supermodular minimization when the objective has bounded total curvature, a parameter $c \in [0,1]$ introduced by Conforti and Cornuéjols that measures how far a monotone set function is from being linear. Their guarantees improve as $c$ decreases. Section 8 of the paper applies this to CSSP: it bounds the curvature of the CSSP objective in terms of the condition number of $A$, so that a well-conditioned matrix has small curvature. This mission formalizes that bound (Lemma 8.2) and the steps of its proof.

Setting

Let $A$ be a real $m \times n$ matrix with columns $c_1, \dots, c_n \in \mathbb{R}^m$, so that $Ax = \sum_i x_i c_i$ for $x \in \mathbb{R}^n$, and let $\|\cdot\|$ be the Euclidean norm. For $S \subseteq [n]$ write $\mathrm{proj}_S(x)$ for the orthogonal projection of $x$ onto $\mathrm{span}\{c_i : i \in S\}$, the point of that subspace nearest to $x$ ($\mathrm{proj}_\emptyset = 0$). The CSSP objective is

$$f^A(S) = \sum_{i=1}^n \|c_i - \mathrm{proj}_S(c_i)\|^2 ,$$

the squared Frobenius distance from $A$ to the nearest matrix whose columns lie in the span of the selected columns.

For a set function $f$ on a finite set $X$, the marginal value of $i$ at $A$ is $f_A(i) = f(A+i) - f(A)$. The function $f$ is monotone decreasing if every $f_A(i) \le 0$. Following (1.1) of the paper, a monotone decreasing $f$ has total curvature at most $c$ when, for every $j \in X$,

$$f_{X-j}(j) \le (1-c)\, f_\emptyset(j),$$

i.e. the last marginal of $j$ is at least a $(1-c)$ fraction of its first marginal in absolute value.

When the columns of $A$ are linearly independent, the condition number is

$$\kappa(A) = \frac{\sup_{\|x\|=1}\|Ax\|}{\inf_{\|x\|=1}\|Ax\|} \in [1, \infty).$$

Formalization targets

Goal: Lemma 8.2 (p. 11)

For $A$ with linearly independent columns,

$$f^A_{[n]-j}(j) \le \frac{1}{\kappa^2(A)}\, f^A_\emptyset(j) \quad \text{for all } j \in [n], \qquad \text{i.e.}\quad c(f^A) \le 1 - \frac{1}{\kappa^2(A)} .$$

Milestones

1. The identity $f^A(S) = \sum_i (\|c_i\|^2 - \|\mathrm{proj}_S(c_i)\|^2)$ (p. 11).
2. Lemma 8.1, monotone part: $S \mapsto \|\mathrm{proj}_S(v)\|^2$ is monotone increasing (p. 11).
3. $f^A$ is monotone decreasing (p. 11).
4. $|f^A_\emptyset(i)| = \sum_j \|\mathrm{proj}_{\{i\}}(c_j)\|^2$ (p. 11).
5. Some unit $x$ has $\|Ax\| \ge \sqrt{|f^A_\emptyset(i)|}$ (p. 12).
6. $|f^A_{[n]-i}(i)| = \|c_i - \mathrm{proj}_{[n]-i}(c_i)\|^2$ (p. 12).
7. Some unit $x'$ has $\|Ax'\| \le \sqrt{|f^A_{[n]-i}(i)|}$ (p. 12).
8. $\kappa^2(A) \ge f^A_\emptyset(i)/f^A_{[n]-i}(i)$ for each $i$, under independence (p. 12).

Items 1–7 hold for arbitrary columns; only item 8 and the goal assume independence.

Significance

Lemma 8.2 turns a spectral quantity into a combinatorial one. With it, any curvature-dependent guarantee for monotone decreasing set functions can be stated for CSSP in terms of $\kappa(A)$. The paper's Lemma 8.3 shows the bound is asymptotically tight, so $1 - \Theta(1/\kappa^2)$ is the right order of the curvature.

There is a caveat that the formalization makes explicit. The paper deduces from Lemma 8.1 that $f^A$ is supermodular, and on that basis quotes its supermodular-minimization algorithm for CSSP with factor $1 + (\kappa^2 - 1)e^{-1} + \kappa^2 O(\epsilon)$ (p. 12). Both the submodularity half of Lemma 8.1 and the supermodularity of $f^A$ are false as printed. With $c_1 = e_1$, $c_2 = e_1 + e_2$ and $v = e_2$ in $\mathbb{R}^2$, the function $F(S) = \|\mathrm{proj}_S(v)\|^2$ has $F(\{1,2\}) + F(\emptyset) = 1 > \tfrac12 = F(\{1\}) + F(\{2\})$. With columns $e_1, e_1+e_2, 2e_2$, $f^A(\{1\}) + f^A(\{2\}) = \tfrac{15}{2} > 7 = f^A(\{1,2\}) + f^A(\emptyset)$. A small perturbation gives a counterexample with independent columns. This is the "suppressor variable" phenomenon of regression. The curvature bound itself does not depend on supermodularity and is true, so it is the goal here. The approximation factor for CSSP is not posed, because it rests on the false claim and is stated with $O(\epsilon)$.

As far as is known, none of these statements has a machine-checked proof. The results are proved on paper, except the two false claims above. The mission adds a formal proof of the curvature–conditioning bound, and makes explicit which parts of §8 survive.

Difficulty

The individual steps are short linear algebra, but they rest on infrastructure that has to be assembled: orthogonal projections onto spans of finite column families, the Pythagorean identities that relate $f^A$ to projections, and the variational characterization of $\sup$ and $\inf$ of $\|Ax\|$ over the unit sphere.

The paper's own route to the large unit vector (milestone 5) does not work as printed. It sets $x_j \propto \|\mathrm{proj}_{\{i\}}(c_j)\|$ and uses $c_i \cdot c_j = \|c_i\|\,\|\mathrm{proj}_{\{i\}}(c_j)\|$, which is false when $c_i \cdot c_j < 0$, so the signs of the inner products have to be accounted for. In the final step, the ratio $f^A_\emptyset(i)/f^A_{[n]-i}(i)$ is a quotient of two negative numbers. Comparing it with $\kappa^2$ requires the supremum to be finite and the infimum to be strictly positive, which is exactly what independence of the columns provides.

Formalization scope

• Columns are a family c : Fin n → EuclideanSpace ℝ (Fin m), and $Ax$ is ∑ j, x j • c j. Subsets of $[n]$ are Finset (Fin n). $[n] - i$ is Finset.univ.erase i and $S + i$ is insert i S.
• $\mathrm{proj}_S$ is Mathlib's Submodule.starProjection onto Submodule.span ℝ (c '' S), which is finite-dimensional, so the argmin of the page exists and is unique.
• $\kappa(A)$ is a real quotient of ⨆ and ⨅ of $\|Ax\|$ over the unit-sphere subtype of $\mathbb{R}^n$. It is used only under LinearIndependent ℝ c, the page's "non-singular", where both are finite and positive. The page's $\kappa = \infty$ for dependent columns is not encoded.
• Monotonicity is the marginal form of §2.1. Curvature is the product form above, which agrees with the ratio of (1.1) whenever $f_\emptyset(j) \ne 0$ and avoids Lean's $x/0 = 0$.
• For $n = 0$ the goal has no $j$ and is empty, as on the page.
• A trivializing reading is ruled out by the encoding: had $\kappa$ evaluated to a junk value $0$, the goal would reduce to monotonicity. Under independence the supremum and infimum are genuine, so the bound $1 - 1/\kappa^2$ is the paper's.
• Not in scope: Lemma 8.3 (stated with $O(\cdot)$), the CSSP approximation factor, and the false submodularity and supermodularity claims.

Useful infrastructure that would be reusable beyond this mission includes Pythagoras for starProjection, monotonicity of projection norms under inclusion of subspaces, and extremal characterizations of $\sup_{\|x\|=1}\|Ax\|$ and $\inf_{\|x\|=1}\|Ax\|$. Contributions of any of the milestones, or of these general lemmas, are welcome.

Selected references

• M. Sviridenko, J. Vondrák, J. Ward, Optimal approximation for submodular and supermodular optimization with bounded curvature, SODA 2015 (version of Oct. 9, 2014); journal version Math. Oper. Res. 42(4), 2017. https://doi.org/10.1287/moor.2016.0842
• M. Conforti, G. Cornuéjols, Submodular set functions, matroids and the greedy algorithm: tight worst-case bounds and some generalizations of the Rado–Edmonds theorem, Discrete Appl. Math. 7(3), 1984. https://doi.org/10.1016/0166-218X(84)90003-9
• C. Boutsidis, P. Drineas, M. Magdon-Ismail, Near-optimal column-based matrix reconstruction, SIAM J. Comput. 43(2), 2014. https://doi.org/10.1137/12086755X
• A. Das, D. Kempe, Submodular meets spectral: greedy algorithms for subset selection, sparse approximation and dictionary selection, ICML 2011. https://arxiv.org/abs/1102.3975

Motivation

Tensor decompositions are a standard tool in signal processing, chemometrics, telecommunications and blind source separation. The canonical polyadic decomposition (CPD) writes a tensor as a sum of rank-one tensors; the rank-$(L_r,L_r,1)$ block term decomposition (BTD) replaces each rank-one term by the outer product of a rank-$L_r$ matrix and a vector. Sorber, Van Barel and De Lathauwer (SIAM J. Optim. 23(2), 2013) introduced a common generalization, the (rank-$L_r\circ$ rank-1) BTD, and derived optimization algorithms that compute all three decompositions of complex tensors. Their quasi-Newton and nonlinear conjugate gradient methods are built on complex Taylor models in $z$ and $\bar z$ (Sorber, Van Barel and De Lathauwer, SIAM J. Optim. 22(3), 2012), and each iteration needs the complex cogradient of the least-squares objective. This mission formalizes the closed form of that cogradient, Theorem 4.4 of the paper. The algorithms were released in the MATLAB toolbox Tensorlab.

Setting

Let $\mathcal T\in\mathbb C^{I_1\times\cdots\times I_N}$ be a tensor of order $N=P+Q$. Fix block sizes $L_1,\dots,L_R$ and put $R'=\sum_r L_r$. The unknowns are complex matrices $A^{(p)}\in\mathbb C^{I_p\times R'}$ for $1\le p\le P$ and $C^{(q)}\in\mathbb C^{I_{P+q}\times R}$ for $1\le q\le Q$, collected in the vector

$$z=(\operatorname{vec}A^{(1)},\dots,\operatorname{vec}A^{(P)},\operatorname{vec}C^{(1)},\dots,\operatorname{vec}C^{(Q)}).$$

Let $E=\operatorname{diag}(1_{1\times L_1},\dots,1_{1\times L_R})$ be the $R\times R'$ block-diagonal matrix of row vectors of ones. The (rank-$L_r\circ$ rank-1) BTD is read as a structured CPD with $N$ factor matrices

$$A^{(1)},\dots,A^{(P)},\qquad A^{(P+q)}=C^{(q)}E\quad(1\le q\le Q).$$

Its residual tensor is $\mathcal F_{\mathrm{BTD}}=\sum_{r'=1}^{R'}a^{(1)}_{r'}\circ\cdots\circ a^{(N)}_{r'}-\mathcal T$, and the objective is the nonlinear least-squares function

$$f_{\mathrm{BTD}}=\tfrac12\|\mathcal F_{\mathrm{BTD}}\|^2 ,$$

with the Frobenius norm induced by $\langle\mathcal T,\mathcal U\rangle=\sum\overline{t_{i_1\cdots i_N}}\,u_{i_1\cdots i_N}$.

For a set $\sigma$ of modes, $V^\sigma=\bigodot_{n\notin\sigma}A^{(n)}$ is the Khatri–Rao product (column-wise Kronecker product) of the factor matrices outside $\sigma$, and $W^\sigma=\ast_{n\notin\sigma}A^{(n)\mathrm H}A^{(n)}$ is the Hadamard product of their Gramians. The mode-$n$ unfolding $T_{(n)}$ arranges the mode-$n$ fibres of $\mathcal T$ as columns, ordered like the rows of $V^{\{n\}}$.

The objective is real-valued, hence not holomorphic in $z$. The cogradient $\partial f/\partial z$ is the vector of Wirtinger derivatives

$$\frac{\partial f}{\partial z_\beta}=\frac12\Big(\frac{\partial f}{\partial x_\beta}-\mathrm i\frac{\partial f}{\partial y_\beta}\Big),\qquad z_\beta=x_\beta+\mathrm i y_\beta,$$

the derivative with respect to $z_\beta$ with $\bar z_\beta$ held constant. Overbars denote entrywise complex conjugation.

Formalization targets

Goal: Theorem 4.4

For every tensor $\mathcal T$ and every $z$, for $1\le p\le P$ and $1\le q\le Q$,

$$2\,\frac{\partial f_{\mathrm{BTD}}}{\partial A^{(p)}}=\overline{A}^{(p)}\,W^{\{p\}}-\overline{T}_{(p)}\,V^{\{p\}},$$ $$2\,\frac{\partial f_{\mathrm{BTD}}}{\partial C^{(q)}}=\overline{C}^{(q)}\,E\,W^{\{P+q\}}\,E^{\mathrm T}-\overline{T}_{(P+q)}\,V^{\{P+q\}}\,E^{\mathrm T}.$$

The identity (4.4), which assembles $\partial f/\partial z$ from these blocks, holds by construction in the formalization.

Milestones

1. Proposition 4.2: for any permutation $p$, $(A^{(p_1)}\odot\cdots\odot A^{(p_N)})^{\mathrm H}(A^{(p_1)}\odot\cdots\odot A^{(p_N)})=\ast_n A^{(n)\mathrm H}A^{(n)}$.
2. Corollary 4.3: $(V^\sigma)^{\mathrm H}V^\sigma=W^\sigma$.
3. The splitting $f=\frac12(f^{(1)}-f^{(2)}-\overline{f^{(2)}}+f^{(3)})$ with $f^{(1)}=\|A^{(n)}V^{\{n\}\mathrm T}\|^2$, $f^{(2)}=\langle T_{(n)},A^{(n)}V^{\{n\}\mathrm T}\rangle$, $f^{(3)}=\|T_{(n)}\|^2$.
4. $\partial f^{(1)}/\partial A^{(n)}=\overline{A}^{(n)}W^{\{n\}}$.
5. $\partial f^{(2)}/\partial A^{(n)}=\overline{T}_{(n)}V^{\{n\}}$ and $\partial\overline{f^{(2)}}/\partial A^{(n)}=\partial f^{(3)}/\partial A^{(n)}=0$.
6. Proposition 4.1: $\operatorname{vec}(C^{(q)}E)^{\mathrm T}=\operatorname{vec}(C^{(q)})^{\mathrm T}(E\otimes\mathbb I)$ and the chain rule $\partial/\partial C^{(q)}=\partial/\partial A^{(P+q)}\cdot E^{\mathrm T}$.

Significance

Theorem 4.4 gives the gradient information every first-order and quasi-Newton method for the CPD, the rank-$(L_r,L_r,1)$ BTD and their generalization needs. It replaces derivatives through the large Khatri–Rao product $V^{\{n\}}$ by small $R'\times R'$ Hadamard products of Gramians, and leaves one tensor–Khatri–Rao product per mode as the dominant cost. Table 4.1 of the paper specializes it to third-order CPDs and BTDs, and the same structure underlies the alternating least squares updates of §4.1.

The paper proves the result. As far as a search of the Prove2Me catalogue and Mathlib shows, it has not been machine-checked. Mathlib has the Kronecker product and Matrix.vec, but no Khatri–Rao product, tensor unfoldings or Wirtinger calculus on matrix-valued unknowns. A formal proof would supply a checked derivation of the gradient formula used in tensor software, together with reusable lemmas about Khatri–Rao products and Wirtinger derivatives of real quadratic forms.

Difficulty

Each step is elementary, but the bookkeeping is easy to get wrong. The cogradient of a real function is not a complex derivative. Conjugates appear exactly where the conventions put them: with $\partial/\partial\bar z$ instead of $\partial/\partial z$, or with the inner product conjugated on the other argument, every overbar in (4.5) moves. The proof also differentiates with respect to $A^{(P+q)}$, which is not an unknown. That step is justified only by treating the structured CPD as an unstructured one and pulling back through $A^{(P+q)}=C^{(q)}E$ with a chain rule for Wirtinger derivatives, which needs real differentiability of the outer function. Finally, the Khatri–Rao identity $V^{\mathrm H}V=W$ requires reindexing a product over a permuted list of modes into a product over all modes.

Formalization scope

Modes are Fin P ⊕ Fin Q. The paper's mode $p$ is Sum.inl (p-1) and its mode $P+q$ is Sum.inr (q-1). The $R'$ columns are Σ r : Fin R, Fin (L r) in lexicographic order, which is the order of (3.5), and all indices are 0-based. vec is column-major (Mathlib's Matrix.vec). The rows of $V^\sigma$ and the columns of $T_{(n)}$ are indexed by tuples of the remaining indices rather than by integers; both use the same order, so the products in (4.5) are unchanged. The residual is defined entrywise and $f_{\mathrm{BTD}}$ is a real number cast to $\mathbb C$. The cogradient is defined coordinatewise from Mathlib's real deriv along the directions $e_\beta$ and $\mathrm i e_\beta$. The functions differentiated in the goal are real polynomials in the real and imaginary parts, so deriv is a genuine derivative there; a junk value of $0$ would make (4.5) false, not vacuous. Milestones 3–5 are stated for the unstructured CPD, in which every mode, including $n>P$, has a free factor matrix. Proposition 4.1's chain rule assumes the outer function is real-differentiable at the point. That is the only hypothesis in the mission.

The paper's qualifiers "$a_r^{(n)}$ nonzero" (3.1) and "rank-$L_r$" (3.5) describe what a decomposition is, and §4 never uses them. They are not imposed: the statements hold for arbitrary complex factor matrices and tensors, and $P$, $Q$, $R$, $L_r$, $I_n$ may be $0$. The goal must not be replaced by its (4.4) part alone, which holds by construction, nor by a version with $\partial/\partial\bar z$ or with the Hermitian transpose $A^{\mathrm H}$ in place of the entrywise conjugate $\overline A$.

A complete development needs Khatri–Rao Gramian identities, derivatives of real-valued polynomial maps along real lines in $\mathbb C^n$, and the Wirtinger chain rule for a $\mathbb C$-linear inner map. All of these can be reused outside this mission. Contributions are welcome on any milestone, and on general Wirtinger-calculus lemmas (linearity, products, conjugation, chain rule).

Selected references

• L. Sorber, M. Van Barel, L. De Lathauwer, Optimization-Based Algorithms for Tensor Decompositions: Canonical Polyadic Decomposition, Decomposition in Rank-$(L_r,L_r,1)$ Terms, and a New Generalization, SIAM J. Optim. 23(2), 695–720, 2013. https://doi.org/10.1137/120868323
• L. Sorber, M. Van Barel, L. De Lathauwer, Unconstrained Optimization of Real Functions in Complex Variables, SIAM J. Optim. 22(3), 879–898, 2012. https://doi.org/10.1137/110832124
• T. G. Kolda, B. W. Bader, Tensor Decompositions and Applications, SIAM Review 51(3), 455–500, 2009. https://doi.org/10.1137/07070111X
• L. De Lathauwer, Decompositions of a Higher-Order Tensor in Block Terms — Part II: Definitions and Uniqueness, SIAM J. Matrix Anal. Appl. 30(3), 1033–1066, 2008. https://doi.org/10.1137/070690729

Motivation

For a $k$-regular graph $G=(V,E)$ on $n$ vertices, the Expander Mixing Lemma (Alon–Chung 1988; see the survey of Hoory, Linial and Wigderson 2006) says that for any vertex sets $A,B$, the number of edges between them is close to its expected value: $\bigl||E(A,B)|-k|A||B|/n\bigr|\le\rho\sqrt{|A||B|}$, where $\rho$ is the largest absolute value of a nontrivial eigenvalue of $kI-\Delta^+$. The deviation, called the discrepancy of $A$ and $B$, measures how random-like the graph is. It is one of the basic tools of spectral graph theory and of the theory of expanders.

Higher-dimensional expansion asks the same questions for simplicial complexes, where the role of edges is played by $d$-dimensional cells. Parzanchevski, Rosenthal and Tessler (Combinatorica 35, 2015) prove, for complexes with a complete skeleton, a Cheeger inequality and a Mixing Lemma: the number of $d$-cells with one vertex in each of $d+1$ disjoint vertex sets is controlled by the spectrum of the upper Laplacian. The Mixing Lemma is the result of this mission. In the same paper it yields a lower bound on Gromov's geometric overlap (Gromov 2010) of complexes with a small spectral spread.

Setting

Let $V=\{0,1,\dots,n-1\}$ and $d\ge1$. A $d$-dimensional complex with a complete skeleton $X$ contains every set of at most $d$ vertices as a cell, together with a set $X^d$ of $d$-cells, each a set of $d+1$ vertices. Cells are oriented by the increasing order of their vertices.

A $(d-1)$-form is a real function on the $d$-element subsets of $V$; the space $\Omega^{d-1}$ carries the inner product $\langle f,g\rangle=\sum_\sigma f(\sigma)g(\sigma)$. Similarly $\Omega^{d-2}$ consists of functions on $(d-1)$-subsets ($\Omega^{-1}=\mathbb R^{\{\varnothing\}}$ when $d=1$) and $\Omega^d$ of functions on $d$-cells. Writing $\sigma_i$ for the $i$-th vertex of $\sigma$, the co-boundary operators are $(\partial^*_j f)(\sigma)=\sum_{i}(-1)^i f(\sigma\setminus\sigma_i)$ (for $j=d$ only on $\sigma\in X^d$), and the boundary operators are $(\partial_j g)(\sigma)=\sum_{v\sim\sigma} g(v\sigma)$, the sum over vertices $v\notin\sigma$ with $\{v\}\cup\sigma$ a cell, where $v\sigma$ is the oriented cell $[v,\sigma_0,\sigma_1,\dots]$.

The upper, lower and full Laplacians on $\Omega^{d-1}$ are

$$\Delta^+=\partial_d\partial^*_d,\qquad \Delta^-=\partial^*_{d-1}\partial_{d-1},\qquad \Delta=\Delta^++\Delta^-.$$

The $(d-1)$-cycles are $Z_{d-1}=\ker\partial_{d-1}$, and the exact forms are $B^{d-1}=\operatorname{im}\partial^*_{d-1}$; they are orthogonal complements in $\Omega^{d-1}$. Write $\mathbb P_{Z_{d-1}}$ and $\mathbb P_{B^{d-1}}$ for the orthogonal projections. For vertex sets $A_0,\dots,A_d$, $F(A_0,\dots,A_d)$ is the set of $d$-cells with exactly one vertex in each $A_i$.

For $\alpha\in\mathbb R$ the operator $\alpha I-\Delta^+$ preserves $Z_{d-1}$, and

$$\rho_\alpha=\max\bigl\{|\mu|\;:\;\mu\text{ an eigenvalue of }(\alpha I-\Delta^+)|_{Z_{d-1}}\bigr\}=\bigl\|(\alpha I-\Delta^+)|_{Z_{d-1}}\bigr\|.$$

Formalization targets

Goal: Theorem (1.4), the Mixing Lemma (p. 16)

For every $\alpha\in\mathbb R$ and any pairwise disjoint vertex sets $A_0,\dots,A_d$ (not necessarily a partition),

$$\left|\,|F(A_0,\dots,A_d)|-\frac{\alpha\cdot|A_0|\cdots|A_d|}{n}\right|\le\rho_\alpha\cdot\bigl(|A_0|\cdots|A_d|\bigr)^{\frac d{d+1}}.$$

The parameter $\alpha$ is free; the introduction's Theorem 1.4 (p. 4) is the case $\alpha=k$, the average degree of a $(d-1)$-cell, and is included as a companion item.

Milestones

The milestones follow the proof of §4.3 and the facts of §§2–3 that it uses:

1. the chain-complex identity $\partial_{d-1}\partial_d=0$ and its dual, so $B^{d-1}\subseteq\ker\Delta^+$ (p. 6), and $B^{d-1}=Z_{d-1}^\perp$ (p. 6);
2. Proposition 3.2: $\Delta^+_{\overline X}=nI-\Delta_X$ for the complement complex (3.4), and $\Delta^-=n\,\mathbb P_{B^{d-1}}$ (3.5);
3. for the signed indicator forms $\varphi=\delta_{A_0,A_1,\dots,A_{d-1}}$ and $\psi=\delta_{A_d,A_1,\dots,A_{d-1}}$: the norm identity (4.5), the count $\langle\varphi,(D-\Delta^+)\psi\rangle=|F(A_0,\dots,A_d)|$ (4.6), its $\alpha$-shift (4.7), the splitting (4.8), the main term (4.9), the spectral bound (4.10)–(4.11);
4. the asymmetric bound $\rho_\alpha\sqrt{|A_0||A_d|}\,|A_1|\cdots|A_{d-1}|$ and the geometric-mean step that symmetrizes it.

Significance

For $d=1$ the theorem is the Expander Mixing Lemma of graphs. For $d\ge2$ it shows that a complex whose upper Laplacian has its nontrivial spectrum concentrated near $\alpha$ distributes its $d$-cells evenly across any $d+1$ disjoint vertex sets, with a bound that depends on the set sizes only through $(\prod|A_i|)^{d/(d+1)}$. The paper uses it to bound the geometric overlap of such complexes from below (Corollary 1.6), and to show that random Linial–Meshulam complexes (Linial–Meshulam 2006) have small discrepancy.

The result is proved in the paper. To the knowledge of this mission, neither it nor the underlying signed simplicial Laplacians with their Hodge decomposition on a complete skeleton has been formalized. The formalization produces the boundary and co-boundary operators, the three Laplacians, the identity for the complement complex, and the projection formula for $\Delta^-$, all of which are standard in the spectral theory of complexes and reusable well beyond this result.

Difficulty

Each step of the proof is a finite computation, but the steps are sign-sensitive. The count (4.6) relies on the orientation conventions: $\psi$ replaces only the first block of $\varphi$, and the off-diagonal terms of $\Delta^+$ carry the signs that make $\varphi(\sigma)\psi(\sigma')=+1$ on neighbouring cells. Proposition 3.2 (1) and (3) are identities between explicit matrices indexed by $d$-subsets that hold only because the skeleton is complete. The spectral step needs the orthogonal decomposition $\Omega^{d-1}=B^{d-1}\oplus Z_{d-1}$, invariance of $Z_{d-1}$ under $\Delta^+$, and the spectral theorem for the restriction. A direct transcription of the graph argument does not apply: the operators act on $(d-1)$-forms rather than on functions on vertices, and the expected count $\alpha\prod|A_i|/n$ is not visible from $\Delta^+$ of $X$ alone.

Formalization scope

Vertices are Fin n. A complex is the structure Complex n d holding its $d$-cells, so the complete skeleton is built in; forms are EuclideanSpace ℝ over the subtype of $k$-element subsets of Fin n, oriented by the order of Fin n. All boundary and co-boundary operators are defined by their formulas, not as adjoints of one another. The standing assumption $d\ge1$ is added to every statement. The disjoint sets $A_0,\dots,A_d$ form an arbitrary pairwise-disjoint family, some possibly empty; no partition and no bound $n\ge d+1$ is assumed. The exponent $d/(d+1)$ is a real power of a nonnegative real.

$\rho_\alpha$ is not defined as a maximum or supremum: the goal is stated for every real $\rho$ that bounds $|\mu|$ for all eigenvalues $\mu$ of $\alpha I-\Delta^+$ on $Z_{d-1}$, which is equivalent to the statement for $\rho_\alpha$. Trivializing formalizations are ruled out by construction: forms cannot be supported outside $d$-subsets (such a form would be a spurious cycle), an empty eigenvalue set cannot be read through a default sSup value, the boundary operators carry orientation signs, and disjoint families are not confused with partitions.

Contributions are welcome at every milestone, in particular a general library for signed boundary operators on finite complexes, the identity $\partial\partial=0$, and the Hodge-type decomposition of forms, which are reusable for the Cheeger inequality of the companion mission and for other work on high-dimensional expanders.

Selected references

• O. Parzanchevski, R. Rosenthal, R. J. Tessler, Isoperimetric inequalities in simplicial complexes, Combinatorica 35 (2015); arXiv:1207.0638v3 (2013). https://arxiv.org/abs/1207.0638 , https://doi.org/10.1007/s00493-014-3002-x
• N. Alon, F. R. K. Chung, Explicit construction of linear sized tolerant networks, Discrete Mathematics 72 (1988). https://doi.org/10.1016/0012-365X(88)90189-6
• S. Hoory, N. Linial, A. Wigderson, Expander graphs and their applications, Bulletin of the AMS 43 (2006). https://doi.org/10.1090/S0273-0979-06-01126-8
• N. Linial, R. Meshulam, Homological connectivity of random 2-complexes, Combinatorica 26 (2006). https://doi.org/10.1007/s00493-006-0027-9
• M. Gromov, Singularities, expanders and topology of maps. Part 2: From combinatorics to topology via algebraic isoperimetry, Geometric and Functional Analysis 20 (2010). https://doi.org/10.1007/s00039-010-0073-8

Motivation

Polynomial optimization problems invariant under permutation of the variables appear throughout global optimization, and their semidefinite (SDP) relaxations, the moment–sum-of-squares hierarchy of Lasserre (2001), grow quickly with the number of variables $n$. Riener, Theobald, Jansson Andrén and Lasserre (arXiv:1103.0486v3, Math. Oper. Res. 2013) show how symmetry reduces the size of these relaxations. Their last section treats a model family in which the reduction is extreme: minimizing one power sum subject to fixed values of the lower power sums.

For this family the paper gives two bounds. A lower bound (Theorem 6.6, the subject of a companion mission) comes from an SDP over an $n \times n$ Hankel matrix. The present mission concerns the complementary upper bound, Theorem 6.7, which comes from an SDP in a single real variable whose matrix has size $m$, the number of constraints plus one, independently of $n$. The idea behind it, that candidate optima have at most $m$ non-zero coordinates, is the same phenomenon studied in the degree principle for symmetric polynomials of Timofte (2003) and Riener (2012).

Setting

Fix $n, m, q \in \mathbb{N}$ and real numbers $\gamma_1, \dots, \gamma_{m-1}$. For $x \in \mathbb{R}^n$ write $s_j(x) = \sum_{i=1}^n x_i^j$ for the power sum of order $j$. Problem (6.4) is

$$\mathrm{P}_{nmq}: \quad \min \sum_{i=1}^n x_i^q \quad \text{s.t.} \quad \sum_{i=1}^n x_i^j = \gamma_j, \quad j = 1, \dots, m-1,$$

with standing assumptions $q \ge m$ and $m \le n + 1$, and optimal value $\min \mathrm{P}_{nmq}$ (an infimum, possibly $+\infty$ if the constraints are inconsistent, or $-\infty$).

For a sequence $s = (s_0, s_1, \dots)$ the Hankel matrix of order $k$ is $H_k(s) = (s_{i+j-2})_{1 \le i,j \le k}$. For a monic real polynomial

$$p = X^m + \sum_{j=0}^{m-1} p_j X^j$$

with complex roots $x_1, \dots, x_m$ (with multiplicity), the Newton sums are $s_k = \sum_{j=1}^m x_j^k$, so $s_0 = m$. Newton's identities express the $s_k$ through the coefficients and back; in particular $s_1, \dots, s_{m-1}$ determine $p_{m-1}, \dots, p_1$, so that, once $s_j = \gamma_j$ for $j < m$ is imposed, the only free coefficient is $p_0$, and every later Newton sum is a known polynomial $s_j = Q_j(p_0)$, $j \ge m$. For $q \le 2m - 2$ the paper defines

$$\mathrm{U}_{nmq} = \min_{p_0} \{\, Q_q(p_0) : H_m(s) \succeq 0 \,\}, \tag{6.7}$$

where $s_0 = m$, $s_j = \gamma_j$ for $1 \le j \le m-1$, and $s_j = Q_j(p_0)$ for $m \le j \le 2m-2$. The constraint $H_m(s) \succeq 0$ is equivalent to the real-rootedness of $p$ (Hermite's criterion).

Formalization targets

Goal: Theorem 6.7

For $m \le n$ and $m \le q \le 2m-2$,

$$\min \mathrm{P}_{nmq} \le \mathrm{U}_{nmq},$$

and if $\mathrm{P}_{nmq}$ has an optimal solution $x^* \in \mathbb{R}^n$ with at most $m$ non-zero entries, then

$$\min \mathrm{P}_{nmq} = \mathrm{U}_{nmq}.$$

Both halves are part of the goal; the equality is what makes (6.7) an exact convex reformulation of (6.4) in the sparse case.

Milestones

In attack order:

1. Newton's identities for a monic polynomial of degree $m$, in the two forms printed on p. 27 ($k \ge m$ and $1 \le k < m$).
2. Given $s_1, \dots, s_{m-1}$ and $p_0$, there is exactly one monic polynomial of degree $m$ with those Newton sums and constant coefficient $p_0$.
3. $Q_j$ is affine in $p_0$ whenever $m \le j \le 2m-1$.
4. Hermite's criterion: all roots of $p$ are real if and only if $H_m(s) \succeq 0$.
5. From the proof of Theorem 6.7: if $p$ is feasible for (6.7) and real-rooted, its roots padded with zeros form a feasible point of $\mathrm{P}_{nmq}$ with objective $s_q$.

Significance

The result. Theorem 6.7 bounds a nonconvex polynomial optimization problem in $n$ variables from above by a semidefinite program in one variable with an $m \times m$ matrix. Combined with Theorem 6.6 it brackets $\min \mathrm{P}_{nmq}$ between two SDP values, and when an $m$-sparse optimum exists, it identifies the optimal value exactly by a problem whose size does not grow with $n$. The milestones are classical facts about univariate polynomials (Newton's identities, Hermite's criterion) that are used far beyond this paper, in real algebraic geometry, the theory of hyperbolic and real-rooted polynomials, and moment problems.

Formalizing it. The inequality has a short printed proof; the equality case has none, so the formalization supplies a complete argument for a stated but unproved claim. Hermite's real-root criterion is, as far as the platform's library search shows, not available in Mathlib, and a formal proof of it is a reusable contribution in its own right. Newton's identities exist in Mathlib for multivariate symmetric polynomials; their form for the coefficients and roots of a monic univariate polynomial, as printed here, is a separate statement. No part of Theorem 6.7 is known to have a machine-checked proof.

Difficulty

The obvious argument for the inequality, sending a feasible $p_0$ to the vector of roots of $p$, needs the roots to be real, and real-rootedness is encoded only through positive semidefiniteness of a Hankel matrix of Newton sums. The step "psd Hankel matrix implies real roots" is Hermite's criterion, a nontrivial classical theorem with no library support. The equality half needs the converse direction of the same criterion together with a construction that the text does not give. Finally, the printed proof starts from "an optimal solution $p_0$ of the SDP (6.7)", which need not exist; the formal statement is about infima and must be proved without assuming attainment.

Formalization scope

Points are Fin n → ℝ; $\gamma$ is a sequence ℕ → ℝ of which only $\gamma_1, \dots, \gamma_{m-1}$ are read. Both $\min \mathrm{P}_{nmq}$ and $\mathrm{U}_{nmq}$ are infima in EReal over their feasible sets, so an infeasible problem has value $+\infty$ and no optimizer is assumed. The Hankel matrix is indexed from $0$: entry $(i,j)$ is $s_{i+j}$, $0 \le i, j < k$, which is the paper's $(s_{i+j-2})_{1 \le i,j \le k}$. Newton sums of a real polynomial are the real part of $\sum_z z^k$ over its complex roots with multiplicity (Polynomial.aroots), which is real for real polynomials.

The feasible set of (6.7) is parametrized by the monic polynomial $p$ itself rather than by $p_0$: $p$ monic of degree $m$, $s_j(p) = \gamma_j$ for $1 \le j \le m-1$, and $H_m(s(p)) \succeq 0$. Milestone 2 shows that this is a bijective change of coordinates $p \leftrightarrow p_0$, and $s_j(p)$ for $j \ge m$ is the paper's $Q_j(p_0)$ by definition. The hypothesis $m \le q$ in the goal is the standing assumption $q \ge m$ of (6.4); with $q \le 2m-2$ it forces $m \ge 2$. In milestone 5, $q \ge 1$ is needed because the padding zeros contribute $0^q$.

Several formalizations would trivialize or falsify the targets and are ruled out: values are never real infima (which are $0$ on empty or unbounded sets); $\mathrm{U}_{nmq}$ is defined by (6.7) with $s_m, \dots, s_{2m-2}$ tied to $p_0$, not by the relaxation (6.6) in which those entries are free and which can have a strictly smaller value; and the equality half is stated in full, not weakened to an inequality.

The needed infrastructure is Vieta's formulas and Newton's identities over Polynomial, Hermite's quadratic form and its signature, and the passage between root multisets and coordinate vectors. All of it is reusable outside this mission; contributions of Hermite's criterion in particular are welcome.

Selected references

• C. Riener, T. Theobald, L. Jansson Andrén, J. B. Lasserre, Exploiting symmetries in SDP-relaxations for polynomial optimization, Math. Oper. Res. 38(1):122–141, 2013. arXiv:1103.0486v3 (the version formalized here): https://arxiv.org/abs/1103.0486, DOI: https://doi.org/10.1287/moor.1120.0558
• J. B. Lasserre, Global optimization with polynomials and the problem of moments, SIAM J. Optim. 11(3):796–817, 2001. https://doi.org/10.1137/S1052623400366802
• V. Timofte, On the positivity of symmetric polynomial functions. Part I: General results, J. Math. Anal. Appl. 284(1):174–190, 2003. https://doi.org/10.1016/S0022-247X(03)00301-9
• C. Riener, On the degree and half-degree principle for symmetric polynomials, J. Pure Appl. Algebra 216(4):850–856, 2012. https://doi.org/10.1016/j.jpaa.2011.08.012
• S. Basu, R. Pollack, M.-F. Roy, Algorithms in Real Algebraic Geometry, 2nd ed., Springer, 2006 (Hermite's quadratic form). https://doi.org/10.1007/3-540-33099-2

Motivation

Many polynomial optimization problems are invariant under permuting the variables: the objective and the constraints are symmetric polynomials. Lasserre's moment–sum-of-squares hierarchy solves such problems by a sequence of semidefinite programs whose size grows quickly with the number of variables and the degree, and symmetry is the main tool for shrinking them. Riener, Theobald, Jansson Andrén and Lasserre (arXiv:1103.0486, Math. Oper. Res. 38(1), 2013, doi:10.1287/moor.1120.0558) develop symmetry reductions of the hierarchy and, in their Section 6, a second route: instead of reducing the moment matrices, describe the orbit space of the symmetry group directly by polynomial inequalities, following Procesi and Schwarz, and optimize over it.

For the symmetric group the orbit space is described by a single Hankel matrix of power sums. This mission formalizes the lower bound that results for a basic family of symmetric problems: minimize one power sum subject to prescribed values of the lower power sums.

Setting

Fix $n \in \mathbb N$ and a point $x = (x_1,\dots,x_n) \in \mathbb R^n$. The power sum of order $j$ is $s_j(x) = \sum_{i=1}^n x_i^{\,j}$; note $s_0(x) = n$.

The problem $\mathrm P_{nmq}$ (6.4). Given $m, q \in \mathbb N$ with $q \ge m$, $m \le n+1$, and prescribed reals $\gamma_1,\dots,\gamma_{m-1}$, minimize $s_q(x)$ over all $x \in \mathbb R^n$ with $s_j(x) = \gamma_j$ for $j = 1,\dots,m-1$. Its optimal value is written $\min \mathrm P_{nmq}$.

The Hankel matrix. For a real sequence $s = (s_0, s_1, \dots)$ and $k \in \mathbb N$, $H_k(s) = (s_{i+j-2})_{1\le i,j\le k}$.

The matrix $J$ (6.1). For a finite group acting on $\mathbb R^n$ with fundamental invariants $\pi_1,\dots,\pi_l$, the matrix $J = (\langle d\pi_i, d\pi_j\rangle)_{1\le i,j\le l}$ has $G$-invariant entries $\sum_k \partial_k\pi_i\,\partial_k\pi_j$. For the symmetric group with $\pi_j = s_j/j$, $J$ is the Hankel matrix $H_n(s)$ of power sums.

The SDP (6.5). For $m \le n+1$ and $m \le q \le 2n-2$,

$$\mathrm L_{nmq} = \min_s \{\, s_q \;:\; H_n(s) \succeq 0,\ s_0 = n,\ s_j = \gamma_j,\ j = 1,\dots,m-1 \,\}.$$

The variable is the sequence $s$; the constraint is a single linear matrix inequality of size $n$.

Formalization targets

Goal: Theorem 6.6

For $m \le n+1$ and $m \le q \le 2n-2$:

(a)

$$\min \mathrm P_{nmq} \ge \mathrm L_{nmq}.$$

(b) If $q = m = 2r$, put $\gamma_0 = n$, $H_r(\gamma) = (\gamma_{i+j-2})_{1\le i,j\le r}$ and $u_r(\gamma)^T = (\gamma_r,\dots,\gamma_{2r-1})$. If $H_r(\gamma) \succ 0$, then

$$\min \mathrm P_{nmq} \ge u_r(\gamma)^T H_r(\gamma)^{-1} u_r(\gamma).$$

The goal is the conjunction of (a) and (b).

Milestones

1. For $\pi_j = s_j / j$, $J = H_n(s)$ (p. 26).
2. The power-sum vector of a feasible point of (6.4) is feasible for (6.5) (proof of (a), p. 26).
3. $H_n(s) \succeq 0 \Rightarrow H_{r+1}(s) \succeq 0$ for $r < n$ (proof of (b), p. 26).
4. The Schur complement bound $u^T H_r(s)^{-1} u \le s_{2r}$ when $H_{r+1}(s) \succeq 0$ and $H_r(s) \succ 0$ (proof of (b), p. 26).

Significance

Part (a) bounds a nonconvex polynomial problem of degree $q$ in $n$ variables by a semidefinite program with one $n\times n$ matrix variable, whose size does not depend on $q$. Part (b) evaluates the bound in closed form in the even case $q = m$, where the SDP reduces to a single scalar constraint on the free moment $s_{2r}$. Together with the upper bound of Theorem 6.7 (a separate mission of this series), it brackets the optimal value of power-sum problems, which appear in moment problems, in the analysis of symmetric inequalities, and as test cases for symmetry-adapted relaxations.

The paper's proof is a few lines, citing the Procesi–Schwarz theorem and Schur's complement. A formal development makes each step explicit: the positive semidefiniteness of the Hankel matrix of power sums of a real point, the principal-submatrix argument, and the Schur complement for a Hankel block. It also fixes a misprint: the bound of (b) is printed as $u_r(\gamma)^T H_r(\gamma)\, u_r(\gamma)$ without the inverse, which is false (for $n = 2$, $m = q = 2$, $\gamma_1 = 1$ the minimum is $1/2$ but the printed bound is $2$). No machine-checked version of this result is known to exist.

Difficulty

The individual steps are classical, and the difficulty is in the bookkeeping. The paper reaches (a) through the Procesi–Schwarz description of the orbit space, a deep result stated for general finite groups; the formal route needs only the positive semidefiniteness of the Hankel matrix of power sums of every real point, which the page leaves implicit. The identity $J = H_n(s)$ requires computing partial derivatives of power sums with the right index shifts. In (b), $H_{r+1}(s)$ must be identified with a $2\times 2$ block matrix whose top-left block is $H_r(\gamma)$ with $\gamma_0 = n$, which needs the prescribed values $\gamma_1,\dots,\gamma_{2r-1}$ and a reindexing of $\{0,\dots,r\}$ as $\{0,\dots,r-1\}\sqcup\{r\}$. The extended-real values must be handled: the problem can be infeasible (value $+\infty$), and for odd $q$ unbounded below.

Formalization scope

Lean namespace SymPolyOpt.PowerSumLB. Conventions:

• A point is x : Fin n → ℝ; $s_j(x)$ is powerSum x j, so powerSum x 0 = n.
• $\gamma \in \mathbb R^{m-1}$ is a sequence γ : ℕ → ℝ of which only $\gamma_1,\dots,\gamma_{m-1}$ are read. The page's "$=\gamma_k$" in (6.4) is read as $\gamma_j$.
• $\min \mathrm P_{nmq}$ and $\mathrm L_{nmq}$ are infima in EReal: $+\infty$ for an infeasible problem, possibly $-\infty$. They are compared only by ≤.
• The SDP variable is a sequence s : ℕ → ℝ; its entries past $s_{2n-2}$ are unconstrained and irrelevant since $q \le 2n-2$.
• $H_k(s)$ is 0-indexed, $H_k(s)_{ij} = s_{i+j}$ for $0 \le i,j < k$, the same matrix as the page's.
• $J$ is formed with the product of partial derivatives as the pairing of (6.1), using Mathlib's MvPolynomial.psum and pderiv.
• $q \le 2n-2$ is natural-number subtraction; for $n \le 1$ it forces $q = 0$.

Disclosed changes relative to the page:

• Theorem 6.6(b) is stated with $H_r(\gamma)^{-1}$, as the paper's proof gives; the printed form without the inverse is false.
• $H_r(\gamma) \succ 0$ is added to (b) and to the Schur milestone, so that the inverse exists.

Trivializing formalizations are ruled out: real-valued infima (which return $0$ on empty or unbounded sets) are not used; $\mathrm L_{nmq}$ is an SDP over arbitrary sequences, not over power-sum vectors of points, so (a) is not definitional; the constraint $s_0 = n$ is kept; Lean's junk inverse of a singular matrix is excluded by the positive-definiteness hypothesis; and the printed false bound is not stated.

Needed infrastructure: Mathlib's Matrix.PosSemidef, Matrix.PosDef, PosSemidef.submatrix and the Schur-complement lemmas for block matrices (fromBlocks). A reusable lemma that Hankel matrices of power sums (and, more generally, of moments of a finitely supported measure) are positive semidefinite would serve other moment-problem missions. Contributions to any milestone are welcome.

Selected references

• C. Riener, T. Theobald, L. Jansson Andrén, J. B. Lasserre, Exploiting symmetries in SDP-relaxations for polynomial optimization, Math. Oper. Res. 38(1):122–141, 2013. arXiv:1103.0486v3, doi:10.1287/moor.1120.0558
• C. Procesi, G. Schwarz, Inequalities defining orbit spaces, Invent. Math. 81:539–554, 1985. doi:10.1007/BF01388587
• J. B. Lasserre, Global optimization with polynomials and the problem of moments, SIAM J. Optim. 11(3):796–817, 2001. doi:10.1137/S1052623400366802

Motivation

Phase recovery (also called phase retrieval) asks to reconstruct a complex signal from the magnitudes of linear measurements of it, the phases being lost. It arises wherever detectors record intensities only: X-ray crystallography, diffraction imaging, astronomy and optics, where a sensor measures $|\langle a_i, x\rangle|^2$ and not $\langle a_i, x\rangle$. The classical algorithms are alternating projections of Gerchberg–Saxton type (Gerchberg and Saxton, 1972; Fienup, 1982), which are fast but non-convex and can stall at stationary points that are not solutions.

Two convex relaxations change this picture. PhaseLift (Chai, Moscoso and Papanicolaou, 2011; Candès, Strohmer and Voroninski, 2013) lifts the signal $x$ to the rank-one matrix $xx^*$ and relaxes the rank constraint to a trace minimization; for random Gaussian measurements it provably recovers $x$ with high probability (Candès and Li, 2014). Waldspurger, d'Aspremont and Mallat (arXiv:1206.0102v3, Math. Programming 2015) instead separate amplitude and phase and obtain PhaseCut, a semidefinite program with the structure of the Goemans–Williamson MaxCut relaxation, for which MaxCut algorithms (block coordinate descent, rounding) apply. The question this mission formalizes is how the two relaxations compare in the noiseless case.

Timeline. An earlier version of the paper showed that PhaseLift tightness implies PhaseCut tightness; Voroninski (2012, working paper, cited in arXiv:1206.0102v3, §4.3) reversed the argument under additional assumptions. Version 3 of the paper (2013) introduces the refined program PhaseCutMod and proves that, for injective $A$, PhaseCutMod is tight whenever PhaseLift is (Corollary 4.3).

Setting

Fix a measurement matrix $A\in\mathbb C^{n\times p}$ with rows $a_1^*,\dots,a_n^*$, and a vector of amplitudes $b\in\mathbb R^n$. The phase recovery problem (1) is

$$\text{find } x\in\mathbb C^p \text{ such that } |Ax| = b,$$

i.e. $|a_i^*x| = b_i$ for every $i$. Write $\mathbf H_m$ for the $m\times m$ Hermitian matrices and $X\succeq 0$ for positive semidefiniteness. When $A$ is injective, $A^*A$ is invertible and the pseudoinverse is $A^\dagger = (A^*A)^{-1}A^*$; $AA^\dagger$ is the orthogonal projector onto the range of $A$. With $\operatorname{diag}(b)$ the diagonal matrix of the $b_i$, set

$$M = \operatorname{diag}(b)(\mathbf I - AA^\dagger)\operatorname{diag}(b),\qquad B = \operatorname{diag}(b)A^{\dagger *}A^\dagger\operatorname{diag}(b).$$

The three programs are

$$\text{(PhaseLift)}\quad \min_{X\in\mathbf H_p} \operatorname{Tr}(X)\ \ \text{s.t.}\ \operatorname{Tr}(a_ia_i^*X) = b_i^2,\ i=1,\dots,n,\ X\succeq 0;$$ $$\text{(PhaseCut)}\quad \min_{U\in\mathbf H_n} \operatorname{Tr}(UM)\ \ \text{s.t.}\ \operatorname{diag}(U) = 1,\ U\succeq 0;$$ $$\text{(PhaseCutMod)}\quad \min_{U\in\mathbf H_n} \operatorname{Tr}(BU)\ \ \text{s.t.}\ \operatorname{Tr}(MU) = 0,\ \operatorname{diag}(U) = 1,\ U\succeq 0.$$

PhaseCutMod selects, among the optimal points of PhaseCut (which in the noiseless case are exactly its points with $\operatorname{Tr}(MU)=0$), one of minimal $\operatorname{Tr}(BU)$. A relaxation is tight when it has a unique optimal solution and that solution has rank one. The two lifted programs are linked by the map

$$\Phi(X) = \operatorname{diag}(b)^{-1}AXA^*\operatorname{diag}(b)^{-1}.$$

In Lean these objects are PhaseCut.Tight.pinv, Mmat, Bmat, Phi, IsSolvable, phaseLiftFeasible/phaseLiftOptimal, phaseCutFeasible/phaseCutOptimal, phaseCutModFeasible/phaseCutModOptimal and IsTight.

Formalization targets

Goal: Corollary 4.3 (p. 13)

If $A$ is injective, $b_i\neq 0$ for all $i$, and (1) has an exact solution, then

$$\text{PhaseLift is tight}\ \Longrightarrow\ \text{PhaseCutMod is tight}.$$

The statement fixes no constants and no probabilistic model: it is a deterministic comparison of two programs, valid for every $A$ and $b$ meeting the three hypotheses.

Milestones

1. $M\succeq 0$ (§2.2, p. 4).
2. For $U, M\succeq 0$, $\operatorname{Tr}(UM)=0$ implies $UM=0$ (proof of Proposition 4.2, p. 12).
3. $\operatorname{Tr}(MU)=0,\ U\succeq 0$ (15) holds iff $U=\Phi(X)$ for some $X\succeq 0$ (16) (proof of Proposition 4.2, p. 12).
4. Proposition 4.2 (p. 12): $\Phi$ is a bijection from the feasible set of PhaseLift onto the feasible set of PhaseCutMod.
5. $\operatorname{Tr}(X) = \operatorname{Tr}(B\Phi(X))$ (proof of Corollary 4.3, p. 13).
6. $\operatorname{Rank}(X) = \operatorname{Rank}(\Phi(X))$ (proof of Corollary 4.3, p. 13).
7. Companion: when (1) is solvable, every optimal solution of PhaseCut is feasible for PhaseCutMod (§4.3, p. 12).

Significance

The result. Corollary 4.3 says that the MaxCut-type relaxation, once refined to PhaseCutMod, is never worse than PhaseLift at exact recovery: every recovery guarantee proved for PhaseLift (for instance the Gaussian-measurement results of Candès–Strohmer–Voroninski and Candès–Li) transfers to PhaseCutMod. Since PhaseCutMod can be approached by solving PhaseCut with $M$ replaced by $M+\gamma B$ for small $\gamma>0$, the algorithms designed for PhaseCut (greedy and block-coordinate SDP solvers) inherit those guarantees. The bijection of Proposition 4.2 also gives a precise dictionary between the lifted signal space $\mathbf H_p$ and the phase space $\mathbf H_n$.

Formalizing it. The result is proved in the paper; to the best of current knowledge none of it has a machine-checked proof. The mission produces a formal statement of the three programs over complex Hermitian matrices, of the pseudoinverse identities for injective matrices, and of the trace and rank transport under the congruence $\Phi$, which are reusable for other lifting arguments in phase retrieval and low-rank SDP.

Difficulty

The corollary itself is short once Proposition 4.2 and the two transport identities are available; the substance lies in Proposition 4.2. The natural first attempt is to check that $\Phi$ maps feasible points to feasible points, which is direct. Surjectivity is the step that requires an argument: a PhaseCutMod-feasible $U$ satisfies only the scalar constraint $\operatorname{Tr}(MU)=0$, and one must first upgrade it to the matrix equation $MU=0$ (using positive semidefiniteness of both factors) and then to a range condition on $\operatorname{diag}(b)U\operatorname{diag}(b)$ before a preimage can be written down. In Lean, an additional cost is the absence of a pseudoinverse and of trace-zero-product lemmas for complex Hermitian matrices in the form needed, and the bookkeeping between complex traces and real objectives.

Formalization scope

Matrices are Matrix (Fin n) (Fin p) ℂ and Matrix (Fin m) (Fin m) ℂ; $X\succeq 0$ is Mathlib's Matrix.PosSemidef (Hermitian included); rank is Matrix.rank. Conventions committed to:

• Injectivity is Function.Injective A.mulVec. $A^\dagger$ is encoded as (Aᴴ * A)⁻¹ * Aᴴ, Mathlib's inverse returning $0$ on singular matrices; it is the pseudoinverse only for injective $A$, so every statement mentioning $M$, $B$ or $\Phi$ assumes injectivity. This adds injectivity, not present on the page, to the claim $M\succeq 0$ (p. 4) and to the companion claim on p. 12.
• $\operatorname{diag}(b)^{-1}$ has entries $b_i^{-1}$ and is used only under $b_i\neq 0$.
• $\operatorname{Tr}(a_ia_i^*X) = a_i^*Xa_i$ is the diagonal entry $(AXA^*)_{ii}$.
• Objectives are real parts of complex traces (the traces are real on the feasible sets); the PhaseCutMod constraint $\operatorname{Tr}(MU)=0$ is an equation in $\mathbb C$.
• Optimal sets are sets of minimizers over the feasible set; tight means the optimal set equals $\{X\}$ for some $X$ of rank one. Weakening this to "some optimal solution has rank one" would be a different, weaker statement and is excluded, as is applying tightness to the feasible set instead of the optimal set.
• Proposition 4.2's "bijection" is Set.BijOn Φ (PhaseLift feasible) (PhaseCutMod feasible).
• The hypotheses of Proposition 4.2 and Corollary 4.3 are the page's, including the solvability of (1).

Infrastructure that a complete development needs, and that is reusable beyond this mission: the pseudoinverse identities $A^\dagger A = \mathbf I$ and $AA^\dagger A = A$ for injective complex $A$, the fact that $\mathbf I - AA^\dagger$ is a Hermitian idempotent, the lemma $\operatorname{Tr}(UM)=0\Rightarrow UM=0$ for complex positive semidefinite matrices, and rank invariance under congruence by injective matrices. Proofs of any milestone are welcome, as are general Mathlib-style lemmas on these points.

Selected references

• I. Waldspurger, A. d'Aspremont, S. Mallat, Phase recovery, MaxCut and complex semidefinite programming, Mathematical Programming 149 (2015); preprint arXiv:1206.0102v3 (2013). https://arxiv.org/abs/1206.0102
• E. J. Candès, T. Strohmer, V. Voroninski, PhaseLift: Exact and stable signal recovery from magnitude measurements via convex programming, Comm. Pure Appl. Math. 66 (2013). https://arxiv.org/abs/1109.4499
• E. J. Candès, X. Li, Solving quadratic equations via PhaseLift when there are about as many equations as unknowns, Found. Comput. Math. 14 (2014). https://arxiv.org/abs/1208.6247
• A. Chai, M. Moscoso, G. Papanicolaou, Array imaging using intensity-only measurements, Inverse Problems 27 (2011). https://doi.org/10.1088/0266-5611/27/1/015005
• V. Voroninski, A comparison between the PhaseLift and PhaseCut algorithms, working paper, 2012 (as cited in arXiv:1206.0102v3, p. 28). https://arxiv.org/abs/1206.0102
• M. X. Goemans, D. P. Williamson, Improved approximation algorithms for maximum cut and satisfiability problems using semidefinite programming, J. ACM 42 (1995). https://doi.org/10.1145/227683.227684
• J. R. Fienup, Phase retrieval algorithms: a comparison, Applied Optics 21 (1982). https://doi.org/10.1364/AO.21.002758