Convergence Analysis of Nonconvex ADMM for Rigid Registration

TL;DR

This paper analyzes the convergence of a nonconvex ADMM algorithm for rigid registration, proving conditions under which it converges to the global optimum and explaining its stability properties.

Contribution

It provides the first convergence analysis of REG-ADMM for rigid registration, including explicit parameter choices and stability conditions based on noise levels.

Findings

REG-ADMM fixed points are stationary points of REG-SDP

Explicit ADMM parameter $ ho$ guarantees convergence to global optimum with clean data

High noise levels can cause REG-ADMM instability regardless of initialization

Abstract

We consider the problem of rigid registration, where we wish to jointly register multiple point sets via rigid transforms. This arises in applications such as sensor network localization, multiview registration, and protein structure determination. The least-squares estimator for this problem can be reduced to a rank-constrained semidefinite program (REG-SDP). It was recently shown that by formally applying the alternating direction method of multipliers (ADMM), we can derive an iterative solver (REG-ADMM) for REG-SDP, wherein each subproblem admits a simple closed-form solution. The empirical success of REG-ADMM has been demonstrated for multiview registration. However, its convergence does not follow from the existing literature on nonconvex ADMM. In this work, we study the convergence of REG-ADMM and our main findings are as follows. We prove that any fixed point of REG-ADMM is a stationary (KKT) point of REG-SDP. Moreover, for clean measurements, we give an explicit formula for the ADMM parameter $\rho$, for which REG-ADMM is guaranteed to converge to the global optimum (with arbitrary initialization). If the noise is low, we can still show that the iterates converge to the global optimum, provided they are initialized sufficiently close to the optimum. On the other hand, if the noise is high, we explain why REG-ADMM becomes unstable if $\rho$ is less than some threshold, irrespective of the initialization. We present simulation results to support our theoretical predictions. The novelty of our analysis lies in the fact that we exploit the notion of tightness of convex relaxation to arrive at our convergence results.

Equations

$${\bar{\text{x}}}_k={\bar{\mathbf{O}}}_i{\text{x}}_{k,i}+{\bar{\text{t}}}_i.$$

$$\underset{({\mathbf{O}}_i),({\text{t}}_i),({\text{z}}_k)}{\min }\sum _{i=1}^M\sum _{k\in P_i}\Vert {\text{z}}_k-\left({\mathbf{O}}_i{\text{x}}_{k,i}+{\text{t}}_i\right){\Vert }^2.$$

$$\underset{({\mathbf{O}}_i)}{\min }\left[\underset{({\text{t}}_i),({\text{z}}_k)}{\min }\sum _{i=1}^M\sum _{k\in P_i}\Vert {\text{z}}_k-({\mathbf{O}}_i{\text{x}}_{k,i}+{\text{t}}_i){\Vert }^2\right].$$

$$\underset{{\mathbf{O}}_1,\dots ,{\mathbf{O}}_M\in \mathbb{O}(d)}{\min }\sum _{i,j=1}^M\mathrm{Tr}\left({\mathbf{C}}_{ij}{\mathbf{O}}_j^{\top }{\mathbf{O}}_i\right).$$

$$\sum _{i=1}^M\sum _{k\in P_i}\Vert {\text{z}}_k^{\ast }-({\mathbf{O}}_i^{\ast }{\text{x}}_{k,i}+{\text{t}}_i^{\ast }){\Vert }^2$$

$$=\mathrm{Tr}\left({\mathbf{C}}_{ij}{\mathbf{O}}_j^{\ast \top }{\mathbf{O}}_i^{\ast }\right).$$

$$\sum _{i,j=1}^M\mathrm{Tr}\left({\mathbf{C}}_{ij}{\mathbf{O}}_j^{\top }{\mathbf{O}}_i\right)=\mathrm{Tr}\left(\mathbf{C}{\mathbf{O}}^{\top }\mathbf{O}\right).$$

$$\underset{\mathbf{G}\in {\mathbb{S}}_{+}^{Md}}{\min }$$

$$[\mathbf{G}{]}_{ii}={\mathbf{I}}_d,i\in [1:M],$$

$$\mathrm{rank}(\mathbf{G})=d.$$

$$\underset{\mathbf{G},\mathbf{H}\in {\mathbb{S}}^{Md}}{\min }$$

$$\mathbf{G}\in \Omega ,\mathbf{H}\in \Theta ,$$

$$\mathbf{G}-\mathbf{H}=0,$$

$${\mathcal{L}}_{\rho }(\mathbf{G},\mathbf{H},\mathbf{\Lambda })=\mathrm{Tr}\left(\mathbf{CG}\right)+\mathrm{Tr}(\mathbf{\Lambda }(\mathbf{G}-\mathbf{H}))+\frac{\rho }{2}\Vert \mathbf{G}-\mathbf{H}{\Vert }^2,$$

$${\mathbf{G}}^{k+1}$$

$${\mathbf{H}}^{k+1}$$

$${\mathbf{\Lambda }}^{k+1}$$

$${\mathbf{G}}^{k+1}$$

$${\mathbf{H}}^{k+1}$$

$${\Pi }_{\Omega }(\mathbf{A})=\sum _{i=1}^d\max ({\lambda }_i,0){\mathbf{u}}_i{\mathbf{u}}_i^{\top }.$$

$$\underset{\mathbf{G}\in {\mathbb{S}}_{+}^{Md}}{\min }$$

$$[\mathbf{G}{]}_{ii}={\mathbf{I}}_d,i\in [1:M].$$

$${\mathbf{G}}^{k+1}$$

$${\mathbf{H}}^{k+1}$$

$${\mathbf{\Lambda }}^{k+1}$$

$$\underset{\mathbf{G},\mathbf{H}\in {\mathbb{S}}^{Md}}{\min }$$

$$\mathbf{G}-\mathbf{H}=0,$$

$${\lambda }_i(\mathbf{X})-\Vert \mathbf{Y}{\Vert }_2\le {\lambda }_i(\mathbf{X}+\mathbf{Y})\le {\lambda }_i(\mathbf{X})+\Vert \mathbf{Y}{\Vert }_2.$$

$$\underset{{\mathbf{O}}_1,\cdots ,{\mathbf{O}}_M\in {\mathbb{R}}^{d\times d}}{\min }$$

$${\mathbf{I}}_d-{\mathbf{O}}_i^{\top }{\mathbf{O}}_i=0,i\in [1:M].$$

$$\mathcal{L}(({\mathbf{O}}_i),({\mathbf{\Lambda }}_i))=\sum _{i,j=1}^M\mathrm{Tr}\left([\mathbf{C}{]}_{ij}{\mathbf{O}}_j^{\top }{\mathbf{O}}_i\right)+\sum _{i=1}^M\mathrm{Tr}({\mathbf{\Lambda }}_i\left({\mathbf{I}}_d-{\mathbf{O}}_i^{\top }{\mathbf{O}}_i)\right),$$

$${\mathbf{O}}_i^{\ast }{\mathbf{\Lambda }}_i^{\ast }=-\sum _{j=1}^M{\mathbf{O}}_j^{\ast }[\mathbf{C}{]}_{ji},i\in [1:M].$$

$${\mathbf{\Lambda }}_i^{\ast }=-\sum _{j=1}^M{\mathbf{O}}_i^{\ast \top }{\mathbf{O}}_j^{\ast }[\mathbf{C}{]}_{ji}=-\sum _{j=1}^M[{\mathbf{G}}^{\ast }{]}_{ij}[\mathbf{C}{]}_{ji}=-[{\mathbf{G}}^{\ast }\mathbf{C}{]}_{ii}.$$

$$\displaystyle\underset{\mathbf{G}}{\min}\$$

$$s.t.$$

$$\mathbf{G}\in {\mathbb{S}}_{+}^{Md}.$$

$$[{\mathbf{\Lambda }}^{\ast }{]}_{ii}=-[\mathbf{C}{\mathbf{G}}^{\ast }{]}_{ii},i\in [1:M].$$

$${\mathbf{H}}^{k+1}={\Pi }_{\Theta }\left({\mathbf{G}}^{k+1}\right).$$

Full text

††journal:

1 Introduction

In rigid registration, we want to infer the global coordinates of a set of points, given the coordinates of overlapping subsets of these points in different local coordinate systems Chaudhury et al. (2015b). The local coordinate systems are related to each other by rigid transforms, which are otherwise unknown. This problem arises in applications such as sensor network localization, multiview registration, protein structure determination, and manifold learning Sanyal et al. (2017b); Cucuringu et al. (2012a); Zhang et al. (2010); Krishnan et al. (2005); Sharp et al. (2004); Cucuringu et al. (2012b); Fang and Toh (2013); Zhang and Zha (2004), where we wish to reconstruct an underlying global structure based on observations of multiple local sub-structures. For instance, consider an adhoc wireless network consisting of geographically distributed sensor nodes with limited radio range. To make sense of the data collected from the sensors, one usually requires the positions of the individual sensors. The positions can be found simply by attaching a GPS with each sensor, but this is often not feasible due to cost, power, and weight considerations. On the other hand, we can estimate (using time-of-arrival) the distance between sensors that are within the radio range of each other Mao et al. (2007). The problem of estimating sensor locations from the available inter-sensor distances is referred to as sensor network localization (SNL) Mao et al. (2007); Shang et al. (2004). Efficient methods for accurately localizing small-to-moderate sized networks have been proposed over the years Soares et al. (2015); Simonetto and Leus (2014); Wang et al. (2008); Biswas et al. (2006). However, these methods typically cannot be used to localize large networks. To address this, scalable divide-and-conquer approaches for SNL have been proposed Zhang et al. (2010); Cucuringu et al. (2012a); Chaudhury et al. (2015a); Sanyal et al. (2017b), where the large network is first subdivided into smaller subnetworks which can be efficiently and accurately localized (pictured in Fig. 1(a)). Each subnetwork (called patch) is localized independent of other subnetworks. Thus, the coordinates returned for a patch will in general be an arbitrarily rotated, flipped, and translated version of the ground-truth coordinates (Fig. 1(b)). Consequently, the network is divided into multiple patches, where each patch can be regarded as constituting a local coordinate system which is related to the global coordinate system by an unknown rigid transform. We now want to assign coordinates to all the nodes in a global coordinate system based on these patch-specific local coordinates.

