机器学习:从零开始学习梯度下降

机器学习算法与Python实战

共 10646字,需浏览 22分钟

 ·

2021-03-14 14:08

作者:SETHNEHA 翻译:王可汗 校对:陈丹



梯度下降是一个需要理解的重要算法,因为它是机器学习和深度学习中使用的许多更先进算法的基础。因此,掌握梯度下降的内部工作原理对任何计划进一步探索机器学习算法的人来说都是非常有益的。

最好的学习方法是实践,因此在本文中,我将一步步介绍梯度下降过程是如何工作的,而不使用像scikit-learn这样的ML库。在日常工作中,使用这些库当然会更快、更简洁,但在学习过程中,我发现手工实现对这个特定算法来说是非常宝贵的。

梯度下降的目标

梯度下降的目标是最小化模型预测与原始数据间的误差。在本文的背景下,我们将着眼于二次多项式模型,也称为二次方程:


将二阶多项式画出来后, 看起来是这样的:



多项式回归

我们这里专门观察多项式回归, 即自变量x与因变量y之间的关系被建模为x的n阶多项式。简单地说,我们的二次多项式的系数a、b、c将被估计、评价和修改,直到我们可以将线准确地拟合到输入的x数据上。梯度下降是这个过程中的优化步骤,它修改和优化这些系数的值。
 
现在,我们将看看如何创建和绘制这样的曲线,并建立一个初始模型来拟合这个数据,然后我们将使用梯度下降来优化和改进它。如果我们能得到一个可以准确描述数据的模型,希望它能够准确预测另一组x值的y值。
 
我们可以开始为二次多项式方程(𝑎𝑥²+𝑏𝑥+𝑐)选择系数,并应用到我们将尝试建模的数据:

coeffs = [2,-5, 4]

这将是我们希望我们的预测模型尽可能接近的真值模型的系数[5] 。接下来,我们需要一个二次多项式的评估函数,在给定一组系数和给定输入𝑥的情况下,返回相应的𝑦。

def eval_2nd_degree(coeffs, x):    """    Function to return the outputof evaluating a second degree polynomial,    given a specific x value.
Args: coeffs: List containingthe coefficients a,b, and c for the polynomial. x: The input x value tothe polynomial.
Returns: y: The correspondingoutput y value for the second degree polynomial.
""" a = (coeffs[0]*(x*x)) b = coeffs[1]*x c = coeffs[2] y = a+b+c return y

当x=3时,我们可以看看它的作用:

coeffs = [2, -5, 4]x=3eval_2nd_degree(coeffs, x)
7

创建数据和基础模型

定义一些x数据(输入),我们希望预测y(输出):

import numpy as npimport matplotlib.pyplot as plt
hundred_xs=np.random.uniform(-10,10,100)print(hundred_xs)
x_y_pairs = []for x in hundred_xs: y =eval_2nd_degree(coeffs, x) x_y_pairs.append((x,y))
xs = []ys = []for a,b in x_y_pairs: xs.append(a) ys.append(b)
plt.figure(figsize=(20,10))plt.plot(xs, ys, 'g+')plt.title('Original data')plt.show()


这很好,但是我们可以通过让事情变得更真实来改进它。你可以添加噪音或“抖动”的值,使他们可以类似于现实世界的数据:

defeval_2nd_degree_jitter(coeffs, x, j):    """    Function to return the noisy output ofevaluating a second degree polynomial,    given a specific x value. Output values canbe within [y−j,y+j].
Args: coeffs: List containing thecoefficients a,b, and c for the polynomial. x: The input x value to the polynomial. j: Jitter parameter, to introduce noiseto output y.
Returns: y: The corresponding jittered output yvalue for the second degree polynomial.
""" a = (coeffs[0]*(x*x)) b = coeffs[1]*x c = coeffs[2] y = a+b+c print(y)
interval = [y-j, y+j] interval_min = interval[0] interval_max = interval[1] print(f"Should get value in the range{interval_min} - {interval_max}") jit_val = random.random() *interval_max # Generate a randomnumber in range 0 to interval max
while interval_min > jit_val: # While the random jittervalue is less than the interval min, jit_val = random.random() *interval_max # it is not in the rightrange. Re-roll the generator until it # give a number greater than the interval min.
return jit_val

