2018-12 文献阅读：Computing the Best Approximation ... Doubly Nonnegative Cone(Sun2018)——半光滑牛顿法SNCG

文献：Computing the Best Approximation Over the Intersection of a Polyhedral Set and the Doubly Nonnegative Cone(2018)

---

作者：Ying Cui，Defeng Sun，Kim-Chuan Toh
阅读原因：机器学习讨论班

0.摘要

本文介绍了一种高效的算法——imABCD算法，来计算在线性等式、不等式（多面体集）和双重非负锥（半正定矩阵且元素非负）的交集上，给定矩阵的最佳逼近问题。此算法并不像BCD算法一样，而是用2-block CD来求解4-block 无约束非光滑对偶问题。算法的关键是用Semismooth Newton法求解子问题（因为最初的4 blocks被合并为2 blocks，因此没有封闭形式）。imABCD算法的迭代复杂度是$O(1/k^2)$。

1. Preliminaries：定义和定理

Def 1（锥cone和凸锥convex cone）

设 $K\subset {\mathbb{R}}^n$，若 $\forall x\in K,\lambda >0$，都有 $\lambda x\in K$，则称K为锥。这样的集合是由原点发射的半直线的并集。称同时是凸集的锥为凸锥。

• 性质：任意多个凸锥（凸集）的交仍是凸锥（凸集）

Def 2（极锥与对偶锥）

设锥K，则
极锥为 K o : = { s ∈ R n : ⟨ s , x ⟩ ≤ 0 , ∀ x ∈ K } K^o:=\{s\in\mathbb{R}^n:\langle s,x\rangle\leq0,\forall x\in K\} Ko:={s∈Rn:⟨s,x⟩≤0,∀x∈K}
对偶锥为 K ∗ : = { s ∈ R n : ⟨ s , x ⟩ ≥ 0 , ∀ x ∈ K } K^*:=\{s\in\mathbb{R}^n:\langle s,x\rangle\geq0,\forall x\in K\} K∗:={s∈Rn:⟨s,x⟩≥0,∀x∈K}

• $K^o=-K^{\ast }$

Def 3（凸函数convex function）

1. 设$f(x):{\mathbb{R}}^n\to \overline{\mathbb{R}}$为广义实函数，定义域为$D_f$。若f(x)的上图像
   e p i f = { ( x , μ ) : x ∈ D f , μ ∈ R , μ ≥ f ( x ) } epif=\{(x,\mu):x\in D_f,\mu\in\mathbb{R},\mu\geq f(x)\} epif={(x,μ):x∈Df​,μ∈R,μ≥f(x)}为凸集，则称函数f(x)为凸函数。称f(x)的有效域为
   d o m f = { x ∣ f ( x ) < + ∞ } , domf=\{x|f(x)<+\infty\}, domf={x∣f(x)<+∞},则实质上
   e p i f = { ( x , μ ) : x ∈ d o m f , μ ∈ R , μ ≥ f ( x ) } epif=\{(x,\mu):x\in domf,\mu\in\mathbb{R},\mu\geq f(x)\} epif={(x,μ):x∈domf,μ∈R,μ≥f(x)}

