数值优化（Numerical Optimization）学习系列-信赖域方法

信赖域方法和线搜索类似都是迭代方法，与其不同的是，每次迭代时，在一个选定的可信赖区域内，选择当前迭代点的近似模型 mkm_k ，然后计算最优步长；如果步长不合适，可以对区域进行缩放。该小结主要介绍：

信赖域方法的基本形式

求解信赖域的基础方法

信赖域方法的收敛性和收敛速度

信赖域方法的扩展

信赖域方法的基本形式

在信赖域方法中，可信赖的区域（Region）的选择很重要，一般都会根据上一步结果进行动态变化；如果上一步太大则缩小，否则进行扩大。

在模型最优化问题中，选择TR方法比LS方法能够较快的收敛，例如


在该例子中，在非凸函数F中，当前步骤TR方法要优于LS。

信赖域方法有几个参数需要选择：

1. 近似模型 mkm_k

2. 可信赖区域 Δk\Delta_k

3. 求解参数 pkp_k

基本形式

在本节中模型选择为二次近似模型，采用函数二阶泰勒展开，即f(xk+p)=fk+gTkp+12pT∇2fkpf(x_k+p)=f_k + g_k^Tp + \frac 12 p^T \nabla^2 f_k p一般情况下会用 BkB_k 去近似Hessian矩阵，即mk=fk+gTk+12pTBkpm_k = f_k + g_k^T + \frac 12 p^T B_k p其中 BkB_k为对称矩阵。

信赖域的基本形式为：min mk(p)=fk+gTk+12pTBkps.t ||p||≤Δkmin \ m_k(p)=f_k + g_k^T + \frac 12 p^T B_k p \quad s.t \ ||p|| \le \Delta_k

该问题为关于p的带约束的最优化问题，参数p被限制在一个球形区域内。如果BkB_k选择为Hessian，则为TR的牛顿方法。

如果||B−1kgk||≤Δk||B_k^{-1}g_k|| \le \Delta_k则 pk=−B−1kgkp_k=-B_k^{-1}g_k为完全步(Full Step)，即球形约束没有作用。

Δk\Delta_k的选择

参数Δk\Delta_k的选择一般会根据上一步的结果进行调整，定义ρk=fk−f(xk+pk)mk(0)−mk(pk)\rho_k=\frac{f_k-f(x_k+p_k)}{m_k(0)-m_k(p_k)}其中分子表示函数实际减小的值；分母表示近似模型减少的值。分析 ρk\rho_k:

1. 如果 ρk\rho_k 小于0，一般情况下分母不可能小于0，因为目标函数求解的是最小值；此时说明分子小于0，即下一个目标点比上一步大，此时需要舍弃。

2. 如果 ρk\rho_k 大于0并且接近1，说明模型和实际的预期比较相符合，此时可以考虑扩大Δk\Delta_k

3. 如果ρk\rho_k大于0但是明显小于1，此时可以不用调整

4. 如果ρk\rho_k大于0，但是接近0，说明模型变化范围比较大，但是实际改变比较小，此时应该收缩或者减少Δk\Delta_k

具体算法如下：

子问题的最优解

为简化形式，将信赖域问题的子问题表示为：min m(p)=f+gTp+12pTBps.t ||p||≤Δmin \ m(p)=f + g^Tp + \frac 12 p^T B p \quad s.t \ ||p|| \le \Delta

该问题为标准的带不等式约束的二次优化问题，可以根据KKT条件（后面会深入介绍）得到该问题的最优解

定理

如果向量 p∗p^* 为子问题的最优解，当且仅当满足 p∗p^* 为可行解，并且存在标量 λ\lambda 满足

(B+λI)=−g(a)λ(Δ−||p∗||)=0(b)(B+λI)正定(c)\begin{align}
& (B+\lambda \mathbf I)=-g \quad (a)\\
& \lambda(\Delta-||p^*||)=0 \quad (b)\\
& (B+\lambda \mathbf I) 正定\quad (c)
\end{align}

其中条件（b）为补充条件，即要么λ为0要么Δ=||p∗||\lambda为0要么\Delta=||p^*||。

从下图中可以看出最优解和参数λ\lambda的关系



当参数Δ=Δ1\Delta=\Delta^1时，最优解为p3p^3此时相当于没有约束，此时λ=0\lambda=0

当参数Δ=Δ1或者Δ2\Delta=\Delta^1或者\Delta^2时，最优解被球形约束限制，此时满足Δ=||p∗||\Delta=||p^*||，根据上面条件（a）有λp∗=−Bp∗−g\lambda p^*=-Bp^*-g

