作者：SETHNEHA

翻译：王可汗

校对：陈丹

本文约3200字，建议阅读15分钟

本文为大家介绍多项式回归中的梯度下降算法。

关键词：梯度下降、多项式回归、模型优化

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

最好的学习方法是实践，因此在本文中，我将一步步介绍梯度下降过程是如何工作的，而不使用像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]*xc = coeffs[2]y = a+b+creturn y

当x=3时，我们可以看看它的作用：

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

创建数据和基础模型

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

import numpy as np
import matplotlib.pyplot as plthundred_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]*xc = coeffs[2]y = a+b+cprint(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 maxwhile 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

测试一下：

7
Should get value in the range 3 - 11
6.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 labelsy_bar: prediction labelsReturns: Calculated MSE loss."""return sum((ys - y_bar)*(ys - y_bar)) /len(ys)initial_model_loss = loss_mse(ys, y_bar)initial_model_loss47922.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 datapointsoutputs_y: actual y output pointsReturns: 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**2partial_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.0001
a_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.9704266522693037
Now 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 datapointsoutputs_y: actual y output pointslr: learning rateReturns: 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_ab_new = coeffs[1] - lr * gradient_bc_new = coeffs[2] - lr * gradient_cnew_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 rateReturns: Tuple containing (updated_model_loss,new_model_coeffs, new_y_bar predictions, saved loss updates)"""losses = []rand_coeffs_to_test = rand_coeffsfor 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

编辑：黄继彦

校对：杨学俊

译者简介

王可汗，清华大学机械工程系直博生在读。曾经有着物理专业的知识背景，研究生期间对数据科学产生浓厚兴趣，对机器学习AI充满好奇。期待着在科研道路上，人工智能与机械工程、计算物理碰撞出别样的火花。希望结交朋友分享更多数据科学的故事，用数据科学的思维看待世界。

翻译组招募信息

工作内容：需要一颗细致的心，将选取好的外文文章翻译成流畅的中文。如果你是数据科学/统计学/计算机类的留学生，或在海外从事相关工作，或对自己外语水平有信心的朋友欢迎加入翻译小组。

你能得到：定期的翻译培训提高志愿者的翻译水平，提高对于数据科学前沿的认知，海外的朋友可以和国内技术应用发展保持联系，THU数据派产学研的背景为志愿者带来好的发展机遇。

其他福利：来自于名企的数据科学工作者，北大清华以及海外等名校学生他们都将成为你在翻译小组的伙伴。

点击文末“阅读原文”加入数据派团队~

转载须知

如需转载，请在开篇显著位置注明作者和出处（转自：数据派ID：DatapiTHU），并在文章结尾放置数据派醒目二维码。有原创标识文章，请发送【文章名称-待授权公众号名称及ID】至联系邮箱，申请白名单授权并按要求编辑。

发布后请将链接反馈至联系邮箱（见下方）。未经许可的转载以及改编者，我们将依法追究其法律责任。

点击“阅读原文”拥抱组织