2. 若凸函数$f(x)$的值域为$(-\infty ,+\infty ]$且不恒为$+\infty$，则称f为正常凸函数（proper convex）。例如，凸集$K$上的示性函数记为${\delta }_K$：
   $${\delta }_K(x)=\{\begin{array}{rlr}0& :& x\in K\\ +\infty & :& x\notin K\end{array}$$由定义和$K$的凸性，${\delta }_K$是正常凸函数。

Def 4（下半连续函数lower semicontinuous，lsc）

称函数 $f(x):S\to \mathbb{R}$是下半连续的，如果对于 { x n } ⊂ S : x n → x \{ x_n\}\sub S:x_n\rightarrow x {xn​}⊂S:xn​→x，都有
$$f(x)\le \underset{n\to \infty }{\lim }f(x_n)$$从定义中可以看出，称f(x)在x点处是下半连续的，如果
f ( x ) = lim inf ⁡ y → x f ( y ) = lim ⁡ ϵ → 0 ( inf ⁡ { f ( y ) ∣ ∣ y − x ∣ ≤ ϵ } ) f(x)=\liminf_{y\rightarrow x}f(y)=\lim_{\epsilon\rightarrow 0}(\inf\{f(y)||y-x|\leq\epsilon\}) f(x)=y→xliminf​f(y)=ϵ→0lim​(inf{f(y)∣∣y−x∣≤ϵ})

• 性质：函数f(x)在$R^n$是下半连续的
  ⇔ ∀ α ∈ R , { x ∣ f ( x ) ≤ α } \Leftrightarrow\forall\alpha\in R,\{x|f(x)\leq\alpha\} ⇔∀α∈R,{x∣f(x)≤α}是闭集
  $\iff epif$是$R^{n+1}$中闭集

Def 5（函数的凸包convex hull，$cof$）

记cof为函数f(x)的凸包，满足
$$epi(cof)=co(epif)$$即“凸包的上图像等于上图像的凸包”。

• 函数f(x)的凸包是可以由f(x)控制的最大的凸函数，若f(x)是凸函数，则它的凸包等于自身。

Def 6（函数的下半连续包$lscf$，凸函数的闭包closed hull，$clf$）

1. cl(epif)为任意函数f的上图像的闭包，则必存在一个函授g，使得
   $$cl(epif)=epi(g)$$则称g为f的下半连续包，记作$lscf$。

2. 记$clf$为$f$的闭包，满足
   $$clf(x)\equiv \{\begin{array}{rlr}lscf(x)& :& 若lscf(x)>-\infty ,\forall x\\ -\infty & :& 若lscf(x)=-\infty ,\exists x\end{array}$$即，对于正常凸函数，闭性与下半连续性等价；而闭的非正常凸函数一定是常函数$+\infty ,-\infty$。

Theorem 7（$lscf$和$clf$的性质）

设 $f:X\to [-\infty ,+\infty ]$为凸函数，则 $lscf$和 $clf$都是凸的。进一步的，若 $f$在一些点取值有限，则 $lscf$和 $clf$都是正常凸的，且有 $clf=lscf$。否则 $lscf$具有以下形式：
$$clf(x)=\{\begin{array}{rlr}-\infty & ,& x\in cl(domf)\\ +\infty & ,& x\notin cl(domf)\end{array}$$在后一种情形中，若 $cl(domf)\ne \varnothing$，则 $lscf(x)=+\infty ,\forall x\notin cl(domf)$，而根据定义， $clf(x)=-\infty$，故 $lscf\not\equiv clf$。而其余情形， $lscf\equiv clf$

Def 8（共轭函数）

设 $f:X\to [-\infty ,+\infty ]$，则f的共轭 $f^{\ast }:\Omega \to [-\infty ,+\infty ]$定义为
f ∗ ( ω ) = sup ⁡ x ∈ X { ⟨ x , ω ⟩ − f ( x ) } f^*(\omega)=\sup_{x\in X}\{\langle x,\omega\rangle-f(x)\} f∗(ω)=x∈Xsup​{⟨x,ω⟩−f(x)}

Theorem 9

对于任意的f，f ${}^{\ast }$是 $\Omega$上的闭凸函数。

Def 10（投影）

1. 定义
   设K是有界闭凸集，${\Pi }_K:Y\to K$是Y到K上的投影，则${\Pi }_K(y)={\mathrm{a}\mathrm{r}\mathrm{g}\mathrm{m}\mathrm{i}\mathrm{n}}_{x\in K}||x-y|{|}^2$（回顾：泛函分析最佳逼近定理）

2. Moreau分解
   设K为闭凸锥，则对于$x,x_1,x_2\in R^n$，有以下性质等价：
   （1）$x=x_1+x_2，x_1\in K，x_2\in K^o$，且$⟨x_1,x_2⟩=0$
   （2）$x_1={\Pi }_K(x),x_2={\Pi }_{K^o}(x)$

• 性质与定理¹：定义
  $$\begin{array}{rl}{d}_K：& R^n\to R\\ & x\mapsto \underset{y\in K}{\min }\frac{1}{2}||x-y|{|}^2\end{array}$$则由Moreau分解，$d_K(x)=\frac{1}{2}||x-{\Pi }_K(x)|{|}^2=\frac{1}{2}||{\Pi }_{K^o}(x)|{|}^2$，且$d_K$是一个可微的凸函数，其梯度为$\nabla d_K(x)={\Pi }_{K^o}(x)$

3. 半正定矩阵$S_{+}^n$上的投影
   设$C\in S^n$，则C有正交对角化：$C=Pdiag({\lambda }_1,...,{\lambda }_n)P^T$。记
   $$C_{+}=P\left[\begin{array}{ccc}\max ({\lambda }_1,0)& & \\ & \ddots & \\ & & \max ({\lambda }_n,0)\end{array}\right]P^T$$
   $$C_{-}=P\left[\begin{array}{ccc}\min ({\lambda }_1,0)& & \\ & \ddots & \\ & & \min ({\lambda }_n,0)\end{array}\right]P^T$$有如下结论：${\Pi }_{S_{+}^n}(C)=C_{+}，{\Pi }_{S_{-}^n}(C)=C_{-}，$其中$(S_{+}^n{)}^o=S_{-}^n$

Def 11（近端映射proximal mapping）

定义函数 $\psi :X\to R$在 $x$点处的近端映射为
P r o x ψ ( x ) ≜ a r g m i n ω ∈ X { ψ ( ω ) + 1 2 ∣ ∣ ω − x ∣ ∣ 2 } \mathrm{Prox}_{\psi}(x)\triangleq\mathrm{argmin}_{\omega\in X}\{\psi(\omega)+\frac{1}{2}||\omega-x||^2\} Proxψ​(x)≜argminω∈X​{ψ(ω)+21​∣∣ω−x∣∣2}

Theorem 12（Rademacher定理）

设 $E\subset {\mathbb{R}}^n$是开集， $F:E\to {\mathbb{R}}^n$是Lipschitz连续函数，则 $F$在 $E$上几乎处处可微。

Def 13（次微分与广义Jacobi矩阵）

设多值函数 $F:{\mathbb{R}}^n\to {\mathbb{R}}^m$是（局部）Lipschitz连续函数，根据Rademacher定理， $F$几乎处处可微。设
D F ≜ { x ∈ R n ∣ F 在 x 处 可 微 } D_F\triangleq\{x\in\mathbb{R}^n|F在x处可微\} DF​≜{x∈Rn∣F在x处可微}则 ${\Omega }_F\triangleq {\mathbb{R}}^n\backslash D_F$是零测度集。
令 $JF(x)$表示 $F$在 $x\in D_f$处的 $m\times n$ Jacobian矩阵。 $F$在 $x\in {\mathbb{R}}^n$处的B(Bouligand)-subdifferential定义为
∂ B F ( x ) ≜ { V ∈ R m × n : V = lim ⁡ x k → x x k ∈ D f J F ( x k ) } \partial_BF(x)\triangleq\{V\in\mathbb{R}^{m\times n}:V=\lim_{{x^k\to x}\atop{x^k\in D_f}} JF(x^k)\} ∂B​F(x)≜{V∈Rm×n:V=xk∈Df​xk→x​lim​JF(xk)}即 $JF(x^k)$的聚点集。从而，Clarke意义下的广义Jacobian是 ${\partial }_{B}F(x)$的凸包，即
$$\partial F(x)=\mathrm{c}\mathrm{o}{\partial }_{B}F(x)$$且 $\partial F(x)$是非空凸紧的、上半连续的（upper-semicontinuous）集合。（非空说明定义中的极限存在性）

• 【例】$f(x)=|x|$在$x=0$处不可微，在其他点处可微，则$$\underset{\genfrac{}{}{0ex}{}{x^k\to 0}{x^k\in D_f}}{\lim }JF(x^k)=-1或1$$故 ∂ B f = { − 1 , 1 } \partial_Bf=\{-1,1\} ∂B​f={−1,1}，从而$${\partial }_f=co{\partial }_{B}f=[-1,1]$$

• 【注】特别的，若函数$f$是凸的，则在$f$的二阶可导点x处，${\nabla }^{2}f(x)$是半正定矩阵。进一步的，若$\nabla f$Lipschitz连续，则对于$\forall x\in D_f$，${\nabla }^{2}f(x)$是半正定矩阵，从而${\partial }^{2}f(x)$中的元素都是半正定矩阵（利用矩阵的连续性，半正定矩阵的极限还是半正定的）

• 【注】对于连续可微函数，
  $$\underset{\genfrac{}{}{0ex}{}{x^k\to x}{x^k\in D_f}}{\lim }JF(x^k)=JF(x)$$即${\partial }_{B}F(x)$只包含x处的Jacobian矩阵这一个点，此时次微分就是微分（Jacobian）。

Def 14（半光滑函数semismooth）

设 $\varphi$是 局部Lipschitz连续函数。称 $\varphi$在x处semismooth，如果满足
（1） $\varphi$在x处方向可微（任意方向上方向导数存在）
（2）对于 $\forall V\in \partial \varphi (x+h)$，有
$$\varphi (x+h)-\varphi (x)-Vh=o(\Vert h\Vert ),h\to 0$$进一步的，若 $\varphi$在x处semismooth，且 $\forall V\in \partial \varphi (x+h)$，有
$$\varphi (x+h)-\varphi (x)-Vh=O(\parallel h{\parallel }^2)$$则称 $\varphi$在x处strongly semismooth。

• 理解：给定方向h，则在这个方向有上述Taylor展开式，其中$V\in \partial \varphi (x+h)$就体现了方向，即$V$来自于 “后方” 【这不由得让我想起了Nesterov gradient方法】。因为并不是所有的$V\in \partial \varphi (x)$都有上述Taylor展开，也就是可以理解为 “依方向Taylor展开” 。例如$f(x)=|x|$在$x=0$处，左侧右侧都可以分别做Taylor展开，但右侧要选取$\partial f(0)=\partial f(0+\delta )=1$，而左侧要选取$\partial f(0)=\partial f(0-\delta )=-1$，才能实现分别做Taylor展开。而如果一个函数在x点处可微，则当$h\to 0$时，$\partial \varphi (x+h)\to \nabla \varphi (x)$，即为Taylor展开式。

• 【注】最典型的（strongly）semismooth的例子是分片光滑函数，例如$y=|x|、y={\Pi }_{S_{+}^n}(x)$等。一般涉及min、max的情形，除了最优点之外都是光滑的，就是分片光滑函数。

• 性质：（strongly）semismooth函数的复合仍是（strongly）semismooth的。

2. 研究问题

凸优化问题（1）及对偶问题DP（3）

• 本文主要的目的是求解以下凸优化问题（1）：
  $$\begin{array}{rl}& \underset{x\in X}{\min }f(x)+\varphi (x)\\ & s.t.A(x)=b,g(x)\in C,x\in K\end{array}$$其中，X、Y、Z是有限维欧式空间，$f:X\to R$是光滑强凸函数，$\varphi :X\to (-\infty ,+\infty ]$是闭的正常凸函数，$A:X\to Y$是线性算子，$b\in Y$是给定数据，$C\subset Z，K\subset X$是两个闭凸锥，$g:X\to Z$是一个光滑的$C$-凸 map，满足以下性质：
  $$g(\alpha x_1+(1-\alpha )x_2)-[\alpha g(x_1)+(1-\alpha )g(x_2)]\in C,\forall x_1,x_2\in g^{-1}(C),\forall \alpha \in (0,1)$$关于g函数，例如不等式约束就包含其中。

• 问题（1）的Dual Problem有如下形式（3）：
  
  这是一个无约束convex nonsmooth problem，其中h是一个可微的凸函数，其梯度Lipschitz连续，${\varphi }_1,...,{\varphi }_4$是proper closed convex functions。问题（3）不仅是（1）的对偶问题，很多优化问题本身就具有（3）的形式，比如robust PCA²，也可以参考PCA 与 Robust PCA区别

• 可以用imABCD算法求解问题（3），迭代复杂度$O(1/k^2)$，将在Section2中详细介绍。

最佳逼近问题（2）

研究问题（1）的兴趣来源于以下最佳逼近问题（2）：
$$\begin{array}{rl}& \underset{X\subset S^n}{\min }\frac{1}{2}||X-G|{|}^2\\ & s.t.AX=b,BX\ge d,X\subset S_{+}^n,X\ge 0\end{array}$$其中，$A:S^n\to R^{m_E}$，$B:S^n\to R^{m_I}$是两个线性算子（从矩阵空间映射到向量空间），$G\in S^n,b\in R^{m_E},d\in R^{m_I}$是给定数据，$S^n,S_{+}^n$分别表示对称矩阵锥和半正定矩阵（凸）锥。

最小二乘协方差调整问题LSCAP

问题（2）是来源于最小二乘协方差矩阵调整问题（Least Squares Covariance Adjustment Problem，LSCAP）³：
$$\begin{array}{rl}& \underset{X\subset S^n}{\min }\frac{1}{2}||X-G|{|}^2\\ & s.t.X\subset S_{+}^n\\ & tr(A_{i}X)=b_i,i=1,...,p\\ & tr(C_{j}X)\le d_j,i=j,...,m\end{array}$$其中，$G,A_i,C_j\in S^n$是给定数据。
LSCAP的几何解释： 给定一个$G\in S^n$，希望（在Frobenius范数意义下）计算G到一般的多面体集交集上的投影，这些多面体集包括线性等式、线性不等式、半正定矩阵锥。
LSCAP的概率解释： 在理论上，协方差矩阵是对称半正定矩阵。但是在实际应用中，根据样本数据计算出的协方差矩阵一般只能保证是对称矩阵，所以就需要对其进行调整（adjustment），即计算它到半正定矩阵锥上的投影，同时还可能满足一些约束条件（例如，某两个变量之间的协方差已知或满足不等式）。这个过程本质上就是求解上述优化问题，得到的最优解是一个有效的协方差矩阵，可以保证后续的分析是有意义的。

3 问题求解与算法

求解对偶问题（3）

将4-blocks变量$({\omega }_1,{\omega }_2,{\omega }_3,{\omega }_4)$分为2-blocks，形式如下（5）：

因为$▽h$全局Lipschitz连续，则存在两个自伴半正定线性算子$\mathcal{Q}、\hat{\mathcal{Q}}:\mathbb{U}\to \mathbb{V}$使得

下面是$\mathcal{Q}、\hat{\mathcal{Q}}$的分解以及假设条件：

下面利用imABCD算法求解对偶问题（3）

【Algorithm】imABCD算法（$O(1/k^2)$）

【注】当${\epsilon }_k\equiv 0$时，此算法退化到mABCD算法。BCD算法可参考【机器学习】坐标下降法（Coordinate descent）和坐标下降法(Coordinate descent)

• 注意！ 算法中求解的是$\hat{h}$而不是$h$，这是因为$h$的形式可能很复杂，性质也未知，不方便求解。而$\hat{h}(\omega ,{\omega }^{\prime })$是一个凸二次函数，是$h$在${\omega }^{\prime }$处的二次近似，文献中给出证明，此算法最终可以收敛到最优解，且迭代复杂度为$O(1/k^2)$。而一般的gradient-type算法都是$O(1/k)$，之所以能加速，是使用了Nesterov加速技巧⁴，即算法中的Step 2.这一加速技巧非常常用，在很多论文中都有提到。
• 在某点处做函数的二次近似的想法很常用，因为凸二次问题求解简单且有很多好的方法，同时又容易保证算法的收敛性（Amazing！）。比如：

1. Newton法：$f(x)\approx f(x^k)+\nabla f(x^k{)}^T(x-x^k)+\frac{1}{2}||x-x^k|{|}_{H(x^k)}^2\triangleq \hat{f}(x;x^k)$
2. imABCD法：$h(\omega )\le h({\omega }^{\prime })+⟨\nabla h({\omega }^{\prime }),(\omega -{\omega }^{\prime })⟩+\frac{1}{2}||\omega -{\omega }^{\prime }|{|}_{\hat{\mathcal{Q}}}^2\triangleq \hat{h}(\omega ;{\omega }^{\prime })$
3. $g$梯度Lipschitz连续$\Rightarrow g(x)\le g(x^{\prime })+⟨\nabla g(x^{\prime }),(x-x^{\prime })⟩+\frac{L}{2}||x-x^{\prime }|{|}^2\triangleq \hat{g}(x;x^{\prime })$

求解最佳逼近问题（2）

回顾最佳逼近问题（2）：

问题（2）的对偶问题为（19）：

为了使用imABCD算法，将(y,S)分为一块，(z,Z)分为一块，这样可以使问题不会退化且精度更高。下面用Newton-type法求解imABCD算法中的两个子问题（7）

第一个子问题的分析：

其中，式(21)和$\nabla \xi$的形式需要借助投影算子的性质得到。而且有如下重要性质：
（1）${\Pi }_{S_{+}^n}(\cdot )$是$S^n$上的处处strongly semismooth函数，而且是1-Lipschitz连续函数，因此几乎处处可微。
（2）strongly semismooth函数的复合仍是strongly semismooth的。

接下来，就是大名鼎鼎的半光滑Newton法出场了！（听说孙老师的看家本领就是它）

【Algorithm】Semismooth Newton-CG算法（SNCG）

此方法的收敛性在其他论文中已经得到证明，主要参考(H. Qi AND D.F. Sun, 2006)⁵和(X.Y. Zhao, D.F. Sun, AND K.-C. Toh, 2010)⁶

• Semismooth Newton-CG法

1. 适用条件：
   \quad 我们知道，经典的Newton法要求函数二阶可导，故Step 1中的广义Hessian阵（一个集合）退化为Hessian阵，而Step 1变为用CG法求解：我们知道，经典的Newton法要求函数二阶可导，故Step 1中的广义Hessian阵（一个集合）退化为Hessian阵，而Step 1变为用CG法求解：
   $$H(x^k)d+\nabla f(x^k)=0$$而当二阶不可微时，Newton法便不再适用。
   \quad 也就是说，在求解一个凸优化问题时，如果函数一阶可微但二阶不可微，但是$\nabla f$是semismooth的，就可以使用Semismooth Newton-CG法进行求解。
   \quad 另一种情形（如子问题二）：若$F$是semismooth的，则在求解非线性方程组$F(x)=0$也可以使用SNCG法。但此时，因为没有优化问题的背景，所以收敛性可能无法保证，一般需要结合加速近端梯度法（accelerated proximal gradient，APG）使用。在求解子问题二时会进一步讨论。
2. 优点：
   （1）大量数值实验表明，SNCG算法的收敛速度很快，很多问题只需几步即可迭代到最优解。从理论分析可知，一般SNCG是超线性收敛的（superlinear），而如果$\nabla f$是strongly semismooth的，则二阶收敛。
   （2）在Step 1中，从（广义hessian）中任意选取V，算法都是有效的。
3. CG法：可参考共轭梯度法的推导与完整算法和梯度下降法和共轭梯度法有何异同？
   （1）用于求解线性方程组$AX=b$，其中$A$是正定矩阵。
   （2）具体操作时，尽量将矩阵A变稀疏或将A分解，计算复杂度可由O($n^2$)下降为O(C$n$)，C为与具体问题相关的常数。

• 在求解第一个子问题时，此算法中有以下几点值得注意⁶：
  （1）精确表示${\partial }^2\xi$是困难的，因此定义
  $${\hat{\partial }}^2\xi (y)\triangleq \mathcal{A}{\Pi }_{S_{+}^n}(G_1+{\mathcal{A}}^{\ast }{y}^j){\mathcal{A}}^{\ast }$$且有性质：$\forall d,{\partial }^2\xi (y)d\subset {\hat{\partial }}^2\xi (y)d$。这意味着，若所有${\hat{\partial }}^2\xi (y)$中的元素都正定，那么${\partial }^2\xi (y)$中的所有元素一定都是正定的。（由命题3.6⁵，在最优解$y^{\ast }$处，$\forall V\in {\partial }^2\xi (y^{\ast })$都正定，故在$y^{\ast }$的一个邻域中，所有$V\in {\partial }^2\xi (y)$都正定。另外，在文献⁶Algorithm 2中，使用的是$V_j+{\epsilon }_{j}I$，也是保证其正定性）
  （2）CG法的终止条件为：$||\mathcal{A}{V}^j{\mathcal{A}}^{\ast }{d}^j+\nabla \xi (y^j)||\le \delta$，$\delta$的形式根据具体问题的收敛性证明的需要而定。
  （3）$\forall j\ge 0$，$d^j$都是$\xi$的下降方向，因此Step 2中的line search是well defined的
  （4）$\exists j_0$正整数，使得$\forall j\ge j_0$都有
  $$\xi (y^j+d^j)\le \xi (y^j)+\mu ⟨\nabla \xi (y^j),d^j⟩$$即$y^{j+1}=y^j+d^j$

第二个子问题的分析：

• 添加近端项

  • 需要注意，近端项$\frac{c}{2}||z-z_0|{|}^2$的作用是保证目标函数关于$z$是强凸的。添加之前，目标函数关于$z$求梯度为：
    $$\mathcal{B}({\mathcal{B}}^{\ast }z+Z+G_2)-d$$而$\mathcal{B}{\mathcal{B}}^{\ast }$是半正定矩阵（有限维Euclid空间中，线性算子作用在向量上可看作矩阵），这可能导致无法解除$z$。而添加近端项之后，梯度变为：
    $$\mathcal{B}({\mathcal{B}}^{\ast }z+Z+G_2)-d+c(z-z_0)$$此时$z$的系数矩阵为：$\mathcal{B}{\mathcal{B}}^{\ast }+cI\succ 0$
  • 在imABCD算法中，$\hat{h}$是$h$的一个上界，所以再扩大仍然是上界，这意味着，添加近端项后算法仍是有效的。

• 与第一个子问题（21）不同，问题（23）是一个有约束的SC${}^1$问题（即，目标函数连续可微，梯度semismooth）。根据文献⁷，我们也可以使用SNCG算法，但每步迭代都需要求解一个严格凸二次规划问题，比较麻烦。

• 另一个思路是将问题（23）看做一个更一般的无约束非光滑凸问题：
  
  具体来说，问题（23）的等价形式为：
  $$\underset{z}{\min }\frac{1}{2}||{\Pi }_{\ge 0}({\mathcal{B}}^{\ast }z+G_2)|{|}^2-⟨d,z⟩+\frac{c}{2}||z-z_0|{|}^2+{\delta }_{\ge 0}(z)$$其中，$\psi (z)={\delta }_{\ge 0}(z)$凸非光滑（${\mathbb{R}}^{n^{+}}$是凸锥，更具体说是凸多面体集），剩下的强凸光滑函数部分对应为$\zeta (z)$。此强凸问题可以用APG算法求解⁴，是全局收敛的，但收敛速度只有线性的⁸。可以参考①近端梯度下降算法(Proximal Gradient Algorithm)、②近端梯度下降proximal gradient descent、③对近端梯度算法（Proximal Gradient Method）的理解和④近端梯度法(Proximal Gradient Method, PG)。具体算法过程为：

【Algorithm】Accelerated Proximal Gradient算法（APG）

1. 主要原理：
   问题（24）的最优解是下面非线性方程组的根：
   
   推导方法有三种：
   （1）上面四篇参考文章的第二篇，核心为：
   $$\begin{array}{rl}\hat{\zeta }(x;x_k)& \triangleq \zeta (x_k)+⟨\nabla \zeta (x_k),(x-x_k)⟩+\frac{L}{2}||x-x_k|{|}^2\\ & =\frac{L}{2}\parallel x-(x_k-\frac{1}{L}\nabla \zeta (x_k)){\parallel }^2+C\end{array}$$其中C是与x无关的常数。用$\hat{\zeta }$近似$\zeta$，则有：
   x ∗ = a r g m i n x { ζ ^ ( x ) + ψ ( x ) } = a r g m i n x { L 2 ∥ x − ( x k − 1 L ∇ ζ ( x k ) ) ∥ 2 + ψ ( x ) } = P r o x ψ / L ( x k − 1 L ∇ ζ ( x k ) ) \begin{aligned} x^*&=\mathrm{argmin}_x\{\widehat{\zeta}(x)+\psi(x)\}\\ &=\mathrm{argmin}_x\{\frac{L}{2}\parallel x-(x_k-\frac{1}{L}\nabla \zeta(x_k))\parallel^2+\psi(x)\}\\ &=\mathrm{Prox}_{\psi/L}(x_k-\frac{1}{L}\nabla \zeta(x_k)) \end{aligned} x∗​=argminx​{ζ ​(x)+ψ(x)}=argminx​{2L​∥x−(xk​−L1​∇ζ(xk​))∥2+ψ(x)}=Proxψ/L​(xk​−L1​∇ζ(xk​))​
   （2）上面四篇参考文章的第三篇，核心为：
   
   （3）原问题$\min \zeta (x)+\psi (x)$，则最优解$x^{\ast }$满足：
   $$0\in \nabla \zeta (x^{\ast })+\partial \psi (x^{\ast })$$转化后的问题为$F(x)\triangleq x-{\mathrm{P}\mathrm{r}\mathrm{o}\mathrm{x}}_{\psi }(x-\nabla \zeta (x))=0$，它的根满足：
   x ∗ = a r g m i n ω { ψ ( ω ) + 1 2 ∥ ω − x ∗ + ∇ ζ ( x ∗ ) ∥ 2 } x^*=\mathrm{argmin}_{\omega}\{\psi(\omega)+\frac{1}{2}\parallel\omega-x^*+\nabla\zeta(x^*)\parallel^2\} x∗=argminω​{ψ(ω)+21​∥ω−x∗+∇ζ(x∗)∥2}则也满足类似的关系：
   $$0\in \partial \psi (x^{\ast })+x^{\ast }-x^{\ast }+\nabla \zeta (x^{\ast })$$因此，两个问题等价。
2. 这意味着，若用APG算法求解问题（24）得到最优解，等价于用它去寻找$F(x)=0$的根，迭代过程是一样的！（为后面APG-SNCG算法做铺垫）

进一步分析

• 值得注意的是，有如下性质存在：
  （1）任意凸函数都是semismooth的。
  （2）投影算子${\Pi }_K(\cdot )$，若$K$是凸集，则${\Pi }_K(\cdot )$是凸函数。所以，${\Pi }_{S_{+}^n}(\cdot )$、${\Pi }_{\ge 0}(\cdot )$都是semismooth的。
  （3）近端算子Prox的性质：当$\psi (x)={\delta }_{\ge 0}(x)$时，${\mathrm{P}\mathrm{r}\mathrm{o}\mathrm{x}}_{c\psi }(x)={\Pi }_{\ge 0}(x)$是凸函数。一般的，Prox函数的性质依赖于$\psi$，例如，当$\psi =||\cdot |{|}_1$时，${\mathrm{P}\mathrm{r}\mathrm{o}\mathrm{x}}_{||\cdot |{|}_1}(x)$为软阈值操作；当$\psi =||\cdot |{|}_2$时，${\mathrm{P}\mathrm{r}\mathrm{o}\mathrm{x}}_{||\cdot |{|}_2}(x)$为semismooth函数。
  （4）一般的，在涉及max、min时，都是分片光滑的，也就是说，${\Pi }_{S_{+}^n}(\cdot )$、${\Pi }_{\ge 0}(\cdot )$都是分片光滑函数，从而是strongly semismooth的。
• 因此，再具体求解子问题（23）时，$F$是semismooth的，因此可以利用SNCG法和APG算法的杂合体进行求解。

【Algorithm】APG-SNCG算法

用流程图总结如下：

• 与子问题（21）不同，$F(x)=0$是一个新问题，并没有优化问题的背景，也就是说，$F$已经不是某个优化问题的梯度，所以没有对应的目标函数，Step1所求出的$d^j$也就不是哪个目标函数的下降方向！而且。同时，Step2中对m设置了上界m${}_0$，以使步长不会过小。因此，Step2中的不等式不一定满足（从而可能转到step3）。
  \quad 另一方面，文献4命题3.6证明了$\partial F(y^{\ast })$中所有元素都是正定矩阵，而文献中$F(y)\triangleq \mathcal{A}{\Pi }_{S_{+}^n}(G+{\mathcal{A}}^{\ast }{y}^j)$与子问题二中的$F$不同（与子问题一$\nabla \xi$相同）。因此，Step1中的$V^j$不一定是正定矩阵，CG法对应的不等式（26）不一定满足（从而可能转到step3）。
• Remark.
  \quad 1. 由$\zeta$是强凸的，Step3 APG是此算法全局收敛性的保障！
  \quad 2. 在不等式（26）成立下，由Step1得到的$d^j$以概率1是$||F(x)||$在$x^j$处的下降方向。

将子问题二分解为smaller decoupled problems（略）

此部分核心在于，将子问题二拆分为多个更容易求解的问题，最终形式为：

4. 数值实验（略）

此部分只列出几种算法

【Algorithm】Accelerated Block Coordinate Gradient Descent算法（ABCGD）

【Algorithm】A four-block inexact enhanced Randomized ABCD算法（eRABCD）

---

1. A dual approach to semidefinite least-squares problems（Jerome Malick，2004） ↩︎

2. Robust principal component analysis: Exact recovery of corrupted low-rank matrices via convex optimization（J. Wright等，2009） ↩︎

3. Least-squares covariance matrix adjustment（S. Boyd AND X. Lin，2005） ↩︎

4. A method of solving a convex programming problem with convergence rate$O(1/k^2)$（Nesterov，1983） ↩︎ ↩︎

5. A quadratically convergent newton method for computing the nearest correlation matrix（H. Qi AND D.F. Sun，2006） ↩︎ ↩︎

6. A Newton-CG augmented Lagrangian method for semidefinite programming（X.Y. Zhao, D.F. Sun, AND K.-C. Toh，2010） ↩︎ ↩︎ ↩︎

7. A Globally Convergent Newton Method for Convex SC${}^1$ Minimization Problems（J.-S. Pang AND L. Qi，1995） ↩︎

8. Convergence rates of inexact proximal-gradient methods for convex optimization（M. Schmidt, N.L. Roux, AND F.R. Bach，2011） ↩︎