Rigid registration also comes up in multiview registration, where the objective is to reconstruct a $3$D model of an object based on partial overlapping scans of the object. Here, the scans can be seen as patches, which are to be registered in a global reference frame via rotations and translations. Similar situation arises in protein conformation, where we are required to determine the $3$D structure of a protein (or other macromolecule) from overlapping fragments Cucuringu et al. (2012b); Fang and Toh (2013).

1.1 Problem Statement

We now formally describe the rigid registration problem and set up the rank constrained semidefinite program (REG-SDP) that is the focus of our analysis. We first clarify our notation: $[m:n]$ denotes the set of integers $\{m,\ldots,n\}$; $\mathbb{O}(d)$ denotes the set of $d\times d$ orthogonal matrices; $\mathbb{S}^{m}$ denotes the set of symmetric matrices of size $m\times m$; $\mathbb{S}^{m}_{+}$ denotes the set of symmetric positive semidefinite matrices of size $m\times m$; $\mathbf{I}_{d}$ denotes the $d\times d$ identity matrix. We will regard a matrix $\mathbf{A}\in\mathbb{R}^{Md\times Md}$ as a block matrix composed of $M\times M$ blocks, where each block is of size $d\times d$; we will refer to the $(i,j)$-th block as $[A]_{ij}$, where $1\leq i,j\leq M$. We also recall that a rigid transform in $\mathbb{R}^{d}$ is a composition of an orthogonal transform (rotation, reflection) and a translation, and is denoted by $(\mathbf{O},\textbf{t})$, where $\mathbf{O}\in\mathbb{O}(d)$ denotes the orthogonal transform, and $\textbf{t}\in\mathbb{R}^{d}$ denotes the translation. Specifically, a rigid transform $(\mathbf{O},\textbf{t})$ maps a point x to the point $\mathbf{O}\textbf{x}+\textbf{t}$.

Suppose a network consists of $N$ nodes in $\mathbb{R}^{d}$, which we label using $\mathcal{S}=[1:N]$. The network is sub-divided into $M$ patches (subsets of $\mathcal{S}$), labeled $P_{1},\cdots,P_{M}$. Let $\bar{\textbf{x}}_{1},\ldots,\bar{\textbf{x}}_{N}\in\mathbb{R}^{d}$ be the true coordinates of the $N$ nodes in some global coordinate system. We associate with each patch a local coordinate system: If node $k$ belongs to $P_{i}$, let $\textbf{x}_{k,i}\in\mathbb{R}^{d}$ be the local coordinate of node $k$ in patch $P_{i}$. In other words, if $(\bar{\mathbf{O}}_{i},\bar{\textbf{t}}_{i})$ is the (unknown) rigid transform (defined with respect to the global coordinate system) associated with patch $P_{i}$, then we have

[TABLE]

In general, the local coordinate measurements are noisy and we cannot expect (1) to hold exactly. Instead, we resort to least-squares fitting to solve the registration problem. More precisely, given the local coordinate measurements $(\textbf{x}_{k,i})$, we consider the following least-squares minimization to estimate the rigid transforms and the global coordinates Chaudhury et al. (2015b):

[TABLE]

The fundamental difficulty of REG stems from the fact that variables $\mathbf{O}_{i}$ are constrained to be in $\mathbb{O}(d)$, which is a nonconvex and disconnected set Absil et al. (2009). For instance, $\mathbb{O}(1)=\{-1,1\}$, and $\mathbb{O}(2)$ is topologically equivalent to the union of two disjoint circles in $\mathbb{R}^{2}$. The disconnectedness makes it difficult to apply local optimization methods; in particular, if one initializes the method on the wrong component, one cannot hope to get close to the global optimum.

To isolate the difficulty in REG, we will now reformulate it in a form that involves only the orthogonal transforms. To this end, note that REG can be rewritten as

[TABLE]

It was observed in Chaudhury et al. (2015b) that if we fix the orthogonal transforms $(\mathbf{O}_{i})$, the minimization problem inside the square brackets in (2) is a quadratic convex optimization in $(\textbf{t}_{i})$ and $(\textbf{z}_{k})$, with a closed-form optimum which is linear in $(\mathbf{O}_{i})$. In other words, if we knew the optimal orthogonal transforms for REG, the optimal values of $(\textbf{t}_{i})$ and $(\textbf{z}_{k})$ could be computed using a simple linear transform. Based on this observation (we refer the reader to Chaudhury et al. (2015b) for further details), the variables $(\textbf{t}_{i})$ and $(\textbf{z}_{k})$ can be eliminated, and REG is reduced to the following optimization involving only the orthogonal transforms:

[TABLE]

Here, $\mathbf{C}_{ij}\in\mathbb{R}^{d\times d}$, and $\mathbf{C}_{ji}=\mathbf{C}_{ij}^{\top}$. More specifically, O-REG is equivalent to REG in the following sense: If $(\mathbf{O}_{i}^{\ast})$, $(\textbf{t}_{i}^{\ast})$, $(\textbf{z}_{k}^{\ast})$ are globally optimal for REG, then

[TABLE]

Note that the objective in O-REG can be regarded as a quadratic form in $(\mathbf{O}_{i})$. That is, defining $\mathbf{O}=\left(\mathbf{O}_{1}\cdots\mathbf{O}_{M}\right)\in\mathbb{R}^{d\times Md}$ and $\mathbf{C}=\left((\mathbf{C}_{ij})_{1\leq i,j\leq M}\right)\in\mathbb{R}^{Md\times Md}$, we have

[TABLE]

Here, matrix $\mathbf{C}$ is symmetric, which follows from the fact that $\mathbf{C}_{ji}=\mathbf{C}_{ij}^{\top}$. Furthermore, $\mathbf{C}$ is positive semidefinite Chaudhury et al. (2015b); this is because the objective in O-REG is simply a rewriting of the objective in REG, which is always nonnegative, being the sum of norms.

We will now formulate O-REG as a rank-constrained semidefinite program. To this end, define Gram matrix $\mathbf{G}=\mathbf{O}^{\top}\mathbf{O}\in\mathbb{S}^{Md}_{+}$. Note that the $(i,j)$-th block of $\mathbf{G}$ is given by $[\mathbf{G}]_{ij}=\mathbf{O}_{i}^{\top}\mathbf{O}_{j}$. In particular, the diagonal blocks of $\mathbf{G}$ are the identity matrix, since $[\mathbf{G}]_{ii}=\mathbf{O}_{i}^{\top}\mathbf{O}_{i}=\mathbf{I}_{d}$. Following this, O-REG can be shown Chaudhury et al. (2015b) to be equivalent to

[TABLE]

REG-SDP is a standard semidefinite program Alizadeh et al. (1997), but with an additional rank constraint, which makes it nonconvex. An apparent drawback in reformulating O-REG as REG-SDP is that the number of variables of optimization have increased from $d\times Md$ to $Md\times Md$. On the plus side, though, the domain of optimization in REG-SDP is connected, which means that we can expect local optimization methods to converge to global optimum of REG-SDP. Moreover, the nonconvexity of the problem is isolated in the rank constraint (i.e. dropping the rank constraint from REG-SDP makes the problem convex), a fact we will exploit in our analysis.

1.2 Algorithm

In this paper, we theoretically analyze the iterative algorithm proposed in Miraj Ahmed and Chaudhury (2017) for solving the nonconvex program REG-SDP. This algorithm is based on ADMM (alternating direction method of multipliers), a method traditionally employed for convex programs, where it has strong convergence guarantees Boyd et al. (2011). Recently, however, there have been a spate of works where applying ADMM to nonconvex programs have yielded impressive results Boyd et al. (2011); Chartrand and Wohlberg (2013); Diamond et al. (2017); Hong et al. (2016); Miraj Ahmed and Chaudhury (2017). But the theory of convergence for nonconvex ADMM has not yet caught up with its convex counterpart. In particular, as we shall see when we discuss related works in Section 1.6, existing results on nonconvex ADMM do not apply to our algorithm.

To see how we can apply ADMM to REG-SDP, we first reformulate REG-SDP as

[TABLE]

where, $\Omega=\{\mathbf{X}\in\mathbb{S}^{Md}_{+}:\mathrm{rank}(\mathbf{X})\leq d\}$, and $\Theta=\{\mathbf{X}\in\mathbb{S}^{Md}:[\mathbf{X}]_{ii}=\mathbf{I}_{d}\}$. The equivalence of 5 to REG-SDP follows from the fact that the rank of any $\mathbf{H}\in\Theta$ is at least $d$, since its diagonal blocks are $\mathbf{I}_{d}$. Moreover, if $\mathbf{H}$ is feasible for (5), the rank of $\mathbf{H}$ cannot be more than $d$, since $\mathbf{H}=\mathbf{G}$ and $\mathbf{G}\in\Omega$. The splitting of the constraint into $\Omega$ and $\Theta$, with an additional linear (consensus) constraint, makes the problem algorithmically amenable to ADMM Boyd et al. (2011). For some fixed $\rho>0$, the augmented Lagrangian for (5) is defined to be

