【Python】淺談梯度下降與實作(中):初階的變形者們

在上一章討論到 BGD 的理論,並花費許多時間實作其功能。
緊接著,我們再來了解其他梯度下降的方法與其優缺點。

========================================

注意事項:
   1. 為了順利顯示數學式,請先安裝 google 套件:
TeX All the Things
   2. 本文中有些許高中數學,不會太硬,請放心咀嚼。

========================================

一、隨機梯度下降(stochastic gradient descent, SGD)

使用 BGD 後,人們發現了許多問題,其中一個就是 — 太慢了!。

當資料量逐漸增多,每次迭代都必須要計算一次所有資料。
若資料從 10 筆增加至 10 萬筆的時候, BGD 的緩慢就會一覽無遺。

先複習一下 BGD 的概念:反復迭代 (1) 式,以得到最佳結果:

$$w_k^{(j+1)}=w_k^{(j)}-\alpha\nabla Loss(w_k) ... (1)$$

$$j 為迭代次數,\alpha 為學習率,k為第k個權重$$

把 (1) 式的梯度部分$\nabla Loss(w_k)$代換一下,
用$W$代表所有權重,$X$表示各權重對應的變量,Loss 使用MSE:

$$W^{(j+1)}=W^{(j)}-\alpha\times\frac{2}{n}\sum_{i=1}^{n}{X^i\left(W^{(j)}X^i-Y_t^i\right)} ... (2)$$

$$i 為資料編號,j 為迭代次數,n為資料筆數,\alpha為學習率$$

看完 BGD 後,讓我們回到 SGD。

SGD 的概念是:「每次隨機取一筆資料來算,反復計算以得到最佳結果。」
即省去 (2) 式中的 $\sum_{i=1}^{n}$ 這個運算,而是改成:

$$W^{(j+1)}=W^{(j)}-\alpha\times 2X^i\left(W^{(j)}X^i-Y_t^i\right) ... (3)$$

$$i 為資料編號,j 為迭代次數,\alpha 為學習率$$

若您不太理解 (3) 式所代表的意義,那麼夏恩舉個例子:
若您有一份含有 10 萬筆的資料,使用 BGD 迭代 10 次計算量為 10萬*10次;使用 SGD 則是 1*10次。

兩者相比,SGD 的確是快多了。

接著來看一下 SGD 在誤差空間中的訓練效能表現。
夏恩先引用上一篇文章中的部分程式碼:

# -*- coding: utf-8 -*-
# 引用部分上一章之程式碼
#
# Shayne, 2019.10.09

import numpy as np
import matplotlib.pyplot as plt

# 給定隨機種子,使每次執行結果保持一致
np.random.seed(1)
    
def getdata(n):
    # n為產生資料量
    x = np.arange(-5, 5.1, 10/(n-1))
    # 給定一個固定的參數,再加上隨機變動值作為雜訊,其變動值介於 +-10 之間
    y = 3*x + 2 + (np.random.rand(len(x))-0.5)*20
    return x, y

def plot_error(x, y):
    a = np.arange(-10, 16, 1)
    b = np.arange(-10, 16, 1)
    mesh = np.meshgrid(a, b)

    sqr_err = 0
    for xs, ys in zip(x, y):
        sqr_err += ((mesh[0]*xs + mesh[1]) - ys) ** 2
    loss = sqr_err/len(x)
    
    plt.contour(mesh[0], mesh[1], loss, 20, cmap=plt.cm.jet)
    plt.xlabel('a')
    plt.ylabel('b')
    plt.axis('scaled')
    plt.title('function loss')

# 實作 BGD 
class my_BGD:    
    def __init__(self, a, b, x, y, alpha):
        self.a = a
        self.b = b
        self.x = x
        self.y = y
        self.alpha = alpha
        
        self.a_old = a
        self.b_old = b
        
        self.loss = None;
    
    # Loss function
    def mse(self):
        sqr_err = ((self.a*self.x + self.b) - self.y) ** 2
        return np.mean(sqr_err)
    
    def gradient(self):
        grad_a = 2 * np.mean((self.a*self.x + self.b - self.y) * (self.x))
        grad_b = 2 * np.mean((self.a*self.x + self.b - self.y) * (1))
        return grad_a, grad_b

    def update(self):
        # 計算梯度
        grad_a, grad_b = self.gradient()
        
        # 梯度更新
        self.a_old = self.a
        self.b_old = self.b
        self.a = self.a - self.alpha * grad_a
        self.b = self.b - self.alpha * grad_b
        self.loss = self.mse();

