八、梯度下降法和拟牛顿法

1、梯度
2、梯度上升和梯度下降
3、梯度下降算法详解
- 3.1 直观解释
- 3.2 梯度下降相关概念
- 3.3 梯度下降的矩阵描述
- 3.4 梯度下降的算法调优
4、梯度下降法大家族
5、梯度下降法和其他无约束优化算法的比较
- July解释
6、拟牛顿法
- 6.1 牛顿法
- 6.2 拟牛顿法

1、梯度

梯度下降经典总结

微积分里边，对多元函数求偏导数，把求得的各个参数的偏导数以向量的形式写出来就是梯度。

对函数f(x,y)f(x,y)f(x,y)分别求x，y的偏导数，求得的梯度向量就是(∂f∂x,∂f∂y)T(∂f∂x,∂f∂y)T (\frac{\partial f}{\partial x},\frac{\partial f}{\partial y }) ^T

梯度向量的意义：就是函数变化增加最快的地方，也就是函数在某点上，沿着其梯度向量的方向就是函数增加最快的地方。或者说沿着梯度向量的方向，更容易找到函数的最大值，换句话说沿着梯度向量相反的方向，梯度减小的最快，更容易找到函数的最小值。

2、梯度上升和梯度下降

机器学习算法中，在最小化损失函数的时候，可以通过梯度下降法来一步步迭代求解，得到最小化的损失函数和模型参数。相反，如果需要求解损失函数的最大值，这时就需要用梯度上升法来迭代了。

梯度上升和梯度下降可以通过添加负号来转换。

3、梯度下降算法详解

3.1 直观解释

首先来看看梯度下降的一个直观的解释。比如我们在一座大山上的某处位置，由于我们不知道怎么下山，于是决定走一步算一步，也就是在每走到一个位置的时候，求解当前位置的梯度，沿着梯度的负方向，也就是当前最陡峭的位置向下走一步，然后继续求解当前位置梯度，向这一步所在位置沿着最陡峭最易下山的位置走一步。这样一步步的走下去，一直走到觉得我们已经到了山脚。当然这样走下去，有可能我们不能走到山脚，而是到了某一个局部的山峰低处。

从上面的解释可以看出，梯度下降不一定能够找到全局的最优解，有可能是一个局部最优解。当然，如果损失函数是凸函数，梯度下降法得到的解就一定是全局最优解。

3.2 梯度下降相关概念

在详细了解梯度下降的算法之前，我们先看看相关的一些概念。

步长（Learning rate）：步长决定了在梯度下降迭代的过程中，每一步沿梯度负方向前进的长度。用上面下山的例子，步长就是在当前这一步所在位置沿着最陡峭最易下山的位置走的那一步的长度。
特征（feature）：指的是样本中输入部分，比如2个单特征的样本（x(0),y(0)）,（x(1),y(1)）（x(0),y(0)）,（x(1),y(1)）（x^{(0)},y^{(0)}）,（x^{(1)},y^{(1)}），则第一个样本特征为x(0)x(0)x^{(0)}，第一个样本输出为y(0)y(0)y^{(0)}。
假设函数（hypothesis function）：在监督学习中，为了拟合输入样本，而使用的假设函数，记为hθ(x)hθ(x)h_θ(x)。比如对于单个特征的m个样本（x(i),y(i)）(i=1,2,...m),（x(i),y(i)）(i=1,2,...m),（x^{(i)},y^{(i)}）(i=1,2,...m),可以采用拟合函数如下： hθ(x)=θ0+θ1xhθ(x)=θ0+θ1xh_θ(x)=θ_0+θ_1x。
损失函数（loss function）：为了评估模型拟合的好坏，通常用损失函数来度量拟合的程度。损失函数极小化，意味着拟合程度最好，对应的模型参数即为最优参数。在线性回归中，损失函数通常为样本输出和假设函数的差取平方。比如对于m个样本（xi,yi）(i=1,2,...m)（xi,yi）(i=1,2,...m)（xi,yi）(i=1,2,...m),采用线性回归，损失函数为：

J(θ0,θ1)=∑i=1m(hθ(xi)−yi)2J(θ0,θ1)=∑i=1m(hθ(xi)−yi)2J(θ0,θ1)=∑i=1m(hθ(xi)−yi)2J(θ0,θ1)=∑i=1m(hθ(xi)−yi)2