[TABLE]

where, $\bm{\Lambda}\in\mathbb{S}^{Md}$ is the dual variable for the constraint $\mathbf{G}-\mathbf{H}=\mathbf{0}$ Boyd et al. (2011). The ADMM algorithm in Miraj Ahmed and Chaudhury (2017), initialized with some $\mathbf{H}^{0}$ and $\bm{\Lambda}^{0}$, involves the following updates for $k\geq 0$:

[TABLE]

As observed in Miraj Ahmed and Chaudhury (2017), the first two sub-problems in REG-ADMM ($\mathbf{G}$-update and $\mathbf{H}$-update steps) can be expressed as matrix projections onto $\Theta$ and $\Omega$ respectively. Namely,

[TABLE]

where, $\Pi_{\Omega}(\cdot)$ and $\Pi_{\Theta}(\cdot)$ denotes projection onto $\Omega$ and $\Theta$. Importantly, these projections are computationally efficient. More specifically, projection on $\Theta$ can be computed simply by setting the diagonal blocks of the input matrix to $\mathbf{I}_{d}$. On the other hand, projection on $\Omega$ can be computed efficiently by computing its top $d$ eigenvalues and retaining the ones that are positive Miraj Ahmed and Chaudhury (2017). That is, if $\mathbf{A}=\lambda_{1}\bm{u}_{1}\bm{u}_{1}^{\top}+\cdots+\lambda_{Md}\bm{u}_{Md}\bm{u}_{Md}^{\top}$, where $\lambda_{1}\geq\cdots\geq\lambda_{Md}$ are the eigenvalues, and $\bm{u}_{1},\ldots,\bm{u}_{Md}$ are the corresponding (orthonormal) eigenvectors, then

[TABLE]

Empirically, we find that REG-ADMM has good convergence properties. In fact, it was shown in Miraj Ahmed and Chaudhury (2017) that its performance is comparable to state-of-the-art methods for the registration of three-dimensional multiview scans. Our aim in this paper is to study theoretical convergence of REG-ADMM.

1.3 Convex Relaxation of REG-SDP

In our analysis of REG-ADMM, we will be leveraging convergence properties of a closely related ADMM algorithm that solves a relaxed (convex) version of REG-SDP. Observe that the nonconvex nature of the rigid registration problem REG is isolated in the rank constraint in REG-SDP. As a result, by simply dropping the rank constraint, we obtain the following convex program, which we call the convex relaxation of REG-SDP:

[TABLE]

An ADMM algorithm to solve C-SDP was proposed in Sanyal et al. (2017a). Initialized with some $\mathbf{H}^{0}$ and $\bm{\Lambda}^{0}$, this involves the following updates for $k\geq 0$:

[TABLE]

The expression for the augmented Lagrangian $\mathcal{L}_{\rho}$ in C-ADMM is identical to that in (6). Note that C-ADMM is similar to REG-ADMM, with the only difference being the $\mathbf{G}$-update step. Namely, the minimization is over the nonconvex set $\Omega$ in REG-SDP, while it is over the closed convex set $\mathbb{S}^{Md}_{+}$ in C-SDP. This difference turns out to be crucial for convergence: C-ADMM converges to the global optimum of C-SDP with arbitrary initializations Sanyal et al. (2017a).

While C-ADMM enjoys strong convergence guarantees, it has the typical drawbacks of a convex relaxation. First, the rank of the global optimum of C-SDP is not guaranteed to be $d$, i.e., it might not even be feasible for REG-SDP. If the rank is greater than $d$, we have to “round” the solution of C-SDP to a rank-$d$ matrix Chaudhury et al. (2015b), which will generally be suboptimal for REG-SDP. Second, because the $\mathbf{G}$-update requires us to optimize over the entire PSD cone $\mathbb{S}^{Md}_{+}$, C-ADMM requires the full eigendecomposition of an $Md\times Md$ matrix at every iteration Sanyal et al. (2017a). REG-ADMM overcomes these drawbacks in the following manner. First, REG-ADMM works directly with the nonconvex problem REG-SDP, and thus obviates the need for any rounding. Second, REG-ADMM requires only the top $d$ eigenvectors in the $\mathbf{G}$-update step, resulting in appreciable computational savings. However, these benefits come at a cost—it is usually difficult to derive theoretical guarantees for nonconvex optimization. As far as we know, convergence guarantees for REG-ADMM do not follow from existing results on nonconvex ADMM (cf. Section 1.6 for further details).

1.4 Numerical Experiments

To understand the challenges involved in the convergence analysis of REG-ADMM, we look at some simulation results for the sensor network localization problem. We consider a two-dimensional network with ten nodes. There are three patches, where each patch contains all the ten nodes. The local coordinate system for each patch is obtained by arbitrarily rotating the global coordinate system. Moreover, we perturb the local coordinates using iid Gaussian noise. Our goal is to estimate the patch rotations up to a global transform. We set up REG-SDP for this problem and solve it using REG-ADMM. Fig. 2 shows the dependence of REG-ADMM on $\rho$ at different noise levels, where by “noise level” we mean the variance of the Gaussian noise. Notice that even when the local coordinates measurements are clean, REG-ADMM may get stuck in a local minimum depending on $\rho$. This dependence of the limit point on $\rho$ is observed both at low and high noise levels. This is in contrast with C-ADMM, where the iterates converge to a global minimum for any positive $\rho$ Sanyal et al. (2017a). Also observe that when the noise is relatively large (Fig. 2(c)), the iterates of the algorithm may oscillate without converging if $\rho$ is small. Such non-attenuating oscillations are not observed when the noise is low (see Fig. 2(b)).

1.5 Contribution

The foregoing simulation results provide the main motivation for this paper, namely, we wish to theoretically justify the observed behavior of REG-ADMM in different noise regimes. We first clarify what we mean by “noise”. As is clear from our discussion in Section 1.1, the data matrix $\mathbf{C}$ in REG-SDP ultimately depends on the local coordinate measurements $\textbf{x}_{k,i}$, $i\in[1:M]$, $k\in P_{i}$. If these measurements are exact, we say that $\mathbf{C}$ is clean; otherwise, we say that $\mathbf{C}$ is noisy. Furthermore, we make a distinction between low and high noise. In this context, we bring in the notion of tightness of the convex relaxation C-SDP. Recall that C-SDP is derived by dropping the rank constraint. Let $\mathbf{G}^{\ast}$ denote a global optimum of C-SDP. If $\mathrm{rank}(\mathbf{G}^{\ast})=d$, then clearly $\mathbf{G}^{\ast}$ is global optimum of REG-SDP as well. That is, global minimizer of the convex relaxation is a global minimizer for the original nonconvex program too. In this case, we say that the relaxation is tight Bandeira et al. (2017). Empirically, we notice the well-known phenomena of phase transition for convex relaxations (e.g. see Mixon (2015); Donoho et al. (2011)), where below a certain noise threshold, the relaxation remains tight (see Fig. 3). If the noise level is below (resp. above) this threshold, we say that the noise is low (resp. high). An informal account of our main findings is as follows:

We prove that any fixed point of REG-ADMM is a stationary (KKT) point of O-REG (Theorem 2.1). 2. 2.

We rigorously establish the existence of low noise regime. That is, we show that there is a noise threshold below which the convex relaxation C-SDP is tight (Theorem 2.2). 3. 3.

At low noise, we show that the REG-ADMM iterates converge to the global optimum, provided they are initialized sufficiently close to the optimum (Theorem 2.3). 4. 4.

If the data matrix is clean, then, for arbitrary initialization, we compute $\rho$ for which REG-ADMM converges to the global optimum of REG-SDP (Theorem 2.4). 5. 5.

At high noise, we give a duality-based explanation of why the iterates exhibit non-attenuating oscillations when $\rho$ is small, and why no such oscillations are observed at low noise.

The novelty of our analysis lies in the fact that we exploit the phenomenon of tightness of convex relaxation to prove convergence of the nonconvex ADMM. We contrast this approach with existing works on nonconvex ADMM (which rely on assumptions that do not apply to REG-ADMM) in the following subsection (Section 1.6).

1.6 Related Work

The rank-restricted subset $\Omega$ of the PSD cone is nonconvex, which implies that standard convergence result for ADMM Boyd et al. (2011) does not directly apply to REG-ADMM. However, we do leverage the convergence of convex ADMM for analyzing the convergence of REG-ADMM when the noise is low. A phase transition phenomena similar to the one cited above has been analyzed in Bandeira et al. (2017), albeit in the context of phase synchronization. Our proof of existence of the phase transition for C-SDP is in the spirit of this analysis.

The theoretical convergence of ADMM for nonconvex problems has been studied in Hong et al. (2016); Chang et al. (2016); Wang et al. (2015). However, a crucial working assumption common to these results does not hold in our case. More precisely, observe that we can rewrite REG-SDP as

[TABLE]

where $\iota_{\Gamma}$ is the indicator function associated with a feasible set $\Gamma$ Boyd et al. (2011), namely, $\iota_{\Gamma}(\mathbf{Y})=0$ if $\mathbf{Y}\in\Gamma$, and $\iota_{\Gamma}(\mathbf{Y})=\infty$ otherwise. Notice that, because of the indicator functions, the objective function in (8) is non-differentiable in both $\mathbf{G}$ and $\mathbf{H}$. This violates a regularity assumption common in existing analyses of nonconvex ADMM, namely, that the objective must be smooth in at least one variable. In these works, convergence results are obtained by proving a monotonic decrease in the augmented Lagrangian. This requires: (i) bounding successive difference in dual variables by successive difference in primal variable, which is where the assumption of smoothness is used; (ii) requiring that the parameter $\rho$ is above a certain threshold. In particular, it is not clear whether this thresholding of the value of $\rho$ is fundamental to convergence, or just an artifact of the analysis.

