Python 實現簡單的梯度下降法
機器學習算法常常可以歸結為求解一個最優化問題，而梯度下降法就是求解最優化問題的一個方法。
梯度下降法（gradient descent）或最速下降法（steepest decent），是求解無約束最優化問題的一種最常用的方法。
梯度下降法實現簡單，是一種迭代算法，每一步會求解目標函數的梯度向量。
本文分為理論和 Python 代碼實踐，希望實現簡單的梯度下降法，相關代碼已放在 GitHub 中。
理論
問題定義
那么什么是目標函數，在機器學習中這常常是一個損失函數。不管怎么稱呼，它就是一個函數 $f(x)$，而梯度下降法的目的就是獲取這個函數的極小值。
下面給出一個較為正式的問題定義。
假設 $f(x)$ 是 $R^n$ 上具有一階連續偏導數的函數。需要求解的無約束最優化問題是：
$$\underset{x\in R^n}{min}f(x)$$
即需要求出目標函數 $f(x)$ 的極小點 $x^*$。
算法思想和推導
要理解梯度下降法，首先要理解梯度和負梯度的概念。
梯度是從 n 維推廣出來的概念，類似於斜率。梯度的本意是一個向量，表示某一函數在該點處的方向導數沿着該方向取得最大值，即函數在該點處沿着該方向（此梯度的方向）變化最快，變化率最大（為該梯度的模）。具體定義和公式可以參考百度定義。
舉個例子再體會一下梯度是表示方向的一個向量：
對於函數 $f(x_1,x_2)=2x_1^3-x_2^2$ 來說，它的梯度就是 $g(x_1,x_2)=[6x_1^2,-2x_2]$。對於給定點 $[x_1, x_2]$ 的附近處，它在 $[6x_1^2,-2x_2]$ 方向變化率最大，而其負梯度方向就是 $[-6x_1^2,2x_2]$。例如，在點 $[2, 3]$ 附近處，它的負梯度方向就是 $[-24, -6]$。在此處，點 $[2, 3]$ 向這個方向移動，會使得 $f(x_1,x_2)=2x_1^3-x_2^2$ 值減小的速率最快。反之，如果點 $[2, 3]$ 向梯度方向 $[24, 6]$ 移動，會使得 $f(x_1,x_2)=2x_1^3-x_2^2$ 值增加的速率最快。
理解了梯度之后，其實就可以很容易推導出梯度下降法的算法過程了。
梯度下降法的思想，就是選取適當的初值 $x_{0}$，不斷迭代更新 $x$ 的值，極小化目標函數，最終收斂。
由於負梯度方向是使函數值下降最快的方向，因此梯度下降在每一步采用負梯度方向更新 $x$ 的值，最終達到函數值最小。
可以看出，梯度下降法采用的是貪心的思想。
根據一階泰勒展開，當 $x$ 趨近於 $x_k$ 時：
$$f(x)\approx f(x_k)+g_{k}(x-x_k)$$
這里，$g_k=g(x_k)=\bigtriangledown f(x_k)$ 是 $f(x)$ 在 $x_k$ 的梯度。
我們假設設定了一個初始值 $x_0$，現在需要確定一個 $x_1$，代入上式可得：
$$f(x_1)\approx f(x_0)+g_{0}(x_1-x_0)$$
假設 $x_1$ 和 $x_0$ 之間的距離一定時，為了讓 $f(x_1)$ 最小（貪心策略），應該有：
$$g_{0}(x_1-x_0)=\left | g_{0} \right | \left | x_1-x_0 \right |cos\theta =-\left | g_{0} \right | \left | x_1-x_0 \right |$$
也就是需要讓 $x_1-x_0$ 和梯度 $g_{0}$ 的夾角 $\theta$ 為 180°，使得 $cos\theta =-1$。換言之，$x_1-x_0$ 和梯度 $g_{0}$ 方向相反。
由於 $x_1-x_0=-\frac{g_0}{\left | g_0 \right |}\left | x_1-x_0 \right |$，那么可以得到：
$$x_1=x_0-\frac{g_0}{\left | g_0 \right |}\left | x_1-x_0 \right |=x_0-g_0\lambda_0$$
其中 $\lambda_0=\frac{\left | x_1-x_0 \right |}{\left | g_0 \right |}$ 定義為學習率，它實際上步長除以梯度的模。因此當學習率一定時，步長其實是一直變化的。當梯度較大時，步長也較大；而當梯度較小時，步長也較小。這往往是我們希望的性質，因為當接近於局部最優解時，梯度變得較小，這時往往也需要步長變得更小，以利於找到局部最優解。
同理，我們可以得到 $x_2=x_1-g_1\lambda_1$ ，依次類推，有：
$$x_{k+1}=x_k-g_k\lambda_k$$
其中，學習率 $\lambda_k$ 要足夠小，使得：
- 滿足泰勒公式所需要的精度。
- 能夠很好地捕捉到極小值。
這是一個顯式表達式，可以不斷求出 $x_{k+1}$，當滿足收斂條件時（如梯度足夠小或者 $x_{k+1}$ 更新變化量足夠小），退出迭代，此時 $f(x_{k+1})$ 就是一個求解出來的最小函數值。
至此完成了梯度下降法邏輯上的推導。
Python 代碼實現
理論已經足夠多了，接下來敲一敲實在的代碼吧。
一維問題
假設我們需要求解的目標函數是：
$$f(x)=x^2+1$$
顯然一眼就知道它的最小值是 $x=0$ 處，但是這里我們需要用梯度下降法的 Python 代碼來實現。
1 #!/usr/bin/env python 2 # -*- coding: utf-8 -*- 3 """ 4 一維問題的梯度下降法示例 5 """ 6 7 8 def func_1d(x): 9 """ 10 目標函數 11 :param x: 自變量，標量 12 :return: 因變量，標量 13 """ 14 return x ** 2 + 1 15 16 17 def grad_1d(x): 18 """ 19 目標函數的梯度 20 :param x: 自變量，標量 21 :return: 因變量，標量 22 """ 23 return x * 2 24 25 26 def gradient_descent_1d(grad, cur_x=0.1, learning_rate=0.01, precision=0.0001, max_iters=10000): 27 """ 28 一維問題的梯度下降法 29 :param grad: 目標函數的梯度 30 :param cur_x: 當前 x 值，通過參數可以提供初始值 31 :param learning_rate: 學習率，也相當於設置的步長 32 :param precision: 設置收斂精度 33 :param max_iters: 最大迭代次數 34 :return: 局部最小值 x* 35 """ 36 for i in range(max_iters): 37 grad_cur = grad(cur_x) 38 if abs(grad_cur) < precision: 39 break # 當梯度趨近為 0 時，視為收斂 40 cur_x = cur_x - grad_cur * learning_rate 41 print("第", i, "次迭代：x 值為 ", cur_x) 42 43 print("局部最小值 x =", cur_x) 44 return cur_x 45 46 47 if __name__ == '__main__': 48 gradient_descent_1d(grad_1d, cur_x=10, learning_rate=0.2, precision=0.000001, max_iters=10000)
其輸出結果如下：
第 0 次迭代：x 值為 6.0
第 1 次迭代：x 值為 3.5999999999999996 第 2 次迭代：x 值為 2.1599999999999997 第 3 次迭代：x 值為 1.2959999999999998 第 4 次迭代：x 值為 0.7775999999999998 第 5 次迭代：x 值為 0.46655999999999986 第 6 次迭代：x 值為 0.2799359999999999 第 7 次迭代：x 值為 0.16796159999999993 第 8 次迭代：x 值為 0.10077695999999996 第 9 次迭代：x 值為 0.06046617599999997 第 10 次迭代：x 值為 0.036279705599999976 第 11 次迭代：x 值為 0.021767823359999987 第 12 次迭代：x 值為 0.013060694015999992 第 13 次迭代：x 值為 0.007836416409599995 第 14 次迭代：x 值為 0.004701849845759997 第 15 次迭代：x 值為 0.002821109907455998 第 16 次迭代：x 值為 0.0016926659444735988 第 17 次迭代：x 值為 0.0010155995666841593 第 18 次迭代：x 值為 0.0006093597400104956 第 19 次迭代：x 值為 0.0003656158440062973 第 20 次迭代：x 值為 0.0002193695064037784 第 21 次迭代：x 值為 0.00013162170384226703 第 22 次迭代：x 值為 7.897302230536021e-05 第 23 次迭代：x 值為 4.7383813383216124e-05 第 24 次迭代：x 值為 2.8430288029929674e-05 第 25 次迭代：x 值為 1.7058172817957805e-05 第 26 次迭代：x 值為 1.0234903690774682e-05 第 27 次迭代：x 值為 6.1409422144648085e-06 第 28 次迭代：x 值為 3.684565328678885e-06 第 29 次迭代：x 值為 2.210739197207331e-06 第 30 次迭代：x 值為 1.3264435183243986e-06 第 31 次迭代：x 值為 7.958661109946391e-07 第 32 次迭代：x 值為 4.775196665967835e-07 局部最小值 x = 4.775196665967835e-07
二維問題
接下來推廣到二維，目標函數設為：
$$f(x,y) = -e^{-(x^2 + y^2)}$$
該函數在 $[0, 0]$ 處有最小值。
1 #!/usr/bin/env python 2 # -*- coding: utf-8 -*- 3 """ 4 二維問題的梯度下降法示例 5 """ 6 import math 7 import numpy as np 8 9 10 def func_2d(x): 11 """ 12 目標函數 13 :param x: 自變量，二維向量 14 :return: 因變量，標量 15 """ 16 return - math.exp(-(x[0] ** 2 + x[1] ** 2)) 17 18 19 def grad_2d(x): 20 """ 21 目標函數的梯度 22 :param x: 自變量，二維向量 23 :return: 因變量，二維向量 24 """ 25 deriv0 = 2 * x[0] * math.exp(-(x[0] ** 2 + x[1] ** 2)) 26 deriv1 = 2 * x[1] * math.exp(-(x[0] ** 2 + x[1] ** 2)) 27 return np.array([deriv0, deriv1]) 28 29 30 def gradient_descent_2d(grad, cur_x=np.array([0.1, 0.1]), learning_rate=0.01, precision=0.0001, max_iters=10000): 31 """ 32 二維問題的梯度下降法 33 :param grad: 目標函數的梯度 34 :param cur_x: 當前 x 值，通過參數可以提供初始值 35 :param learning_rate: 學習率，也相當於設置的步長 36 :param precision: 設置收斂精度 37 :param max_iters: 最大迭代次數 38 :return: 局部最小值 x* 39 """ 40 print(f"{cur_x} 作為初始值開始迭代...") 41 for i in range(max_iters): 42 grad_cur = grad(cur_x) 43 if np.linalg.norm(grad_cur, ord=2) < precision: 44 break # 當梯度趨近為 0 時，視為收斂 45 cur_x = cur_x - grad_cur * learning_rate 46 print("第", i, "次迭代：x 值為 ", cur_x) 47 48 print("局部最小值 x =", cur_x) 49 return cur_x 50 51 52 if __name__ == '__main__': 53 gradient_descent_2d(grad_2d, cur_x=np.array([1, -1]), learning_rate=0.2, precision=0.000001, max_iters=10000)
$x_0$ 的初始值設為 $[1,-1]$ ，運行后的結果如下：
[ 1 -1] 作為初始值開始迭代...
第 0 次迭代：x 值為 [ 0.94586589 -0.94586589] 第 1 次迭代：x 值為 [ 0.88265443 -0.88265443] 第 2 次迭代：x 值為 [ 0.80832661 -0.80832661] 第 3 次迭代：x 值為 [ 0.72080448 -0.72080448] 第 4 次迭代：x 值為 [ 0.61880589 -0.61880589] 第 5 次迭代：x 值為 [ 0.50372222 -0.50372222] 第 6 次迭代：x 值為 [ 0.3824228 -0.3824228] 第 7 次迭代：x 值為 [ 0.26824673 -0.26824673] 第 8 次迭代：x 值為 [ 0.17532999 -0.17532999] 第 9 次迭代：x 值為 [ 0.10937992 -0.10937992] 第 10 次迭代：x 值為 [ 0.06666242 -0.06666242] 第 11 次迭代：x 值為 [ 0.04023339 -0.04023339] 第 12 次迭代：x 值為 [ 0.02419205 -0.02419205] 第 13 次迭代：x 值為 [ 0.01452655 -0.01452655] 第 14 次迭代：x 值為 [ 0.00871838 -0.00871838] 第 15 次迭代：x 值為 [ 0.00523156 -0.00523156] 第 16 次迭代：x 值為 [ 0.00313905 -0.00313905] 第 17 次迭代：x 值為 [ 0.00188346 -0.00188346] 第 18 次迭代：x 值為 [ 0.00113008 -0.00113008] 第 19 次迭代：x 值為 [ 0.00067805 -0.00067805] 第 20 次迭代：x 值為 [ 0.00040683 -0.00040683] 第 21 次迭代：x 值為 [ 0.0002441 -0.0002441] 第 22 次迭代：x 值為 [ 0.00014646 -0.00014646] 第 23 次迭代：x 值為 [ 8.78751305e-05 -8.78751305e-05] 第 24 次迭代：x 值為 [ 5.27250788e-05 -5.27250788e-05] 第 25 次迭代：x 值為 [ 3.16350474e-05 -3.16350474e-05] 第 26 次迭代：x 值為 [ 1.89810285e-05 -1.89810285e-05] 第 27 次迭代：x 值為 [ 1.13886171e-05 -1.13886171e-05] 第 28 次迭代：x 值為 [ 6.83317026e-06 -6.83317026e-06] 第 29 次迭代：x 值為 [ 4.09990215e-06 -4.09990215e-06] 第 30 次迭代：x 值為 [ 2.45994129e-06 -2.45994129e-06] 第 31 次迭代：x 值為 [ 1.47596478e-06 -1.47596478e-06] 第 32 次迭代：x 值為 [ 8.85578865e-07 -8.85578865e-07] 第 33 次迭代：x 值為 [ 5.31347319e-07 -5.31347319e-07] 第 34 次迭代：x 值為 [ 3.18808392e-07 -3.18808392e-07] 局部最小值 x = [ 3.18808392e-07 -3.18808392e-07]
我們再試着以初始值 $[3,-3]$ 處開始尋找最小值，即：
gradient_descent_2d(grad_2d, cur_x=np.array([3, -3]), learning_rate=0.2, precision=0.000001, max_iters=10000)
結果可能出乎人意料：
[ 3 -3] 作為初始值開始迭代...
局部最小值 x = [ 3 -3]
梯度下降法沒有找到真正的極小值點！
如果仔細觀察目標函數的圖像，以及梯度下降法的算法原理，你就很容易發現問題所在了。在 $[3, -3]$ 處的梯度就幾乎為 0 了！
print(grad_2d(np.array([3, -3])))
[ 9.13798785e-08 -9.13798785e-08]
由於“梯度過小”，梯度下降法可能無法確定前進的方向了。即使人為增加收斂條件中的精度，也會由於梯度過小，導致迭代中前進的步長距離過小，循環時間過長。
梯度下降法的局限性
梯度下降法實現簡單，原理也易於理解，但它有自身的局限性，因此有了后面很多算法對它的改進。
對於梯度過小的情況，梯度下降法可能難以求解。
此外，梯度下降法適合求解只有一個局部最優解的目標函數，對於存在多個局部最優解的目標函數，一般情況下梯度下降法不保證得到全局最優解（由於凸函數有個性質是只存在一個局部最優解，所有也有文獻的提法是：當目標函數是凸函數時，梯度下降法的解才是全局最優解）。
由於泰勒公式的展開是近似公式，要求迭代步長要足夠小，因此梯度下降法的收斂速度並非很快的。
總結
以上是對用 Python 實現簡單梯度下降法的思考與總結，整個代碼示例參見 GitHub，有何建議和問題請留下您的反饋，謝謝！
原文作者：雨先生
原文鏈接：https://www.cnblogs.com/noluye/p/11108513.html
許可協議：知識共享署名-非商業性使用 4.0 國際許可協議