测试一下:

7Should get value in the range 3 - 116.233537936801398

这个更新的函数将接受二阶多项式的输入和为这个输入添加噪声的抖动值,以给我们一个更现实的输出,而不仅仅是一个完美的曲线:

x_y_pairs = []for x in hundred_xs:    y  =eval_2nd_degree_jitter(coeffs, x, j)    x_y_pairs.append((x,y))
xs = []ys = []for a,b in x_y_pairs: xs.append(a) ys.append(b)
plt.figure(figsize=(20,10))plt.plot(xs, ys, 'g+')plt.title('Original data')plt.show()


当我们建立我们的预测模型,并通过梯度下降来优化它时,我们希望能够得到尽可能接近这些值的结果。
 
建模第一步:尝试一个随机模型

建模的第一步是为二次多项式(𝑦=𝑎𝑥²+𝑏𝑥+𝑐)生成和存储随机系数。这将是我们的初始模型,它很可能不那么精确,我们的目标是改进它,直到它与数据足够吻合。

rand_coeffs=(random.randrange(-10,10),random.randrange(-10,10),random.randrange(-10,10))rand_coeffs
(7, 6, 3)

通过计算输入值的预测输出值来检查这个模型的准确性:

y_bar =eval_2nd_degree(rand_coeffs, hundred_xs)
plt.figure(figsize=(20,10))plt.plot(xs, ys, 'g+', label ='original')plt.plot(xs, y_bar, 'ro',label='prediction')plt.title('Original data vsfirst prediction')plt.legend(loc="lowerright")plt.show()


从上面的图中可以明显看出,这个带有随机系数的新模型并不完全符合我们的数据。为了定量描述这个模型有多不正确,我们计算了模型的均方误差。这是实际输出和预测输出之差的平方和的平均值:


def loss_mse(ys, y_bar):    """    Calculates MSE loss.
Args: ys: training data labels y_bar: prediction labels
Returns: Calculated MSE loss. """ return sum((ys - y_bar)*(ys - y_bar)) /len(ys)
initial_model_loss = loss_mse(ys, y_bar)
initial_model_loss
47922.39790821987

相当大的数。现在让我们看看我们能否通过梯度下降优化模型来改善这个相当高的损失值。
 
梯度下降和损失减少

我们希望改进我们的模型。因此,我们想要改变它的系数a, b和c,以减少误差。因此我们需要知道每个系数是如何影响误差的。这是通过计算损失函数对每个单独系数的偏导数来实现的。

在这个案例中,我们使用MSE作为我们的损失函数——这即是我们希望计算偏导数的函数:

 
我们模型的输出预测如下:


因此,损失可以重新表述为:

 
在这个特定的例子中,我们损失函数的偏导数如下:


  • 如果你计算每个导数的值,你会得到每个系数的梯度。
  • 这些值给出了损失函数相对于每个特定系数的斜率。
  • 它们表明你应该增加还是减少它来减少损失,以及这样做的安全程度。
 
给定系数𝑎,𝑏和𝑐,计算的梯度𝑔𝑎,𝑔𝑏和𝑔𝑐和学习率𝑙𝑟,通常会更新系数,更新的值定义如下:


一旦您将新模型应用到数据中,您的损失应该会减少。
 
减少损失

我们需要一个梯度计算函数,在给定一个二次多项式的系数,以及一组输入𝑥和一组相应的实际输出𝑦,将返回每个系数的梯度。

defcalc_gradient_2nd_poly(rand_coeffs, hundred_xs, ys):    """    calculates the gradient for a second degreepolynomial.
Args: coeffs: a,b and c, for a 2nd degreepolynomial [ y = ax^2 + bx + c ] inputs_x: x input datapoints outputs_y: actual y output points
Returns: Calculated gradients for the 2nddegree polynomial, as a tuple of its parts for a,b,c respectively.
"""
a_s = [] b_s = [] c_s = []