再來是修改計算梯度的方式,以符合 SGD 的要求,其主要差異為:

1. 把原本 BGD 的 np.mean 拿掉。

2. 新增隨機抽取資料的過程。

完整的 SGD 類別設計,除了 gradient 和 BGD 不同,其他都一樣。
這樣的話,這裡就直接繼承 my_BGD 然後再修改 gradient 函數,實作如下:

# 實作 SGD
class my_SGD(my_BGD):    
    def __init__(self, a, b, x, y, alpha):
        super().__init__(a, b, x, y, alpha)
    
    def gradient(self):
        # 隨機取一筆資料進行計算
        idx = np.random.randint(len(x))
        
        grad_a = 2 * (self.a*self.x[idx] + self.b - self.y[idx]) * (self.x[idx])
        grad_b = 2 * (self.a*self.x[idx] + self.b - self.y[idx]) * (1)
        return grad_a, grad_b

最後補上迭代用的主程式,就能看見 SGD 的最佳化過程:

# 隨機產生一組資料
x, y = getdata(51)

# 繪製誤差空間底圖
plot_error(x, y)

# 初始設定
alpha = 0.05

# 從 -9 開始,能看得更明顯
a = -9; b = -9

# 類別 SGD 初始化
mlclass = my_SGD(a, b, x, y, alpha) # <=== 後續的主程式,只要改這句話就好!

plt.plot(a, b, 'ro-')
plt.title('Initial, loss='+'{:.2f}'.format(mlclass.mse())+'\na='+
          '{:.2f}'.format(a)+', b='+'{:.2f}'.format(b))

# 開始迭代
for i in range(1, 11):
    mlclass.update()
    print('iter='+str(i)+', loss='+'{:.2f}'.format(mlclass.loss))
    
    plt.plot((mlclass.a_old, mlclass.a), (mlclass.b_old, mlclass.b), 'ro-')
    plt.title('iter='+str(i)+', loss='+'{:.2f}'.format(mlclass.loss)+'\na='+
              '{:.2f}'.format(mlclass.a)+', b='+'{:.2f}'.format(mlclass.b))

觀察 SGD 的迭代過程,其優點和缺點一樣明顯:計算速度快,但很不穩定。

備註:本例透過 GIF 圖只能看到「不穩定」,無法看出迭代快慢。
           請讀者自行將數據量增大,就能觀察到迭代速度的差異。

二、小批量梯度下降(mini-batch gradient descent, MBGD)

考慮到只取出一筆資料風險太高;取全部資料又太慢。
折衷的方法就是取出一點點資料,例如一次隨機取出 64 筆。

剛提到的 SGD 的數學式,如 (3) 式:

$$W^{(j+1)}=W^{(j)}-\alpha\times 2X^i\left(W^{(j)}X^i-Y_t^i\right) ... (3)$$

$$i 為資料編號,j 為迭代次數,\alpha為學習率$$

到了 MBGD 我們又要把剛才拿掉的 $\sum_{i=1}^{n}$ 再加回來,但要限制取出資料的數量。
把 (3) 式稍作修改,成為 (4) 式,如下:

$$W^{(j+1)}=W^{(j)}-\alpha\times\frac{2}{m}\sum_{i=1}^{m}{X^i\left(W^{(j)}X^i-Y_t^i\right)} ... (4)$$

$$i 為資料編號,j 為迭代次數,m為資料筆數,\alpha為學習率$$

有別於 BGD 中的資料筆數是指「所有」的資料量,(4) 式中的 m 為我們所指定的資料筆數,例如$m=64$。
MBGD 目前在深度學習中是主流的算法,後續的其他演化多是基於 MBGD 的小批量概念。

這裡也可以來點題外話,在建構深度學習模型時常用的幾個名詞,一起介紹一下:

假設資料有 100 筆,用 MBGD 求最佳解,一次隨機取出 10 筆。
此時,batch size 就等於10;每次取出 10 筆來計算的動作,就是一個 iteration;
經過 10 次計算之後,就可以遍歷所有資料,稱為完成一個 epoch。

另外,MBGD 也有隨機取出資料的概念,因此 SGD 可說是 MBGD 的一個特例。
所以大部分的時候,當我們提到 SGD 時,其實都是指 MBGD。

好的,言歸正傳,讓我們來完成 MBGD 的實作。

因為都是使用 np.mean 來取平均,MBGD 仍然可以繼承之前 BGD 的內容,
但需要修改的地方主要有兩個:

1. 代入 batch_size 參數,限制每次迭代的資料數量,並修改計算損失與梯度的函數。

2. 需要新增更新批次資料的方法(update_batch),每次迭代後,都要換上一批新的資料。

夏恩知道在 tensorflow 裡面,批次更新的這一段的操作十分複雜。
但本恩也沒打算要實作 tensorflow... (可能會寫到斷氣)

總之,意思到了就好。
以下實作每次取出 N 筆資料進行迭代,以觀察迭代效果。

# 請注意:
# 本類別不針對資料跑完後做"任何處理",即跑到最後若超出陣列範圍,loss 會變成 nan。
#
# Lazy Shayne, 2019.10.09

class my_MBGD(my_BGD):
    def __init__(self, a, b, x, y, alpha, batch_size):
        super().__init__(a, b, x, y, alpha)
        
        self.idx = 0
        self.batch_size = batch_size
        
        # 使用 np.random.permutation 給定資料取出順序
        self.suffle_idx = np.random.permutation(len(x))       
        self.update_batch();

    # 更新批次
    def update_batch(self):
        #每次更新時,採滾動的方式依次取出 N 筆資料
        idx = self.suffle_idx[self.idx:self.idx+self.batch_size];
        self.idx += self.batch_size

        self.x_batch = self.x[idx]
        self.y_batch = self.y[idx]
        
    # Loss function
    def mse(self):       
        sqr_err = ((self.a*self.x_batch + self.b) - self.y_batch) ** 2
        return np.mean(sqr_err)
    
    def gradient(self):        
        grad_a = 2 * np.mean((self.a*self.x_batch + self.b - self.y_batch) * (self.x_batch))
        grad_b = 2 * np.mean((self.a*self.x_batch + self.b - self.y_batch) * (1)) 
        self.update_batch();
        return grad_a, grad_b

主程式的部分,僅需抽換部分程式碼,
在這邊修改的程式主要為新增批量大小之變數。

# 參數設定
alpha = 0.05
batch_size = 5

# 類別 MBGD 初始化
mlclass = my_MBGD(a, b, x, y, alpha, batch_size)

改完之後,查看一下迭代結果:
從下圖看來 MBGD 的迭代結果比 SGD 穩定多了。

這邊要特別說明的地方是:
實際上我們在訓練的過程中,跑完所有資料時,並非像是夏恩這樣不負責任讓它變成 NaN!

而是取決於我們給定的 epoch 值,重頭再跑過幾輪;或設定無限迴圈,直到某個停止條件才結束;
甚至還可以設定每次重頭來的時候,都要再次打亂所有資料等,非常多種訓練手法。

詳細運作情形,可以回頭看一下夏恩的 TF 筆記系列,在那邊有詳細的說明。

三、動量(momentum)

延續上一個梯度下降法,MBGD。
如同剛才所講,現在夏恩要改用 SGD 來稱呼 MBGD了
把 SGD 加上一個動量項來抑制雜訊,故此方法又稱為 SGDM (SGD with Momentum)。

不論是 BGD 還是 SGD,在迭代的過程中,都容易收斂在局部最小值。
因此人們在 SGD 中加入動量的概念,以擺脫迭代停滯的問題,如 (5) 式。

$$W^{(j+1)}=W^{(j)}-\alpha\nabla L(W^{(j)})+\gamma(W^{(j)}-W^{(j-1)})... (5)$$

夏恩也看過有人把 (5) 式寫成這樣,移項代換,意思是一樣的:

$$V^{(j)}=\gamma V^{(j-1)}+\alpha\nabla L(W^{(j)})$$

$$W^{(j+1)}=W^{(j)}-V^{(j)}... (5.2)$$

$$其中,\nabla L(W^{(j)})=\frac{2}{m}\sum_{i=1}^{m}{X^i\left(W^{(j)}X^i-Y_t^i\right)}$$

$$i 為資料編號,j 為迭代次數,m為資料筆數,\alpha為學習率,\gamma為動量參數$$

直接看 (5.2) 式可能會略感困惑,但只要寫開來看就明白,以下$G^{(n)}$表迭代時的梯度:

第一次迭代:$W^1=W^0-\alpha G^0$
   對照 (5) 式與 (5.2) 式,$V^0=\alpha G^0$

第二次迭代:$W^2=W^1-\alpha G^1+\gamma(W^1-W^0)$
   代換上式的$(W^1-W^0)$,可得$W^2=W^1-\alpha G^1+\gamma(-\alpha G^0)$

   再對照一次,$V^1=\alpha G^1+\gamma V^0$
   代換$V^0$後,得$V^1=\alpha G^1+\gamma\alpha G^0$

第三次迭代:$W^3=W^2-\alpha G^2+\gamma(W^2-W^1)$
   代換上式的$(W^2-W^1)$,可得$W^3=W^2-\alpha G^2+\gamma(-\alpha G^1+\gamma(-\alpha G^0))$

   最後一次對照,$V^2=\alpha G^2+\gamma V^1$
   代換$V^1$後,得$V^2=\alpha G^2+\gamma\alpha G^1+(\gamma)^2\alpha G^0$

除了煩人的數學式之外,夏恩對於繪圖也是略懂:

參考資料:深度学习优化函数详解(4)-- momentum 动量法

一般的 SGD 走的是紅色路線,初始化$W^0$,然後計算梯度走到$W^1$,再到$W^2$。
圖中的綠線,就是我們剛剛算的動量,所以把動量加上去,就可以得到橘色的路徑。
SGD 和 SGDM 在$W^0$走到$W^1$的過程是一樣的,之後則有動量的差異而產生分歧。

不管數學式子怎麼寫,重點在於「本次」與「上次」的權重值的差異:$W^{(j)}-W^{(j-1)}$。
算出差異值後,再乘上係數$\gamma$,一般是給 0.9,也可以改成 0.85 或 0.95,看個人需求。

但請不要讓$\gamma$大於1,$\gamma$應介於$0~1$之間,是一個「衰退係數」。
從 (5.2) 式中,可看到每次的迭代都會不斷地乘上$\gamma$,
當$\gamma$介於$0~1$之間時,表示距離愈遙遠的迭代結果,影響性就會愈低。

以下則是夏恩對於動量這個名詞的初步理解:

權重值的變化,可以理解為權重的「位移量」,位移時變率即為速度。
係數$\gamma$所扮演的角色為「質量」,$「質量」\times「速度」$等於「動量」。
當質量愈大或速度愈大,動量也就愈大,而給定動量的影響在於:

當前後更新的權重同向,則會加速更新步伐,差異愈大則步伐愈大;
當前後更新的權重反向,則會抑制更新步伐,差異愈大則步伐愈小。

講這麼多,終於到了實作的時候了。
好消息是 SGDM 的實作遠比理解數學式容易...

首先繼承剛才的 MBGD 類別,
主要是需新增一個動量項,及覆寫梯度下降 update 函數,程式如下:

class my_SGDM(my_MBGD):
    def __init__(self, a, b, x, y, alpha, batch_size, gamma):
        super().__init__(a, b, x, y, alpha, batch_size)
        
        # Momentum parameter
        self.gamma = gamma
        
        # 動量累加項
        self.sum_Ma = 0
        self.sum_Mb = 0
          
    def update(self):
        self.a_old = self.a 
        self.b_old = self.b
        
        # 計算梯度
        grad_a, grad_b = self.gradient()

        # 動量累加
        self.sum_Ma = self.gamma * self.sum_Ma + self.alpha * grad_a
        self.sum_Mb = self.gamma * self.sum_Mb + self.alpha * grad_b
        
        # 梯度更新
        self.a -= self.sum_Ma
        self.b -= self.sum_Mb 

        self.loss = self.mse();

