数值最优化—无约束问题信赖域算法

一、参考

《数值最优化算法与理论》

二、线性搜索与信赖域算法区别

数值最优化—无约束问题最速下降法和Newton法 和 数值最优化—无约束问题共轭梯度法 中介绍的求解无约束问题的几类算法的共同点是基于方向和步长确定一个下降方向$d^{(k)}$，然后从$x^{(k)}$出发，沿方向$d^{(k)}$进行线性搜索确定步长${\alpha }_k$，得到下一个迭代点$x^{(k+1)}=x^{(k)}+{\alpha }_k{d}^{(k)}$。称这类算法为线性搜索算法。

下面介绍求解无约束问题：
$$minf(x),(x\in R^n)$$
的另一类算法—信赖域算法。信赖域算法的基本思想是：在当前迭代点$x^{(k)}$的附近用一个简单函数近似目标函数$f$。用该近似函数在$x^{(k)}$的某个邻域内的极小值点作为下一个迭代点。与线性搜索型算法比较，信赖域算法在每次迭代同时确定搜索方向和步长。

三、信赖域算法

由于信赖域方法用$f$的某个近似函数在$x^{(k)}$邻域的极小值点作为下一个迭代点，因此，设计信赖域算法需要考虑如下相关问题：

1. 目标函数$f$的（简单）近似形式。
2. 点$x^{(k)}$的邻域（称为信赖域）大小的确定。
3. 函数值序列的下降性检测。
4. 近似问题（称为信赖域子问题）的求解。

考虑到一般非线性函数在任何点的邻域内都可以用二次函数近似，而且二次函数的极小值问题相对容易求解。因此，在信赖域算法中，通常我们采用二次函数作为目标函数$f$的近似，用该二次函数在$x^{(k)}$某邻域内的极小值点作为下一个迭代点，即$x^{(k+1)}$为下面问题的解：
$$minf(x^{(k)})+\nabla f(x^{(k)}{)}^T(x-x^{(k)})+\frac{1}{2}(x-x^{(k)}{)}^T{B}_k(x-x^{(k)})$$$$s.t.||x-x^{(k)}||\le {\Delta }_k$$
其中，$B_k\in R^{n\times n}$是$f$在$x$处的Hessian矩阵或其近似，参数${\Delta }_k>0$控制$x^{(k)}$的邻域大小。若令$d=x-x^{(k)}$，则上面的问题可改写为如下等价的二次函数极小值问题：
$$minf(x^{(k)})+\nabla f(x^{(k)}{)}^{T}d+\frac{1}{2}{d}^T{B}_{k}d\stackrel{△}{=}{q}_k(d)$$$$s.t.||d||\le {\Delta }_k$$
设$d^{(k)}$是问题的解。则下一个迭代点为$x^{(k+1)}=x^{(k)}+d^{(k)}$。问题中的可行域：
D = { d ∈ R n   ∣   ∣ ∣ d ∣ ∣ ≤ Δ k } D = \{d \in R^n \ | \ ||d|| \le \Delta _k \} D={d∈Rn ∣ ∣∣d∣∣≤Δk​}
称为信赖域，参数$\Delta >0$称为信赖域半径，等价的二次函数极小值问题称为信赖域子问题，其中范数可任意选取。本章采用Euclid范数。

信赖域子问题是仅有一个不等式约束的二次函数极小值问题。

由函数的二阶Taylor展开可知，当${\Delta }_k$充分小时，二次函数$q$在信赖域$D$中是$f$的一个很好的近似。另一方面，当$f$本身是一个二次函数或近似与二次函数时，$q$可能在某一个较大范围内也是$f$的一个很好的近似。如何确定信赖域半径是信赖域算法的一个重要环节。一方面，我们希望在$D$中$q$与$f$的近似程度好，另一方面希望信赖域半径尽可能大。为了合理确定信赖域半径，我们定义$\Delta F_k$为$f$在第$k$步的实际下降量。即：
$$\Delta f_k=f(x^{(k)})-f(x^{(k)}+d^{(k)})$$
其中，$d^{(k)}$是信赖域子问题的解。令${\Delta }_{q_k}$为对应的预测下降量，即：
$${\Delta }_{q_k}=f(x^{(k)})-q_k(d^{(k)})$$
注意到$r_k$从某种程度反映了二次函数$q_k(d^{(k)})$与目标函数$f(x^{(k)}+d^{(k)})$的近似程度。若$r_k$接近于1，则可认为二次函数$q_k(d^{(k)})$在信赖域$D$上与目标函数$f(x^{(k)}+d^{(k)})$的近似程度很好。反之，若$r_k$离1较远，我们认为$q_k(d^{(k)})$在信赖域$D$上与目标函数$f(x^{(k)}+d^{(k)})$的近似程度不好。基于上述观察，我们可以用$r_k$与1的近似程度作为对信赖域半径是否合适的准则。即，可通过如下方式调整信赖域半径。