如果能够找到这样的p满足这些条件就能找到最优解。

基于柯西点（Cauchy-Point）的算法

在实际求解过程中不一定找到最优解，而是找到一个充分下降的点满足全局收敛即可。柯西点就是满足该条件的点（pckp_k^c）

柯西点（Caychy-Point）算法

求解步骤如下

1. 计算原子问题的线性蜕化问题，寻找向量pskp_k^s满足psk=arg min fk+gTkps.t. ||p||≤Δkp_k^s=arg \ min \ f_k+g_k^T p \quad s.t. \ ||p|| \le \Delta_k

2. 寻找标量τk>0\tau_k > 0满足 min mk(τpsk)min \ m_k(\tau p_k^s) 同时满足信赖域的约束，即τk=arg min mk(τpsk)s.t. ||τpsk||≤Δk\tau_k = arg \ min \ m_k(\tau p_k^s) \quad s.t. \ ||\tau p_k^s|| \le \Delta_k

3. pck=τkpskp_k^c = \tau_k p_k^s

分别解释如下：

1. 计算pskp_k^s 从上述步骤1中可以看出求解步骤1有必使解，psk=−Δk||gk||gkp_k^s = - \frac {\Delta_k}{||g_k||} g_k。两种思路，一是退化后的函数为线性函数，而且是递减的，只要取下界即可。二是利用KKT条件也可以推出。

2. 将求解得到的pskp_k^s代入子问题可以得到，min m(p)=f+gTp+12pTBps.t ||p||≤Δmin m(τpsk)=f−gTτΔk||gk||gk+12τ2Δk||gk||gTkBkΔk||gk||gks.t. ||τpsk||≤Δkmin \ m(p)=f + g^Tp + \frac 12 p^T B p \quad s.t \ ||p|| \le \Delta \\ min \ m(\tau p_k^s)=f - g^T \tau \frac {\Delta_k}{||g_k||} g_k + \frac 12 \tau ^2\frac {\Delta_k}{||g_k||} g_k^T B_k \frac {\Delta_k}{||g_k||} g_k \\ s.t. \ ||\tau p_k^s|| \le \Delta_k是一个关于τ\tau 的二次函数，如果

1）gTkBkgk≤0g_k^T B_k g_k \le 0，是一个关于τ\tau的递减函数，直接得到τ=1\tau = 1

2)如果gTkBkgk>0g_k^T B_k g_k > 0 求导可以得到τk=||gk||3ΔkgTkBkgk,同时τ需要满足τ\e1\tau_k = \frac {||g_k||^3}{\Delta_k g_k^T B_k g_k},同时\tau需要满足 \tau \e 1

柯西点很容易计算，但是如果只利用柯西点，相当于只利用了梯度方向，相当于线搜索的扩展而已，即收敛速度为线性收敛。

Dogleg算法

该方法适用于当Bk正定时B_k正定时，寻找半径Δ\Delta和最优解p∗(Δ)p^*(\Delta)之间的关系。

在之前定义了完全步(Full Step)，即||B−1kgk||≤Δk||B_k^{-1}g_k|| \le \Delta_k则 pBk=−B−1kgkp_k^B=-B_k^{-1}g_k

该算法的思路为，p∗(Δ)p^*(\Delta)表示不同Δ\Delta条件下的最优解，

1）p∗(Δ)=pB||pB||≤Δp^*(\Delta) = p^B \quad ||p^B|| \le \Delta，否则

2）p(τ)=⎧⎩⎨τpu,pu+(τ−1)(pB−pu),0<τ<11<τ<2p(\tau) =
\begin{cases}
\tau p^u, & 0 < \tau < 1 \\[2ex]
p^u + (\tau-1)(p^B-p^u), & 1<\tau < 2
\end{cases} 其中pup^u定义为pu=−gTggTBggp^u=-\frac{g^Tg}{g^TBg} g，沿着梯度方向的最优解。

此为Dogleg方法，即寻找一个折线即两个线段的交点，即一条线段是沿着负梯度方向的最优值；二是沿着pU到pBp^U 到 p^B方向。在此折线上寻找最优解。

由定理可以证明||p(τ)||是一个关于τ的递增函数，而m(p(τ))是一个关于τ的递减函数||p(\tau)||是一个关于\tau的递增函数，而m(p(\tau))是一个关于\tau的递减函数，又是一个线性函数，因此只要计算p(τ)和||p||=Δp(\tau)和||p||=\Delta的交点即可。从下图可以清晰看到

二维子空间最优化

