ICode9

精准搜索请尝试: 精确搜索
首页 > 其他分享> 文章详细

牛顿法

2021-10-26 02:00:41  阅读:248  来源: 互联网

标签:begin end equation 牛顿 delta lambda


牛顿法和拟牛顿法

最早接触牛顿法是本科时候学《数值计算方法》,现在重新捡起来。

牛顿法

考虑无约束最优化问题

\[\begin{equation} \min \limits _{x \in \mathbb{R} ^ n} f(x) \end{equation} \]

其中, \(x^*\) 为目标函数的极小点。

假设 \(f(x)\) 具有二阶连续偏导数,若第 \(k\) 次迭代值为\(x^(k)\),则可将 \(f(x)\) 在 \(x^{(k)}\) 附近进行二阶泰勒展开:

\[\begin{equation} f(x) = f(x^{(k)}) + g^T_k(x-x^{(k)})+1/2(x-x^{(k)})^T H(x^{(k)})(x-x^{(k)}) \end{equation} \]

这里, \(g_k = g(x^{(k)})=\nabla f(x^{(k)})\) 是 \(f(x)\) 的梯度向量在点 \(x^{(x)}\) 的值, \(H(x^{(k)})\) 是 \(f(x)\) 的黑塞矩阵(hessian matrix):

\[\begin{equation} H(x) = \begin{bmatrix} \frac{\partial^2 f}{\partial x_i \partial x_j} \end{bmatrix} _{n \times n} \end{equation} \]

在点 \(x^{(k)}\) 的值。函数 \(f(x)\)有极值的必要条件是在极值点处一阶导数为 0,即梯度向量为 0。当 \(H(x^{(k)})\) 是正定矩阵时,函数 \(f(x)\) 的极值为极小值。
牛顿法利用极值点的必要条件:

\[\begin{equation} \nabla f(x) = 0 \end{equation} \]

每次迭代中从 \(x^{(k)}\) 开始,求目标函数的极小点,作为 \(k+1\) 次迭代值 \(x^{(k+1)}\) 。具体地, 假设 \(x^{(k+1)}\) 满足:

\[\begin{equation} \nabla f(x^{(k+1)})=0 \end{equation} \]

由式 \((2)\) 有(经过移项变换)

\[\begin{equation} \nabla f(x)=g_{k}+H_{k}(x-x^{(k)}) \end{equation} \]

其中 \(H_{k}=H(x^{(k)})\) 。这样, 式 (5) 成为

\[\begin{equation} g_{k}+H_{k}(x^{(k+1)}-x^{(k)})=0 \end{equation} \]

因此,

\[\begin{equation} x^{(k+1)}=x^{(k)}-H_{k}^{-1} g_{k} \end{equation} \]

或者

\[\begin{equation} x^{(k+1)}=x^{(k)}+p_{k} \end{equation} \]

其中,

\[\begin{equation} H_{k} p_{k}=-g_{k} \end{equation} \]

公式 \((8)\) 作为迭代公式即牛顿法。

算法:
输入:目标函数 \(f(x)\), 梯度 \(g(x)=\nabla f(x)\), 黑塞矩阵 \(H(x)\), 精度要求 \(\varepsilon\);
输出:\(f(x)\) 的极小点 \(x^{*}\) 。

  • (1) 取初始点 \(x^{(0)}\), 置 \(k=0\) 。
  • (2) 计算 \(g_{k}=g(x^{(k)})\) 。
  • (3) 若 \(||g_{k}||<\varepsilon\), 则停止计算, 得近似解 \(x^{*}=x^{(k)}\) 。
  • (4) 计算 \(H_{k}=H(x^{(k)})\), 并求 \(p_{k}\)

\[H_{k} p_{k}=-g_{k} \]

  • (5) 置 \(x^{(k+1)}=x^{(k)}+p_{k}\) 。
  • (6) 置 \(k=k+1\), 转 \((2)\) 。
    步骤 (4) 求 \(p_{k}, p_{k}=-H_{k}^{-1} g_{k}\), 要求 \(H_{k}^{-1}\), 计算比较复杂, 所以有其他改进的方法。

拟牛顿法

针对上述步骤(4)中计算复杂度的问题,考虑使用一个n阶矩阵 \(G_k=G(x_{(k)})\) 来代替 \(H^{-1}_k\) 。
在牛顿 \((6)\) 中取 \(x=x^{(k+1)}\), 即得

\[\begin{equation} g_{k+1}-g_{k}=H_{k}(x^{(k+1)}-x^{(k)}) \end{equation} \]

记 \(y_{k}=g_{k+1}-g_{k}, \delta_{k}=x^{(k+1)}-x^{(k)}\), 则

\[\begin{equation} y_{k}=H_{k} \delta_{k} \ \text{或} \ H_{k}^{-1} y_{k}=\delta_{k} \end{equation} \]

称为拟牛顿条件

如果 \(H_{k}\) 是正定的 \((H_{k}^{-1}\) 也是正定的),那么可以保证牛顿法搜索方向 \(p_{k}\) 是下降 方向。这是因为搜索方向是 \(p_{k}=-H_{k}^{-1} g_{k}\),由式 \((B.8)\) 有

\[\begin{equation} x=x^{(k)}+\lambda p_{k}=x^{(k)}-\lambda H_{k}^{-1} g_{k} \end{equation} \]

所以 \(f(x)\) 在 \(x^{(k)}\) 的泰勒展开式 \((B .2)\) 可以近似写成:

\[\begin{equation} f(x)=f(x^{(k)})-\lambda g_{k}^{\mathrm{T}} H_{k}^{-1} g_{k} \end{equation} \]