给定常数$\eta ,{\eta }_1,{\eta }_2\in (0,1)$满足$\eta <{\eta }_1<{\eta }_2$（一般地，$\eta$接近于或等于0，${\eta }_2$接近于1，如${\eta }_2=\frac{3}{4}$）。

① 若$r_k\ge {\eta }_2$，我们可认为$q$在$D$中是$f$的一个很好的近似，或者说得到一个成功的迭代点$x^{(k+1)}=x^{(k)}+d^{(k)}$。此时，$q$有可能在一个更大一点的区域也是$f$的一个很好的近似。因此，我们可将信赖域半径扩大，即令${\Delta }_{k+1}>{\Delta }_k$。

② 若${\eta }_1<r_k<{\eta }_2$，即$q$在$D$中是$f$的一个好的近似，或者说得到一个好的迭代点$x^{(k+1)}=x^{(k)}+d^{(k)}$。此时可以保持信赖域半径不变，即令${\Delta }_{k+1}={\Delta }_k$（为了下一次得到一个更成功的点，也可减少信赖域半径，即令${\Delta }_{k+1}<{\Delta }_k$）。

③ 若$r_k\le \eta$，即$q$在$D$中是$f$的一个不好的近似，或者说得到一个不成功的迭代点$x^{(k+1)}=x^{(k)}+d^{(k)}$。说明信赖域半径过大，此时需要减少信赖域半径，即令${\Delta }_{k+1}<{\Delta }_k$。

在上面的基础上，给出求解无约束问题的信赖域算法如下：

0. （初始化）取初始点$x^{(0)}\in R^n，\bar{\Delta }>0,{\Delta }_0\in (0,\bar{\Delta })，\eta \in [0,\frac{1}{4})$，精度$ϵ>0$。令$k=0$。
1. （收敛性检测）若$||\nabla f(x^{(k)})||\le ϵ$，则算法终止。得问题的解$x^{(k)}$。否则，转2。
2. （子问题求解）解信赖域子问题得解$d^{(k)}$。
3. （信赖域修正）计算$r_k$。
   若$r_k>\frac{3}{4}$，则令 Δ k + 1 = m i n { 2 Δ k , Δ ˉ } \Delta _{k+1} = min\{2\Delta _k,\bar \Delta \} Δk+1​=min{2Δk​,Δˉ}。
   若$\frac{1}{4}<r_k<\frac{3}{4}$，则令${\Delta }_{k+1}={\Delta }_k$。
   若$\eta <r_k<\frac{1}{4}$，则令${\Delta }_{k+1}=\frac{1}{2}{\Delta }_k$。
4. （可接受检测）若$r_k\le \eta$，令$x^{(k+1)}=x^{(k)}，k=k+1$，转2；否则令$x^{(k+1)}=x^{(k)}+d^{(k)}，k=k+1$。转1。

注：算法中的常数1/4,3/4,1/2是根据经验选取的。实际计算时，可根据问题对它们进行调整。

四、信赖域子问题的求解

