python 中的数值梯度下降
Numeric Gradient Descent in python
我写这段代码是为了得到向量函数的梯度下降
。
其中:f是函数,X0是起点,eta是步长。
它基本上由两部分组成,第一部分获取函数的梯度,第二部分迭代 x,减去梯度。
问题是您通常无法收敛某些函数,例如:
如果我们取
,梯度下降不会收敛到[20,25]
我需要更改或添加什么?
def descenso_grad(f,X0,eta):
def grad(f,X):
import numpy as np
def partial(g,k,X):
h=1e-9
Y=np.copy(X)
X[k-1]=X[k-1]+h
dp=(g(X)-g(Y))/h
return dp
grd=[]
for i in np.arange(0,len(X)):
ai=partial(f,i+1,X)
grd.append(ai)
return grd
#iterations
i=0
while True:
i=i+1
X0=X0-eta*np.array(grad(f,X0))
if np.linalg.norm(grad(f,X0))<10e-8 or i>400: break
return X0
你的梯度下降实现是一个很好的基本实现,但你的梯度有时会振荡和爆炸。首先我们应该明确你的梯度下降并不总是发散。对于eta
和X0
的某些组合,它实际上收敛了。
但首先让我建议对代码进行一些修改:
import numpy as np
语句应该在文件的顶部,而不是在函数中。一般来说,任何导入语句都应该在代码的开头,这样它们只执行一次
- 最好不要写嵌套函数,而是分开写:可以把
partial
函数写在grad
函数外面,grad
函数写在descendo_grad
函数。更利于调试
- 我强烈建议传递学习率(
eta
)、步数(steps
)和容差(在代码中设置为10e-8,或1e)等参数-7) 作为 descendo_grad
函数的参数。这样您就可以比较它们对结果的影响。
无论如何,这是我将在这个答案中使用的代码的实现:
import numpy as np
def partial(g, k, X):
h = 1e-9
Y = np.copy(X)
X[k - 1] = X[k - 1] + h
dp = (g(X) - g(Y)) / h
return dp
def grad(f, X):
grd = []
for i in np.arange(0, len(X)):
ai = partial(f, i + 1, X)
grd.append(ai)
return grd
def descenso_grad(f,X0,eta, steps, tolerance=1e-7):
#iterations
i=0
while True:
i=i+1
X0=X0-eta*np.array(grad(f,X0))
if np.linalg.norm(grad(f,X0))<tolerance or i>steps: break
return X0
def f(X):
return (X[0]-20)**4 + (X[1]-25)**4
现在,关于收敛。我说过你的实现并不总是有分歧。确实:
X0 = [2, 30]
eta = 0.001
steps = 400
xmin = descenso_grad(f, X0, eta, steps)
print(xmin)
将打印 [20.55359068 25.55258024]
但是:
X0 = [2, 0]
eta = 0.001
steps = 400
xmin = descenso_grad(f, X0, eta, steps)
print(xmin)
实际上会偏离 [ 2.42462695e+01 -3.54879793e+10]
1) 发生了什么
您的渐变实际上是围绕 y 轴振荡的。让我们计算 f
在 X0 = [2, 0]
:
处的梯度
print(grad(f, X0))
我们得到 grad(f, X0) = [-23328.00067961216, -62500.01024454831]
,这是相当高但方向正确。
现在让我们计算梯度下降的下一步:
eta = 0.001
X1=X0-eta*np.array(grad(f,X0))
print(X1)
我们得到 X1 = [25.32800068 62.50001025]
。我们可以看到,在 x 轴上,我们实际上更接近最小值,但在 y 轴上,梯度下降跳到了最小值的另一侧,甚至离它更远。实际上,X0[1]
距离最小值 (X0[1] - Xmin[1] = 25
) 在其左侧 的距离为 25,而 X0[1]
现在的距离为 65-25 = 40
但在它的右边*。由于f
画出的曲线绕y轴呈简单的U形,因此X1
中f
取的值会比之前高(为了简化,我们忽略了x 坐标)。
如果我们看一下接下来的步骤,我们可以清楚地看到围绕最小值的爆炸性振荡:
X0 = [2, 0]
eta = 0.001
steps = 10
#record values taken by X[1] during gradient descent
curve_y = [X0[1]]
i = 0
while True:
i = i + 1
X0 = X0 - eta * np.array(grad(f, X0))
curve_y.append(X0[1])
if np.linalg.norm(grad(f, X0)) < 10e-8 or i > steps: break
print(curve_y)
我们得到 [0, 62.50001024554831, -148.43710232226067, 20719.6258707022, -35487979280.37413]
。我们可以看到 X1 在围绕它振荡时离最小值越来越远。
为了说明这一点,我们假设 x 轴上的值是固定的,并且只看 y 轴上发生的情况。图片以黑色显示了在梯度下降的每一步中获取的函数值的振荡(合成数据仅用于说明目的)。由于更新值太大,梯度下降使我们在每一步都离最小值更远:
请注意,我们作为示例给出的梯度下降只有 5 步,而我们编程了 10 步。这是因为当函数取的值太大时,python 无法成功区分 f(X[1])
和 f(X[1]+h)
,因此它计算出一个等于零的梯度:
x = 24 # for the example
y = -35487979280.37413
z = f([x, y+h]) - f([x, y])
print(z)
我们得到 0.0
。这个问题是关于计算机的计算精度的问题,但我们稍后再谈。
因此,这些振荡是由于以下因素的组合造成的:
- 关于 y 轴的局部梯度值非常高
eta
的值太大,无法补偿更新中的爆炸梯度。
如果这是真的,如果我们使用较小的学习率,我们可能会收敛。让我们检查一下:
X0 = [2, 0]
# divide eta by 100
eta = 0.0001
steps = 400
xmin = descenso_grad(f, X0, eta, steps)
print(xmin)
我们会得到[18.25061287 23.24796497]
。 我们可能需要更多步骤,但这次我们正在收敛!!
2) 如何避免?
A) 在您的特定情况下
由于函数形状简单,没有局部最小值或鞍点,我们可以通过简单地剪切梯度值来避免这个问题。这意味着我们为梯度范数定义了一个最大值:
def grad_clipped(f, X, clip):
grd = []
for i in np.arange(0, len(X)):
ai = partial(f, i + 1, X)
if ai<0:
ai = max(ai, -1*clip)
else:
ai = min(ai, clip)
grd.append(ai)
return grd
def descenso_grad_clipped(f,X0,eta, steps, clip=100, tolerance=10e-8):
#iterations
i=0
while True:
i=i+1
X0=X0-eta*np.array(grad_clipped(f,X0, clip))
if np.linalg.norm(grad_clipped(f,X0, clip))<tolerance or i>steps: break
return X0
让我们用发散的例子来测试它:
X0 = [2, 0]
eta = 0.001
steps = 400
clip=100
xmin = descenso_grad_clipped(f, X0, eta, clip, steps)
print(xmin)
这次我们的目标是一致的:[19.31583901 24.20307188]
。请注意,这可能会减慢该过程,因为梯度下降将采取较小的步骤。这里我们可以通过增加步数来更接近真实的最小值。
请注意,此技术还解决了我们在函数值过高时遇到的数值计算问题。
B) 一般
一般来说,梯度下降算法试图避免很多警告(梯度爆炸或消失、鞍点、局部最小值...)。 Adam、RMSprop、Adagrad 等反向传播算法试图避免这些警告。
我不打算深入探讨细节,因为这值得一整篇文章,但是您可以使用以下两个资源(我建议按给定的顺序阅读它们)来加深您对该主题的理解:
- A good article on towardsdatascience.com 解释梯度下降的基础知识及其最常见的缺陷
- An overview of gradient descent algorithms
我写这段代码是为了得到向量函数的梯度下降。
其中:f是函数,X0是起点,eta是步长。
它基本上由两部分组成,第一部分获取函数的梯度,第二部分迭代 x,减去梯度。
问题是您通常无法收敛某些函数,例如:
如果我们取,梯度下降不会收敛到[20,25]
我需要更改或添加什么?
def descenso_grad(f,X0,eta):
def grad(f,X):
import numpy as np
def partial(g,k,X):
h=1e-9
Y=np.copy(X)
X[k-1]=X[k-1]+h
dp=(g(X)-g(Y))/h
return dp
grd=[]
for i in np.arange(0,len(X)):
ai=partial(f,i+1,X)
grd.append(ai)
return grd
#iterations
i=0
while True:
i=i+1
X0=X0-eta*np.array(grad(f,X0))
if np.linalg.norm(grad(f,X0))<10e-8 or i>400: break
return X0
你的梯度下降实现是一个很好的基本实现,但你的梯度有时会振荡和爆炸。首先我们应该明确你的梯度下降并不总是发散。对于eta
和X0
的某些组合,它实际上收敛了。
但首先让我建议对代码进行一些修改:
import numpy as np
语句应该在文件的顶部,而不是在函数中。一般来说,任何导入语句都应该在代码的开头,这样它们只执行一次- 最好不要写嵌套函数,而是分开写:可以把
partial
函数写在grad
函数外面,grad
函数写在descendo_grad
函数。更利于调试 - 我强烈建议传递学习率(
eta
)、步数(steps
)和容差(在代码中设置为10e-8,或1e)等参数-7) 作为descendo_grad
函数的参数。这样您就可以比较它们对结果的影响。
无论如何,这是我将在这个答案中使用的代码的实现:
import numpy as np
def partial(g, k, X):
h = 1e-9
Y = np.copy(X)
X[k - 1] = X[k - 1] + h
dp = (g(X) - g(Y)) / h
return dp
def grad(f, X):
grd = []
for i in np.arange(0, len(X)):
ai = partial(f, i + 1, X)
grd.append(ai)
return grd
def descenso_grad(f,X0,eta, steps, tolerance=1e-7):
#iterations
i=0
while True:
i=i+1
X0=X0-eta*np.array(grad(f,X0))
if np.linalg.norm(grad(f,X0))<tolerance or i>steps: break
return X0
def f(X):
return (X[0]-20)**4 + (X[1]-25)**4
现在,关于收敛。我说过你的实现并不总是有分歧。确实:
X0 = [2, 30]
eta = 0.001
steps = 400
xmin = descenso_grad(f, X0, eta, steps)
print(xmin)
将打印 [20.55359068 25.55258024]
但是:
X0 = [2, 0]
eta = 0.001
steps = 400
xmin = descenso_grad(f, X0, eta, steps)
print(xmin)
实际上会偏离 [ 2.42462695e+01 -3.54879793e+10]
1) 发生了什么
您的渐变实际上是围绕 y 轴振荡的。让我们计算 f
在 X0 = [2, 0]
:
print(grad(f, X0))
我们得到 grad(f, X0) = [-23328.00067961216, -62500.01024454831]
,这是相当高但方向正确。
现在让我们计算梯度下降的下一步:
eta = 0.001
X1=X0-eta*np.array(grad(f,X0))
print(X1)
我们得到 X1 = [25.32800068 62.50001025]
。我们可以看到,在 x 轴上,我们实际上更接近最小值,但在 y 轴上,梯度下降跳到了最小值的另一侧,甚至离它更远。实际上,X0[1]
距离最小值 (X0[1] - Xmin[1] = 25
) 在其左侧 的距离为 25,而 X0[1]
现在的距离为 65-25 = 40
但在它的右边*。由于f
画出的曲线绕y轴呈简单的U形,因此X1
中f
取的值会比之前高(为了简化,我们忽略了x 坐标)。
如果我们看一下接下来的步骤,我们可以清楚地看到围绕最小值的爆炸性振荡:
X0 = [2, 0]
eta = 0.001
steps = 10
#record values taken by X[1] during gradient descent
curve_y = [X0[1]]
i = 0
while True:
i = i + 1
X0 = X0 - eta * np.array(grad(f, X0))
curve_y.append(X0[1])
if np.linalg.norm(grad(f, X0)) < 10e-8 or i > steps: break
print(curve_y)
我们得到 [0, 62.50001024554831, -148.43710232226067, 20719.6258707022, -35487979280.37413]
。我们可以看到 X1 在围绕它振荡时离最小值越来越远。
为了说明这一点,我们假设 x 轴上的值是固定的,并且只看 y 轴上发生的情况。图片以黑色显示了在梯度下降的每一步中获取的函数值的振荡(合成数据仅用于说明目的)。由于更新值太大,梯度下降使我们在每一步都离最小值更远:
请注意,我们作为示例给出的梯度下降只有 5 步,而我们编程了 10 步。这是因为当函数取的值太大时,python 无法成功区分 f(X[1])
和 f(X[1]+h)
,因此它计算出一个等于零的梯度:
x = 24 # for the example
y = -35487979280.37413
z = f([x, y+h]) - f([x, y])
print(z)
我们得到 0.0
。这个问题是关于计算机的计算精度的问题,但我们稍后再谈。
因此,这些振荡是由于以下因素的组合造成的:
- 关于 y 轴的局部梯度值非常高
eta
的值太大,无法补偿更新中的爆炸梯度。
如果这是真的,如果我们使用较小的学习率,我们可能会收敛。让我们检查一下:
X0 = [2, 0]
# divide eta by 100
eta = 0.0001
steps = 400
xmin = descenso_grad(f, X0, eta, steps)
print(xmin)
我们会得到[18.25061287 23.24796497]
。 我们可能需要更多步骤,但这次我们正在收敛!!
2) 如何避免?
A) 在您的特定情况下
由于函数形状简单,没有局部最小值或鞍点,我们可以通过简单地剪切梯度值来避免这个问题。这意味着我们为梯度范数定义了一个最大值:
def grad_clipped(f, X, clip):
grd = []
for i in np.arange(0, len(X)):
ai = partial(f, i + 1, X)
if ai<0:
ai = max(ai, -1*clip)
else:
ai = min(ai, clip)
grd.append(ai)
return grd
def descenso_grad_clipped(f,X0,eta, steps, clip=100, tolerance=10e-8):
#iterations
i=0
while True:
i=i+1
X0=X0-eta*np.array(grad_clipped(f,X0, clip))
if np.linalg.norm(grad_clipped(f,X0, clip))<tolerance or i>steps: break
return X0
让我们用发散的例子来测试它:
X0 = [2, 0]
eta = 0.001
steps = 400
clip=100
xmin = descenso_grad_clipped(f, X0, eta, clip, steps)
print(xmin)
这次我们的目标是一致的:[19.31583901 24.20307188]
。请注意,这可能会减慢该过程,因为梯度下降将采取较小的步骤。这里我们可以通过增加步数来更接近真实的最小值。
请注意,此技术还解决了我们在函数值过高时遇到的数值计算问题。
B) 一般
一般来说,梯度下降算法试图避免很多警告(梯度爆炸或消失、鞍点、局部最小值...)。 Adam、RMSprop、Adagrad 等反向传播算法试图避免这些警告。
我不打算深入探讨细节,因为这值得一整篇文章,但是您可以使用以下两个资源(我建议按给定的顺序阅读它们)来加深您对该主题的理解:
- A good article on towardsdatascience.com 解释梯度下降的基础知识及其最常见的缺陷
- An overview of gradient descent algorithms