1 最优化概论

(1) 最优化的目标

最优化问题指的是找出实数函数的极大值或极小值,该函数称为目标函数。由于定位\(f(x)\)的极大值与找出\(-f(x)\)的极小值等价,在推导计算方式时仅考虑最小化问题就足够了。极少的优化问题,比如最小二乘法,可以给出封闭的解析解(由正规方程得到)。然而,大多数优化问题,只能给出数值解,需要通过数值迭代算法一步一步地得到。

(2) 有约束和无约束优化

一些优化问题在要求目标函数最小化的同时还要求满足一些等式或者不等式的约束。比如SVM模型的求解就是有约束优化问题,需要用到非线性规划中的拉格朗日乘子和KKT条件。这里我们仅介绍无约束优化,有约束优化放在后面的章节讲解。

(3) 线性和非线性规划

线性函数是指目标函数和约束都为线性的优化问题,非线性规划是指目标函数和约束有一个为非线性的优化问题。线性规划一般在运筹学(经济模型、图论网络流等)中有重要运用,而非线性规划在机器学习中有着重要的运用。我们把主要目光放在非线性规划。

(4) 凸优化和非凸优化

按照斯坦福 Boyd 教授(编写凸优化圣经《convex optimization》的那位)的观点,优化问题的分水岭不是线性和非线性,而是凸和非凸。这句话侧面说明了凸优化做为一种特殊的优化问题,显得非常重要,尤其是在机器学习领域。那么为什么凸优化会如此重要呢?首先我们拉看什么是凸函数。凸函数的定义如下:
\(Ω\)为凸集,如果对任意的\(x_1,x_2 \in Ω\)以及每一个\(α(0\leqslant \alpha \leqslant 1)\),有 \(f(αx_1+(1-\alpha)x_2)<=αf(x_1) + (1-α)f(x_2)\),则称定义在凸集上的函数\(f\)是凸的(convex)。
\(Ω\)为凸集,如果对每一个\(α(0<α<1)\)以及\(x_1,x_2 \in Ω\)\(x_1\neq x_2\),有\(f(αx_1+(1-α)x_2)<αf(x_1) + (1-α)f(x_2)\)则称\(f\)是严格凸的(strictly convex)。
以下展示了几个凸函数的图像例子,从几何角度看没如果图形中两点的连线处处都不在图形的下方,则函数是凸的。或者做为二维空间中的函数,如果函数的图形是碗状的,这个函数就是凸的。

数值优化:一阶和二阶优化算法(Pytorch实现)-LMLPHP

那么凸函数有什么神奇的性质值得我们为之兴奋呢?我们有定理:\(f\)是至少含有一个内点的凸集\(Ω\)上的凸函数,当前仅当\(f\)的Hessian矩阵\(\bm{H}\)在整个\(Ω\)上是半正定的。
此处Hessian矩阵正是函数的曲率概念在\(\R^n\)上的推广,凸函数在每个方向上都有正(至少是非负)的曲率。如果一个函数的Hessian矩阵在一个小区域内是半正定的,则称该函数是局部凸的;如果Hessian矩阵在这个区域内是正定的(但不妨碍我说\(\bm{H}\)在整个\(Ω\)上是半正定的,细品),则称这个函数是严格局部凸的(locally strictly convex)。
以下插入一下函数极值点的必要和充分条件的介绍:

函数极值点的必要和充分条件