J(θ_0,θ_1)=∑i=1m(h_θ(x_i)−y_i)2J(θ_0,θ_1)=∑i=1m(h_θ(x_i)−y_i)^2

其中xixix_i表示第i个样本特征，yiyiy_i表示第iii个样本对应的输出，hθ(xi)" role="presentation" style="position: relative;">hθ(xi)hθ(xi)h_θ(xi)为假设函数。

3.3 梯度下降的矩阵描述

3.4 梯度下降的算法调优

在使用梯度下降时，需要进行调优。哪些地方需要调优呢？

算法的步长选择。在前面的算法描述中，我提到取步长为1，但是实际上取值取决于数据样本，可以多取一些值，从大到小，分别运行算法，看看迭代效果，如果损失函数在变小，说明取值有效，否则要增大步长。前面说了。步长太大，会导致迭代过快，甚至有可能错过最优解。步长太小，迭代速度太慢，很长时间算法都不能结束。所以算法的步长需要多次运行后才能得到一个较为优的值。
算法参数的初始值选择。初始值不同，获得的最小值也有可能不同，因此梯度下降求得的只是局部最小值；当然如果损失函数是凸函数则一定是最优解。由于有局部最优解的风险，需要多次用不同初始值运行算法，关键损失函数的最小值，选择损失函数最小化的初值。
归一化。由于样本不同特征的取值范围不一样，可能导致迭代很慢，为了减少特征取值的影响，可以对特征数据归一化，也就是对于每个特征x，求出它的期望x¯¯¯x¯和标准差std(x)，然后转化为：

x−xstd(x)x−xstd(x)

\frac{x−x}{std(x)}

这样特征的新期望为0，新方差为1，迭代次数可以大大加快。

4、梯度下降法大家族

普通的梯度下降算法在更新回归系数时要遍历整个数据集，是一种批处理方法，这样训练数据特别忙庞大时，可能出现如下问题：

1）收敛过程可能非常慢；

2）如果误差曲面上有多个局极小值，那么不能保证这个过程会找到全局最小值。

为了解决上面的问题，实际中我们应用的是梯度下降的一种变体被称为随机梯度下降。

5、梯度下降法和其他无约束优化算法的比较

在机器学习中的无约束优化算法中，除了梯度下降以外，还有最小二乘法、牛顿法和拟牛顿法。

梯度下降法：

需要选择步长
迭代求解
当样本量很大的时候，求解速度快

最小二乘法：

不需要选择步长
计算解析解
如果样本量不算很大，且存在解析解，最小二乘法比起梯度下降法要有优势，计算速度很快。但是如果样本量很大，用最小二乘法由于需要求一个超级大的逆矩阵，这时就很难或者很慢才能求解解析解了，

拟牛顿法 / 牛顿法：

二阶收敛，收敛速度快
同样也是迭代求解，但是梯度下降法是梯度求解，而牛顿法是用二阶的Hessian矩阵的逆矩阵或伪矩阵求解
牛顿法收敛更快，但是每次迭代的时间比梯度下降法长

最小二乘：使误差「所谓误差，当然是观察值与实际真实值的差量」平方和达到最小以寻求估计值的方法，就叫做最小二乘法，用最小二乘法得到的估计，叫做最小二乘估计。

勒让德在论文中对最小二乘法的优良性做了几点说明：

最小二乘使得误差平方和最小，并在各个方程的误差之间建立了一种平衡，从而防止某一个极端误差取得支配地位
计算中只要求偏导后求解线性方程组，计算过程明确便捷
最小二乘可以导出算术平均值作为估计值

July解释

下面是一个典型的机器学习的过程，首先给出一个输入数据，我们的算法会通过一系列的过程得到一个估计的函数，这个函数有能力对没有见过的新数据给出一个新的估计，也被称为构建一个模型。

我们用X1，X2..XnX1，X2..XnX_1，X_2..X_n 去描述feature里面的分量，比如X1=房间的面积，X2=房间的朝向等等，我们可以做出一个估计函数：

θ在这儿称为参数，在这儿的意思是调整feature中每个分量的影响力，就是到底是房屋的面积更重要还是房屋的地段更重要。为了如果我们令X0 = 1，就可以用向量的方式来表示了：

我们程序也需要一个机制去评估我们θ是否比较好，所以说需要对我们做出的h函数进行评估，一般这个进行评估的函数称为损失函数（loss function），描述h函数不好的程度，在下面，我们称这个函数为J函数

在这儿我们可以做出下面的一个损失函数：