在主程式的部分,則是先降低學習率,並多給一個動量參數 M:

# 參數設定
alpha = 0.01
batch_size = 5
gamma = 0.8

# 類別 SGDM 初始化
mlclass = my_SGDM(a, b, x, y, alpha, batch_size, gamma)

從上圖中可以看到在梯度下降時,有一個加速 >> 減速 >> 轉向 >> 再加速的過程,
經過愈陡峭的地方,就會使得下次的更新步伐大增,反之則是給予抑制。

看起來是不是很有動量的感覺呢!

四、NAG

NAG 是 Nesterov accelerated gradient 的縮寫。

針對既有的動量法 SGDM,數學家 Nesterov 再提出新的改進方式。
其核心想法為,每次梯度更新時,「再多往前走一些!」。

延續上述的 (5.2) 式,在其基礎上新增一點東西為 (6) 式,如下:

$$V^{(j)}=\gamma V^{(j-1)}+\alpha\nabla L(W^{(j)}-\gamma V^{(j-1)})$$

$$W^{(j+1)}=W^{(j)}-V^{(j)}... (6)$$

原本是計算$W^{(j)}$的梯度,現在改成算$W^{(j)}-\gamma V^{(j-1)}$。
NAG 大部分的東西都和 SGDM 相同,我看這次就不要再推公式了,因此直接畫張圖吧!

參考資料:深度學習優化函數詳解(5)-- Nesterov accelerated gradient (NAG)

NAG 和 SGDM 在$W^0$走到$W^1$的過程也是一樣的。

然後在計算$W^2$時,不去計算$W^1$的梯度,而是先算$W^1$加上$\gamma(W^1-W^0)$的值,
把$W^1$挪移至新的$(W^1)'$,此時才計算$(W^1)'$的梯度值。
最後才是將$(W^1)'$的梯度值加上$\gamma(W^1-W^0)$,得到上圖中$W^2_{NAG}$的位置。

「多往前走一些」就是 NAG 最大的特色。

把數學式拿出來仔細看看,$\nabla L(W^{(j)}-\gamma V^{(j-1)})$。
其實就是多考慮到二階導數 — 梯度的變化率,也就是俗稱的「加速度」。
因此,NAG 在迭代的過程中會比單純的動量更能提早反應變化與轉向,進而提升收斂的速度。

夏恩的 NAG 實作如下,先繼承 SGDM 的類別,然後僅需覆寫計算梯度的函數即可:

class my_NAG(my_SGDM):
    def __init__(self, a, b, x, y, alpha, batch_size, gamma):
        super().__init__(a, b, x, y, alpha, batch_size, gamma)
    
    def gradient(self):
        
        # 往前多走一點
        a = self.a - self.sum_Ma
        b = self.b - self.sum_Mb
            
        grad_a = 2 * np.mean((a*self.x_batch + b - self.y_batch) * (self.x_batch))
        grad_b = 2 * np.mean((a*self.x_batch + b - self.y_batch) * (1)) 
        self.update_batch();
        return grad_a, grad_b

主程式也和之前相同,更換宣告的類別即可:

# 參數設定
alpha = 0.01
batch_size = 5
gamma = 0.8

# 類別 NAG 初始化
mlclass = my_NAG(a, b, x, y, alpha, batch_size, gamma)

以上四種初階的變形,主要是為了解決 BGD 的緩慢與梯度消失的情況。
但問題解決之後,緊跟著又是另外一個問題 — 學習率太難設定了!

權重們各自對應的項次不同,卻都使用相同的學習率。
最明顯的影響就是在迭代時,高次項權重的更新步伐大,低次項卻緩慢的蠕動,整個收斂過程沒有效率。

此外,不同的資料與深度學習模型的特性差異甚大,
工程師們總是花費大量的時間,只為了找到合適的學習率。

為了讓模型能在訓練過程中可以自行調整學習率,
之後又演化出多種自適應(Adaptive)的梯度下降法。

有關自適應梯度下降法的介紹,請接著下一篇:

【Python】淺談梯度下降與實作(下):猙獰的變形者們