y_bars = eval_2nd_degree(rand_coeffs,hundred_xs)
for x,y,y_bar in list(zip(hundred_xs, ys,y_bars)): # take tuple of (xdatapoint, actual y label, predicted y label) x_squared = x**2 partial_a = x_squared * (y - y_bar) a_s.append(partial_a) partial_b = x * (y-y_bar) b_s.append(partial_b) partial_c = (y-y_bar) c_s.append(partial_c)
num = [i for i in y_bars] n = len(num)
gradient_a = (-2 / n) * sum(a_s) gradient_b = (-2 / n) * sum(b_s) gradient_c = (-2 / n) * sum(c_s) return(gradient_a, gradient_b,gradient_c) # return calculatedgradients as a a tuple of its 3 parts

我们现在要:

  • 使用上面的函数来计算我们表现不佳的随机模型的梯度。
  • 相应调整模型系数。
  • 验证模型的损失现在更小了——梯度下降起作用了!
 
让我们设定一个实验的初始学习速率。这应该保持在很小的范围内,以避免错过全局最小值,但也不能小到要花费无穷长的[8] 时间或陷入局部最小值。lr = 0.0001是一个很好的起点。

calc_grad= calc_gradient_2nd_poly(rand_coeffs, hundred_xs, ys)
lr =0.0001a_new= rand_coeffs[0] - lr * calc_grad[0]b_new= rand_coeffs[1] - lr * calc_grad[1]c_new= rand_coeffs[2] - lr * calc_grad[2]
new_model_coeffs= (a_new, b_new, c_new)print(f"Newmodel coeffs: {new_model_coeffs}")print("")
#updatewith these new coeffs:new_y_bar= eval_2nd_degree(new_model_coeffs, hundred_xs)updated_model_loss= loss_mse(ys, new_y_bar)
print(f"Nowhave smaller model loss: {updated_model_loss} vs {original_model_loss}")
New model coeffs: 5.290395171471687 5.903335222089396 2.9704266522693037Now have smaller model loss: 23402.14716735533 vs 47922.39790821987

通过将训练数据、原始随机模型和更新后的低损失模型一起绘制出来,来可视化这个改进:

plt.figure(figsize=(20,10))plt.plot(xs, ys, 'g+', label ='original model')plt.plot(xs, y_bar, 'ro', label= 'first prediction')plt.plot(xs, new_y_bar, 'b.',label = 'updated prediction')plt.title('Original model vs1st prediction vs updated prediction with lower loss')plt.legend(loc="lower right")plt.show()


多轮迭代梯度下降

我们几乎准备好了。最后一步是在多个轮(周期或迭代)中迭代地执行梯度下降。我们希望每一轮迭代都能看到在降低损失和更好地对原始数据进行模型拟合方面的改进。
 
让我们对上面的calc_gradient_2nd_poly函数进行改进,使其更适用于梯度下降迭代过程:

defcalc_gradient_2nd_poly_for_GD(coeffs, inputs_x, outputs_y, lr):    """    calculates the gradient for a second degreepolynomial.
Args: coeffs: a,b and c, for a 2nd degreepolynomial [ y = ax^2 + bx + c ] inputs_x: x input datapoints outputs_y: actual y output points lr: learning rate
Returns: Calculated gradients for the 2nddegree polynomial, as a tuple of its parts for a,b,c respectively.
""" a_s = [] b_s = [] c_s = []
y_bars = eval_2nd_degree(coeffs, inputs_x)
for x,y,y_bar in list(zip(inputs_x,outputs_y, y_bars)): # take tuple of(x datapoint, actual y label, predicted y label) x_squared = x**2 partial_a = x_squared * (y - y_bar) a_s.append(partial_a) partial_b = x * (y-y_bar) b_s.append(partial_b) partial_c = (y-y_bar) c_s.append(partial_c)
num = [i for i in y_bars] n = len(num)
gradient_a = (-2 / n) * sum(a_s) gradient_b = (-2 / n) * sum(b_s) gradient_c = (-2 / n) * sum(c_s)