信赖域算法中子问题的求解是算法实现的关键。子问题是一个目标为二次函数的约束优化问题。下面介绍采用Euclid范数时，精确求解与非精确求解子问题的特殊算法。为了方便起见，省略迭代指标，用$x\in R^n$表示当前迭代点，此时信赖域子问题为：
$$\{\begin{array}{l}minf(x)+\nabla f(x{)}^{T}d+\frac{1}{2}{d}^{T}Bd\stackrel{△}{=}{q}_k(d)\\ \\ s.t.||d||\le \Delta \end{array}$$

1. 精确求解方法

不难发现，若$B$正定且$\bar{d}\stackrel{△}{=}-B^{-1}\nabla f(x)$满足$||d||\le \Delta$，即无约束问题：
$$minf(x)+\nabla f(x{)}^{T}d+\frac{1}{2}{d}^{T}Bd$$
的解是信赖域子问题的可行点，则$\bar{d}$是信赖域子问题的解。

定理：$d^{\ast }$是信赖域子问题的全局最优解当且仅当$d^{\ast }$可行，存在常数${\lambda }^{\ast }\ge 0$满足$B+{\lambda }^{\ast }I$半正定，且有：
$$\{\begin{array}{l}(B+{\lambda }^{\ast }I)d^{\ast }=-\nabla f(x)\\ \\ {\lambda }^{\ast }(\Delta -||d^{\ast }||)=0\end{array}$$
利用此定理构造求解子问题的算法。设矩阵$B+{\lambda }^{\ast }I$正定。若线性方程组：
$$Bd+\nabla f(x)=0$$
的解$\bar{d}$满足$||\bar{d}||\le \Delta$，则$d^{\ast }=\bar{d}$。此情形对应于${\lambda }^{\ast }=0$且$B$正定。否则，必有${\lambda }^{\ast }>0$。此时求解信赖域子问题等价于解如下方程组：
$$\{\begin{array}{l}(B+\lambda I)d=-\nabla f(x)\\ \\ ||d||=\Delta \end{array}$$
取$\lambda >0$充分大使得$B+\lambda I$正定，可通过解如下关于$\lambda$的一元非线性方程：
$${\varphi }_1(\lambda )=||(B+\lambda I{)}^{-1}\nabla f(x)||-\Delta =0$$
得到解${\lambda }^{\ast }$。然后得子问题的解$d^{\ast }=d({\lambda }^{\ast })=-(B+\lambda I{)}^{-1}\nabla f(x)$。

基于此求信赖域子问题的方法称为子问题的精确解法。利用矩阵的对角化分解可知，${\varphi }_1(\lambda )$是非线性程度高的系统。为了简化方程的计算，定义：
$${\varphi }_2(\lambda )=\frac{1}{\Delta }-\frac{1}{||d(\lambda )||}$$
${\varphi }_2(\lambda )$近似为线性方程系统。该方程可应用牛顿迭代法建立其迭代计算式为：
$${\lambda }_{l+1}={\lambda }_l-\frac{{\varphi }_2({\lambda }_l)}{{\varphi }_2^{\prime }({\lambda }_l)}={\lambda }_l+(\frac{||d_l||}{||q_l||}{)}^2(\frac{||d_l||-\Delta }{\Delta })$$
根据如上分析可建立信赖域子问题精确求解的计算步骤如下：

0. 给定${\lambda }_0>0,\Delta >0$。令$l=0$。
1. 若${\lambda }_l$是问题的解，则解线性问题得解$d^{(l)}$。否则，转2。
2. 作Cholesky（乔列斯基）分解$B+{\lambda }^{(l)}I=R^{T}R$。解方程组：
   $$R^{T}R{d}_l=-\nabla f(x)，R^T{q}_l=d_l$$
   得解$d_l,q_l$。
3. $${\lambda }_{l+1}={\lambda }_l+(\frac{||d_l||}{||q_l||}{)}^2(\frac{||d_l||-\Delta }{\Delta })$$
4. 令$l=l+1$。转1。

2. 折线方法（Dogleg Method 狗腿法）

上面介绍的信赖域子问题精确求解计算量较大。而且当$B+{\lambda }^{\ast }I$非正定时，子问题的求解较为复杂。

我们考虑非精确求解子问题，下面我们介绍非精确求解信赖域子问题的折线方法。有信赖域子问题知其解$d^{\ast }$是信赖域半径$\Delta$的函数，记为$d^{\ast }(\Delta )$。几何上为一条曲线，称其为最优解曲线。如下图：

折线法的思想是用一条折线代替精确解曲线$d^{\ast }(\Delta )$。为了构造合适的折线，首先分析子问题解的特性。若$B$正定且$||-B^{-1}\nabla f(x)||\le \Delta$，则精确解为：
$$d^{\ast }(\Delta )=-b^{-1}\nabla f(x)\stackrel{△}{=}{d}^B$$
否则，我们可作如下分析。当$\Delta$很小时，信赖域子问题中目标函数的二次项的作用不大。因此，可用函数$f(x)$的一次（线性）近似，此时约束在$||d||\le \Delta$下，子问题的近似解为$d^{\ast }(\Delta )\approx -\Delta \frac{\nabla f(x)}{||\nabla f(x)||}$。这显示可采用最速下降方向作为解曲线的近似。基于最速下降的考虑和省略信赖域约束，我们取子问题解的形式为$d=-\tau \nabla f(x)$。通过目标函数极小化可确定子问题沿最速下降方向的解为：
$$d^U=-\frac{\nabla f(x{)}^{T}f(x)}{\nabla f(x{)}^{T}B\nabla f(x)}\nabla f(x)$$
另一方面，为了得到收敛性能更好的搜索方向，我们可以选择牛顿方向$d^B$。由两方向$d^U$和$d^B$可分段构造折线搜索方向（如上图）。记由上图方式构造出的折线为$\tilde{d}(\tau )$，数学上的定义为：
$$\tilde{d}(\tau )=\{\begin{array}{rlr}& \tau d^U,& 0\le \tau \le 1\\ & & \\ & d^U+(\tau -1)(d^B-d^U),& 1\le \tau \le 2\end{array}$$
上式说明当$\tau$较小时方向选用最速下降方向$d^U$，否则取$d^U$和$d^B$的组合方向。