而我们对于任意一个无约束优化问题,函数的最值是要满足一阶必要条件和二阶必要条件的。
一阶必要条件:\(Ω\)\(\R^n\)的一个子集并且\(f\)\(Ω\)上的函数。如果\(\bm{x}^*\)\(f\)\(Ω\)上的相对极小点,那么对\(\bm{x}^*\)点处的任意一个可行的方向\(\bm{d}∈\R^n\),有\(∇f(\bm{x}^*)\bm{d}>=0\)。一个非常重要的特殊情形发生在\(\bm{x}^*\)\(Ω\)内部时(\(\bm{x}^*\)\(Ω\)的内点,\(Ω=\R^n\)就对应这种情形)。在这种情况下,从\(\bm{x^*}\)发散出去的每个方向都是可行方向,因此对所有的\(\bm{d}∈\R^n\),都有\(\nabla f(\bm{x}^{*})\bm{d}>=0\),这就意味着\(\nabla f(\bm{x}^*)=0\)
二阶必要条件:\(Ω\)\(\R^n\)的一个子集并且\(f\)\(Ω\)上的函数。如果\(\bm{x}^*\)\(f\)\(Ω\)上的相对极小点,那么对\(\bm{x}^*\)处的任意一个可行方向\(\bm{d}∈Ω\),有:
\(∇f(\bm{x}^*)\bm{d}>=0\)
② 如果\(∇f(\bm{x}^*)\bm{d}=0\),那么\(d^T∇^2f(\bm{x}^*)\bm{d}>=0\)

同样的,我们在无约束情形下,设\(\bm{x}^*\)是集合\(Ω\)的内点。并且设\(\bm{x}^*\)是函数\(f\)\(Ω\)上的一个内点,那么:
\(\nabla f(\bm{x}^*)=0\)
② 对所有\(\bm{d}\)\(\bm{d}^T∇^2f(\bm{x}^*)\bm{d}\geqslant0\)(这个条件等价于说明Hessian矩阵\(\bm{H}\)是半正定的)

二阶充分条件:稍微加强一下二阶必要条件的条件②,我们就能得到\(\bm{x^*}\)是相对极小点的条件:
\(∇f(\bm{x}^*)=0\)
\(\bm{H}(\bm{x}^*)\)正定。
那么\(\bm{x}^*\)\(f\)的一个严格相对极小点(因为严格正定,不存在Hessain矩阵\(\bm{H}\)特征值为0的困扰)


而上面说了,如果Hessian矩阵在这个区域内是正定的,则称这个函数是严格局部凸的(locally strictly convex),故我们看出上面说的的二阶充分条件要求在每个点\(\bm{x}^*\)处的函数是严格局部凸的。
推广之,设\(f\)是定义在凸集\(Ω\)上的凸函数,那么使函数达到极小值的点集\(Γ\)是凸集,并且\(f\)的相对极小点也是全局极小点。

这下大家应该知道凸函数的好处了,凸函数没有“坑坑洼洼”,相对极小就是全局极小,这样找到相对极小点就可以收功,便于设计出高效的优化算法,如我们求解SVM中的SMO算法(SVM是个凸优化问题)。而有很多“坑坑洼洼”的函数想要找到全局极小点是NP-hard问题,只能采用遗传算法、退火算法这类启发式算法进行求解。我们在深度学习中的大多数函数(可以把带激活层的神经网络当成一个嵌套的函数)是非凸的,不过我们找到这类函数的全局最小值意义不大,一般我们找到局部极小拟合程度就足够好,从而可以解决我们的问题了。因此,在神经网络中我们一般不采用启发式算法来优化,多是采用随机梯度下降、拟牛顿法、动量法等“更正统”的优化算法来找到局部最优解以近似全局最优解。

2 不使用导数的无约束优化——Fibonacci 搜索(也称黄金分割搜索)

(1) 线搜索算法

在一条已知的直线(只有一个变量)上确定极小点的过程,被称为线搜索(line search)。对于一般不能解析地求极小值的非线性函数,这一过程实际上是采用一些巧妙的沿直线搜索的方法来实现的。这些线搜索技巧实际上就是求解一维极小化问题的方法,因为高维问题最终是通过进行一系列逐次线搜索来求解的,所以这些线搜索方法是非线性规划算法的基石。

(2) 黄金分割搜索