a_new = coeffs[0] - lr * gradient_a b_new = coeffs[1] - lr * gradient_b c_new = coeffs[2] - lr * gradient_c
new_model_coeffs = (a_new, b_new, c_new)
#update with these new coeffs: new_y_bar = eval_2nd_degree(new_model_coeffs,inputs_x)
updated_model_loss = loss_mse(outputs_y,new_y_bar) return updated_model_loss,new_model_coeffs, new_y_bar

这将作为gradient_descent函数的一部分被调用:

def gradient_descent(epochs,lr):    """    Perform gradient descent for a seconddegree polynomial.
Args: epochs: number of iterations to performof finding new coefficients and updatingt loss. lr: specified learning rate
Returns: Tuple containing (updated_model_loss,new_model_coeffs, new_y_bar predictions, saved loss updates)
""" losses = [] rand_coeffs_to_test = rand_coeffs for i in range(epochs): loss =calc_gradient_2nd_poly_for_GD(rand_coeffs_to_test, hundred_xs, ys, lr) rand_coeffs_to_test = loss[1] losses.append(loss[0]) print(losses) return loss[0], loss[1], loss[2],losses #(updated_model_loss,new_model_coeffs, new_y_bar, saved loss updates)

最后,让我们训练1500轮,看看我们的模型是否学到了什么:

GD = gradient_descent(1500,0.0001)
plt.figure(figsize=(20,10))plt.plot(xs, ys, 'g+', label ='original')plt.plot(xs, GD[2], 'b.', label= 'final_prediction')plt.title('Original vs Finalprediction after Gradient Descent')plt.legend(loc="lowerright")plt.show()

 
经过完整的训练后,这个训练过的模型显示出了巨大的改进。我们可以通过检查它的最终预测系数a、b和c来进一步检验:

print(f"FinalCoefficients predicted: {GD[1]}")print(f"OriginalCoefficients: {coeffs}")
Final Coefficients predicted: (2.0133237089326155, -4.9936501002139275, 3.1596042252126195)Original Coefficients: [2, -5, 4]

不太远!这是对初始随机模型的一大改进。通过观察训练时减少损失的情况,我们可以获得更深入的见解:

plt.figure(figsize=(20,10))plt.plot(GD[3], 'b-', label ='loss')plt.title('Loss over 1500iterations')plt.legend(loc="lowerright")plt.xlabel('Iterations')plt.ylabel('MSE')plt.show()


我们观察到为了得到更精确的系数,模型损耗接近于零。我们还可以看到,在大约400轮之后,损失并没有明显的改善——绝对不需要1500轮[11] 。另一种策略是在训练步骤中添加某种条件,当达到某个最小损失阈值时停止训练。这将防止过度训练和潜在的模型过拟合。
 
我希望你喜欢这篇关于多项式回归的梯度下降的文章。某些概念乍一看可能会让人望而生畏,但随着时间的推移,如果我们坚持足够长的时间,我们会逐渐熟悉一个问题的“具体细节”。我发现这个练习对我来说确实是这样的,我觉得这是一个很值得学习的经验。

原文标题:
PolynomialRegression — Gradient Descent from Scratch
原文链接:
https://towardsdatascience.com/polynomial-regression-gradient-descent-from-scratch-279db2936fe9?source=collection_home

也可以加一下老胡的微信
围观朋友圈~~~


推荐阅读

(点击标题可跳转阅读)

麻省理工学院计算机课程【中文版】
【清华大学王东老师】现代机器学习技术导论.pdf
机器学习中令你事半功倍的pipeline处理机制
机器学习避坑指南:训练集/测试集分布一致性检查
机器学习深度研究:特征选择中几个重要的统计学概念

老铁,三连支持一下,好吗?↓↓↓

浏览 65
点赞
评论
收藏
分享

手机扫一扫分享

分享
举报
评论
图片
表情
推荐
点赞
评论
收藏
分享

手机扫一扫分享

分享
举报