【数值优化】局部搜索之牛顿法

除了前面说的梯度下降法,牛顿法也是机器学习中用的比较多的一种优化算法。

牛顿法求方程解

牛顿法又称为牛顿-拉弗森方法(Newton-Raphson method),单变量下又称为切线法。
在牛顿迭代法中,先任选一点 x0 作为根的近似估值,过点 (x0, f(x0)) 做曲线的切线,那么切线的斜率即为 f’(x0)。然后取切线与x轴的交点 x1作为新的估值点,以此类推,我们可以找到一点 xn,使得 |f(xn)| 趋于0,也就是说如果我们将上述操作进行下去的话,那么我们离多项式的根就越来越近了。
$$X_{n+1} = X_{n}-\frac{函数值}{函数一阶偏导值}$$

牛顿法应用于最优化

牛顿法求解非线性函数的最优值点。那牛顿法和极值求解有关系?看起来牛顿法只能求零点啊? 一阶导零点不就是函数的极值或者驻点?
将求f(x) = 0的解修改为f’(x) = 0得到牛顿法求极值的迭代公式

$X_{n+1} = X_{n}-\frac{函数一阶偏导值}{函数二阶偏导值}$

对于高维函数,用牛顿法求极值也是这个形式,只不过这里的y’和y’’都变成了向量和矩阵。
y’就变成了对所有参数的偏导数组成的向量

y’’写成矩阵的形式,称为Hessian矩阵

对比标准的牛顿法和梯度法
从本质上去看,牛顿法是二阶收敛,梯度下降是一阶收敛,所以牛顿法就更快。如果更通俗地说的话,比如你想找一条最短的路径走到一个盆地的最底部,梯度下降法每次只从你当前所处位置选一个坡度最大的方向走一步,牛顿法在选择方向时,不仅会考虑坡度是否够大,还会考虑你走了一步之后,坡度是否会变得更大。所以,可以说牛顿法比梯度下降法看得更远一点,能更快地走到最底部。(牛顿法目光更加长远,所以少走弯路;相对而言,梯度下降法只考虑了局部的最优,没有全局思想。)
从几何上说,牛顿法就是用一个二次曲面去拟合你当前所处位置的局部曲面,而梯度下降法是用一个平面去拟合当前的局部曲面,通常情况下,二次曲面的拟合会比平面更好,所以牛顿法选择的下降路径会更符合真实的最优下降路径。

牛顿迭代法的优缺点

1、优点

  • 二阶收敛,收敛速度快;
  • 精确度高

2、缺点

  • 作为一种迭代算法,每一步都需要求解目标函数的Hessian矩阵的逆矩阵,计算比较复杂。二阶方法实践中对高维数据不可行。
  • 可能发生被零除错误。当函数在它的零点附近,导函数的绝对值非常小时,运算会出现被零除错误。
  • 可能出现死循环。当函数在它的零点有拐点时,可能会使迭代陷入死循环。
  • 定步长迭代。改进是阻尼牛顿法。

牛顿法的代码实现

牛顿法的速度相当快,而且能高度逼近最优值。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
# -*- coding: utf-8 -*-
'''
使用Newton法求解函数极值
'''
import numpy as np
import matplotlib.pyplot as plt

#========================================================
# 给定目标函数
#========================================================

X1=np.arange(-1.5,1.5+0.05,0.05)
X2=np.arange(-3.5,2+0.05,0.05)
[x1,x2]=np.meshgrid(X1,X2)
f = x1**2+x2**2+3*x1+4*x2-26; # 给定的函数
plt.contour(x1, x2, f, 20) # 画出函数的20条轮廓线

#========================================================
# 函数一阶导与二阶导
#========================================================

def f1(x):
'''
一阶导
'''
return np.array([2*x[0]+3,2*x[1]+4])

def f2(x):
'''
二阶导,Hessian矩阵
'''
return np.array([[2,0],[0,2]])

#========================================================
# 牛顿法
#========================================================

def newton(x0):
N = 1000 # 迭代次数
i = 1

W = np.zeros((2, N)) # 存储中间
W[:,0] = x0
x = x0

delta = 1 # 斜率

while i<N and delta>0.1:
p = np.dot(np.linalg.inv(f2(x)), f1(x)) # 使用inv函数计算Hessian矩阵的逆矩阵然后与一阶导向量点积
x = x - p
delta = sum((x-x0))
print('第'+str(i)+'次迭代结果: %s', x)

# 存储中间结果
W[:, i] = x
x0 = x

i=i+1
W = W[:, 0:i] # 记录迭代点
return W

#========================================================
# 主程序
#========================================================
if __name__ == '__main__':
x0 = np.array([1,1])
W=newton(x0)

plt.plot(W[0,:],W[1,:],'g*',W[0,:],W[1,:]) # 画出迭代点收敛的轨迹
plt.show()
0%