换言之，我们把对x(i)的估计值与真实值y(i)差的平方和作为损失函数，前面乘上的1/2是为了在求导的时候，这个系数就不见了。

如何调整θ以使得J(θ)取得最小值有很多方法，其中有最小二乘法(min square)，是一种完全是数学描述的方法，另外一种就是梯度下降法。

梯度下降法的算法流程如下：

首先对θ赋值，这个值可以是随机的，也可以让θ是一个全零的向量。
改变θ的值，使得J(θ)按梯度下降的方向进行减少。

为了描述的更清楚，给出下面的图：

这是一个表示参数θ与误差函数J(θ)J(θ)J(θ)的关系图，红色的部分是表示J(θ)J(θ)J(θ)有着比较高的取值，我们需要的是，能够让J(θ)的值尽量的低，也就是达到深蓝色的部分。θ0θ0θ_0，θ1θ1θ_1表示θθθ向量的两个维度。

在上面提到梯度下降法的第一步是给θ给一个初值，假设随机给的初值是在图上的十字点。

然后我们将θ按照梯度下降的方向进行调整，就会使得J(θ)" role="presentation" style="position: relative;">J(θ)J(θ)J(θ)往更低的方向进行变化，如下图所示，算法的结束将是在θθθ下降到无法继续下降为止。

当然，可能梯度下降的最终点并非是全局最小点，即也可能是一个局部最小点，如下图所示：

上面这张图就是描述的一个局部最小点，这是我们重新选择了一个初始点得到的，看来我们这个算法将会在很大的程度上被初始点的选择影响而陷入局部最小点。

下面我将用一个例子描述一下梯度减少的过程，对于我们的函数J(θ)求偏导J：

下面是更新的过程，也就是θi会向着梯度最小的方向进行减少。θi表示更新之前的值，-后面的部分表示按梯度方向减少的量，α表示步长，也就是每次按照梯度减少的方向变化多少。

一个很重要的地方值得注意的是，梯度是有方向的，对于一个向量θ，每一维分量θi都可以求出一个梯度的方向，我们就可以找到一个整体的方向，在变化的时候，我们就朝着下降最多的方向进行变化就可以达到一个最小点，不管它是局部的还是全局的。

用更简单的数学语言进行描述步骤2）是这样的：

梯度下降法找到的一定是下降最快的方向么？机器学习 ML基础中
梯度下降法并不是下降最快的方向，它只是目标函数在当前的点的切平面（当然高维问题不能叫平面）上下降最快的方向。在practical implementation中，牛顿方向（考虑海森矩阵）才一般被认为是下降最快的方向，可以达到superlinear的收敛速度。梯度下降类的算法的收敛速度一般是linear甚至sublinear的（在某些带复杂约束的问题）。

6、拟牛顿法

6.1 牛顿法

假设x∗" role="presentation" style="position: relative;">x∗x∗x^*为函数f(x)f(x)f(x)的根，那么有f(x∗)=0f(x∗)=0f(x^*)=0，在xkxkx_k处进行一阶展开，有：

f(x)=f(xk)+f′(xk)(x−xk)f(x)=f(xk)+f′(xk)(x−xk)

f(x)=f(x_k)+f'(x_k)(x-x_k)

假定xk+1xk+1x_{k+1}为该方程的根，则有：

f(xk+1)=f(xk)+f′(xk)(xk+1−xk)=0f(xk+1)=f(xk)+f′(xk)(xk+1−xk)=0

f(x_{k+1})=f(x_k)+f'(x_k)(x_{k+1}-x_k)=0

那么就有：

xk+1=xk−f(xk)f′(xk)xk+1=xk−f(xk)f′(xk)