We do not make such smoothness assumptions in our analysis. We can afford to do this since we are analyzing a special class of problems, as opposed to the more general setups in Hong et al. (2016); Chang et al. (2016); Wang et al. (2015). Instead of showing a monotonic decrease in the augmented Lagrangian, our analysis relies on the phenomenon of tightness of convex relaxation. This provides more insights into the convergence behavior of the algorithm. For instance, our explanation in Section 2.6 shows that the instability of the algorithm (in the high noise regime) for low values of $\rho$ is fundamental, while suggesting why this instability is not observed in the low noise regime.

1.7 Organization

The rest of the paper is organized as follows. The main results are discussed in Section 2, which consists of six subsections. We conclude with a summary of our findings in Section 3. Appendix contains proofs of some subsidiary results that we use in our analysis.

2 Convergence Analysis

In this section we state and prove the main results of the paper. We set up notations and state some linear algebraic preliminaries in Section 2.1. We then prove some results on duality theory in Section 2.2 that we will use in subsequent analysis. In Section 2.3, we prove that any fixed point of REG-ADMM is a KKT point of O-REG. Sections 2.4 and 2.5 deal with low noise regime: in Section 2.4 we show the existence of, and explicitly compute, the noise threshold below which C-SDP is tight, thus theoretically establishing the existence of the low noise regime; in Section 2.5 we analyze convergence of REG-ADMM in this regime. Finally, in Section 2.6, we give a duality-based explanation of the instability of REG-ADMM in high noise regime for small values of $\rho$.

2.1 Notations and Preliminaries

For a matrix $\mathbf{X}$, $\lVert\cdot\rVert$ and $\lVert\cdot\rVert_{2}$ denote the Frobenius and spectral norms; the latter is simply the largest singular value $\sigma_{\mathrm{max}}(\mathbf{X})$. We note that $\lVert\mathbf{X}\rVert_{2}\leq\lVert\mathbf{X}\rVert$. $\mathrm{Tr}(\mathbf{A})$ denotes the trace of $\mathbf{A}$; and \big{\langle}\mathbf{X},\mathbf{Y}\big{\rangle}=\mathrm{Tr}(\mathbf{X}\mathbf{Y}) is the inner product between two symmetric matrices $\mathbf{X}$ and $\mathbf{Y}$. For a symmetric matrix $\mathbf{X}$, $\lambda_{i}(\mathbf{X})$ denotes the $i$-th eigenvalue, where we assume the eigenvalues $(\lambda_{i})$ to be arranged in a nondecreasing order.

Given a matrix $\mathbf{X}$, we refer to the subspace spanned by its columns as the rangespace of $\mathbf{X}$; $\mathrm{rank}(\mathbf{X})$ is the dimension of the rangespace of $\mathbf{X}$. The subspace spanned by the vectors $\mathbf{v}$ such that $\mathbf{X}\mathbf{v}=\mathbf{0}$ is referred to as the nullspace of $\mathbf{X}$; $\mathrm{nullity}(\mathbf{X})$ is the dimension of nullspace of $\mathbf{X}$.

Following are some of the linear algebraic facts that we shall use in our analysis:

Suppose $\mathbf{X},\mathbf{Y}\in\mathbb{R}^{n\times n}$, such that $\mathbf{X}\mathbf{Y}=\mathbf{0}$. If $\mathrm{rank}(\mathbf{Y})=d$, then $\mathrm{nullity}(\mathbf{X})\geq d$. Similarly, if $\mathrm{nullity}(\mathbf{X})=d$, then $\mathrm{rank}(\mathbf{Y})\leq d$. 2. 2.

Suppose $\mathbf{X},\mathbf{Y}\in\mathbb{S}^{n}_{+}$. Then, $\mathrm{Tr}(\mathbf{X}\mathbf{Y})=\mathbf{0}\iff\mathbf{X}\mathbf{Y}=\mathbf{0}$. 3. 3.

(Weyl’s Theorem Bhatia (2013)) Suppose $\mathbf{X},\mathbf{Y}\in\mathbb{S}^{n}$. Then for $i\in[1:n]$,

[TABLE] 4. 4.

Suppose $\mathbf{X},\mathbf{Y}\in\mathbb{S}^{n}$, such that $\mathbf{X}\mathbf{Y}=\mathbf{0}$. Then $\mathbf{X}$ and $\mathbf{Y}$ are simultaneously diagonalizable.

We now note in the following proposition some properties of positive semidefinite matrices with identity matrix on the diagonal blocks. The proof is deferred to the Appendix.

Proposition 1

Suppose $\mathbf{G}\in\mathbb{S}^{Md}_{+}$ and $[\mathbf{G}]_{ii}=\mathbf{I}_{d}$. Then

1. $(a)$

$\mathrm{rank}(\mathbf{G})\geq d$; 2. $(b)$

$\sigma_{\mathrm{max}}\left([\mathbf{G}]_{ij}\right)\leq 1,\>\>i,j\in[1:M]$; 3. $(c)$

If $\mathrm{rank}(\mathbf{G})=d$, then any nonzero eigenvalue of $\mathbf{G}$ is $M$.

2.2 Duality

We start by discussing some results on Lagrange duality, which plays an important role throughout our analysis. We use duality (i) for establishing KKT (Karush-Kuhn-Tucker) conditions at fixed points of REG-ADMM, (ii) for establishing the existence of a low noise regime, (iii) for establishing convergence in low noise regime, and (iv) for explaining the phenomenon of oscillation observed at high noise. We first derive the KKT conditions for O-REG by writing O-REG as a nonlinear program Bertsekas (1999):

[TABLE]

The Lagrangian for (9) is

[TABLE]

where the symmetric matrix $\bm{\Lambda}_{i}\in\mathbb{R}^{d\times d}$ is the multiplier for the $i$-th equality constraint in (9). The KKT point of (9) has the following characterization. (Recall that $\mathbf{G}$ is the Gram matrix of the $\mathbf{O}_{i}$’s, whose $(i,j)$-th block is $[\mathbf{G}]_{ij}=\mathbf{O}_{i}^{\top}\mathbf{O}_{j}$.)

Lemma 1

The variables $\mathbf{O}^{\ast}_{1},\ldots,\mathbf{O}^{\ast}_{M}\in\mathbb{R}^{d\times d}$ constitute a KKT point of (9) if and only if, for $i\in[1:M]$,

•

$(a)$* $[\mathbf{G}^{\ast}]_{ii}=\mathbf{I}_{d}$, and*

•

$(b)$* $[\mathbf{C}\mathbf{G}^{\ast}]_{ii}=[\mathbf{G}^{\ast}\mathbf{C}]_{ii}$.*

In this case, we say that $\mathbf{G}^{\ast}$ is a KKT point of O-REG.

Proof

For a minimization problem with equality constraints, KKT conditions amount to primal feasibility, and stationarity of Lagrangian (10) with respect to the primal variables Bertsekas (1999) (i.e., the partial derivative of (10) with respect to $\mathbf{O}_{i}$ should vanish). Primal feasibility gives us condition (a). Vanishing of partial derivative of (10) with respect to $\mathbf{O}_{i}$ gives us

[TABLE]

Left multiplying by $\mathbf{O}_{i}^{\ast\top}$ and using the primal feasibility condition that $\mathbf{O}_{i}^{\ast\top}\mathbf{O}_{i}^{\ast}=\mathbf{I}_{d}$, we obtain

[TABLE]

Also, note that $\bm{\Lambda}_{i}^{\ast\top}=-[\mathbf{C}\mathbf{G}^{\ast}]_{ii}$. Since $\bm{\Lambda}_{i}^{\ast}$ is symmetric, condition (b) follows immediately. Conversely, given conditions (a) and (b) on $\mathbf{G}^{\ast}$, it is not difficult to see that $\bm{\Lambda}_{i}$’s defined as $\bm{\Lambda}_{i}=-[\mathbf{C}\mathbf{G}^{\ast}]_{ii}$, and $\mathbf{O}_{i}$’s defined such that $\mathbf{O}_{j}^{\top}\mathbf{O}_{k}=[\mathbf{G}^{\ast}]_{jk}\ \forall j,k\in[1:M]$, together satisfy the KKT conditions. ∎

It is not difficult to check that the gradients of the constraints in (9) are linearly independent. As a result, the KKT conditions necessarily hold at any local minimum of (9) Bertsekas (1999). We now give the KKT conditions for the convex program C-SDP.

Lemma 2

$\mathbf{G}^{\ast}$* is a KKT point of C-SDP if there exists a block diagonal matrix $\bm{\Lambda}^{\ast}\in\mathbb{S}^{Md}$ such that*

•

$(a)$* $\mathbf{G}^{\ast}\in\mathbb{S}^{Md}_{+}$, $[\mathbf{G}^{\ast}]_{ii}=\mathbf{I}_{d}$,*

•

$(b)$* $\mathbf{C}+\bm{\Lambda}^{\ast}\in\mathbb{S}_{+}^{Md}$, and*

•

$(c)$* $(\mathbf{C}+\bm{\Lambda}^{\ast})\>\mathbf{G}^{\ast}=\mathbf{0}$.*

(We call $\bm{\Lambda}^{\ast}$ the dual variable corresponding to $\mathbf{G}^{\ast}$.)

Proof

