非线性优化

非线性优化

非线性最小二乘

  最简单最小二乘问题:
$$\underset{x}{\min }\frac{1}{2}\parallel f\left(x\right){\parallel }_2^2$$
  这里自变量$x\in {\mathbb{R}}^n$，$f$是任意非线性函数，设$f\left(x\right)$有$m$维：$f\left(x\right)\in {\mathbb{R}}^m$。如果$f$是个数学形式上很简单的函数，令目标函数的导数为零，然后求解$x$的最优值，就和求二元函数的极值一样：
$$\frac{df}{dx}=0$$
  解此方程，就得到了导数为零处的极值。可能是极大值、极小值或者鞍点处的值，只要逐个比较它们的函数值大小即可。
  对于不方便直接求解的最小二乘问题，可以用迭代的方式，从一个初始值出发，不断更新当前的最优化变量，是目标函数下降。

1. 给定某个初始值$x_0$
2. 对于第$k$次迭代，寻找一个增量$\Delta x_k$，使得$\parallel \left(x_k+\Delta x_k\right){\parallel }_2^2$达到极小值。
3. 若$\Delta x_k$足够小，则停止。
4. 否则，令$x_{k+1}=x_k+\Delta x_k$，返回第2步

一阶和二阶梯度法

  求解增量最直观的方式是将目标函数在$x$附近进行泰勒展开：
$$\parallel f\left(x+\Delta x\right){\parallel }_2^2\approx \parallel f\left(x\right){\parallel }_2^2+J\left(x\right)\Delta x+\frac{1}{2}{\Delta x}^{T}H\Delta x$$
  $J$是$\parallel f\left(x\right){\parallel }_2^2$关于$x$的导数（雅克比矩阵），而$H$则是二阶导数（海塞矩阵）。可以保留泰勒展开的一阶或二阶项，对应的求解方法则为一阶梯度或二阶梯度法。
  如果保留一阶梯度，那么增量的解就为：
$${\Delta x}^{\ast }=-J^T\left(x\right)$$
  直观意义非常简单，只要沿着反向梯度方向前进即可。
  如果
$${\Delta x}^{\ast }=-\lambda J^T\left(x\right)$$
  称为最速下降法。
  如果保留二阶梯度信息，那么增量方程为：
$${\Delta x}^{\ast }=argmin\parallel f\left(x\right){\parallel }_2^2+J\left(x\right)+\frac{1}{2}{\Delta x}^{T}H\Delta x$$
  对右侧等式关于$\Delta x$的导数并令它为零，就得到了增量的解：
$$H\Delta x=-J^T$$
  该方法又称为牛顿法。

高斯牛顿法

  高斯牛顿的思想是将$f\left(x\right)$进行一阶的泰勒展开：
$$f\left(x+\Delta x\right)\approx f\left(x\right)+J\left(x\right)\Delta x$$
  $J\left(x\right)$为$f\left(x\right)$关于$x$的导数，实际上是一个$m\times n$矩阵，也是一个雅可比矩阵。根据前面的思路，当前的目标是寻找一个下降矢量$\Delta x$，使得$\parallel f\left(x+\Delta x\right){\parallel }^2$达到最小。为了求$\Delta x$，需要解一个线性的最小二乘问题：
$${\Delta x}^{\ast }=arg\underset{\Delta x}{\min }\frac{1}{2}\parallel f\left(x\right)+J\left(x\right)\Delta x{\parallel }^2$$
  先展开目标函数的平方项：
$$\parallel f\left(x\right)+J\left(x\right)\Delta x{\parallel }^2={\left(f\left(x\right)+J\left(x\right)\Delta x\right)}^T\left(f\left(x\right)+J\left(x\right)\Delta x\right)$$
$$=\parallel f\left(x\right){\parallel }_2^2+2{f\left(x\right)}^{T}J\left(x\right)\Delta x+{\Delta x}^T{J\left(x\right)}^{T}J\left(x\right)\Delta x$$
  求上式关于$\Delta x$的导数，并令其为零：
$$2{J\left(x\right)}^{T}f\left(x\right)+2{J\left(x\right)}^{T}J\left(x\right)\Delta x=0$$
$$\Rightarrow {J\left(x\right)}^{T}J\left(x\right)\Delta x=-{J\left(x\right)}^{T}f\left(x\right)$$
  要求解的变量是$\Delta x$，因此这是一个线性方程组，称为增量方程也可以称为高斯牛顿方程（Gauss Newton equation）或者正规方程（Normal equation）。
  对比牛顿法用$J^{T}J$作为牛顿法中的二阶Hessian矩阵的近似，从而省略了计算$H$的过程。求解增量方程是整个优化问题的核心所在。高斯牛顿法的算法步骤可以写成：

