前言
牛顿法和拟牛顿法是两种常用的优化方法,可以用来求解函数的根以及最优化。
牛顿法
考虑无约束优化问题
min
x
∈
R
n
f
(
x
)
\min_{x\in R^n} f(x)
x∈Rnminf(x)
x
∗
x^*
x∗为目标函数的极小点。
假设f(x)具有二阶连续偏导数,若第k次迭代值为
x
(
k
)
x^{(k)}
x(k),则可将f(x)在
x
(
k
)
x^{(k)}
x(k)附近进行二阶泰勒展开:
f
(
x
)
=
f
(
x
(
k
)
)
+
g
k
T
(
x
−
x
(
k
)
)
+
1
2
(
x
−
x
(
k
)
)
T
H
(
x
(
k
)
)
(
x
−
x
(
k
)
)
f(x)=f(x^{(k)})+g_k^T(x-x^{(k)})+\frac{1}{2}(x-x^{(k)})^TH(x^{(k)})(x-x^{(k)})
f(x)=f(x(k))+gkT(x−x(k))+21(x−x(k))TH(x(k))(x−x(k))
其中,
g
k
=
g
(
x
(
x
)
)
=
∇
f
(
x
(
k
)
)
g_k=g(x^{(x)})=\nabla f(x^{(k)})
gk=g(x(x))=∇f(x(k))是
f
(
x
)
f(x)
f(x)的梯度向量在点
x
(
k
)
x^{(k)}
x(k)的值,
H
(
x
(
k
)
)
H(x^{(k)})
H(x(k))是
f
(
x
)
f(x)
f(x)的黑塞矩阵
H
(
x
)
=
[
∂
2
f
∂
x
i
∂
x
j
]
n
×
n
H(x)=[\frac{\partial^2f}{\partial x_i\partial x_j}]_{n\times n}
H(x)=[∂xi∂xj∂2f]n×n
在点
x
(
k
)
x^{(k)}
x(k)的值。函数
f
(
x
)
f(x)
f(x)有极值的必要条件是在极值处一阶导数为0,即梯度向量为0,当
H
(
x
(
k
)
)
H(x^{(k)})
H(x(k))是正定矩阵时,函数
f
(
x
)
f(x)
f(x)的极值为极小值。
牛顿法利用极小点的必要条件
∇
f
(
x
)
=
0
\nabla f(x)=0
∇f(x)=0
每次迭代中从点
x
k
x^{k}
xk开始,求目标函数的极小点,作为第k+1次迭代值
x
(
k
+
1
)
x^{(k+1)}
x(k+1)。具体地,假设
x
(
k
+
1
)
x^{(k+1)}
x(k+1)满足:
∇
f
(
x
(
k
+
1
)
)
=
0
\nabla f(x^{(k+1)})=0
∇f(x(k+1))=0
由之前的泰勒展开式可以得到:
∇
f
(
x
)
=
g
k
+
H
k
(
x
−
x
(
k
)
)
\nabla f(x)=g_k+H_k(x-x^{(k)})
∇f(x)=gk+Hk(x−x(k))
(这就是数值对向量求导的应用,我之前在矩阵求导里头总结过)。
其中,
H
k
=
H
(
x
(
k
)
)
H_k=H({x^{(k)}})
Hk=H(x(k)).
于是,
g
k
+
H
k
(
x
(
k
+
1
)
−
x
(
k
)
)
=
0
g_k+H_k(x^{(k+1)}-x^{(k)})=0
gk+Hk(x(k+1)−x(k))=0
因此,
x
(
k
+
1
)
=
x
(
k
)
−
H
k
−
1
g
k
x^{(k+1)}=x^{(k)}-H^{-1}_k g_k
x(k+1)=x(k)−Hk−1gk
或者
x
(
k
+
1
)
=
x
(
k
)
+
p
k
x^{(k+1)}=x^{(k)}+p_k
x(k+1)=x(k)+pk
其中,
H
k
p
k
=
−
g
k
H_kp_k=-g_k
Hkpk=−gk
用上式做迭代公式的算法就是牛顿法。
算法的流程如下:
输入:目标函数
f
(
x
)
f(x)
f(x),梯度
g
(
x
)
=
∇
f
(
x
)
g(x)=\nabla f(x)
g(x)=∇f(x),黑塞矩阵
H
(
x
)
H(x)
H(x),精度要求
ϵ
\epsilon
ϵ;
输出:
f
(
x
)
f(x)
f(x)的极小点
- 取初始值 x ( 0 ) x^{(0)} x(0),置k=0
- 计算 g k = g ( x ( k ) ) g_k=g(x^{(k)}) gk=g(x(k))
- 若 ∣ ∣ g k ∣ ∣ < ϵ ||g_k||\lt \epsilon ∣∣gk∣∣<ϵ,则停止计算,得近似解 x ∗ = x ( k ) x^*=x^{(k)} x∗=x(k)
- 计算
H
k
=
H
(
x
(
k
)
)
H_k=H(x^{(k)})
Hk=H(x(k)),并求
p
k
p_k
pk
H k p k = − g k H_kp_k=-g_k Hkpk=−gk - 置 x k + 1 = x ( k ) + p k x^{{k+1}}=x^{(k)}+p_k xk+1=x(k)+pk
- 置 k = k + 1 k=k+1 k=k+1,转 2 。
优缺点
- 牛顿法的优点是收敛速度快,因为求解过程考虑了二阶导数,而传统的梯度下降法只考虑了一阶导数。
- 缺点是计算复杂度高,因为要求解黑塞矩阵的逆,所以复杂度较高。机器学习中的有关矩阵的操作一般复杂度都会比较高,所以一般都会采用迭代的方法来求解近似解。
- 其次,牛顿法要求函数必须是凸的,否则上述求解过程不会收敛。
拟牛顿法
针对黑塞矩阵的逆的求解复杂问题,人们又提出了拟牛顿法,思路就是利用一个n阶矩阵
G
k
=
G
(
x
(
k
)
)
G_k=G(x^{(k)})
Gk=G(x(k))来近似替代
H
k
−
1
=
H
−
1
(
x
(
k
)
)
H^{-1}_k=H^{-1}(x^{(k)})
Hk−1=H−1(x(k))
考虑下式:
∇
f
(
x
)
=
g
k
+
H
k
(
x
−
x
(
x
)
)
\nabla f(x)=g_k+H_k(x-x^{(x)})
∇f(x)=gk+Hk(x−x(x))
其中,
g
k
=
∇
f
(
x
(
x
)
)
g_k=\nabla f(x^{(x)})
gk=∇f(x(x)),将上式中的
x
x
x替换为
x
(
k
)
x^{(k)}
x(k),可得
g
k
+
1
−
g
k
=
H
k
(
x
(
k
+
1
)
−
x
(
k
)
)
g_{k+1}-g_k=H_k(x^{(k+1)}-x^{(k)})
gk+1−gk=Hk(x(k+1)−x(k))
记
y
k
=
g
k
+
1
−
g
k
y_k=g_{k+1}-g_k
yk=gk+1−gk,
δ
k
=
x
(
k
+
1
)
−
x
(
k
)
\delta_k=x^{(k+1)}-x^{(k)}
δk=x(k+1)−x(k),则
y
k
=
H
k
δ
k
y_k=H_k\delta_k
yk=Hkδk
或
H
k
−
1
y
k
=
δ
k
H_k^{-1}y_k=\delta_k
Hk−1yk=δk
上式被称为拟牛顿条件
为了减少计算,我们不直接计算 H k − 1 H_k^{-1} Hk−1,选取 H k − 1 H_k^{-1} Hk−1的近似有两种选择,于是就有了两种不同的算法。
-
DFP算法
在DFP算法中,我们用 G k G_k Gk作为 H k − 1 H_k^{-1} Hk−1的近似, G k G_k Gk的迭代式如下:
G k + 1 = G k + P k + Q k G_{k+1}=G_k+P_k+Q_k Gk+1=Gk+Pk+Qk
于是有
G k + 1 y k = G k y k + P k y k + Q k y k G_{k+1}y_k=G_ky_k+P_ky_k+Q_ky_k Gk+1yk=Gkyk+Pkyk+Qkyk
其中, P k P_k Pk和 Q k Q_k Qk的满足以下条件:
P k y k = δ k P_ky_k=\delta_k Pkyk=δk
Q k y k = − G k y k Q_ky_k=-G_ky_k Qkyk=−Gkyk
为了满足以上条件,取
P k = δ k δ k T δ k T y k P_k=\frac{\delta_k\delta_k^T}{\delta_k^Ty_k} Pk=δkTykδkδkT
Q k = − G k y k y k T G k y k T G k y k Q_k=-\frac{G_ky_ky_k^TG_k}{y_k^TG_ky_k} Qk=−ykTGkykGkykykTGk
于是, G k G_k Gk的迭代式如下:
G k + 1 = G k + δ k δ k T δ k T y k − G k y k y k T G k y k T G k y k G_{k+1}=G_k+\frac{\delta_k\delta_k^T}{\delta_k^Ty_k}-\frac{G_ky_ky_k^TG_k}{y_k^TG_ky_k} Gk+1=Gk+δkTykδkδkT−ykTGkykGkykykTGk
(这种思想的确很厉害,一般人难以想到) -
BFGS算法
BFGS的思想与DFP相似,都是采用近似的方法减少计算复杂度,不过近似的思路不太一样。
BFGS采用 B k B_k Bk近似 H H H,拟牛顿条件如下:
B k + 1 δ k = y k B_{k+1}\delta_k=y_k Bk+1δk=yk
迭代式如下:
B k + 1 = B k + P k + Q k B_{k+1}=B_k+P_k+Q_k Bk+1=Bk+Pk+Qk
B k + 1 δ k = B k δ k + P k δ k + Q k δ k B_{k+1}\delta_k=B_k\delta_k+P_k\delta_k+Q_k\delta_k Bk+1δk=Bkδk+Pkδk+Qkδk
P k P_k Pk和 Q k Q_k Qk满足如下条件:
P k δ k = y k P_k\delta_k=y_k Pkδk=yk
Q k δ k = − B k δ k Q_k\delta_k=-B_k\delta_k Qkδk=−Bkδk
于是,最终的迭代式如下:
B k + 1 = B k + y k y k T y k T δ k − B k δ k δ k T B k δ k T B k δ k B_{k+1}=B_k+\frac{y_ky_k^T}{y_k^T\delta_k}-\frac{B_k\delta_k\delta_k^TB_k}{\delta_k^TB_k\delta_k} Bk+1=Bk+ykTδkykykT−δkTBkδkBkδkδkTBk
总结
其实牛顿法和拟牛顿法的核心思想还是近似,尤其是对于导数和矩阵的近似,这在其他的机器学习问题里多有体现,以后遇到相关的问题再讨论。关于相关的代码等以后有时间了再写,对于牛顿法和拟牛顿法我目前也没有完全弄懂,等以后理解更深刻了再重写这篇博文。
参考资料
- 统计机器学习
- 牛顿法 知乎