因 \(H_{k}^{-1}\) 正定, 故有 \(g_{k}^{\mathrm{T}} H_{k}^{-1} g_{k}>0\) 。当 \(\lambda\) 为一个充分小的正数时, 总有 \(f(x)<f(x^{(k)})\), 也就是说 \(p_{k}\) 是下降方向。
拟牛顿法将 \(G_{k}\) 作为 \(H_{k}^{-1}\) 的近似, 要求矩阵 \(G_{k}\) 满足同样的条件。首先, 每次迭 代矩阵 \(G_{k}\) 是正定的。同时, \(G_{k}\) 满足下面的拟牛顿条件:

\[\begin{equation} G_{k+1}\ y_{k}=\delta_{k} \end{equation} \]

根据拟牛顿条件选择 \(G_k\) 作为 \(H^{-1}_k\) 的近似,或选择 \(B_k\) 作为 \(H_k\) 的近似的算法为拟牛顿法。不同的合理的选择有不同的实现方法,接下来只介绍BFGS(Broyden-Fletcher-Goldfarb-Shanno),一个比较流行的算法。

假设每一步的迭代过程 \(B_{k+1}\) 是由 \(B_k\) 加上两个附加项构成的,即:

\[\begin{equation} B_{k+1}=B_k + P_k + Q_k \end{equation} \]

且满足拟牛顿条件:

\[\begin{equation} B_{k+1} \delta _k =B_k \delta _k + P_k \delta _k + Q_k \delta _k \end{equation} \]

为了满足上述拟牛顿条件,可以令:

\[\begin{equation} \begin{aligned} P_{k} \delta_{k}&=y_{k} \\ Q_{k} \delta_{k}&=-B_{k} \delta_{k} \end{aligned} \end{equation} \]

找出适合条件的 \(P_{k}\) 和 \(Q_{k}\), 得到 BFGS 算法矩阵 \(B_{k+1}\) 的迭代公式:

\[\begin{equation} B_{k+1}=B_{k}+\frac{y_{k} y_{k}^{\mathrm{T}}}{y_{k}^{\mathrm{T}} \delta_{k}}-\frac{B_{k} \delta_{k} \delta_{k}^{\mathrm{T}} B_{k}}{\delta_{k}^{\mathrm{T}} B_{k} \delta_{k}} \end{equation} \]

可以证明, 如果初始矩阵 \(B_{0}\) 是正定的, 则迭代过程中的每个矩阵 \(B_{k}\) 都是正定的。

算法:
输入:目标函数 \(f(x), g(x)=\nabla f(x)\), 精度要求 \(\varepsilon\);
输出: \(f(x)\) 的极小点 \(x^{*}\) 。

  • (1) 选定初始点 \(x^{(0)}\), 取 \(B_{0}\) 为正定对称矩阵, 置 \(k=0\) 。
  • (2) 计算 \(g_{k}=g\left(x^{(k)}\right)\) 。若 \(||g_{k}||<\varepsilon\), 则停止计算, 得近似解 \(x^{*}=x^{(k)}\); 否则转 (3)。
  • (3) 由 \(B_{k} p_{k}=-g_{k}\) 求出 \(p_{k}\) 。
  • (4) 一维捜索:求 \(\lambda_{k}\) 使得

\[f\left(x^{(k)}+\lambda_{k} p_{k}\right)=\min _{\lambda \geqslant 0} f\left(x^{(k)}+\lambda p_{k}\right) \]

  • (5) 置 \(x^{(k+1)}=x^{(k)}+\lambda_{k} p_{k}\) 。
  • (6) 计算 \(g_{k+1}=g\left(x^{(k+1)}\right)\), 若 \(||g_{k+1}||<\varepsilon\), 则停止计算, 得近似解 \(x^{*}=x^{(k+1)}\); 否则, 按式 (B.30) 算出 \(B_{k+1}\) 。
  • (7) 置 \(k=k+1\), 转 (3)。
# newton algorithm written in python
def f(x):
    return math.exp(x) + x**2 + x*3 + 5.0

def df(x):
    return math.exp(x) + 2*x + 3.0

def ddf(x):
    return math.exp(x) + 2.0
    
def newton_optimize(x_init=0, eps=0.01):
    epoch = 1
    x_0 = x_init
    x_1 = x_0 - df(x_0)/ddf(x_0)
    while abs(x_1-x_0) > eps:
        x_0 = x_1
        x_1 = x_0 - df(x_0)/ddf(x_0)
        epoch+=1
    return x_1, epoch

if __name__ == "__main__":
    x_opt, epoch = newton_optimize(0, 1e-4)
    print(epoch, x_opt, f(x_opt))

标签:begin,end,equation,牛顿,delta,lambda
来源: https://www.cnblogs.com/clearhanhui/p/15464034.html

本站声明: 1. iCode9 技术分享网(下文简称本站)提供的所有内容,仅供技术学习、探讨和分享;
2. 关于本站的所有留言、评论、转载及引用,纯属内容发起人的个人观点,与本站观点和立场无关;
3. 关于本站的所有言论和文字,纯属内容发起人的个人观点,与本站观点和立场无关;
4. 本站文章均是网友提供,不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属;如您发现该文章侵犯了您的权益,可联系我们第一时间进行删除;
5. 本站为非盈利性的个人网站,所有内容不会用来进行牟利,也不会利用任何形式的广告来间接获益,纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。

专注分享技术,共同学习,共同进步。侵权联系[81616952@qq.com]

Copyright (C)ICode9.com, All Rights Reserved.

ICode9版权所有