精华内容
下载资源
问答
  • 事实上这个方法(思想)在实际中是被用于求解线性方程的,当然单纯的梯度下降形式并没有被直接采用,被广泛用于求解对称、正系数矩阵方程组的方法是,基于梯度下降原理实现的共轭梯度法,这一方法在大...

    共轭梯度法(1)--的基本原理

    之前已经搞明白了,梯度下降法的基本原理,当然解释的调度是从求函数极值的角度出发的,事实上从这个角度来理解,个人感觉是一个最为直接的理解角度,其完完全全是建立在多变量函数的微分系统中的。事实上这个方法(思想)在实际中是被用于求解线性方程的,当然单纯的梯度下降形式并没有被直接采用,被广泛用于求解对称、正系数矩阵方程组的方法是,基于梯度下降原理实现的共轭梯度法,这一方法在大型稀疏线性方程组的求解中,有极高的效率。但是,之前对于梯度下降思想的描述,好像认为这个方法是用于求函数极值的方法,似乎与线性方程组求解问题无关。因此,对于对称、正定系数矩阵的二次型的理解是很必要的,它能让我们明白矩阵系统中线性方程组与函数极值的关系,能够解释为什么一个求函数极值的思想可以被用于求解线性方程组系统。

    $Ax = b$对应的二次型

    这里就不严格去给出严格的数学表示了,因为对于非数学专业的我们来说,好像没啥必要,知晓核心的思想更为重要。注意本次针对的系数矩阵$A$都是对称、正定的,线性方程组系统$Ax = b$对应的是二次型对于自变量导数为0时的表达式,不严谨地可将其看作是函数吧,即二次型为

    $$

    f(x) = \frac{1}{2} x^{T} A x - b^{T} x + c \tag{1}

    $$

    对于$Ax = b$来说,$c$是一个0向量。我们知道求$f(x)$极小值的方法最直接的一种就是利用其导数为0,获得极值点从而获得极小值,因此,求导得到

    $$

    f'(x) = \frac{1}{2}A^T x + \frac{1}{2}A x - b = 0 \tag{2}

    $$

    因此,求解线性方程组$Ax = b$就等价于求解$f'(x) = 0$这个表达式,那么一阶导数为0这个表达换一种说法,不就是去求二次型$f(x)$的函数的极值嘛!这也就解释了为什么一个求解函数极值的方法却可以应用到线性方程组的求解当中来。

    为了更好的理解线性方程组$Ax = b$的求解与函数极值之间的关系,这里举一个简单的,可以被可视化的例子,若方程组为:

    $$

    \begin{bmatrix}

    2 & 2\\

    2 & 5

    \end{bmatrix}

    \begin{bmatrix}

    x\\

    y

    \end{bmatrix}

    =

    \begin{bmatrix}

    6\\

    3

    \end{bmatrix} \tag{3}

    $$

    这样一个对称正定的方程组,为了待求的未知数为了更符合函数表达,这里用了$x, y$来表示,依据二次型公式,得到

    $$

    f(x, y) = \frac{1}{2}

    \begin{bmatrix}

    x\\

    y

    \end{bmatrix}

    \begin{bmatrix}

    2 & 2\\

    2 & 5

    \end{bmatrix}

    \begin{bmatrix}

    x & y

    \end{bmatrix} -

    \begin{bmatrix}

    6 & 3

    \end{bmatrix}

    \begin{bmatrix}

    x\\

    y

    \end{bmatrix}

    $$

    这样一个$2 \times 2$的线性方程组对应的二次型就是一个二元函数,为了得到表达式,这里利用sympy库进行符号运算。import numpy as np

    from sympy import *

    # 先定义符号变量x, y

    x, y = symbols('x y')

    A = np.array([[2, 2], [2, 5]])

    b = np.array([6, 3])

    # 自变量向量这里用x1表示

    x1 = np.array([x, y])

    f_xy = 1/2 * x1.T @ A @ x1 - b.T @ x1

    init_printing(use_unicode = True) # 更美观的显示

    pprint(f_xy.expand())2 2

    1.0⋅x + 2.0⋅x⋅y - 6⋅x + 2.5⋅y - 3⋅y

    利用函数图像来直观地可视化函数极值情况。import matplotlib.pyplot as plt

    import numpy as np

    from mpl_toolkits.mplot3d import Axes3D

    x = np.linspace(-10, 10, 100)

    y = np.linspace(-10, 10, 100)

    xx, yy = np.meshgrid(x, y)

    f_xy = xx * (xx + yy) - 6 * xx + yy * (xx + 2.5 * yy) - 3 * yy

    # 这里用等值线云图与伪色彩图来可视化,也顺便对比下两者的区别

    fig, ax = plt.subplots(1, 2)

    cf1 = ax[0].pcolormesh(xx, yy, f_xy, cmap='RdBu_r')

    cs1 = ax[0].contour(x, y, f_xy, colors='gray', levels=15)

    ax[0].clabel(cs1, inline=True, colors='k')

    fig.colorbar(cf1, orientation='vertical')

    ax[0].set_title('function by pcolormesh')

    cf2 = ax[1].contourf(x, y, f_xy, cmap='RdBu_r')

    cs2 = ax[1].contour(x, y, f_xy, colors='gray', levels=15)

    ax[1].clabel(cs2, inline=True, colors='k')

    ax[1].set_title('function by contourf')

    fig.savefig('function2.jpg', dpi=600)

    plt.show()

    # 三维的函数图像

    fig1 = plt.figure()

    ax = Axes3D(fig1)

    cd = ax.plot_surface(xx, yy, f_xy, cmap='RdBu_r')

    ax.set_title('the function of 3D image')

    fig1.colorbar(cd)

    fig1.savefig('image3d.jpg', dpi=600)

    plt.show()

    通过函数图像的可视化,可以看到线性方程组本身是对称正定的,那么对应的二次型的函数图像则必定都位于在$z= 0 $平面以上,函数均是大于0的,且从函数的等值线及三维图像中可以看到是存在极小值的。明白了这个关系之后,对于梯度下降法或者共轭梯度法应用于方程组(系数矩阵为对称正定)的求解就不足为奇了。

    Gradient descent method求解方程组

    共轭梯度法的实现全来自于或者说是基于梯度下降法,因此,理解梯度下降法求解线性方程系统很有必要。在之前的梯度下降法的原理中,已经明白了这个方法对于函数极值的求解,是通过沿着梯度方向负向进行的,每一次的计算实质上获得是自变量的增量向量,而现在面对方程组时,由于之前的二次型的内容已经得出了线性方程组与多元函数极值之间的关系,即$f'(x_{i}) = Ax_{i} - b = 0$。那么梯度下降法可以被这样描述与构建:每一次产生的近似值$x_{i}$向量,近似解自然会有一个残差向量(residual vector)为$r_{i}$, 下表$i$表示第几步的迭代。

    $$

    r_{i} = b - A x_{i}\\

    r_{i} = - f'(x_{i}) \tag{4}

    $$

    式(4)中的$r_{i}$不就是梯度向量嘛,表示了点下一步搜索或者移动的方向,当然要想确定自变量的增量大小,还需要知道步长,这一点在梯度下降原理中已经给出了解释;在这里步长设定为$\alpha_{i}$,所以,产生的新的点$x_{i + 1}$为:

    $$

    x_{i + 1} = x_{i} + \alpha_{i} r_{i} \tag{5}

    $$

    到这里好像似乎与之前的原理中阐述并无大的区别,实际上现在问题在于如何确定步长$\alpha_i$,而在之前的介绍中步长是一个固定的,对于迭代计算来说,固定的步长的计算效率显然是不能够令人满意的。对于函数的极小值的确定来说,在这里每一步获得的新的函数值为$f(x_{i+1})$,对于整个求解过来说,只要这个更新的函数值达到最小值那么任务就算是完成了,依据式(5),函数值可以看做是关于步长的函数,为了获得极小值,通过导数为0来完成,即

    $$

    \frac{\mathrm{d}f'(x_{i + 1})}{\mathrm{d} \alpha_{i}} = 0 \tag{6}

    $$

    结合式(5)来看,这是一个复合函数,利用链式求导法则,得到

    $$

    式(6) = f'(x_{i+1}) ^{\mathrm{T}} \cdot \frac{\mathrm{d}x_{i + 1}}{\mathrm{d}\alpha_{i}} = 0 \tag{7}

    $$

    这个式子中,$f'(x)$不就是梯度向量嘛,依据式(4)给出的定义,在这里它也称作残差向量$-r_{i + 1}$,$x_{i + 1}$求导的结果可以根据式(5)为$r_{i}$,所以,要保证函数取得极值,就要使得前后两步的梯度向量(也就是残差向量$r$)保证正交!这里写为内积形式:

    $$

    r_{i + 1} ^ {\mathrm{T}} \cdot r_{i} = 0 \tag{8}

    $$

    那么依据这个关系,带入残差向量的表达,得到:

    $$

    (b- A x_{i + 1}) ^ {\mathrm{T}} \cdot r_{i} = 0\\

    (b - A (x + \alpha_{i} r_{i}))^{\mathrm{T}} \cdot r_{i} = 0

    \tag{9}

    $$

    通过这个式子就可以得到步长$\alpha_{i}$的计算公式:

    $$

    \alpha_{i} = \frac{r_{i} ^{\mathrm{T}} r_{i}}{r_{i} ^{\mathrm{T}} A r_{i}} \tag{10}

    $$

    到这里,所有量的表达都得到了解决,那么对于线性方程组系统来说,梯度下降法的迭代过程为:

    $$

    repeat:\\

    r_{i} = b - A x_{i}\\

    \alpha_{i} = \frac{r_{i}^{\mathrm{T}} r_{i}}{r_{i} ^{\mathrm{T}} A r_{i}}\\

    x_{i + 1} = x_{i} + \alpha_{i} r_{i}\\

    until: x_{i}的增量幅度满足设定精度\\

    end

    \tag{11}

    $$

    迭代计算的代码如下:import numpy as np

    # 每一步的残差向量r0都会被以行向量的形式存放在r矩阵中

    def gd_method(A, b, x0, tol, max_time):

    r0 = b - np.dot(A, x0)

    number = 0

    r = np.zeros([max_time, b.shape[0]])

    x = np.zeros([max_time, b.shape[0]])

    while np.linalg.norm(r0, np.inf) > tol:

    r0 = b - np.dot(A, x0)

    r[number, :] = r0

    alpha = np.dot(r0.T, r0) / np.dot(np.dot(r0,A), r0)

    x0 = x0 + alpha * r0

    x[number, :] = x0

    number += 1

    if number > max_time:

    print('warning:the result is not stable!')

    return

    return x, number, r

    这里对一个$3\times3$大小的对称正定矩阵进行验算,梯度下降算法迭代了67次得到结果,$[0.99983945,0.99976565, 1.99978575]$,而这个方程组的精确解为$[1,1,2]$。可以迭代结果还是比较可信与可靠的。A = np.array([[4, -2, -1], [-2, 4, -2], [-1, -2, 3]])

    b = np.array([0, -2, 3])

    x0 = np.array([1, 1, 1])

    x_gd = gd_method(A,b, x0, 0.0001,500)

    print(x_gd[0][x_gd[1]-1,:])

    print(x_gd[1])

    print(x_gd[2])[0.99983945 0.99976565 1.99978575]

    67

    [[-1. -2. 3. ]

    [-0.39130435 0.43478261 0.15942029]

    [ 0.09201888 0.02436289 0.15942029]

    ...

    [ 0. 0. 0. ]

    [ 0. 0. 0. ]

    [ 0. 0. 0. ]]x = np.linalg.solve(A,b)

    print(x)[1. 1. 2.]

    依据上述的分析以及这个算例的计算,实际上对于线性方程组迭代的迭代计算会有一个新的认识。首先,从最为原始的迭代的角度来看,通过每一迭代步中都会产生一个残差向量,似乎在方程组的求解角度来看,残差向量仅仅也就表示一个误差吧,好像也并没有很直观的一个认识。但是从梯度下降的算法角度来看,线性方程组迭代步中的残差向量实际上就是二次型(函数)的负梯度向量,那么也就是说,这个向量$r$不断地在修正点移动的方向,使得移动方向不断的向函数的极小值点处靠近!由于gd_method函数中会记录每一步的参差向量,因为这个是一个3个未知数的方程组,因此残差向量$r$本身就是一个$R^{3}$空间中的向量,同时返回得到的$x$矩阵中存放着每一步的近似解,这样进行可视化,可以看出梯度下降中前后的参差向量确实是正交的。因为三元函数无法进行可视化,不是很容易看出梯度下降的搜寻过程。可能二元函数的化对于理解更方便吧。xs = x_gd[0][0:x_gd[1]-1,:]

    plt.style.use('ggplot')

    fig = plt.figure()

    ax = Axes3D(fig)

    ax.plot(xs[:,0], xs[:,1], xs[:,2], marker='o')

    ax.set_title('point in every step')

    ax.set_xlabel('x')

    ax.set_ylabel('y')

    ax.set_zlabel('z')

    fig.savefig('gdimage.jpg', dpi=600)

    plt.show()

    在上述的分析中,我们会发现实际好像似乎梯度下降这个迭代算法的效率并不是很高,很简单的方程组都要迭代67次,显然这样的效率去求解大型的稀疏方程组是不能令人满意的。这个方法之所以没有被广泛地应用,主要原因是因为每一步的残差向量也就是梯度向量,都必须与相邻步的残差向量保持正交,如果为了便于理解,以二维空间为例,就是每一次的搜索(移动)方向都是相互垂直的,那么间隔的残差向量必定平行,也就说从初始值向精确解移动是呈折线的方式,一个搜索方向会被重复了好多次,这样使得迭代效率并不高。但是梯度下降法实现的思想确是很有意义的,为函数极值的确定与线性方程组的求解之间建立了联系。真正有意义的方法是在这个方法基础上实现的共轭梯度法,它对搜索方向进行了共轭处理,在理解了梯度下降的原理后,再对共轭梯度法进行学习可能更好一点。

    本站文章如非特别说明,均为原创。未经本人同意,请勿转载!转载请务必注明出处!

    展开全文
  • 本文开始介绍一种全新的优化方法——共轭梯度法。最初,共轭梯度法是用来求解线性方程 的,称为线性共轭梯度法。后来,有人把这种方法扩展到了非线性优化问题中,称为非线性共轭梯度法。本文从基本的线性共轭梯度讲...

    本文开始介绍一种全新的优化方法——共轭梯度法。

    最初,共轭梯度法是用来求解线性方程

    的,称为线性共轭梯度法。后来,有人把这种方法扩展到了非线性优化问题中,称为非线性共轭梯度法。本文从基本的线性共轭梯度讲起,介绍它的思想和算法流程。

    问题描述

    共轭梯度法并不适用于任意线性方程,它要求系数矩阵对称且正定。我们会在后面逐渐提到为什么会有这样苛刻的要求。

    比较巧的是,在非线性优化问题中,所有的二阶方法的每次迭代都需要求解一个形如

    这样的方程。比如,前面提到的牛顿系列方法,需要求解
    ,SLAM中常用的高斯牛顿法,需要求解
    ,这些系数矩阵
    通常都满足对称正定的要求(如果不正定,则想办法让它们正定,比如加上一个单位阵)。

    之前,我们都默认直接对系数矩阵求逆,得到

    ,而没有考虑求逆的性能代价。即便是通过矩阵分解的方式求逆,计算量仍然很大。线性共轭梯度法的提出,让我们能够通过迭代的方式求
    ,而不必对系数矩阵求逆。

    具体如何实现呢?我们先额外构造一个二次函数

    如果对该函数求导,并令导数为零,我们得到

    你会发现,使得该函数导数为0的

    恰好是线性方程
    的解。换句话说,
    的极小值点(
    是正定矩阵,意味着
    有极小值)恰好是线性方程
    的解。于是,我们可以把求解线性方程的问题转化为求二次函数极值的问题。

    寻找思路

    如何求二次函数极值?可以用之前讲过的最速下降法、牛顿法之类的方法吗?

    显然,牛顿法我们不用考虑了,因为刚才提到,求解线性方程其实是牛顿法的每次迭代都要面对的问题,所以它们不是一个层次上的问题。最速下降法倒是可以考虑,我们只需要求解当前位置的导数

    ,然后在该方向上确定一个最优的步长,不断迭代即可。回顾一下线搜索(二)那篇文章中给出的最速下降法的收敛路径,如下图所示。

    143d4d930faf58b6b978c27dce8b0264.png

    对于一个简单的二次函数,最速下降法竟然走出这么蜿蜒曲折的路径,这要用在求解线性方程上,恐怕比直接矩阵求逆还要慢。

    我们仍然以上图所示的二维目标函数为例,最好的结果是走出下面这条路径。

    e6fd0854332d42f8d5524c12c685a226.png

    我们这样解读这条新的路径:整个路径只有两段,即两次迭代,对应着目标函数只有两维。第一次迭代从

    ,沿着坐标轴
    的方向前进,达到了目标函数在该方向上可能达到的最优值。第二次迭代从
    ,沿着坐标轴
    的方向前进,达到了目标函数在该方向上可能达到的最优值。两次迭代后就达到了极小值点。

    重新审视这一路线,可以发现,每次迭代的前进方向是互相正交的,因为都是沿着坐标轴前进的。但最速下降法的搜索路径也是互相正交的,为什么就不行呢?如果我们把目标函数稍微旋转一下,就会发现问题的关键。

    da7a480c13ab1a75e891f4055114444a.png

    这是另一个二次函数,它与前面的二次函数很不一样,椭圆的轴不与坐标轴平行。(这里补充一点,之所以我们画出的二次函数的图形都是椭圆,是因为我们只考察系数矩阵为对称阵的情况。)此时,仍然沿着坐标轴搜索,这就是不折不扣的最速下降法了。这说明,能不能在两次迭代内达到极小值,与是否沿着坐标轴搜索没什么关系,关键在于椭圆的轴是否与搜索方向一致。在第一种情况中,椭圆的轴与坐标轴平行,搜索方向也和坐标轴平行,所以两次迭代即可。第二种情况中,椭圆的轴不与坐标轴平行,搜索方向依旧与坐标轴平行,所以两次迭代是不够的。这给了我们一些提示,如果我们先算出椭圆长短轴的方向,然后沿着这些方向去搜索,就能在两个迭代后达到极小值。

    从上面的二维函数推广到N维函数,想象N维空间中的椭圆,我们仍然可以计算出每一维的轴向,然后沿着这些轴搜索,就能确保在N次迭代后达到极小值。那如何求椭圆的轴呢?很简单,椭圆的轴其实就是系数矩阵的特征向量,对

    做特征值分解就可以了。说到这里,其实有点怪怪的,特征值分解的计算量恐怕不小吧,为了避免求逆却整出来个特征值分解,岂不是舍近求远。

    那么,有没有不求特征值,也能在N次迭代达到极小值点的方法呢?比如,能不能走出下面这样的路径。

    edb40559a3b2c8a66967a6e779fc771b.png

    第一次迭代前进的方向并不与椭圆的轴平行,在这种情况下,第二次迭代却能一次达到最优点。异想天开?当然不是,共轭梯度法就是要走出这样的路线。

    我们用数学工具探究一下这种路径需要满足的条件。对于二次函数中的任意一个点来说,该点与函数极值点的连线正是牛顿法中每次迭代的线搜索方向。回顾线搜索与信赖域中的式(2.14),套用这里的系数矩阵,可以得到图中的向量

    。如果假设起始点
    的搜索方向是最速下降方向,即
    ,于是有

    其中,最后两个梯度的乘积等于0是因为最速下降法的每次迭代方向都是垂直的。我们从二维例子中推出了一个看起来有些意义的结果,这个结果表明,连续的两次迭代的前进方法关于系数矩阵

    正交。也就是说,其中的一个方向经过
    的变换后与另一个方向正交。

    接下来的问题是如何把这一结论推广到N维目标函数中去,我们先不关注证明,直接看一下结论。

    共轭梯度法

    对于一组向量

    ,如果任意两个向量间满足

    则称这组向量与对称正定矩阵

    共轭。

    可见,刚才我们发现的正是共轭的性质,只不过,推广到高维空间后,任意两个向量之间都要关于

    正交,不仅仅是相邻的。

    共轭向量的好处正如我们前面发现的,依次沿着每个向量优化,最终可以在N个迭代后达到极小值。书中有严格的证明,但我们无意在此着墨,我认为结合上面的例子便可以理解其中的思想。

    现在的问题是,给定任意一个对称正定矩阵,如何求取一组共轭向量?

    显然,

    的特征向量就是一组共轭向量,因为
    对于
    的任意两个特征向量
    都成立,满足式(5.4),这就从公式上验证了我们前面的发现。不过,之前也说了,实际中是不能用特征向量的,因为计算太复杂。(请注意,这里我们用到了
    这一性质,只有对称矩阵有此性质,呼应了文章开头提到的对矩阵
    的特殊要求。这一性质贯穿了整个共轭梯度法的证明过程。)

    实际中,需要一种不依赖于矩阵分解,可以在每次迭代时以较小的代价计算出

    的方法。回顾上文的最后一张图,当时我们假设起始点
    的搜索方向是最速下降方向。这其实不是随便假设的,这样做会使得后面的计算变得简单。接下来,我们就从
    开始,推导出后续每次迭代所需的

    满足式(5.4)的共轭向量有无数个,而我们只需要取其中的一个。在

    的时刻,我们已知
    ,注意到这两个向量一定是线性无关的,因为
    是该方向上使得目标函数最小的步长,所以在该方向上目标函数不可能继续下降。那么就可以在这两个向量构成的平面上生成新的
    ,令
    ,将
    带入式(5.4),即可求得
    ,从而得到

    在此基础上,结合数学归纳法,我们可以证明,令

    其中

    即可使得生成的所有

    满足式(5.4)的要求,从而得到一组共轭向量。有了每次迭代的方向
    ,就可以计算出
    在该方向上的最优步长。计算方式很简单,我们用
    替换
    中的
    ,得到
    ,将其对
    求导,令导数等于0,即可得到最优步长

    于是,整个共轭梯度法的流程整理如下

    其中,用

    替代了
    。从
    开始,
    ,依次迭代下去。

    实际中,为了降低计算量,一般会对式(5.13)稍做改动,得到下面的算法流程

    这一版本的共轭梯度算法涉及的矩阵乘法和加法会更少一些。

    收敛速度

    优化算法的收敛速度是个值得关心的问题。前面提到了,共轭梯度法一定会在N次迭代内收敛到极值。但在N次迭代的过程中,收敛速度怎样呢?

    我们用线搜索(二):收敛性和收敛速度中提到过的向量的矩阵范数来衡量共轭梯度法的收敛速度。比如,评估

    随着
    的增加的变化率。这里省略证明过程,直接给出结论

    其中,

    是矩阵
    的第
    个特征值,
    是任意一个关于
    阶多项式。这个等式可以理解为,算法在第
    次迭代后,残差的最大值与一个
    阶多项式的值有关,这个
    阶多项式是所有
    阶多项式里面结果最小的那个,但每个多项式的输入必须是使其取得最大值的
    。单从这个不等式其实是看不出什么东西的,但却可以推导出一个惊人的结论。

    如果矩阵

    的N个特征值只能取
    个互不相同的值,这
    个值分别是
    。我们定义一个
    阶多项式

    可以发现该多项式满足

    。我们再定义另一个多项式

    此时有,

    是它的一个根。而且,
    是一个
    次多项式。如果我们把这个多项式带入到式(5.32)中,可以得到

    其中,第二行用特定的

    阶多项式
    代替了任意的
    阶多项式
    ,相当于放大了表达式的值。第四行,根据
    的定义,输入任意的
    都等于0。因此,我们最终发现,在第
    次迭代后,
    就已经收敛到了极小值

    这是个非常好的性质,它说明,共轭梯度法在N次迭代后收敛其实是最坏情况下的结果,当特征值存在重叠时,只需要更少的迭代次数就能收敛。

    不过,以上还没有回答我们想要了解的问题,在迭代过程中,算法的收敛速度如何。下面的式子在一定程度上给出了回答(这里仍然省略证明)

    在第

    次迭代后,残差不大于
    。这个结论也是很有价值的,这里的
    是从小到大排列的特征值,
    越小,残差就越小。对于同样的
    次迭代,如果矩阵
    的第
    大的特征值比矩阵
    的第
    大的特征值小,矩阵
    就会收敛得更快。另一方面,也可以理解为,整个收敛过程中的收敛速度是变化的,与特征值的大小分布有关。如果特征值的取值呈聚类状分布,比如有5个较大的特征值,还有几个接近于1的特征值,那么迭代过程很可能像下图这样

    da4a722b25f3b32bcf38863c97a4f75d.png

    在7次迭代后就达到了比较小的值。这是因为除了5个较大的特征值,其它特征值都集中在一个聚类中。根据前面得出的结论,如果

    只有
    个不同大小的特征值,
    次迭代后就能收敛。现在的情况与之类似,虽然聚类中的特征值并非完全一样,但我们仍然可以将其当做近似相等的特征值,所以优化在6次迭代后理论上就会接近于极值。当然,从图中可以看到,实际情况下,第7次迭代才接近极值,这也是合理的。

    如果我们进一步对式(5.34)放缩,可以得到一个更宽松的收敛性质

    在该式中,收敛速度受

    的条件数制约,条件数越小,收敛速度越快。在矩阵的条件数中我们介绍了条件数的定义和计算方式,即
    。这说明,特征值分布越集中,收敛速度越快,分布越分散,收敛速度越慢。极端情况下,当
    时,所有特征值都相等,
    的条件数为1,此时只需要一次迭代即可收敛到极值点。

    既然共轭梯度法的收敛速度如此依赖于特征值分布,我们理应做点什么,以改善那些特征值分布不佳的情况。

    预处理策略

    最容易想到,也是最有效的方法就是,对

    做个变换,让变换后的
    具有更为集中的特征值。比如,令

    针对

    的优化目标函数为

    其中,新优化方程中的系数矩阵为

    。如果我们找到合适的
    ,使得
    具有比较集中的特征值以及较小的条件数,那么对
    的优化就能很快收敛。

    但问题的关键在于如何选择合适的

    。极端情况下,如果令
    ,可以发现,
    恰好是
    的Cholesky分解,即
    。但我们不能这样做,因为对
    做Cholesky分解无异于直接求
    的逆,计算量太大。一种更快的方案是Incomplete Cholesky分解,该分解得到使得
    的稀疏矩阵
    ,该矩阵可保证
    ,从而得到良好的特征值分布。

    考虑预处理之后,式(5.23)转化为如下算法

    开始,
    ,依次迭代下去。不过,预处理虽然改善了收敛速度,却多了求解
    的步骤,会增加一部分计算量。

    总结

    本文介绍了线性共轭梯度法,这是一种以迭代的方式求解线性方程的方法。这种方法的妙处在于,不需要矩阵分解和求逆,只需要上次迭代的步长和当前的梯度,即可计算下次的步长。计算量极低,且不占用额外的内存空间,适用于求解大规模的线性方程。

    参考资料

    Chapter 16 Preconditioning Illinois math course

    上一篇:信赖域(二):确切方法
    下一篇:共轭梯度法(二):非线性共轭梯度
    展开全文
  • 共轭梯度法学习笔记

    千次阅读 2019-06-27 20:07:28
    共轭梯度法的基本思想是把共轭性与最速下降法相结合,利用已知点处的梯度构造一组共轭方向,并沿这组方向进行搜索,求出目标函数的极小点。根据共轭方向的基本性质,这种方法具有二次终止性。 用共轭梯度法求解简单...

    共轭梯度法是介于最速下降法与牛顿法之间的一个方法,它仅需一阶导数信息,但克服了最速下降法收敛慢的缺点,又避免了牛顿法需要存储和计算Hesse矩阵并求逆的缺点。共轭梯度法的基本思想是把共轭性与最速下降法相结合,利用已知点处的梯度构造一组共轭方向,并沿这组方向进行搜索,求出目标函数的极小点。根据共轭方向的基本性质,这种方法具有二次终止性。

    用共轭梯度法求解简单多元函数问题:

    引入泰勒展开式概念:

    f(x)=f\left(x_{0}\right)+f^{\prime}\left(x_{0}\right)\left(x-x_{0}\right)+\frac{f^{(2)}\left(x_{0}\right)}{2 !}\left(x-x_{0}\right)^{2}+\cdots+\frac{f^{(n)}\left(x_{0}\right)}{n !}\left(x-x_{0}\right)^{n}

    则对x_{1}x_{2}求一阶偏导得到在x^{(0)}点的梯度方向g_{0}=\nabla f\left(x^{(0)}\right)。那么$g_{0}$$p_{0}$是什么关系呢?注意这道题是求最小值问题,梯度方向是函数局部上升最快方向,则反方向为局部下降最快方向。

    题目中\beta有以下几种常见形式

    \beta_{j}=\frac{\left\|\nabla f\left(y^{(j+1)}\right)\right\|^{2}}{\left\|\nabla f\left(y^{(j)}\right)\right\|^{2}}

    \beta_{j}=\frac{\boldsymbol{g}_{j+1}^{\mathrm{T}}\left(\boldsymbol{g}_{i+1}-\boldsymbol{g}_{j}\right)}{\boldsymbol{g}_{j}^{\mathrm{T}} \boldsymbol{g}_{j}}

    \beta_{j}=\frac{\boldsymbol{g}_{i+1}^{\mathrm{T}}\left(\boldsymbol{g}_{i+1}-\boldsymbol{g}_{i}\right)}{\boldsymbol{d}^{(j) \mathrm{T}}\left(\boldsymbol{g}_{j+1}-\boldsymbol{g}_{j}\right)}

    \beta_{j}=\frac{\boldsymbol{d}^{(j) T} \nabla^{2} f\left(\boldsymbol{x}^{(j+1)}\right) \boldsymbol{g}_{j+1}}{\boldsymbol{d}^{(j) T} \nabla^{2} f\left(\boldsymbol{x}^{(j+1)}\right) \boldsymbol{d}^{(j)}}

    参考网址:

    1、数学优化入门:梯度下降法、牛顿法、共轭梯度法

    展开全文
  • 码字不易,转载请注明出处http://xlindo.comhttps://zhuanlan.zhihu.com/p/52692238优化中常用两种梯度方法,...据我的理解,应该是用普通梯度法来阐述梯度下降思想,用共轭梯度法作引申。1 问题提出2 梯度法2.1 ...

    码字不易,转载请注明出处

    http://xlindo.com

    https://zhuanlan.zhihu.com/p/52692238

    优化中常用两种梯度方法,普通梯度法(我们老师叫的Einfaches Gradienten-Verfahren)以及共轭梯度法(Conjugate Gradients)。

    据我的理解,应该是用普通梯度法来阐述梯度下降思想,用共轭梯度法作引申。

    1 问题提出

    2 梯度法

    2.1 梯度下降法 simple gradient method

    梯度下降法的思想很简单,每次下降方向选择梯度反方向,这样能保证每次值都在减小。

    但问题是“锯齿”现象,所以效果很差,实际情况中一般不用。

    x_k = 6 # The algorithm starts at x=6

    alpha = 0.01 # step size multiplier

    precision = 0.00001 # end mark

    previous_step_size = 1

    max_iters = 10000 # maximum number of iterations

    iters = 0 #iteration counter

    df = lambda x: 4 * x**3 - 9 * x**2

    while previous_step_size > precision and iters < max_iters:

    prev_x_k = x_k

    g_k = df(prev_x_k)

    d_k = -g_k

    x_k += alpha * d_k

    previous_step_size = abs(x_k - prev_x_k)

    iters+=1

    print("The local minimum occurs at", x_k)

    print("The number of iterations is", iters)

    The local minimum occurs at 2.2499646074278457

    The number of iterations is 70

    2.2 共轭梯度法 conjugate gradient method

    该方法基于梯度下降法。

    由于梯度下降法每次下降的方向为负梯度方向,这并不能保证其是最优的方向。通过共轭方向的计算,保证第二步开始的下降方向在一个圆锥内,这能极大的提高下降的效率。

    x_k = 6 # The algorithm starts at x=6

    alpha = 0.01 # step size multiplier

    precision = 0.00001 # end mark

    previous_step_size = 1

    max_iters = 10000 # maximum number of iterations

    iters = 0 #iteration counter

    beta = 0

    df = lambda x: 4 * x**3 - 9 * x**2

    g_k = df(x_k)

    while previous_step_size > precision and iters < max_iters:

    prev_x_k = x_k

    prev_g_k = g_k

    if 0 == iters:

    d_k = -g_k

    else:

    g_k = df(x_k)

    beta = g_k*g_k/(prev_g_k*prev_g_k)

    d_k = -g_k + beta*d_k

    x_k += alpha * d_k

    previous_step_size = abs(x_k - prev_x_k)

    iters+=1

    print("The local minimum occurs at", x_k)

    print("The number of iterations is", iters)

    The local minimum occurs at 2.2500110335395793

    The number of iterations is 22

    2.3 线搜索 line search

    线搜加入了对步长的计算。

    名为Backtracking Line Search(BLS)的梯度下降法以Armijo条件为方法,在以2.2为基础的搜索方向上选取不越过最优点的最大步长: 1. 定义

    ,及

    2. 从

    开始迭代

    ,计算

    其中Armijo条件为

    # -*- coding: utf-8 -*-

    # min ( )=(x-3)**2

    import matplotlib.pyplot as plt

    def f(x):

    '''The function we want to minimize'''

    return (x-3)**2

    def f_grad(x):

    '''gradient of function f'''

    return (x-3)*2

    x = 6

    y = f(x)

    MAX_ITER = 300

    curve = [y]

    i = 0

    step = 0.1

    previous_step_size = 1.0

    precision = 1e-4

    #下面展示的是我之前用的方法,看上去貌似还挺合理的,但是很慢

    while previous_step_size > precision and i < MAX_ITER:

    prev_x = x

    gradient = f_grad(x)

    x = x - gradient * step

    new_y = f(x)

    if new_y > y: #如果出现divergence的迹象,就减小step size

    step *= 0.8

    curve.append(new_y)

    previous_step_size = abs(x - prev_x)

    i += 1

    print("The local minimum occurs at", x)

    print("The number of iterations is", i)

    plt.figure()

    plt.show()

    plt.plot(curve, 'r*-')

    plt.xlabel('iterations')

    plt.ylabel('objective function value')

    #下面展示的是backtracking line search,速度很快

    x = 6

    y = f(x)

    alpha = 0.25

    beta = 0.8

    curve2 = [y]

    i = 0

    previous_step_size = 1.0

    precision = 1e-4

    while previous_step_size > precision and i < MAX_ITER:

    prev_x = x

    gradient = f_grad(x)

    step = 1.0

    while f(x - step * gradient) > f(x) - alpha * step * gradient**2:

    step *= beta

    x = x - step * gradient

    new_y = f(x)

    curve2.append(new_y)

    previous_step_size = abs(x - prev_x)

    i += 1

    print("The local minimum occurs at", x)

    print("The number of iterations is", i)

    plt.plot(curve2, 'bo-')

    plt.legend(['gradient descent', 'BLS'])

    plt.show()

    运行结果如图

    The local minimum occurs at 3.0003987683987354

    The number of iterations is 40

    The local minimum occurs at 3.000008885903001

    The number of iterations is 10

    3 其他方法

    其他方法还有如Riccati与Collocation,以后有时间再单独总结。

    展开全文
  • 最优化方法-共轭梯度法

    千次阅读 2019-03-31 17:16:23
    最优化方法-共轭梯度法 1.简介 共轭梯度法最初由Hesteness和Stiefel于1952年为求解线性方程组而提出的。其基本思想是把共轭性与最速下降方法相结合,利用已知点处的梯度构造一组共轭方向,并沿这组方向进行搜素,求...
  • 共轭梯度法关键是要找正交向量寻找方向,去不断逼近解。 其本质是最小二乘解的思想 最小二乘解 其中A系数矩阵是确定的,Ax是永远都取不到向量 b的,取得到那就是不用最小二乘解 我要求AX和b最小的距离,就是...
  • 提出了一种分析频谱的新方法,其主要思想是采用共轭梯度法训练傅里叶基神经网络权值,根据权值获得信号的幅度谱和相位谱,并给出了基于Matlab语言的频谱分析应用实例。仿真结果表明,与FFT相比,该方法具有计算精度...
  • 牛顿主要解决非线性优化问题,其收敛速度比梯度下降速度更快,其需要解决的问题可以描述为:对于目标函数f(x),在无约束条件下求他的最小值minf(x)。 牛顿的主要思想:在现有的极小值估计值的附近对f(x...
  • 最优化方法——最速下降法,阻尼牛顿法,共轭梯度法 目录 最优化方法——最速下降法,阻尼牛顿法,共轭梯度法 1、不精确一维搜素 1.1 Wolfe-Powell 准则 2、不精确一维搜索算法计算步骤 3、最速下降法 3.1 ...
  • 实验的题目和要求一、所属课程名称:最优化方法二、实验日期:2010年5月10日~2010年5月15日三、实验目的掌握最速下降法,牛顿法和共轭梯度法的算法思想,并能上机编程实现相应的算法。二、实验要求用MATLAB实现最速...
  • 牛顿(Newton's method)是一种在实数域和复数域上近似求解方程的方法,,它使用函数f(x)的泰勒级数的前面几项来寻找方程f(y)=0的根。 牛顿最初由艾萨克·牛顿在《Method of Fluxions》(1671年完成,... 基本思想3....
  • 1. 梯度下降(Gradient Descent) ...梯度下降的优化思想是用当前位置负梯度方向作为搜索方向,因为该方向为当前位置的最快下降方向,所以也被称为是”最速下降“。最速下降越接近目标值,步长越小,前进...
  • 梯度:有时候也称之为斜度,也就是一个曲面沿着给定方向的...梯度下降的优化思想是用当前位置负梯度方向作为搜索方向,因为该方向为当前位置的最快下降方向,所以也被称为是”最速下降“。最速下降越接近目标值.
  • 梯度下降的优化思想是用当前位置负梯度方向作为搜索方向,因为该方向为当前位置的最快下降方向,所以也被称为是”最速下降“。最速下降越接近目标值,步长越小,前进越慢。梯度下降的搜索迭代示意图如下图所...
  • (一)共轭梯度算法

    万次阅读 2019-06-19 16:37:26
    共轭梯度算法是干什么? 共轭梯度算法是一种迭代算法,在一次次的对待中最终...共轭梯度法是属于最小化类的迭代方法。为了求解Ax = b这样的一次函数,可以先构造一个二次齐函数 f(x)=12xTAx−bTxf(x) = \frac{1}{...
  • 学习共轭梯度算法,最好是在了解了梯度下降和(拟)牛顿之后,再学习,省事而且三种方法对比记忆,效果会更好,废话不多说,进入正题。 链接抛上: 梯度下降算法原理及其在线性回归中的应用 最优化算法之牛顿和...
  • 本文首发于公众微信号-AI研究订阅号,来源...接下来介绍的共轭方向是介于最速下降和Newton之间的一种方法,它克服了最速下降的锯齿现象,从而提高了收敛速度;它的迭代公式也比较简单,不必计算目标函数的...
  • 常见的: 1.梯度下降:全批度下降,随机梯度下降(SGD),小批量梯度下降(batch SGD) 2.牛顿法:优化函数的二阶导数信息,海森矩阵求解困难,还有海森...4.共轭梯度法共轭梯度法是介于梯度下降法(最速下降法...
  • 1 无约束优化算法思路 2 无约束优化的几种方法 2.1 最速下降法(梯度法) ...注意:共轭方向法仅仅是一种思想,其思想描述如下,由此衍生的共轭梯度法 看了好多人的讲解,可脑袋里还是没...
  • [基础知识点]PCG算法以及代码解析

    千次阅读 2020-05-08 23:19:07
    共轭梯度法思想是,选择一个优化方向后,本次选择的步长能够将这个方向的误差更新完,在以后的优化更新过程中不再需要朝这个方向更新了。由于每次将一个方向优化到了极小,后面的优化过程将不再影响之前优化方向上...
  • 当分析一个问题时,直接求解该问题可能会非常的复杂 一种很常见的思想是通过已有的信息,逐步将问题逼近目标答案 渐近的思想在很多算法中都可以见到 ... 比如梯度下降算法、牛顿算法以及共轭梯度法 5、级...
  • 第四章 无约束最优化方法 优化算法中心思想4.1 梯度法(Steepest Descent Method)4.2 牛顿法4.2.1 基本牛顿法4.2.2 阻尼牛顿法4.3 变尺度法(拟牛顿法)(Quasi-Newton's Method)4.4 共轭梯度法(Conjugate Gradient ...
  • 预处理的基本思想是通过修改病态...预处理共轭梯度共轭梯度法的收敛速度正比于 ,因此对于条件数比较差的矩阵收敛速度会很慢,我们需要选择矩阵 使得 ,这样就具有更快的收敛速度了。为了保证预处理后的矩阵 是对...
  • 常见的最优化方法有梯度下降法、牛顿法和拟牛顿法、共轭梯度法等等。1. 梯度下降法(Gradient Descent)梯度下降法实现简单,当目标函数是凸函数时,梯度下降法的解是全局解。一般情况下,其解不保证是全局最优解,...
  • 几种常见的优化算法

    千次阅读 2019-04-24 23:28:55
    目录 神经网络优化最重要的思想: 1. 梯度下降法(Gradient Descent)(一阶) ...3. 共轭梯度法(Conjugate Gradient) 4. 启发式优化方法 5.待补充 神经网络优化最重要的思想: 梯度的反向传...
  • 由于求解本问题等价于一个严格凸泛函的极小化问题,所以本文中所提到的多重网格法均采用以下三种非线性无约束最优化方法:Polack-Ribiere共轭梯度法,Hooke-Jeeves模式搜索法及不含线搜索的SSC梯度法作为非线性磨光...
  • 常用优化算法

    2018-12-18 10:03:47
    常见的最优化方法有梯度下降法、牛顿法和拟牛顿法、共轭梯度法等等。 1. 梯度下降法(Gradient Descent) 梯度下降法是最常用的一种优化算法。其核心思想是:在当前位置寻找梯度下降最快的方向,来逐渐逼...
  • 最优化方法——BGFS变尺度算法

    千次阅读 2019-04-18 22:07:29
    目录 ...最优化方法——最速下降法,阻尼牛顿法,共轭梯度法 1、不精确一维搜素 1.1 Wolfe-Powell 准则 2、不精确一维搜索算法计算步骤 1、BGFS基本思想 2、BGFS变尺度算法的计算公式 ...
  • 《最优化方法》复习题一、 简述题1、怎样判断一个函数是否为凸函数.(例如: 判断函数是否为凸函数)2、写出几种迭代的收敛条件.3、熟练...简述共轭梯度法的基本思想.写出Goldstein、Wolfe非精确一维线性搜索的公式。5...
  • 主要内容包括 (精确或非精确) 线搜索技术, 最速下降法与 (修正) 牛顿法, 共轭梯度法, 拟牛顿法, 信赖域方法, 非线性最小二乘问题的解 法, 约束优化问题的最优性条件, 罚函数法, 可行方向法, 二次规划问题的解法, ...

空空如也

空空如也

1 2 3
收藏数 42
精华内容 16
关键字:

共轭梯度法思想