求解线搜索问题的一个最普遍的方法是本节所描述的斐波那契搜索方法。一旦解的范围已知,黄金分割搜索是一种有效找出单变量函数\(f(x)\)的最小值的方式。我们假设\(f\)是一个单峰函数,在区间\([a,b]\)上具有相对极小。选择区间内的两点\(x_1\)\(x_2\),使得\(a<x_1<x_2<b\)。我们将使用新的更小的区间替换原始区间。根据以下法则该区间可以继续括住极小值。如果\(f(x_1)\leqslant f(x_2)\),则在下一步中保持区间\([a, x_2]\)。如果\(f(x_1)>f(x_2)\),则保持\([x_1, b]\)。如下图所示。

数值优化:一阶和二阶优化算法(Pytorch实现)-LMLPHP

不过,我们如何将\(x_1\)\(x_2\)放置在区间\([a,b]\)上呢?我们在选择\(x_1\)\(x2\)时有两个标准:
(a) 关于区间保持对称(由于我们不知道极小在区间的哪一侧)
(b) 选择\(x_1\)\(x_2\)使得不管在下一步中使用哪种选择,\(x_1\)\(x_2\)都是下一步中的某个采样点。为了简化讨论,我们以\([a,b]=[0, 1]\)为例子,可以推广到其他区间。即要求\(x_1 = 1 - x_1\)(关于区间中心对称),\(x_1 = x_2^2\)。如下图所示,如果新区间为\([0, x_2]\),标准(b)保证原始的\(x_1\)将会在下个区间中变为\(x_2\),因而仅仅需要进行依次函数求值,即\(f(x_1g)\),(这里\(g\)\(x_2\)的初始值)同样,如果新的区间为\([x_1, 1]\),则\(x_2\)变为新的"\(x_1\)"。这种重用函数求值的能力意味着在第一步后,每步仅需要目标函数的单次求值。每轮迭代演示如下:
数值优化:一阶和二阶优化算法(Pytorch实现)-LMLPHP
根据上图所示,我们需要选择黄金分割搜索的比例,即\(x_2\)所放置的位置。 旧区间和新区间的比例为\(1/g = (1+ \sqrt{5})/2\),即黄金分割。这样,每轮放置的\(x_1 = 1-g=(1+ \sqrt{5})/2,x_2 = g = ( \sqrt{5} − 1)/2=0.618\),如下图所示:
数值优化:一阶和二阶优化算法(Pytorch实现)-LMLPHP
故黄金分割算法如下:

import numpy as np
import math
def gss(f, a, b, k):
    g = (math.sqrt(5)-1)/2
    # 计算x1和x2
    x1 = a + (1-g)*(b-a)
    x2 = a + g*(b-a)
    f1, f2 = f(x1), f(x2)
    for i in range(k):
        if f1 < f2 :
            # 依次更新b, x2, x1
            b = x2
            x2 = x1
            # 这里代码设计的很巧妙,b是已经更新后的新b
            x1 = a + (1-g)*(b-a)
            f2 = f1
            f1 = f(x1)
        else:
            a = x1
            x1 = x2
            x2 = a + g*(b-a)
            f1 = f2
            f2 = f(x2)
    y = (a+b)/2
    return(a, b), y
if __name__ == '__main__':
    a, b = 0, 1
    k = 15
    (a,b), y = gss(lambda x: x**6-11*x**3+17*x**2-7*x+1, a, b, k)
    print("(%.4f, %.4f)"%(a, b), y)

算法的运行结果如下:

(0.2834, 0.2841) 0.28375198388070366

可以看到函数\(f(x)=x^6-11x^3+17x^2-7x+1\)在区间\([0,1]\)上的最小值在\(0.2834\)\(0.2841\)之间,可以近似为\(0.28375\)

3 使用一阶导数的无约束优化——梯度下降法

\(f\)是多元函数,\(\bm{x}^{(t)}\)\(\bm{x}^{(t+1)}\)都是向量。梯度下降法的迭代式为:

\[\bm{x}^{(t+1)} = \bm{x}^{(t)}-η∇f(\bm{x}^{(t)})\]
10-18 04:49