凸优化学习笔记：等式约束凸优化问题的Newton方法、不可行初始点的Newton方法

基础预备：

无约束优化：线搜索最速下降

无约束优化：修正阻尼牛顿法

无约束优化：Hessian-Free Optimization 拟牛顿类算法（BFGS，L-BFGS）

凸优化学习笔记：Lagrange函数、对偶函数、对偶问题、KKT条件

约束优化：约束优化的三种序列无约束优化方法

等式约束凸优化问题

考虑下述等式约束凸优化问题的求解方法：
$$\begin{array}{ll}\min & f(x)\\ \text{ s.t. }& Ax=b\end{array}\text{\hspace{1em}(1)}$$
其中$f:{\mathbb{R}}^n\to \mathbb{R}$是二次连续可微凸函数，$A\in {\mathbb{R}}^{p\times n}$，$rank(A)=p<n$。根据KKT条件可知，$x^{\ast }\in domf$是优化问题（1）最优解的充要条件是，存在${\nu }^{\ast }\in {\mathbb{R}}^p$满足：
$$A{x}^{\ast }=b,\nabla f\left(x^{\ast }\right)+A^T{\nu }^{\ast }=0\text{\hspace{1em}(2)}$$
因此，求解等式约束优化问题（1）等价于求解一个含有$n+p$个变量$x^{\ast },{\nu }^{\ast }$的$n+p$个方程的问题。一般第一组方程$A{x}^{\ast }=b$称为原可行方程，是线性的，而第二组方程$\nabla f\left(x^{\ast }\right)+A^T{\nu }^{\ast }=0$称为对偶可行方程，通常是非线性的，很难求出其解析解。将$\nabla f\left(x^{\ast }\right)$线性化就是等式约束Newton方法的基本思想。

当然，任何等式约束优化问题都可以通过消除等式约束转化为等价的无约束优化问题，然后利用之前讲到的各种无约束优化问题求解方法去解；或者是利用无约束优化方法求解对偶问题，然后根据对偶解得到原问题（1）的解。但是，消除方法（或对偶方法）通常会破坏问题的结构，比如稀疏性，这样可能会导致忽视掉原问题的一些结构优势，造成耗时。

等式约束凸二次规划

考虑等式约束凸二次规划问题：
$$\begin{array}{ll}\min & f(x)=(1/2)x^{T}Px+q^{T}x+r\\ \text{ s.t. }& Ax=b\end{array}\text{\hspace{1em}(3)}$$
其中，$P\in {\mathbf{S}}_{+}^n,A\in {\mathbb{R}}^{p\times n}$。该问题是扩展Newton方法处理等式约束问题的基础。此时最优性条件可以写为：
$$A{x}^{\ast }=b,P{x}^{\ast }+q+A^T{\nu }^{\ast }=0\text{\hspace{1em}(4)}$$
用矩阵表示为：
$$\left[\begin{array}{cc}P& A^T\\ A& 0\end{array}\right]\left[\begin{array}{l}{x}^{\star }\\ {\nu }^{\star }\end{array}\right]=\left[\begin{array}{c}-q\\ b\end{array}\right]\text{\hspace{1em}(5)}$$
这个$n+p$个线性方程组称为等式约束优化问题（3）的KKT系统，系统矩阵称为KKT矩阵。如果KKT矩阵非奇异，存在唯一最优的原对偶对 $(x^{\ast },{\nu }^{\ast })$；如果KKT矩阵奇异，但是有解，任何解都构成最优对偶对 $(x^{\ast },{\nu }^{\ast })$；如果KKT矩阵无解，那么二次优化问题或者无解或者无下界。已证明，$P$正定是KKT矩阵非奇异的充分条件。

等式约束的Newton方法

该方法和无约束优化的Newton方法几乎一样，初始点必须可行，即满足$x\in domf且Ax=b$，根据等式约束的需要对Newton方向的定义进行适当修改，即要保证Newton方向$\Delta x_{nt}$是可行方向，满足$A\Delta x_{nt}=0$。

Newton方向：基于二阶近似的定义

前面指出，一般情况下（2）是非线性方程，很难有解析解。为了导出等式约束问题（1）在**可行点$x$**处的Newton方向$\Delta x_{nt}$，我们对目标函数$f$在$x$附近做二阶Taylor展开：
$$\begin{array}{ll}\text{ min }& \hat{f}(x+v)=f(x)+\nabla f(x{)}^{T}v+(1/2)v^T{\nabla }^{2}f(x)v\\ \text{ s.t. }& A(x+v)=b\end{array}\text{\hspace{1em}(6)}$$
该问题变量为$v$，且是一个等式约束凸二次规划问题，可以用解析解求解。我们假定该系统的KKT矩阵非奇异，$v=\Delta x_{nt}$，在此基础上定义$x$处的Newton方向$\Delta x_{nt}$为问题（6）的解。这样一来，通过求解问题（6）得到$\Delta x_{nt}$，这是一个向量，被称为Newton方向，将其加到$x$上就得到了以$f$的二阶近似$\hat{f}(x+v)$为目标函数的对应问题的最优解。