KKT conditions for the convex program C-SDP are just the conditions of primal feasibility, dual feasibility, and complementary slackness Boyd and Vandenberghe (2004). Condition (a) is the primal feasibility condition, condition (b) is the dual feasibility condition, and condition (c) is the complementary slackness condition that follows from strong duality Alizadeh et al. (1997). Now, to see why $\bm{\Lambda}^{\ast}$ is block diagonal, we write C-SDP as a standard semidefinite program Alizadeh et al. (1997)

[TABLE]

Here $\mathbf{A}_{k}\in\mathbb{S}^{Md}$, and $\mathrm{Tr}(\mathbf{A}_{k}\mathbf{G})=\mathbf{b}_{k},k=[1:m]$ collectively encode the condition that $[\mathbf{G}]_{ii}=\mathbf{I}_{d}$. From Alizadeh et al. (1997), $\bm{\Lambda}^{\ast}$ is of the form $\sum_{k=1}^{m}y_{k}\mathbf{A}_{k}$, for some scalars $y_{k}$, and it is not difficult to see that $\sum_{k=1}^{m}y_{k}\mathbf{A}_{k}$ forms a symmetric block diagonal matrix in this case. ∎

Observe that the identity matrix $\mathbf{I}_{Md}$ is strictly feasible for C-SDP. Thus, by Slater’s condition Boyd and Vandenberghe (2004), the KKT conditions in Lemma 2 are both necessary and sufficient for $\mathbf{G}^{\ast}$ to be a global minimizer of C-SDP Boyd and Vandenberghe (2004). Note that we refer to $\bm{\Lambda}^{\ast}$ as the dual variable corresponding to $\mathbf{G}^{\ast}$: this is because, given a global minimizer $\mathbf{G}^{\ast}$, there is a unique block diagonal matrix $\bm{\Lambda}^{\ast}$ satisying conditions in Lemma 2. This is formalized in the following proposition, which expresses dual optimum in terms of the primal optimum. The proof is deferred to the Appendix.

Proposition 2

Let $\mathbf{G}^{\ast}$ and $\bm{\Lambda}^{\ast}$ be as in Lemma 2. Then

[TABLE]

Finally, we consider the dual variables generated by C-ADMM and REG-ADMM. For both these algorithms, we will assume that the dual initialization $\bm{\Lambda}^{0}$ is block diagonal.

Proposition 3

Let $\bm{\Lambda}^{0}$ be a block diagonal matrix in REG-ADMM (or C-ADMM). Then

1. (i)

$\bm{\Lambda}^{k}$* is block diagonal at every $k\geq 1$;* 2. (ii)

The $\mathbf{H}$-update step in (7) reduces to

[TABLE]

The proof is deferred to the Appendix.

2.3 Fixed Point Analysis

The first result we prove is that any fixed point of REG-ADMM is a stationary (KKT) point of O-REG. Before we do so, we clarify what we mean by a “fixed point”.

Definition 1

Suppose $\mathcal{T}$ is an iterative procedure initialized with $\textbf{z}^{0}$, and working as follows: $\textbf{z}^{k+1}=\mathcal{T}(\textbf{z}^{k})$. Then, $\textbf{z}^{k_{0}}$ is a fixed point of $\mathcal{T}$ if $\textbf{z}^{k_{0}}=\mathcal{T}(\textbf{z}^{k_{0}})$.

In other words, the iterative procedure stops updating the variables when the fixed point is reached. We now characterize fixed point of REG-ADMM.

Proposition 4

Suppose REG-ADMM is initialized with $(\mathbf{H}^{0},\bm{\Lambda}^{0})$, where $\bm{\Lambda}^{0}$ is a block diagonal matrix. Then $(\mathbf{H}^{k_{0}},\bm{\Lambda}^{k_{0}})$ is a fixed point of REG-ADMM if and only if $\mathbf{G}^{k_{0}+1}=\mathbf{H}^{k_{0}}$.

Proof

If $\mathbf{G}^{k_{0}+1}=\mathbf{H}^{k_{0}}$, then clearly $[\mathbf{G}^{k_{0}+1}]_{ii}=\mathbf{I}_{d}$. This means that $\mathbf{H}^{k_{0}+1}=\Pi_{\Theta}\left(\mathbf{G}^{k_{0}+1}\right)=\mathbf{G}^{k_{0}+1}=\mathbf{H}^{k_{0}}$. This, in turn, implies that $\bm{\Lambda}^{k_{0}+1}=\bm{\Lambda}^{k_{0}}+\rho(\mathbf{G}^{k_{0}+1}-\mathbf{H}^{k_{0}+1})=\bm{\Lambda}^{k_{0}}$. Thus, we have that $(\mathbf{H}^{k_{0}+1},\bm{\Lambda}^{k_{0}+1})=(\mathbf{H}^{k_{0}},\bm{\Lambda}^{k_{0}})$.

Conversely, if $(\mathbf{H}^{k_{0}},\bm{\Lambda}^{k_{0}})$ is a fixed point of REG-ADMM, then by definition $(\mathbf{H}^{k_{0}+1},\bm{\Lambda}^{k_{0}+1})=(\mathbf{H}^{k_{0}},\bm{\Lambda}^{k_{0}})$. Now, $\bm{\Lambda}^{k_{0}+1}=\bm{\Lambda}^{k_{0}}$ implies that $\mathbf{G}^{k_{0}+1}-\mathbf{H}^{k_{0}+1}=\mathbf{0}$. Thus, we have that $\mathbf{G}^{k_{0}+1}=\mathbf{H}^{k_{0}+1}=\mathbf{H}^{k_{0}}$. ∎

We now prove the following theorem, which, combined with Proposition 4, implies stationarity of any fixed point of REG-ADMM.

Theorem 2.1

Suppose REG-ADMM is initialized with $(\mathbf{H}^{0},\bm{\Lambda}^{0})$, where $\bm{\Lambda}^{0}$ is a block diagonal matrix. Suppose $\mathbf{G}^{k_{0}+1}=\mathbf{H}^{k_{0}}$. Then $\mathbf{H}^{k_{0}}$ is a KKT point of O-REG.

Proof

Recall that $\mathbf{G}^{k_{0}+1}$ is formed from top-$d$ eigendecomposition of $\mathbf{H}^{k_{0}}-\rho^{-1}\left(\mathbf{C}+\bm{\Lambda}^{k_{0}}\right)$. This implies that the rangespace of $\left(\mathbf{H}^{k_{0}}-\rho^{-1}\left(\mathbf{C}+\bm{\Lambda}^{k_{0}}\right)-\mathbf{G}^{k_{0}+1}\right)$ is orthogonal to rangespace of $\mathbf{G}^{k_{0}+1}$, which gives us

[TABLE]

Or,

[TABLE]

Since $\mathbf{G}^{k_{0}+1}=\mathbf{H}^{k_{0}}$, the left hand side of (15) is $\mathbf{0}$, and we have

[TABLE]

Observe that $[\bm{\Lambda}^{k_{0}}\mathbf{H}^{k_{0}}]_{ii}$ is symmetric (since $\bm{\Lambda}^{k_{0}}$ is block diagonal, and $[\mathbf{H}^{k_{0}}]_{ii}=\mathbf{I}_{d}$). Thus, we have that $[\mathbf{C}\mathbf{H}^{k_{0}}]_{ii}$ is symmetric. Using Lemma 1, we can now conclude that $\mathbf{H}^{k_{0}}$ is a KKT point of O-REG. ∎

2.4 Existence of Low Noise Regime: Tightness of C-SDP

In Section 2.5, we will derive convergence results for REG-ADMM in the low noise regime. Recall from our discussion in Section 1.5, that we say we are in the low noise regime when the convex relaxation C-SDP is tight; that is, when global minimizer of C-SDP is also a global minimizer of REG-SDP. Before we prove convergence results for REG-ADMM in the low noise regime, it is imperative that we theoretically establish the empirical observation (Fig. 3) that there is indeed a noise threshold below which C-SDP is tight. This is the task for this section. Recall that the data matrix $\mathbf{C}$ in REG-SDP ultimately depends on the local coordinate measurements $\textbf{x}_{k,i}$, $i\in[1:M]$, $k\in P_{i}$. If these measurements are exact, we say that $\mathbf{C}$ is clean; otherwise, we say that $\mathbf{C}$ is noisy.

Our analysis for tightness of REG-SDP relies on the assumption that nullity of the clean data matrix is $d$. Before stating the main results of this section, we examine some properties of the clean data matrix that will help demystify this assumption. To do so, we rely on our discussion in Section 1.1. Let $\bar{\textbf{x}}_{1},\ldots,\bar{\textbf{x}}_{N}\in\mathbb{R}^{d}$ be the true coordinates of the $N$ nodes in some global coordinate system. Let the local coordinate system associated with patch $P_{i}$ be related to the global coordinate system via rigid transform $(\bar{\mathbf{O}}_{i},\bar{\textbf{t}}_{i})$, $i\in[1:M]$. It follows from equation (1) that

[TABLE]

Let $\mathbf{C}_{0}$ be the clean data matrix formed from exact local coordinate measurements. Define $\mathbf{G}_{0}$ be defined by

[TABLE]

Then, from (3), and from the fact that REG-SDP is an exact reformulation of O-REG, we get that

[TABLE]

Now, since $\mathbf{C}_{0},\mathbf{G}_{0}\in\mathbb{S}^{Md}_{+}$, we get that

[TABLE]