在DogLeg算法中，可以理解为最优解p∗p^*可以表示为pU和pB的扩展子空间，即p∗=λ1g+λ2B−1gp^U和p^B的扩展子空间，即p^* = \lambda_1g + \lambda_2B^{-1}g，则原子问题可以表示为：

min m(p)=f+gTp+12pTBps.t ||p||≤Δp∈span[g,−B−1g]min \ m(p)=f + g^Tp + \frac 12 p^T B p \quad \\ s.t \ ||p|| \le \Delta \\ p \in span[g, -B^{-1}g]

迭代算法

根据子问题的最优解形式可以得到

1）当λ=0时\lambda=0时需要满足Bp∗=−gB正定||p∗||≠ΔB p^* = -g \\ B 正定 \\ ||p^*|| \ne \Delta

2) 当λ≠0时\lambda\ne0时需要满足(B+λI)p∗=−gB+λI正定||p∗||=Δ(B+\lambda \mathbf I)p^* = -g \\ B+\lambda \mathbf I 正定 \\ ||p^*|| = \Delta即||p∗||=||−(B+λI)−1g||=Δ||p^*||=||-(B+\lambda \mathbf I)^{-1}g||=\Delta

通过定义p(λ)=−(B+λI)−1gp(\lambda)=-(B+\lambda \mathbf I)^{-1}g，寻找特定的λ使得||p(λ)||=Δ\lambda使得||p(\lambda)||=\Delta

相关定义

由于B正定因此可以进行正交分解,B=QΛQT； Λ=diag(λ1,λ2,...,λn)并且λ1≤λ2≤...≤λnB=Q\Lambda Q^T；\ \Lambda=diag(\lambda_1,\lambda_2,...,\lambda_n)并且\lambda_1 \le \lambda_2 \le ... \le \lambda_n

B+λI=Q(Λ+λI)QTp(λ)=Q(Λ+λI)−1QTgp(λ)2=−∑1n(qTjg)2(λj+λ)2B+\lambda \mathbf I = Q (\Lambda + \lambda \mathbf I)Q^T \\ p(\lambda)=Q (\Lambda + \lambda \mathbf I)^{-1}Q^T g \\
p(\lambda)^2=-\sum_1^n{\frac {(q_j^Tg)^2}{(\lambda_j+\lambda)^2}}由于是正交分解因此qiqj=1q_i q_j=1

迭代算法

1）由上述定义p(λ)2=−∑1n(qTjg)2(λj+λ)2p(\lambda)^2=-\sum_1^n{\frac {(q_j^Tg)^2}{(\lambda_j+\lambda)^2}}可以看出

2） 当λ>−λ1时\lambda > - \lambda_1时函数值是一个单调递减函数，特别当limλ→∞(p(λ)2)→0\lim_{\lambda \to \infty}(p(\lambda)^2) \to 0

3)因此我们需要寻找λ∈(−∞,λ1)使得||p||=Δ\lambda \in (-\infty, \lambda_1)使得||p||=\Delta

λ\lambda求解

在qTig≠0q_i^Tg \ne 0

1）通过牛顿迭代法求解ϕ(λ)=||p(λ)||−Δ=0\phi(\lambda)=||p(\lambda)||-\Delta=0

2）近似算法求解ϕ1(λ)=C1λ+λ1+C2phi1(λ)=1Δ−λ+λ3C3\phi_1(\lambda)=\frac {C1}{\lambda+\lambda_1}+C2 \\
phi_1(\lambda)=\frac {1}{\Delta}-\frac{\lambda+\lambda_3}{C3}

当qTig=0q_i^Tg = 0是一个Hard Case，此时可以添加一个误差因子进行近似

p(λ)=Q(Λ+λI)−1QTg+τEp(\lambda)=Q (\Lambda + \lambda \mathbf I)^{-1}Q^T g + \tau \mathbf E

信赖域的其他扩展

poor scaling问题，可以扩展为球形或者椭圆形

可以构造Scaled的算法min m(p)=f+gTp+12pTBps.t ||Dp||≤Δmin \ m(p)=f + g^Tp + \frac 12 p^T B p \quad s.t \ ||Dp|| \le \Delta

矩阵D的构造和∂2f/∂xi\partial ^2f/\partial x_i相关

可以构造通用的柯西点废

其他Norm方法，例如||p||1≤Δ||p||∞≤Δ||p||_1 \le \Delta \\||p||_\infty \le \Delta

总结

通过该节需要了解

1. 信赖域方法和线搜索方法的不同

2. 信赖域方法的基本形式

3. 信赖域方法的柯西点算法、DogLeg算法和最优解迭代算法

4. 信赖域方法收敛