x_{k+1}=x_k-\frac{f(x_k)}{f'(x_k)}

于是就有了一个递归方程，可以通过迭代的方式不断逼近f(x)f(x)f(x)的解

牛顿法的思想：因为函数的二阶导数反映了函数的凹凸性，而二阶导数越大，一阶导数的变化越大，所以在搜索中，利用一阶导/二阶导的形式，可以做些搜索方向的修正。

牛顿法的迭代过程：已知xk处的函数值，和该点的一阶导和二阶导，可以根据这三个信息，可以拟合得到蓝色的线，并且可以求得其极小点，极值点对应xk+dk，这就是一次迭代，就是xk->xk+dk，在xk+dk点，再次通过二次曲线拟合，直至迭代至收敛

牛顿法最优值的步骤如下：

随机选取起始点x0x0x 0
计算目标函数f(x)f(x)f(x) 在该点xkxkx_ k 的一阶导数和海森矩阵；
依据迭代公式xk+1=xk−H−1kf'kxk+1=xk−Hk−1f′kx_{ k+1} =x _k −H ^{−1}_ k f ′ _k 更新xxx 值

如果E(f(xk+1)−f(xk))<ϵE(f(xk+1)−f(xk))<ϵE(f(x _{k+1} )−f(x_ k ))，则收敛返回，否则继续步骤2,3直至收敛

6.2 拟牛顿法

雅克比矩阵：函数对各自变量的一阶导数

Hessian矩阵：函数对自变量的二次微分

当Hessian是正定的(所有特征值都是正的)，则该临界点是局部极小点。
当Hessian是负定的(所有特征值都是负的)，则该临界点是局部极大点。

我们这里讨论的点都是临界点，即f’(x) = 0的点，且函数都是连续可导的。

拟牛顿法的思想：

我们可以看到，当我们的特征特别多的时候，求海森矩阵的逆的运算量是非常大且慢的，这对于在实际应用中是不可忍受的，因此我们想能否用一个矩阵来代替海森矩阵的逆呢，这就是拟牛顿法的基本思路。

求Hessian矩阵的逆会影响计算效率，同时，只有搜索方向满足和负梯度的夹角小于90°即可，所有可以用近似矩阵代替Hessian矩阵，只要其满足正定、容易求逆即可。

我们要选择一个矩阵来代替海森矩阵的逆，那么我们首先要研究一下海森矩阵需要具有什么样的特征才能保证牛顿法成功的应用。通过上面的描述我们将下面的公式称为拟牛顿条件：

因此对于选择的代替矩阵GkGkG_k，需要满足两个条件：

拟牛顿条件，即Gk(f′(xk+1)−f′(xk))=xk+1−xkGk(f′(xk+1)−f′(xk))=xk+1−xkG_k(f'(x^{k+1})-f'(x_k))=x^{k+1}-x^{k}
要保证GkGkG_k为正定矩阵，因为只有保证为正定矩阵，才能保证牛顿法的搜索方向是下降的。

如何得到替代矩阵：

Goldfeld等于提出，可以用正定矩阵Gk+vkIGk+vkIG_k+v_kI来代替Hessian矩阵GkGkG_k，计算修正的牛顿方向。

根据下式：

f′(xk)=f′(xk+1)+f′′(xk+1)(xk−xk+1)f′(xk)=f′(xk+1)+f″(xk+1)(xk−xk+1)

f'(x_k)=f'(x_{k+1})+f''(x_{k+1})(x_k-x_{k+1})

可得：

gk=gk+1+Hk+1(xk−xk+1)gk=gk+1+Hk+1(xk−xk+1)

g_k=g_{k+1}+H_{k+1}(x_k-x_{k+1})

所以有：

gk+1−gk=Hk+1(xk+1−xk)gk+1−gk=Hk+1(xk+1−xk)

g_{k+1}-g_{k}=H_{k+1}(x_{k+1}-x_k)

为了简化下面的符合表达式，令：

yk=gk+1−gkyk=gk+1−gk

y_k=g_{k+1}-g_{k}

sk=xk+1−xksk=xk+1−xk

s_k=x_{k+1}-x_k

Bk+1=f′′(xk+1)Bk+1=f″(xk+1)

B_{k+1}=f''(x_{k+1})

Dk+1=B−1k+1Dk+1=Bk+1−1

D_{k+1}=B_{k+1}^{-1}

于是有：

yk=Bk+1⋅skyk=Bk+1⋅sk

y_k=B_{k+1}\cdot s_k

sk=Dk+1⋅yksk=Dk+1⋅yk

s_k=D_{k+1} \cdot y_k

其中：Bk+1Bk+1B_{k+1}为Hessian矩阵，Dk+1Dk+1D_{k+1}为Hessian矩阵的逆

我们可以看出，Hessian矩阵满足上面两个式子

如何在上面两个条件的基础上构造Hessian矩阵呢：DFP和BFGS

1、DFP
该算法的核心思想在于通过迭代的方法，逐步改善近似海塞矩阵。我们使用D来表示近似海塞矩阵的逆矩阵。

2、BFGS

BFGS算法是用直接逼近海塞矩阵的方式来构造近似海塞矩阵，同样，我们使用迭代的方式来逐步逼近。