From (20) and the fact that $\mathrm{rank}(\mathbf{G}_{0})=d$, it follows that $\mathrm{nullity}(\mathbf{C}_{0})\geq d$. The assumption we make on $\mathbf{C}_{0}$ is that $\mathrm{nullity}(\mathbf{C}_{0})=d$. This is equivalent to the condition that the body graph corresponding to the registration problem is affinely rigid Chaudhury et al. (2015b). Loosely speaking, the registration problem is well-posed if $\mathrm{nullity}(\mathbf{C}_{0})=d$. For example, this is true when each patch contains all the $N$ points (the setup for the simulations in Section 1.4). More generally, $\mathrm{nullity}(\mathbf{C}_{0})=d$ if there exists an ordering of the patches such that $P_{1}$ contains at least $d+1$ points, and $P_{i}$ has at least $d+1$ points in common with $P_{1}\cup P_{2}\cup\cdots\cup P_{i-1}$, for $i\geq 2$ Chaudhury et al. (2015b) (this occurs naturally in applications like multiview registration Miraj Ahmed and Chaudhury (2017)).

We now state a result that we require to prove the main result of this section (Theorem 2.2). This result is also of independent interest in the context of stability of C-SDP. In particular, the lemma quantifies the perturbation in the optimum of C-SDP due to a perturbation in $\mathbf{C}_{0}$. We defer the proof to the Appendix.

Lemma 3

Suppose $\mathrm{nullity}(\mathbf{C}_{0})=d$. Let $\mathbf{C}$ be the corresponding noisy data matrix, and let $\mathbf{W}=\mathbf{C}-\mathbf{C}_{0}$. Let $\mathbf{G}^{\ast}$ be a global optimum for C-SDP corresponding to the noisy data $\mathbf{C}$. Let $\Delta=\mathbf{G}^{\ast}-\mathbf{G}_{0}$, where $\mathbf{G}_{0}$ is the ground-truth solution. Then

[TABLE]

We now state the main result of this section concerning tightness of C-SDP.

Theorem 2.2

Suppose $\mathrm{nullity}(\mathbf{C}_{0})=d$. Let $\mathbf{C}=\mathbf{C}_{0}+\mathbf{W}$ be the noisy data matrix. Then there exists (explicitly computable) $\eta>0$ such that $\mathrm{rank}(\mathbf{G}^{\ast})=d$ (i.e. C-SDP is tight) if $\lVert\mathbf{W}\rVert<\eta$.

Proof

Let $\bm{\Lambda}^{\ast}$ be the dual corresponding to $\mathbf{G}^{\ast}$. From the fact that $(\mathbf{C}+\bm{\Lambda}^{\ast})\mathbf{G}^{\ast}=\mathbf{0}$ (Lemma 2(b)), and that $\mathrm{rank}(\mathbf{G}^{\ast})\geq d$ (Proposition 1(a)), we have that $\mathrm{nullity}(\mathbf{C}+\bm{\Lambda}^{\ast})\geq d$. This also tells us the following: if we prove $\mathrm{nullity}(\mathbf{C}+\bm{\Lambda}^{\ast})=d$, then it necessarily follows that $\mathrm{rank}(\mathbf{G}^{\ast})=d$, implying that C-SDP is tight. Thus, our goal is to prove that there exists $\eta>0$ such that $\mathrm{nullity}(\mathbf{C}+\bm{\Lambda}^{\ast})=d$ if $\lVert\mathbf{W}\rVert<\eta$.

From Proposition 2, we have that

[TABLE]

Here $\mathrm{bd}:\mathbb{R}^{Md\times Md}\rightarrow\mathbb{R}^{Md\times Md}$ is the linear operator that leaves the diagonal blocks untouched, and sets every other element to [math]. To get the last equality, we have added a superfluous term $\mathrm{bd}(\mathbf{C}_{0}\mathbf{G}_{0})$. We can do this because $\mathbf{C}_{0}\mathbf{G}_{0}=\mathbf{0}$, see (20).

Thus, in summary, $\mathbf{C}+\bm{\Lambda}^{\ast}=\mathbf{C}_{0}+\mathbf{A}$, where $\mathbf{A}=\mathbf{W}-\mathrm{bd}(\mathbf{W}\mathbf{G}^{\ast})-\mathrm{bd}(\mathbf{C}_{0}\bm{\Delta})$, and $\bm{\Delta}=\mathbf{G}^{\ast}-\mathbf{G}_{0}$ is the perturbation in the optimum of C-SDP due to pertubation in the data matrix. Now, by Weyl’s theorem,

[TABLE]

We have already established that $\mathrm{nullity}(\mathbf{C}+\bm{\Lambda}^{\ast})\geq d$. This implies that the bottom $d$ eigenvalues of $(\mathbf{C}+\bm{\Lambda}^{\ast})$ are guaranteed to be [math], i.e., $\lambda_{1}(\mathbf{C}+\bm{\Lambda}^{\ast})=\cdots=\lambda_{d}(\mathbf{C}+\bm{\Lambda}^{\ast})=0$. To prove $\mathrm{nullity}(\mathbf{C}+\bm{\Lambda}^{\ast})=d$, we need to show that $\lambda_{d+1}(\mathbf{C}+\bm{\Lambda}^{\ast})>0$. We can do this using (21), if we upperbound $\lVert\mathbf{A}\rVert_{2}$.

[TABLE]

where we have used the results that $\lVert\cdot\rVert_{2}\leq\lVert\cdot\rVert$, that $\lVert\cdot\rVert$ obeys triangle inequality, that $\lVert\mathrm{bd}(\cdot)\rVert\leq\lVert\cdot\rVert$, and that $\lVert\cdot\rVert$ is sub-multiplicative Bhatia (2013). We now need an upper bound on $\lVert\bm{\Delta}\rVert$, which is the magnitude of perturbation in the optimum of C-SDP when the data matrix is perturbed from $\mathbf{C}_{0}$ to $\mathbf{C}=\mathbf{C}_{0}+\mathbf{W}$. This is where we use Lemma 3, which gives us

[TABLE]

Thus, we have

[TABLE]

where, the fact that $\lVert\mathbf{G}^{\ast}\rVert\leq M\sqrt{d}$ follows from Proposition 1(b). The terms in the paranthesis are constant. Thus, by controlling $\lVert\mathbf{W}\rVert$, we can control $\lVert\mathbf{A}\rVert_{2}$. In particular, if

[TABLE]

then $\lambda_{d+1}(\mathbf{C}_{0})-\lVert\mathbf{A}\rVert_{2}>0$. This, from equation (21), implies that $\lambda_{d+1}(\mathbf{C}+\bm{\Lambda}^{\ast})>0$. ∎

When the noise is low enough for the convex relaxation C-SDP to be tight, we say that we are in the low noise regime. Theorem 2.2 establishes the validity of this notion; that is, we do have a non-zero noise level $\eta$, given by (22), below which C-SDP is tight.

2.5 Convergence in Low Noise Regime

We start with a result on C-ADMM, which we will leverage to get results for REG-ADMM in the low noise regime. Suppose $\mathbf{A}^{k}=\mathbf{H}^{k}-\rho^{-1}(\mathbf{C}+\bm{\Lambda}^{k})$, $\mathbf{A}^{*}=\mathbf{H}^{*}-\rho^{-1}(\mathbf{C}+\bm{\Lambda}^{*})$, where $\mathbf{H}^{\ast}$ is a global minimum for the convex relaxation C-SDP, and $\bm{\Lambda}^{\ast}$ is the corresponding optimal dual as in Lemma 2.

Lemma 4

For C-ADMM, $\lVert\mathbf{A}^{k+1}-\mathbf{A}^{*}\rVert^{2}\leq\lVert\mathbf{A}^{k}-\mathbf{A}^{*}\rVert^{2}$, for every $k\geq 0$.

Proof of Lemma 4 is deferred to the Appendix. We now state our result on convergence of REG-ADMM in the low noise regime ($\mathbf{W}$ and $\eta$ are as in Theorem 2.2).

Theorem 2.3

Suppose $\lVert\mathbf{W}\rVert<\eta$, which implies that the convex relaxation is tight, i.e., $\mathrm{rank}(\mathbf{G}^{\ast})=d$, where $\mathbf{G}^{\ast}$ is a global optimum for C-SDP. Let $\bm{\Lambda}^{\ast}$ be the corresponding optimum dual variable as in Lemma 2. Suppose $\mathbf{H}^{0}$ and $\bm{\Lambda}^{0}$ are such that

[TABLE]

Then REG-ADMM converges to the global optimum, that is, $\mathbf{H}^{k}\rightarrow\mathbf{G}^{\ast}$ and $\bm{\Lambda}^{k}\rightarrow\bm{\Lambda}^{\ast}$ as $k\rightarrow\infty$.

Proof

Since $\lVert\mathbf{W}\rVert<\eta$, we have that $\mathrm{nullity}(\mathbf{C}+\bm{\Lambda}^{\ast})=d$ (see proof of Theorem 2.2), which, in particular implies that $\lambda_{d+1}(\mathbf{C}+\bm{\Lambda}^{\ast})>0$. We now deduce the eigenvalues of $\mathbf{A}^{\ast}$. From KKT condition (a) in Lemma 2, we have that

[TABLE]

In particular, this means that $\mathbf{C}+\bm{\Lambda}^{\ast}$, and $\mathbf{G}^{\ast}$ are simultaneously diagonalizable, and since their product is $\mathbf{0}$, their rangespaces are orthogonal. Now, since $\mathrm{rank}(\mathbf{G}^{\ast})=d$, the non-zero eigenvalues of $\mathbf{G}^{\ast}$ are all $M$ (Proposition 1(c)). Moreover, from KKT condition (b) in Lemma 2, we have that $\mathbf{C}+\bm{\Lambda}^{\ast}\in\mathbb{S}^{Md}_{+}$. Putting everything together, we have that

• •