已知$Ax=b$，根据KKT条件，Newton方向$\Delta x_{nt}$可由以下方程确定：
$$\left[\begin{array}{cc}{\nabla }^{2}f(x)& A^T\\ A& 0\end{array}\right]\left[\begin{array}{c}\Delta x_{\mathrm{nt}}\\ w\end{array}\right]=\left[\begin{array}{c}-\nabla f(x)\\ 0\end{array}\right]\text{\hspace{1em}(7)}$$
其中$\omega$是该二次问题的最优对偶变量，$\Delta x_{nt}$只在KKT矩阵非奇异时有定义。如同用无约束优化问题的Newton方法一样，当目标函数$f$是严格二次函数时，Newton修正向量$x+\Delta x_{nt}$是等式约束优化问题的准确最优解，此时$\omega$就是原问题最优的对偶变量。也就是说，当$f$接近二次时，$x+\Delta x_{nt}$可以说是最优解$x^{\ast }$的很好估计，$\omega$就是最优对偶变量${\nu }^{\ast }$的很好估计。

Newton减量：用于设计停止准则

将等式约束优化问题的Newton减量定义为：
$$\lambda (x)={\left(\Delta x_{\mathrm{nt}}^T{\nabla }^{2}f(x)\Delta x_{\mathrm{nt}}\right)}^{1/2}\text{\hspace{1em}(8)}$$
已知$\hat{f}(x+v)=f(x)+\nabla f(x{)}^{T}v+(1/2)v^T{\nabla }^{2}f(x)v$ 表示$f$在$x$处的二阶Taylor近似，则原目标函数$f(x)$和近似模型$\hat{f}(x+v)$之间的差值为：
f(x)−inf⁡{f^(x+v)∣A(x+v)=b}=λ(x)2/2(9)f(x)-\inf \{\hat{f}(x+v) \mid A(x+v)=b\}=\lambda(x)^2 / 2 \tag{9} f(x)−inf{f^​(x+v)∣A(x+v)=b}=λ(x)2/2(9)
这说明，$\lambda (x{)}^2/2$对$x$处的$f(x)-p^{\ast }$给出了基于二阶近似模型的一个估计，而$\lambda (x)$（或者$\lambda (x{)}^2/2$）也可以作为设计停止准则的基础。

Newton方向$\Delta x_{nt}$总是可行下降方向。首先，满足$A\Delta x_{nt}=0$说明该方向是可行方向；其次，已知$f$沿方向$\Delta x_{nt}$的方向导数为
$${\frac{d}{dt}f\left(x+t\Delta x_{\mathrm{nt}}\right)|}_{t=0}=\nabla f(x{)}^T\Delta x_{\mathrm{nt}}=-\lambda (x{)}^2\text{\hspace{1em}(10)}$$
是小于0的，说明其为下降方向。

算法框架

从上述分析可知，Newton方向$\Delta x_{nt}$总是可行下降方向，因此该方法是一种可行下降方法，该方法需要每个$x$处的KKT矩阵可逆。

不可行初始点的Newton方法

根据上文分析已知，上述Newton方法是一个可行下降方法，要求初始点$x\in domf$满足$Ax=b$。下面介绍一种能从不可行初始点开始进行迭代的Newton方法，即$x\in domf$但$Ax\ne b$。

不可行点的Newton方向

从等式约束优化问题（1）的最优性条件（2）出发，假定$x$为当前点且不可行，我们需要基于该点$x$找到一个方向$\Delta x$使得$x+\Delta x$满足（至少近似满足）最优性条件，即$x+\Delta x\approx x^{\ast }$。为找到方向$\Delta x$，用$x+\Delta x$代替$x^{\ast }$，用$\omega$代替${\nu }^{\ast }$，将原函数$f$的梯度进行一阶泰勒展开：
$$\nabla f(x+\Delta x)\approx \nabla f(x)+{\nabla }^{2}f(x)\Delta x\text{\hspace{1em}(11)}$$
并带入（2）中，得到：
$$A(x+\Delta x)=b,\nabla f(x)+{\nabla }^{2}f(x)\Delta x+A^{T}w=0\text{\hspace{1em}(12)}$$
可将其写成关于变量$\Delta x$和$\omega$的线性方程组形式：
$$\left[\begin{array}{cc}{\nabla }^{2}f(x)& A^T\\ A& 0\end{array}\right]\left[\begin{array}{c}\Delta x\\ w\end{array}\right]=-\left[\begin{array}{c}\nabla f(x)\\ Ax-b\end{array}\right]\text{\hspace{1em}(13)}$$
此时的Newton方向$\Delta x$和在可行点$x$处定义的Newton方向不一样，但可以看出，若残差向量$Ax-b=0$，方程（13）就退化成在可行点$x$处定义Newton方向的方程（7）。为避免混淆，我们用符号$\Delta x_{nt}$表示（13）所定义的方向$\Delta x$，称为不可行点的Newton方向。