1. 给定初始值$x_0$。
2. 对于第$k$次迭代，求出当前的雅可比矩阵$J\left(x_k\right)$和误差$f\left(x_k\right)$。
3. 求解增量方程：${J\left(x\right)}^{T}J\left(x\right)\Delta x=-{J\left(x\right)}^{T}f\left(x\right)$
4. 若$\Delta x_k$足够小，则停止。否则，令$x_{k+1}=x_k+\Delta x_k$，返回第2步。

列文伯格—马夸尔特方法

  高斯牛顿法中采用的近似二阶泰勒展开只能在展开点附近有较好的近似效果，所以很自然的想到应该给$\Delta x$添加一个信赖区域（Trust Region），不能让它太大而使得近似不准确。非线性优化中有一系列这类方法，这类方法被称为信赖区域方法（Trust Region Method），在信赖区域里边认为近似是有效的；出了这个区域，近似可能会出问题。
  如何确定这个信赖区域范围呢？一个较好的方法是根据近似模型跟实际函数之间的差异来确定：如果差异小，就让范围尽可能大；如果差异大，就缩小这个近似范围。
$$\rho =\frac{f\left(x+\Delta x\right)-f\left(x\right)}{J\left(x\right)\Delta x}$$
  使用上式来判断泰勒近似是否够好。$\rho$的分子是实际函数下降的值，分母是近似模型下降的值。如果$\rho$接近于1，则近似是好的。如果$\rho$太小，说明实际减小的值远小于近似减小的值，认为近似较差，需要缩小近似范围。反之，如果$\rho$比较大，则说明实际下降的比预计的跟大，可以放大近似范围。步骤如下：

1. 给定初始值$x_0$以及初始优化半径$\mu$。
2. 对于第$k$次迭代，求解：
   $$\underset{{\Delta x}_k}{\min }\frac{1}{2}\parallel f\left(x_k\right)+J\left(x_k\right){\Delta x}_k{\parallel }^2，s.t.\parallel D{\Delta x}_k{\parallel }^2\le \mu$$
3. 计算$\rho$
4. $\rho >\frac{3}{4}，\mu =2\mu$
5. $\rho <\frac{1}{4}，\mu =0.5\mu$
6. 如果$\rho$大于某阈值，则认为近似可行。令$x_{k+1}=x_k+{\Delta x}_k$
7. 判断算法是否收敛。不收敛回第2步，否则结束。

  近似范围扩大的倍数和阈值都是经验值，可以调节。在第2步中，我们把增量限定于一个半径为$\mu$的球中，认为只是在这个球内才是有效的。带上$D$之后这个球可以看成一个椭球。在列文伯格提出的优化方法中，把$D$取成单位阵$I$，相当于直接把$\Delta x$约束在一个球中。随后，马夸尔提出将$D$取成非负数对角阵-----实际中通常用$J^{T}J$的对角元素平方根，使得在梯度小的维度上约束范围更大一些。
  步骤2中是带不等式约束的优化问题，使用拉格朗日乘子将它转化为一个无约束优化问题：
$$\underset{{\Delta x}_k}{\min }\frac{1}{2}\parallel f\left(x_k\right)+J\left(x_k\right){\Delta x}_k{\parallel }^2+\frac{\lambda }{2}\parallel D\Delta x{\parallel }^2$$
   $\lambda$为拉格朗日乘子。类似于高斯牛顿法中的做法，展开后发现问题的核心乃是计算增量的线性方程：
$${J\left(x\right)}^{T}J\left(x\right)\Delta x+\lambda D^{T}D\Delta x=-{J\left(x\right)}^{T}f\left(x\right)$$
$$\Rightarrow \left(H\left(x\right)+\lambda D^{T}D\right)\Delta x=-{J\left(x\right)}^{T}f\left(x\right)$$
  令$g=-{J\left(x\right)}^{T}f\left(x\right)$
$$\Rightarrow \left(H\left(x\right)+\lambda D^{T}D\right)\Delta x=g$$
  根据上式可以得到，当$\lambda$比较大时，$H$占主要地位，说明二次近似模型在该范围内是比较好的，列文伯格—马夸尔特法方法更接近于高斯牛顿法。另一方面，当$\lambda$比较大时，$\lambda I$占据主要地位，列文伯格—马夸尔特法更接近于一阶梯度下降法（即最速下降），这说明附近的二次近似不够好。列文伯格—马夸尔特方法的求解方式，可在一定程度上避免线性方程组的系数矩阵的非奇异和病态问题，提供更稳定，更准确的增量$\Delta x$