$\mathbf{A}^{\ast}$ has $d$ positive eigenvalues, each of them equal to $M$.

• •

$\mathbf{A}^{\ast}$ has $(M-1)d$ negative eigenvalues, which are the negative of the non-zero eigenvalues of $\rho^{-1}(\mathbf{C}+\bm{\Lambda}^{\ast})$.

Let

[TABLE]

By Weyl’s theorem, we deduce that for any $\mathbf{A}\in\mathcal{B}^{\ast}$, the top $d$ eigenvalues of $\mathbf{A}$ would lie in the interval $\left(M-\rho^{-1}\lambda_{d+1}(\mathbf{C}+\bm{\Lambda}^{\ast}),M+\rho^{-1}\lambda_{d+1}(\mathbf{C}+\bm{\Lambda}^{\ast})\right)$, and other eigenvalues of $\mathbf{A}$ would lie in the interval \big{(}-\rho^{-1}\left(\lambda_{Md}(\mathbf{C}+\bm{\Lambda}^{\ast})+\lambda_{d+1}(\mathbf{C}+\bm{\Lambda}^{\ast})\right),0\big{)}. In other words, for any $\mathbf{A}\in\mathcal{B}^{\ast}$, only the top $d$ eigenvalues of $\mathbf{A}$ would be nonnegative. Let

[TABLE]

Suppose $\mathbf{A}^{k}\in\mathcal{B}^{\ast}$. Then, because only the top $d$ eigenvalues of $\mathbf{A}^{k}$ are nonnegative, we have

[TABLE]

That is, projection on the nonconvex set $\Omega$ is same as projection on the convex set $\mathbb{S}^{Md}_{+}$.

Now, from Lemma 4 we infer that, if $\mathbf{A}^{k}\in\mathcal{B}^{\ast}$, then $\mathbf{A}^{k+1}\in\mathcal{B}^{\ast}$. Thus, every subsequent projection on $\Omega$ is equivalent to projection on $\mathbb{S}^{Md}_{+}$. Note that condition (23) in the hypothesis of the theorem implies that $\mathbf{A}^{0}\in\mathcal{B}^{\ast}$. Thus, the iterates generated by REG-ADMM initialized with $\mathbf{H}_{0}$, $\bm{\Lambda}_{0}$ is the same as the iterates generated by C-ADMM initialized with $\mathbf{H}_{0}$, $\bm{\Lambda}_{0}$. Now, since C-ADMM converges to global optimum Sanyal et al. (2017a), we deduce that REG-ADMM converges to global optimum. ∎

In summary, when the conditions in Theorem 2.3 are satisfied, REG-ADMM and C-ADMM generate the same iterates, which allows us to infer the convergence of REG-ADMM from that of C-ADMM. Thus, we reap the convergence benefits of C-ADMM, while retaining the computational advantages of REG-ADMM (see discussion in Section 1.3).

Observe that Theorem 2.3 implies a tradeoff, namely, if we intialize the dual sufficiently close to the optimal, we can be lax with the primal initialization. This principle is brought to the fore in the clean case, where we know that the dual optimum $\bm{\Lambda}^{\ast}=\mathbf{0}$ (this is because, in the clean case, $\mathbf{G}_{0}$ is the primal global optimum, and $\mathbf{C}_{0}\mathbf{G}_{0}=\mathbf{0}$; see (20) and Proposition 2).

Theorem 2.4

Let $\mathrm{nullity}(\mathbf{C}_{0})=d$ and $\bm{\Lambda}^{0}=\mathbf{0}$. Then given any primal initialization $\mathbf{H}^{0}$, REG-ADMM converges to global optimum provided

[TABLE]

Proof

For the clean case, $\bm{\Lambda}^{\ast}=\mathbf{0}$. Thus, with $\bm{\Lambda}^{0}=\mathbf{0}$, the condition in Theorem 2.3 reduces to

[TABLE]

Now,

[TABLE]

Since $\lVert\mathbf{G}^{\ast}\rVert=M\sqrt{d}$, we have

[TABLE]

Thus, the algorithm would converge to the global optimum if

[TABLE]

which proves the result. ∎

In Fig. 2(a), REG-ADMM gets stuck in a local optimum for $\rho=1$, but converges to the global optimum if $\rho$ is set using Theorem 2.4. The result also sheds light on the robustness of the algorithm to primal initializations, as observed in Miraj Ahmed and Chaudhury (2017). Note that we have an upper bound on $\rho$, in contrast to existing results in the literature which prescribe a lower bound. This phenomenon is peculiar to low noise data. In contrast, small values of $\rho$ in the high noise regime may lead to non-attenuating oscillations (Fig. 2(c)).

2.6 Oscillations in High Noise Regime

Before discussing oscillations in the high-noise regime, we will derive a necessary condition that a fixed-point of REG-ADMM must satisfy.

Proposition 5

Suppose $(\mathbf{H}^{\ast},\bm{\Lambda}^{\ast})$ is a fixed-point of REG-ADMM algorithm with parameter $\rho$. Then $(\mathbf{C}+\bm{\Lambda}^{\ast})\mathbf{H}^{\ast}=\mathbf{0}$.

Proof

Suppose we initialize REG-ADMM with $\mathbf{H}^{0}=\mathbf{H}^{\ast}$, $\bm{\Lambda}^{0}=\bm{\Lambda}^{\ast}$. Then, by Proposition 4, $\mathbf{G}^{1}=\mathbf{H}^{0}$. Now, since $\mathbf{G}^{1}$ is formed from top-$d$ eigendecomposition of $\mathbf{H}^{0}-\rho^{-1}\left(\mathbf{C}+\bm{\Lambda}^{0}\right)$ (see (7)), the range space of $\left(\mathbf{H}^{0}-\rho^{-1}\left(\mathbf{C}+\bm{\Lambda}^{0}\right)-\mathbf{G}^{1}\right)$ is orthogonal to range space of $\mathbf{G}^{1}$, which gives us

[TABLE]

Or,

[TABLE]

Since $\mathbf{G}^{1}=\mathbf{H}^{0}\ (=\mathbf{H}^{\ast})$, the right hand side of the preceding equation is $\mathbf{0}$, and we obtain the desired result. ∎

Let us now focus on the high-noise regime. Suppose REG-ADMM is initialized with $(\mathbf{H}^{0},\bm{\Lambda}^{0})$ such that $\mathbf{H}^{0}$ is feasible for REG-SDP, $\bm{\Lambda}^{0}$ is block-diagonal, and $(\mathbf{C}+\bm{\Lambda}^{0})\mathbf{H}^{0}=\mathbf{0}$. In the high-noise regime (i.e., when the convex relaxation is not tight), $\mathbf{C}+\bm{\Lambda}^{0}$ must have at least one negative eigenvalue (say $-\mu^{2}$). Indeed, if $\mathbf{C}+\bm{\Lambda}^{0}$ were positive semidefinite, it would follow from Lemma 2 that the convex relaxation is tight. Now, the $\mathbf{G}$-update is

[TABLE]

From Proposition 1, we know that the $d$ non-zero eigenvalues of $\mathbf{H}^{0}$ are $M$. Now, observe that $-\rho^{-1}(\mathbf{C}+\bm{\Lambda}^{0})$ has a positive eigenvalue $\rho^{-1}\mu^{2}$. If $\rho$ is sufficiently small, we would have $\rho^{-1}\mu^{2}>M$. Then, since $\mathbf{G}^{1}$ is determined by the top $d$ eigenvalues of $\mathbf{H}^{0}-\rho^{-1}(\mathbf{C}+\bm{\Lambda}^{0})$, we would have that the top eigenvalue of $\mathbf{G}^{1}$ is strictly bigger than $M$, making $\mathbf{G}^{1}$ infeasible for REG-SDP (since, by Proposition 1, any $\mathbf{G}$ feasible for REG-SDP necessarily has all non-zero eigenvalues equal to $M$). This would imply that $\mathbf{H}^{1}\neq\mathbf{G}^{1}$, and consequently, that $\bm{\Lambda}^{1}\neq\bm{\Lambda}^{0}$. That is, for small value of $\rho$, we see that REG-ADMM does not stabilize even when $\mathbf{H}^{0}$ and $\bm{\Lambda}^{0}$ satisfy $(\mathbf{C}+\bm{\Lambda}^{0})\mathbf{H}^{0}=\mathbf{0}$, a property that any candidate for a fixed-point of REG-ADMM must satisfy. Put differently, there can be no fixed-point to which REG-ADMM converges if $\rho$ is sufficiently small.

Observe that the argument above depended on the existence of a negative eigenvalue of $\mathbf{C}+\bm{\Lambda}^{0}$. This argument does not hold in the low-noise regime because we can simultaneously have the properties that $(\mathbf{C}+\bm{\Lambda}^{0})\mathbf{H}^{0}=\mathbf{0}$, and that all the eigenvalues of $\mathbf{C}+\bm{\Lambda}^{0}$ are nonnegative (which holds when $\mathbf{H}^{0}$ is global optimum; see condition (b) in Lemma 2). This suggests why the instability is not observed for low values of $\rho$ in the low-noise regime.

3 Conclusion

