文章目录
Momentum
Momentum的迭代公式为:
vt=γvt−1+η∇θJ(θ)θ=θ−vt
其中J(⋅)一般为损失函数。我们知道,一般的梯度下降,是没有γvt−1这一项的,有了这一项之后,θ的更新和前一次更新的路径有关,使得每一次更新的方向不会出现剧烈变化,所以这种方法在函数分布呈梭子状的时候非常有效。
先来看看这个函数利用梯度下降的效果吧。
import matplotlib.pyplot as plt
import numpy as np
"""
z = x^2 + 50 y ^2
2x
100y
"""
partial_x = lambda x: 2 * x
partial_y = lambda y: 100 * y
partial = lambda x: np.array([partial_x(x[0]),
partial_y(x[1])])
f = lambda x: x[0] ** 2 + 50 * x[1] ** 2
class Decent:
def __init__(self, function):
self.__function = function
@property
def function(self):
return self.__function
def __call__(self, x, grad, alpha=0.4, beta=0.7):
t = 1
fx = self.function(x)
dist = - grad @ grad
while True:
dx = x - t * grad
fdx = self.function(dx)
if fdx <= fx + alpha * t * dist:
break
else:
t *= beta
return dx
grad_decent = Decent(f)
x = np.array([30., 15.])
process = []
while True:
grad = partial(x)
if np.sqrt(grad @ grad) < 1e-7:
break
else:
process.append(x)
x = grad_decent(x, grad)
process = np.array(process)
print(len(process))
x = np.linspace(-40, 40, 1000)
y = np.linspace(-20, 20, 500)
fig, ax= plt.subplots()
X, Y = np.meshgrid(x, y)
ax.contour(X, Y, f([X, Y]), colors='black')
ax.plot(process[:, 0], process[:, 1])
plt.show()
怎么说呢,有点震荡?289步1e-7的误差
x = np.array([30., 15.])
process = []
v = 0
gamma = 0.7
eta = 0.016
while True:
grad = partial(x)
v = gamma * v + eta * grad
if np.sqrt(grad @ grad) < 1e-7:
break
else:
process.append(x)
x = x - v
用117步,话说,这个参数是不是难调啊,感觉一般η很小啊。
还有一个很赞的分析,在博客:
路遥知马力-Momentum
Nesterov accelerated gradient
NGD的迭代公式是:
vt=γvt−1+η∇θJ(θ−γvt−1)θ=θ−vt
和上面的区别就是,第t步更新,我们关心的是下一步(一个近似)的梯度,而不是当前点的梯度,我之前以为这是有一个搜索的过程的,但是实际上没有,所以真的是这个式子具有前瞻性?或许真的和上面博客里说的那样,因为后面的部分可以看成一个二阶导的近似。
x = np.array([30., 15.])
process = []
v = 0
gamma = 0.7
eta = 0.013
while True:
grad = partial(x-gamma*v)
v = gamma * v + eta * grad
if np.sqrt(grad @ grad) < 1e-7:
break
else:
process.append(x)
x = x - v
感觉没有momentum好用啊
NESTEROV 的另外一个方法?
在那个overview里面,引用的是
但是里面的方面感觉不是NGD啊,不过的确是一种下降方法,所以讲一下吧。
假设f(x)满足其一阶导函数一致连续的凸函数,比如用以下条件表示:
∥f′(x)−f′(y)∥≤L∥x−y∥,∀x,y∈E
由此可以推得(不晓得这个0.5哪来的,虽然有点像二阶泰勒展式,但是呢,凸函数好像没有这性质吧,去掉0.5就可以直接证出来了,而且这个0.5对证明没有什么大的影响吧,因为只要让L=0.5L就可以了啊):
f(y)≤f(x)+<f′(x),y−x>+0.5L∥y−x∥2(2)
为了解决min{f(x)∣x∈E},且最优解X∗非空的情况,我们可以:
- 首先选择一个点y0∈E,并令
k=0,a0=1,x−1=y0,α−1=∥y0−z∥/∥f′(y0)−f′(z)∥
其中z是E中不同于y0的任意点,且f′(y0)̸=f′(z) - 第k 步:
a) 计算最小的i满足:
f(yk)−f(yk−2−iαk−1f′(yk))≥2−i−1αk−1∥f′(yk)∥2(4)
b) 令
αk=2−iαk−1,xk=yk−αkf′(yk)ak+1=(1+4ak2+1)/2yk+1=xk+(ak−1)(xk−xk−1)/ak+1.
作者证明了这个方法的收敛率是O(1/k2)的。
即在满足上面提到的假设,且利用上面给出的方法所求,可以证明,对于任意的k≥0:
f(xk)−f∗≤C/(k+2)2
其中C=4L∥y0−x∗∥2并且f∗=f(x∗),x∗∈X∗。
还有一些关于收敛步长的分析就不贴了。
证明:
令yk(α)−αf′(yk), 可以得到(通过(2)):
f(yk)−f(yk(α))≥0.5α(2−αL)∥f′(yk)∥2
结果就是, 只要2−iαk−1≤L−1,不等式(4)就成立,也就是说αk≥0.5L−1,∀k≥0, 否则2−iαk−1>L−1。
令pl=(ak−1)(xk−1−xk),则yk+1=xk−pk/ak+1,于是:
pk+1−xk+1=pk−xk+ak+1αk+1f′(yk+1)
于是:
∥pk+1−xk+1+x∗∥2=∥pk−xk+x∗∥2+2(ak+1−1)αk+1<f′(yk+1,pk>+2ak+1αk+1<f′(yk+1,x∗−yk+1>+ak+12αk+12∥f′(yk+1)∥2
利用不等式(4)和f(x)的凸性,可得:
<f′(yk+1),yk+1−x∗>0.5αk+1∥f′(yk+1)∥2≥f(xk+1)−f∗+0.5αk+1∥f′(yk+1)∥2(5)≤f(yk+1)−f(xk+1)≤f(xk)−f(xk+1)−ak+1−1<f′(yk+1,pk>(6)
其中第一个不等式,先利用凸函数得性质:
f∗≥f(yk+1)+<f′(yk+1),x∗−yk+1)
再利用不等式(4):
f(yk+1)−f(xk+1)≥0.5αk+1∥f′(yk+1)∥2
代入这俩个不等式可得:
∥pk+1−xk+1+x∗∥2−∥pk−xk+x∗∥2≤2(ak+1−1)αk+1<f′(yk+1),pk>−2ak+1αk+1(f(xk+1−f∗)+(ak+12−ak+1)αk+12∥f′(yk+1)∥2≤−2ak+1αk+1(f(xk+1)−f∗)+2(ak+12−ak+1)αk+1(f(xk)−f(xk+1))=2αk+1ak2(f(xk)−f∗)−2αk+1ak+12(f(xk+1)−f∗)≤2αkak2(f(xk)−f∗)−2αk+1ak+12(f(xk+1)−f∗)
其中第一个不等式用到了(5), 第二个不等式用到了(6), 等式用到了ak+12−ak+1=ak2,最后一步用到了αk≤αk+1且f(xk)≥f∗。
因此:
2αk+1ak+12(f(xk+1)−f∗)≤2αk+1ak+12(f(xk+1)−f∗)+∥pk+1−xk+1+x∗∥2≤2αkak(f(xk)−f∗)+∥pk−xk+x∗∥2≤2α0a02(f(x0)−f∗)+∥p0−x0+x∗∥2≤∥y0−x∗∥2.
最后一个不等式成立是因为p0=0,x0=y0,左边第一项大于等于0.
又αk≥0.5L−1,ak+1≥ak+0.5≥1+0.5(k+1),所以:
f(xk+1)−f∗≤C/(k+3)2
证毕。
class Decent:
def __init__(self, function, grad):
self.__function = function
self.__grad = grad
self.process = []
@property
def function(self):
return self.__function
@property
def grad(self):
return self.__grad
def __call__(self, y, z):
def find_i(y, alpha):
i = 0
fy = self.function(y)
fdy = self.grad(y)
fdynorm = fdy @ fdy
while True:
temp = self.function(y - 2 ** (-i) * alpha * fdy)
if fy - temp > 2 ** (-i -1) * alpha * fdynorm:
return i, fdy
else:
i += 1
a = 1
x = y
fdy = self.grad(y)
fdz = self.grad(z)
alpha = np.sqrt((y-z) @ (y-z) /
(fdy-fdz) @ (fdy - fdz))
k = 0
while True:
k += 1
self.process.append(x)
i, fdy = find_i(y, alpha)
if np.sqrt(fdy @ fdy) < 1:
print(k)
return x
alpha = 2 ** (-i) * alpha
x_old = np.array(x, dtype=float)
x = y - alpha * fdy
a_old = a
a = (1 + np.sqrt(4 * a ** 2 + 1)) / 2
y = x + (a_old - 1) * (x - x_old) / a
grad_decent = Decent(f, partial)
x = np.array([30., 15.])
z = np.array([200., 10.])
grad_decent(x, z)
process = np.array(grad_decent.process)
x = np.linspace(-40, 40, 1000)
y = np.linspace(-20, 20, 500)
X, Y = np.meshgrid(x, y)
fig, ax = plt.subplots()
ax.contour(X, Y, f([X, Y]), colors="black")
ax.scatter(process[:, 0], process[:, 1])
ax.plot(process[:, 0], process[:, 1])
plt.show()
用了30步就能到达上面的情况,不过呢,如果想让∥f′(x)∥≤1e−7得1000多步,主要是因为会来回振荡。