基于等式约束问题的原对偶方法

原对偶方法是指同时修改原变量$x$和对偶变量$\nu$使KKT条件（近似）满足的方法。换种思路考虑，将$x$和$\nu$都视为变量，可将KKT条件表示成$r(x^{\ast },{\nu }^{\ast })=0$，其中$r:{\mathbb{R}}^n\times {\mathbb{R}}^p\to {\mathbb{R}}^n\times {\mathbb{R}}^p$由下式定义：
$$r(x,\nu )=\left(r_{\text{dual }}(x,\nu ),r_{\text{pri }}(x,\nu )\right)\text{\hspace{1em}(14)}$$
其中，$r_{\text{dual }}(x,\nu )=\nabla f(x)+A^T\nu ,r_{\text{pri }}(x,\nu )=Ax-b$，分别表示对偶残差和原残差。将$r$在当前估计$y$附近进行一阶Taylor展开，得到：
$$r(y+z)\approx \hat{r}(y+z)=r(y)+\mathrm{Dr}(y)z\text{\hspace{1em}(15)}$$
其中，$\mathrm{Dr}(y)\in {\mathbf{R}}^{(n+p)\times (n+p)}$表示函数$r$在$y$处一阶导。根据KKT条件，将原对偶Newton方向$\Delta y_{pd}=(\Delta x_{pd},\Delta {\nu }_{pd})$定义为使$\hat{r}(y+z)=0$的步长$z$，即
$$Dr(y)\Delta y_{\mathrm{pd}}=-r(y)\text{\hspace{1em}(16)}$$
进一步可写成：
$$\left[\begin{array}{cc}{\nabla }^{2}f(x)& A^T\\ A& 0\end{array}\right]\left[\begin{array}{c}\Delta x_{\mathrm{pd}}\\ \Delta {\nu }_{\mathrm{pd}}\end{array}\right]=-\left[\begin{array}{c}{r}_{\text{dual }}\\ r_{\text{pri }}\end{array}\right]=-\left[\begin{array}{c}\nabla f(x)+A^T\nu \\ Ax-b\end{array}\right]\text{\hspace{1em}(17)}$$
令${\nu }^{+}=\nu +\Delta {\nu }_{pd}$，（17）又可写为：
$$\left[\begin{array}{cc}{\nabla }^{2}f(x)& A^T\\ A& 0\end{array}\right]\left[\begin{array}{c}\Delta x_{\mathrm{pd}}\\ {\nu }^{+}\end{array}\right]=-\left[\begin{array}{c}\nabla f(x)\\ Ax-b\end{array}\right]\text{\hspace{1em}(18)}$$
观察（13）（17）和（18）可以发现$\Delta x_{\mathrm{nt}}=\Delta x_{\mathrm{pd}}$，修正后的对偶变量$ w=\nu^{+}=\nu+\Delta \nu_{\mathrm{pd}}$，这说明最初定义的不可行Newton方向的方程（13）计算出的是不可行点Newton方向和修正后的对偶变量，而方程（17）表明不可行点Newton方向和相应的对偶方向是以原对偶残差为等式右边项的方程组的解，对偶变量的当前值对计算对偶方向或者其修正值都不起作用。

停止准则

结合（13）可知$f$沿方向$\Delta x_{nt}$的方向导数为
$$\begin{array}{rl}{\frac{d}{dt}f(x+t\Delta x)|}_{t=0}& =\nabla f(x{)}^T\Delta x\\ & =-\Delta x^T\left({\nabla }^{2}f(x)\Delta x+A^{T}w\right)\\ & =-\Delta x^T{\nabla }^{2}f(x)\Delta x+(Ax-b{)}^{T}w\end{array}\text{\hspace{1em}(19)}$$
可以发现，除非$x$可行，则上式不一定为负，这说明在不可行点$x$处的Newton方向不一定是$f$的下降方向。但我们可以计算残差范数沿$\Delta y_{pd}$的方向导数，即
$${\frac{d}{dt}{\Vert r\left(y+t\Delta y_{\mathrm{pd}}\right)\Vert }_2^2|}_{t=0}=2r(y{)}^T\mathrm{Dr}(y)\Delta y_{\mathrm{pd}}=-2r(y{)}^{T}r(y)=-||r(y)|{|}_2\text{\hspace{1em}(19)}$$
这表明残差范数沿不可行点Newton方向和相应的对偶方向$\Delta y_{pd}=(\Delta x_{pd},\Delta {\nu }_{pd})$下降，可用来设计停止准则。

算法框架

---

参考文献

Convex Optimization – Boyd and Vandenberghe