In this work, we analyzed the convergence behavior of REG-ADMM in different noise regimes. Existing results on nonconvex ADMM do not apply to REG-ADMM as they rely on certain smoothness assumptions that are not satisfied by REG-SDP. We bypassed these assumptions, and exploited the tightness phenomenon of convex relaxation to guide our analysis. We started with a fixed point analysis, where we proved that any fixed point of REG-ADMM converges to a stationary (KKT) point of REG-SDP. To further refine the result, we looked at the behavior of REG-ADMM when the noise is low. In particular, we precisely defined what is meant by “low” noise by establishing tightness of the convex relaxation C-SDP below a certain noise threshold. We then proved that, by initializing the primal and dual variables sufficiently close to the optimum, the iterates of REG-ADMM are guaranteed to converge to the global optimum. By applying this result to the clean case, we showed that given any primal initialization, we can explicitly compute values of $\rho$ for which the algorithm converges to the global optimum. For high noise, we showed that for sufficiently small $\rho$, the iterates generated by REG-ADMM do not stabilize, even if initialization of REG-ADMM satisfies necessary property of a fixed-point.

4 Appendix

4.1 Proof of Proposition 1

$(a)$ Consider the first $d$ columns of $\mathbf{G}$. Since $\mathbf{G}_{11}=\mathbf{I}_{d}$, we conclude that the first $d$ columns of $\mathbf{G}$ are linearly independent.

$(b)$ Let $\textbf{x}_{i}$, $\textbf{x}_{j}\in\mathbb{R}^{d}$, $\lVert\textbf{x}_{i}\rVert=\lVert\textbf{x}_{j}\rVert=1$ . Regard any vector in $\mathbb{R}^{Md}$ as consisting of $M$ vectors in $\mathbb{R}^{d}$ stacked vertically. Construct $\textbf{x}\in\mathbb{R}^{Md}$ with $i$-th block as $\textbf{x}_{i}$, $j$-th block as $-\textbf{x}_{j}$, and every other block as $\mathbf{0}$. Then

[TABLE]

where we have used the fact that $[\mathbf{G}]_{ji}=[\mathbf{G}]_{ij}^{\top}$. Since $\mathbf{G}\in\mathbb{S}^{Md}_{+}$, we get that

[TABLE]

This gives us

[TABLE]

Similarly, by replacing $-\textbf{x}_{j}$ with $\textbf{x}_{j}$ in the $j$-th block of x, we get

[TABLE]

Putting these together, we have

[TABLE]

Unit vectors $\textbf{x}_{i}$, $\textbf{x}_{j}\in\mathbb{R}^{Md}$ were arbitrary, and thus, the result follows.

$(c)$ Consider the spectral decomposition of $\mathbf{G}$,

[TABLE]

where $\alpha_{i}>0$ are the non-zero eigenvalues, and $\mathbf{v}_{i}$ are the corresponding orthonormal eigenvectors. Let

[TABLE]

Regard $\mathbf{B}$ as a block-row, where each element is of size $d\times d$. Notationally, $\mathbf{B}=[\mathbf{B}_{1}\cdots\mathbf{B}_{M}]$, where $\mathbf{B}_{i}\in\mathbb{R}^{d\times d}$. Thus,

[TABLE]

In particular, since $[\mathbf{G}]_{ii}=\mathbf{I}_{d}$, we get that

[TABLE]

Now, suppose $\mathbf{v}_{1}=(\mathbf{v}_{11},\cdots,\mathbf{v}_{1M})$, where $\mathbf{v}_{1j}\in\mathbb{R}^{d}$. Note that $\sqrt{\alpha_{1}}\mathbf{v}_{1j}^{\top}$ forms the first row of $\mathbf{B}_{j}$. From orthogonality of $\mathbf{B}_{j}$, we get that, for every $j$,

[TABLE]

Or,

[TABLE]

Now, since $\lVert\mathbf{v}_{1}\rVert^{2}=1$, we have

[TABLE]

from where we get that $\alpha_{1}=M$. Similarly, $\alpha_{2}=\cdots=\alpha_{d}=M$.

4.2 Proof of Proposition 2

Lemma 2 tells us that $(\mathbf{C}+\bm{\Lambda}^{\ast})\>\mathbf{G}^{\ast}=\mathbf{0}$. That is,

[TABLE]

Lemma 2 also tells us that $\bm{\Lambda}^{\ast}$ is a symmetric block-diagonal matrix, and that $[\mathbf{G}^{\ast}]_{ii}=\mathbf{I}_{d}$. Using these facts, and comparing the diagonal blocks of the left-hand side and the right-hand side, we get

[TABLE]

Thus,

[TABLE]

where, $\mathrm{bd}$ is the linear operator that leaves the diagonal blocks untouched and sets other elements to [math].

4.3 Proof of Proposition 3

We prove the proposition for REG-ADMM by induction. The proof for C-ADMM is exactly the same since the $\mathbf{H}$-update and $\bm{\Lambda}$-update steps are identical for REG-ADMM and C-ADMM. Clearly, the proposition holds for $k=0$. Assume that the proposition holds for $k=k_{0}$. We will show that it then has to hold for $k=k_{0}+1$. Consider the $\mathbf{H}$-update step in (7)

[TABLE]

We know that $\Theta$ is the set of symmetric matrices for which $d\times d$ diagonal blocks are $\mathbf{I}_{d}$. It is clear that projection of a matrix on $\Theta$ is obtained by setting the diagonal blocks of the matrix to $\mathbf{I}_{d}$. In other words, $\Pi_{\Theta}(\cdot)$ affects only the diagonal blocks of its argument. By induction hypothesis, $\bm{\Lambda}^{k_{0}}$ is block diagonal, and thus adding it to $\mathbf{G}^{k}_{0}$ in the $\mathbf{H}$-update step does not affect the projection on $\Theta$, or,

[TABLE]

Now, consider the $\bm{\Lambda}$-update step,

[TABLE]

Since $\Pi_{\Theta}(\cdot)$ affects only the diagonal blocks of its arguments, it is clear that the off-diagonal blocks of $\mathbf{H}^{k_{0}+1}$ and $\mathbf{G}^{k_{0}+1}$ are the same. Thus, $\mathbf{G}^{k_{0}+1}-\mathbf{H}^{k_{0}+1}$ is a block diagonal matrix, which along with the hypothesis that $\bm{\Lambda}^{k_{0}}$ is block diagonal, implies that $\bm{\Lambda}^{k_{0}+1}$ is block diagonal.

4.4 Proof of Lemma 3

Since $\mathrm{rank}(\mathbf{G}_{0})=d$ we can write the spectral decomposition of the ground truth solution as (see Proposition 1(c))

[TABLE]

where $\mathbf{s}_{i},\ i\in[1:d]$ are orthonormal. Similarly, we can write the optimal solution corresponding to the perturbed data $\mathbf{C}$ as

[TABLE]

where $\mathbf{g}_{i},\ i\in[1:d]$ are orthonormal.

Let $\mathbf{P}$ be the orthoprojector on $\mathrm{nullspace}(\mathbf{C}_{0})$. Since $\mathbf{C}_{0}\mathbf{G}_{0}=\mathbf{0}$, $\mathrm{rank}(\mathbf{G}_{0})=d$, and $\mathrm{nullity}(\mathbf{C}_{0})=d$, we deduce that

[TABLE]

Clearly, $\mathbf{R}=\mathbf{I}-\mathbf{P}$ is the orthoprojector on $\mathrm{range}(\mathbf{C}_{0})$, and

[TABLE]

Let $\mathbf{h}_{i}:=\mathbf{R}\mathbf{g}_{i}$. As will be apparent, we want to lowerbound \big{\langle}\mathbf{C}_{0},\mathbf{G}^{\ast}\big{\rangle}. Now,

[TABLE]

Thus, to get a lowerbound on \big{\langle}\mathbf{C}_{0},\mathbf{G}^{\ast}\big{\rangle}, we just have to lowerbound $\sum_{i=1}^{Md}\alpha_{i}\,\lVert\mathbf{h}_{i}\rVert^{2}$. To do this, we first reformulate $\lVert\mathbf{h}_{i}\rVert^{2}$ as,

[TABLE]

where the second equality holds because $\mathbf{R}^{2}=\mathbf{R}$, since $\mathbf{R}$ is an orthoprojector. Thus,

[TABLE]

Now, to upperbound \big{\langle}\mathbf{G}_{0},\mathbf{G}^{\ast}\big{\rangle}, we note that

[TABLE]

where, to get the first inequality, we use the fact that any singular value of $[\mathbf{G}^{\ast}]_{ij}$ is at most $1$ for every $1\leq i,j\leq M$ (Proposition 1). Combining (26) and (27), we get

[TABLE]

Plugging this in (25) gives us,

[TABLE]

Now,

[TABLE]

where the first inequality is Cauchy-Schwarz (with inner product between two symmetric matrices defined as the trace of their product), and the last inequality is due to the fact that the term in first paranthesis is negative by optimality of $\mathbf{G}^{\ast}$, and that the first term in second paranthesis is [math]. So,

[TABLE]

which finally gives us,

[TABLE]

4.5 Proof of Lemma 4

The proof of this lemma essentially follows the convergence proof for convex ADMM presented in Boyd et al. (2011). Note that

[TABLE]

where we have used the fact that, for every $k$,

[TABLE]

because $\mathbf{H}^{k}-\mathbf{H}^{\ast}$ has $\mathbf{0}$ on the diagonal blocks, and $\bm{\Lambda}^{k}-\bm{\Lambda}^{\ast}$ is block diagonal. Further algebraic manipulations give us

[TABLE]

Observe that the first two terms in (29) are always nonnegative. Nonnegativity of the third term in (29) follows from the convex projection property, i.e.,

[TABLE]

To see why the fourth term in (29) is nonnegative, observe that the Lagrangian for the convex program C-SDP is given by

[TABLE]

We have already seen in the proof of Lemma 2 that strong duality holds for C-SDP. Thus, from the saddle point property of the Lagrangian at global optimum Boyd and Vandenberghe (2004), we get

[TABLE]

which gives us

[TABLE]