深度學習優化算法（牛頓法-->梯度下降法-->Nadam）

目錄

• 一、牛頓法與擬牛頓法
  • 1、牛頓法
    • 1.1 原始牛頓法（假設f凸函數且兩階連續可導，Hessian矩陣非奇異）
      • 算法1.1 牛頓法
    • 1.2 阻尼牛頓法
      • 算法1.2 阻尼牛頓法
  • 2、擬牛頓法(如何不用二階導數構造海森矩陣)
    • 2.1 擬牛頓條件（擬牛頓方程，割線條件）
  • 3、擬牛頓法之DFP算法
    • 算法2.1 DFP算法
  • 4、擬牛頓法之BFGS算法
    • 算法2.2 BFGS算法(I)
    • 算法2.2 BFGS算法(II)
  • 5、L-BFGS算法
    • 算法2.4 （$D_k·g_k$的快速算法）
• 二、梯度下降法
  • 2.1 梯度
  • 2.2 梯度下降與梯度上升
  • 2.3 梯度下降法算法詳解
    • 2.3.1 梯度下降法的相關概念
    • 2.3.2 梯度下降法的詳細算法（代數法和矩陣法）
      • 梯度下降法的代數描述
      • 梯度下降法的矩陣方式描述
    • 3.3 梯度下降的算法調優
  • 4.梯度下降法和其他無約束優化算法的比較
• 三、梯度下降法的改進
  • 3.1 優化算法通用框架
  • 3.2 SGD with Momentum
  • 3.3 SGD with Newterov Acceleration
  • 3.3 AdaGrad
  • 3.4 AdaDelta/RMSProp
  • 3.5 Adam
  • 3.6 Nadam
• 四、如何選擇優化算法
  • 4.1 Adam存在的幾個問題
    • 4.1.1 可能不收斂
    • 4.1.2 可能錯過全局最優解
    • 4.1.3 小結
  • 4.2 Adam+SGD
    • 4.2.1 切換算法后用什么樣的學習率
    • 4.2.2 什么時候切算法
  • 4.3 優化算法常用tricks
  • 4.4 為什么深度學習不使用牛頓法或者擬牛頓法

一、牛頓法與擬牛頓法

擬牛頓法(Quasi-Newton Methods)是求解非線性優化問題最有效的方法之一，於20世紀50年代提出。DFP、BFGS和L-BFGS算法都是重要的擬牛頓法。考慮如下無約束的極小化問題\(\underset{x} {min} f(x)\)，其中\({\tt x}=(x_1,x_2,...,x_N)^T \in R^N\)這里假設\(f\)為凸函數，且兩階連續可微，此外記極小問題的解為\(x^*\)

1、牛頓法

1.1 原始牛頓法（假設f凸函數且兩階連續可導，Hessian矩陣非奇異）

為簡單起見，首先考慮\(N=1\)的簡單清醒。牛頓法的基本思想是：在現有極小點估計值的附近對\(f(x)\)做二階泰勒展開，進而找到極小點的下一個估計值，設\(x_k\)為當前的極小點的估計值，則

\(\varphi(x)=f(x_k) + f'(x_k)(x-x_k)+\frac{1}{2}f''(x_k)(x-x_k)^2\)

表示\(f(x)\)在\(x_k\)附近的二階泰勒展開式（略去了關於\(x-x_k\)的高階項），由於求的是最值，由極值必要條件可值

\(\varphi'(x)=0\)， 即 \(f'(x)+f''(x)(x-x_k)=0\)

從而求得\(x = x_k - \frac{f'(x_k)}{f''(x_k)}\)

如是若給定初始值\(x_0\)則可以構造如下的迭代格式\(x_{k+1} = x_k - \frac{f'(x_k)}{f''(x_k)}\)

產生序列\(\{x_k\}\)來逼近\(f(x)\)的極小點。在一定條件下，\(\{x_k\}\)可以收斂到\(f(x)\)的極小點

對於\(N>1\)的情形，二階泰勒展開式可以做推廣，此時

\(\varphi({\tt x})=f({\tt x}_k) +\nabla f({\tt x}_k)({\tt x}-{\tt x}_k)+\frac{1}{2}({\tt x}-{\tt x}_k)^T\nabla^2 f({\tt x}_k)({\tt x}-{\tt x}_k)\)

其中\(\nabla f\)為\(f\)的梯度向量，\(\nabla^2f\)為\(f\)的海森矩陣(Hessian matrix)，其定義分別為：

\(\nabla f=\begin{bmatrix} \partial f/ \partial x_1 \\ \partial f/ \partial x_2 \\ \vdots \\\partial f/ \partial x_N \end{bmatrix}\), \(\nabla^2 f=\begin{bmatrix} \partial^2 f/ \partial x_1^2 & \partial^2 f/ \partial x_1 \partial x_2 & \cdots & \partial^2 f/ \partial x_1 \partial x_N\\ \partial^2 f/ \partial x_2 \partial x_1 & \partial^2 f/ \partial x_2^2 & \cdots & \partial^2 f/ \partial x_2 \partial x_N \\ \vdots & \vdots & \cdots & \vdots \\\partial^2 f/ \partial x_N \partial x_1 & \partial^2 f/ \partial x_N \partial x_2 & \cdots & \partial^2 f/ \partial x_N^2 \end{bmatrix}\)

以下分別將其簡記為\(g\) 和\(H\)，特別的，若\(f\)的混合偏導數可交換次序，即對\(\forall i,j, \partial^2 f/\partial x_i \partial x_j = \partial^2 f/\partial x_j \partial x_i\)，則海森矩陣\(H\)為對稱矩陣。其中\(\nabla f(x_k)\)和\(\nabla^2f(x_k)\)表示\(\tt x\)去值為\(x_k\)后得到的實值向量和矩陣，分別記為\(g_k\)和\(H_k\)

​ 同樣地，由於是求極小點，極值必要條件要求它為\(\varphi(x)\)的駐點，即\(\nabla\varphi(x)=0\)，亦即

\(g_k + H_k({\tt x}-{\tt x}_k)=0\)

進一步，若矩陣\(H_k\)非奇異，則可解得\({\tt x}={\tt x}_k - H_k^{-1}·g_k\)，於是給定初始值\(\tt x_0\)，同樣可以構造出迭代格式：\({\tt x}_{k+1}={\tt x}_k - H_k^-·g_k\)

這就是原始牛頓迭代法，其迭代格式中的搜索方向\(d_k=-H_k^{-1}·g_k\)稱為牛頓方向

算法1.1 牛頓法

1. 給定初值\({\tt x}_0\)和精度閾值\(\epsilon\)，並令\(k:=0\)
2. 計算\(g_k\)和\(H_k\)
3. 若\(||g_k||<\epsilon\)，則停止迭代；否則確定搜索方向\(d_k=-H_k^{-1}·g_k\)
4. 計算新的迭代點\({\tt x}_{k+1} := {\tt x}_k + d_k\)
5. 令\(k:=k+1\)，轉至步驟2

當目標函數是二次函數時，由於二階泰勒展開函數與原目標函數不是接近而是完全相同的二次式，海森矩陣退化成一個常數矩陣，從任一初始值出發，只需一步迭代就可以達到\(f(x)\)的極小點\(x^*\)，因此牛頓法是一種具有二次收斂性的算法。至於非二次函數，若函數的二次型較強，或迭代點已進入極小點的鄰域，其手鏈數獨也是很快的，這是牛頓法的主要優點。

但原始牛頓法由於迭代公式中沒有步長因子，而是定長補償迭代，對於非二次型目標函數有時會使函數值上升，即出現\(f({\tt x}_{k+1}) > f({\tt x}_k)\)的情況，這表明原始牛頓法不能保證函數值穩定的下降，在嚴重的情況下甚至可能造成迭代點\(\{{\tt x}_k\}\)的發散而導致計算失敗

1.2 阻尼牛頓法

為了消除原始牛頓法不能保證函數值穩定的下降的弊病，人們提出了阻尼牛頓法，每次迭代的方向仍采用\(d_k\)，但每次迭代需沿此方向作一維搜索(line search)，尋求最優的步長因子\(\lambda_k\)即

\(\lambda_k = \underset{\lambda \in R}{argmin}f({\tt x}_k + \lambda d_k)\)

算法1.2 阻尼牛頓法

1. 給定初值\({\tt x}_0\)和精度閾值\(\epsilon\)，並令\(k:=0\)
2. 計算\(g_k\)和\(H_k\)
3. 若\(||g_k||<\epsilon\)，則停止迭代；否則確定搜索方向\(d_k=-H_k^{-1}·g_k\)
4. 利用\(\lambda_k = \underset{\lambda \in R}{argmin}f({\tt x}_k + \lambda d_k)\)得到步長\(\lambda_k\)，並令\({\tt x}_{k+1} := {\tt x}_k + \lambda_kd_k\)
5. 令\(k:=k+1\)，轉至步驟2

注1.1: 算法中涉及到\(H_k^{-1}\)的計算，實際應用中通常不直接對\(H_k\)求逆，而是將其轉化為求解線性代數方程組\(H_kd_k=-g_k\)，此時可根據稀疏矩陣\(H_k\)的形態來選擇合適的迭代法，如預條件共軛梯度法(PCG)、代數多重網格法(AMG)等

牛頓法是梯度下降法的進一步發展，梯度法利用目標函數的一階偏導信息，以負梯度方向作為搜索方向，只考慮目標函數在迭代點的局部性質，而牛頓法不僅使用目標函數的一階偏導數，還進一步利用了目標函數的二階偏導數，這樣就考慮了梯度變化的趨勢，因而能更全面的確定合適的搜索方向，以加快收斂，它具有二階收斂數獨。但牛頓法主要存在以下兩個缺點：

1）對目標函數有較嚴格的要求，函數必須具有連續一、二階偏導數，海森矩陣必須正定

2）計算相當復雜，除需要計算梯度外，還需要計算二階偏導數矩陣和它的逆矩陣，計算量和存儲量均很大，且均以維度\(N\)的平方比增加，當N很大時這個問題更加突出。

2、擬牛頓法(如何不用二階導數構造海森矩陣)

如上所述，牛頓法雖然收斂速度快，但是計算過程中需要計算目標函數的二階偏導數，計算量較大。而且有時目標函數的海森矩陣無法保持正定，從而使得牛頓法失效。為了客服這兩個問題，人們提出了擬牛頓法。這個方法的基本思想是：不用二階偏導數二構造出近似海森矩陣的正定矩陣，在“擬牛頓”的條件下優化目標函數。不同的構造方法就產生了不同的擬牛頓法

也有人把擬牛頓法翻譯成准牛頓法，其實是表示類似於牛頓法的意思，因為只是對算法中用來計算搜索方向的海森矩陣作了近似計算。

2.1 擬牛頓條件（擬牛頓方程，割線條件）

符號約定，使用\(B \approx H\)，\(D\approx H^{-1}\)

設經過\(k+1\)次迭代后得到\({\tt x}_{k+1}\)，此時將目標函數\(f(x)\)在\({\tt x}_{k+1}\)附近作泰勒展開，取二階近似，得到

\(f({\tt x})=f({\tt x}_{k+1}) + \nabla f({\tt x}_{k+1})({\tt x} - {\tt x}_{k+1}) + \frac{1}{2}({\tt x} - {\tt x}_{k+1})^T\nabla^2f({\tt x}_{k+1})({\tt x} - {\tt x}_{k+1})\)

兩邊同時作用一個梯度算子\(\nabla\)可得：\(\nabla f({\tt x}) \approx \nabla f({\tt x}_{k+1}) + H_{k+1}({\tt x} - {\tt x}_{k+1})\)

取\({\tt x} = {\tt x_k}\)，並整理可得：\(g_{k+1}-g_k \approx H_{k+1}·({\tt x}_{k+1} - {\tt x}_k)\)

若引入記號\({\tt s}_k = {\tt x}_{k+1} - {\tt x}_k\)，\({\tt y}_k = g_{k+1} - g_k\)

則2.16可以寫成緊湊形式\({\tt y}_k \approx H_{k+1}{\tt s}_k\)或\({\tt s}_k \approx H_{k+1}^{-1}{\tt y}_k\)

這就是所謂的擬牛頓法，它對迭代過程中的海森矩陣\(H_{k+1}\)作約束，因此對\(H_{k+1}\)做近似的\(B_{k+1}\)以及對\(H_{k+1}^{-1}\)做近似的\(D_{k+1}\)可以將上面的式子改寫為\({\tt y}_k = B_{k+1}{\tt s}_k\)或者\({\tt s}_k = D_{k+1}{\tt y}_k （2.21）\)

3、擬牛頓法之DFP算法

DFP算法由Davidon於1959年首先提出，后經Fletcher和Powell加以發展和完善，是最早的擬牛頓法，該算法的核心是：通過迭代的方法，對\(H_{k+1}^{-1}\)做近似，迭代格式為：

\(D_{k+1}=D_k + \Delta D_k, k=0,1,2,... （2.22）\)

其中的\(D_0\)通常取為單位矩陣\(I\)，因此關鍵是每一步的矯正矩陣$ D\Delta_k$如何構造。

注意迭代格式將嵌套在算法1.2中，因此我們才想\(\Delta D_k\)可能於\({\tt s}_k,{\tt y}_k\)和\(D_k\)發生關聯。這么我們采用待定法，即首先將\(\Delta D_k\)待定稱某中形式，然后結合擬牛頓條件（2.21）來進行推導。那將$ \Delta D_k$待定成什么形式呢？這個說起來比較tricky，我們將其待定為：

\(\Delta D_k = \alpha \tt uu^T + \beta vv^T （2.23）\)

其中，\(\alpha, \beta\)是待定系數，\(\tt u,v \in R^N\)為待定向量。從形式上看，這種待定公式至少保證了矩陣\(\Delta D_k\)的對稱性，因為\(\tt uu^T\) 和 $ \tt vv^T$都是對稱矩陣。將(2.23)帶入(2.22)並結合指導條件(2.21)，可得

\({\tt s}_k = D_k{\tt y}_k + \alpha{\tt uu^Ty}_k + \beta{vv^Ty}_k (2.24) \\ = D_k{\tt y}_k + {\tt u}(\alpha{\tt u^Ty}_k) + {\tt v}(\beta{v^Ty}_k)=D_k{\tt y}_k + (\alpha{\tt u^Ty}_k){\tt u} +(\beta{v^Ty}_k){\tt v}\)

其中括號中的\(\alpha{\tt u^Ty}_k\)和\(\beta{v^Ty}_k\)是兩個數，不妨做出簡單賦值，令

\(\alpha{\tt u^Ty}_k=1\)和\(\beta{v^Ty}_k=-1 （2.26）\)

即\(\alpha = \frac{1}{{\tt u^Ty}_k}\)和\(\beta = \frac{-1}{{v^Ty}_k}\)，其中\({\tt u}\)和\({\tt v}\)仍有待確定

將(2.26)帶入(2.25)得到\({\tt u} - {\tt v} = {\tt s}_k - D_k{\tt y}_k\)

要上式成立，不妨直接取\(\tt u=s_k, v=D_ky_k （2.29）\)

再將(2.29)帶入(2.27)便得\(\alpha = \frac{1}{{\tt s_k^Ty}_k}\)和\(\beta = \frac{-1}{(D_k{\tt y}_k)^T{\tt y}_k}= \frac{-1}{{\tt y}_k^T(D_k){\tt y}_k}\)（這里\(D_k\)是對稱的）

這樣我們便可以將\(\Delta D_k\)構造出來了，即\(\Delta D_k=\frac{{\tt s}_k{\tt s}_k^T}{{\tt s}_k^T{\tt y}_k} - \frac{D_k{\tt y}_k{\tt y}_k^TD_k}{{\tt y}_k^TD_k{\tt y}_k}\)

算法2.1 DFP算法

1. 給定初值\({\tt x}_0\)和精度閾值\(\epsilon\)，並令\(D_0:=I,k:=0\)
2. 確定搜索方向\(d_k=-D_k·g_k\)(這里\(d_k\)是定義參數，稱作牛頓方向)
3. 利用\(\lambda_k = \underset{\lambda \in R}{argmin}f({\tt x}_k + \lambda d_k)\)得到步長\(\lambda_k\)，並令\({\tt s}_k = \lambda_kd_k, {\tt x}_{k+1} := {\tt x}_k + {\tt s}_k\)
4. 若\(||g_{k+1}||<\epsilon\)，則算法結束
5. 計算\({\tt y}_k = g_{k+1} - g_k\)
6. 計算\(D_{k+1}=D_k + \Delta D_k = D_k + \frac{{\tt s}_k{\tt s}_k^T}{{\tt s}_k^T{\tt y}_k} - \frac{D_k{\tt y}_k{\tt y}_k^TD_k}{{\tt y}_k^TD_k{\tt y}_k}\)
7. 令\(k:=k+1\)，轉至步驟2

4、擬牛頓法之BFGS算法

BFGS算法是以發明者Broyden，Fletcher，Goldfard和Shanno四個人的名字的首字母命名的，比DFP性能更佳，目前它已成為求解無約束非線性優化問題最常用的方法之一。BFGS算法已有較完善的局部收斂理論，對其全局收斂性的研究也取得了重要成果。BFGS算法中核心公式的推導過程於DFP完全類似，只是互換了\({\tt s}_k,{\tt y}_k\)的位置。

需要注意的是，BFGS算法是直接逼近海森矩陣，即\(B_k \approx H_k\)，仍采用迭代方法，設迭代格式為：

\(B_{k+1}=B_k + \Delta B_k, k=0,1,2,...  (2.32)\)

其中\(B_0\)也常取單位矩陣\(I\)，因此關鍵是每一步的矯正矩陣\(\Delta B_k\)如何構造，同樣將其待定為

\(\Delta B_k = \alpha \tt uu^T + \beta vv^T （2.33）\)

將(2.33)帶入(2.32)並結合指導條件(2.20)可得\({\tt y}_k = B_{k+1}{\tt s}_k = B_k{\tt s}_k + (\alpha{\tt u}^T{\tt s}_k){\tt u}+ (\beta{\tt v}^T{\tt s}_k){\tt v}\)

令\(\alpha{\tt u^Ts}_k=1\)和$\beta{v^Ts}_k=-1 \(，以及\)u={\tt y}_k, {\tt v}=B_k{\tt s}_k\(，可算的即\)\alpha = \frac{1}{{\tt y}_k^T{\tt s}_k}\(和\)\beta = \frac{-1}{{\tt s}_k^TB_k{\tt s}_k}$

從而便可以得到\(\Delta B_k=\frac{{\tt y}_k{\tt y}_k^T}{{\tt y}_k^T{\tt s}_k} - \frac{B_k{\tt s}_k{\tt s}_k^TB_k}{{\tt s}_k^TB_k{\tt s}_k}\)

算法2.2 BFGS算法(I)

1. 給定初值\({\tt x}_0\)和精度閾值\(\epsilon\)，並令\(B_0:=I,k:=0\)
2. 確定搜索方向\(d_k=-B_k^{-1}·g_k\)(這里\(d_k\)是定義參數，稱作牛頓方向)
3. 利用\(\lambda_k = \underset{\lambda \in R}{argmin}f({\tt x}_k + \lambda d_k)\)得到步長\(\lambda_k\)，並令\({\tt s}_k = \lambda_kd_k, {\tt x}_{k+1} := {\tt x}_k + {\tt s}_k\)
4. 若\(||g_{k+1}||<\epsilon\)，則算法結束
5. 計算\({\tt y}_k = g_{k+1} - g_k\)
6. 計算\(B_{k+1}=B_k + \Delta B_k = B_k + \frac{{\tt y}_k{\tt y}_k^T}{{\tt y}_k^T{\tt s}_k} - \frac{B_k{\tt s}_k{\tt s}_k^TB_k}{{\tt s}_k^TB_k{\tt s}_k}\)
7. 令\(k:=k+1\)，轉至步驟2

算法2.2中步2通常是通過求解線性代數方程組\(B_kd_k=-g_k\)來進行，然而更一般的做法是通過步6的遞推關系應用Sherman-Morrison公式，直接給出\(B_{k+1}^{-1}\)和\(B_k^{-1}\)之間的關系式

\(B_{k+1}^{-1}= (I-\frac{s_ky_k^T}{y_k^Ts_k})B_k^{-1}(I-\frac{y_ks_k^T}{y_k^Ts_k}) + \frac{s_ks_k^T}{y_k^Ts_k} \\=B_k^{-1} + (\frac{1}{s_k^ty_k} + \frac{y_k^TB_k{-1}y_k}{(s_k^Ty_k)^2})s_ks_k^T -\frac{1}{s_k^Ty_k}(B_k^{-1}y_ks_k^T+s_ky_k^TB_k^{-1}) (2.39)\)

Sherman-Morrison公式：設\(A \in R^n\)為非奇異方陣，\(\text{u,v} \in R^n\)，若\(1+{\tt v}^TA^{-1}{\tt u} \neq 0\)，則有：

\((A+{\tt uv}^T)^{-1} = A^{-1} - \frac{A^-1{\tt uv}^TA^{-1}}{1+{\tt v}^TA^{-1}{\tt u}}\)

且為了避免出現矩陣求擬的符號，我們統一將\(B_i^{-1}\)用\(D_i\)替換，整個算法不再需要求解線性代數方程組，由矩陣-向量運算就可以完成了

算法2.2 BFGS算法(II)

1. 給定初值\({\tt x}_0\)和精度閾值\(\epsilon\)，並令\(D_0:=I,k:=0\)
2. 確定搜索方向\(d_k=-B_k^{-1}·g_k\)(這里\(d_k\)是定義參數，稱作牛頓方向)
3. 利用\(\lambda_k = \underset{\lambda \in R}{argmin}f({\tt x}_k + \lambda d_k)\)得到步長\(\lambda_k\)，並令\({\tt s}_k = \lambda_kd_k, {\tt x}_{k+1} := {\tt x}_k + {\tt s}_k\)
4. 若\(||g_{k+1}||<\epsilon\)，則算法結束
5. 計算\({\tt y}_k = g_{k+1} - g_k\)
6. 計算\(D_{k+1} = (I-\frac{s_ky_k^T}{y_k^Ts_k})D_k(I-\frac{y_ks_k^T}{y_k^Ts_k}) + \frac{s_ks_k^T}{y_k^Ts_k}\)（與DFP算法的區別點）
7. 令\(k:=k+1\)，轉至步驟2

5、L-BFGS算法

在BFGS算法中，需要用到一個\(N \times N\)的矩陣\(D_k\)，當N很大時，存儲這個矩陣將變得很耗計算機資源，L-BFGS算法不再存儲完整的矩陣\(D_k\)，而是存儲計算過程中的向量序列\(\{s_i\},\{y_i\}\)，需要矩陣\(D_k\)時，利用向量序列\(\{s_i\},\{y_i\}\)的計算來代替。而且向量序列\(\{s_i\},\{y_i\}\)也不是所有的都存，而是固定存最新的m個（參數m時可指定參數），每次計算\(D_k\)時，只利用最新的m個\(\{s_i\},\{y_i\}\)。顯然這樣以來我們將存儲由原來的\(O(N^2)\)降到了\(O(mN)\)

由算法2.3的更新公式\(D_{k+1} = (I-\frac{s_ky_k^T}{y_k^Ts_k})D_k(I-\frac{y_ks_k^T}{y_k^Ts_k}) + \frac{s_ks_k^T}{y_k^Ts_k}\)，可令\(\rho_k=\frac{1}{y_k^ts_k}, V_k=I-\rho_ky_ks_k^T\)，於是更新公式可以改寫為：\(D_{k+1}=V_k^TD_kV_k + \rho_ks_ks_k^T （2.43）\)

如果給定初始矩陣\(D_0\)，則利用（2.43）依次可得

\(D_0=V_0^TD_0V_0 + \rho_0s_0s_0^T\)

\(D_1=V_1^TD_1V_1 + \rho_1s_1s_1^T = V_1^T(V_0^TD_0V_0 + \rho_0s_0s_0^T)V_1 + \rho_1s_1s_1^T\)

...

\(D_{k+1}= (V_k^TV_{k-1}^T...V_1^TV_0^T)D_0(V_0V_1...V_{k-1}V_k) \\+(V_k^TV_{k-1}^T...V_2^Tv_1^T)(\rho_0s_0s_0^T)(V_1V_2...V_{k-1}V_k) \\+...\\+(V_k^TV_{k-1}^T)(\rho_{k-2}s_{k-2}s_{k-2}^T)(V_{k-1}V_k)\\+V_k^T(\rho_{k-1}s_{k-1}s_{k-1}^T)V_k \\+\rho_ks_ks_k^T\)

由上式可見，計算\(D_{k+1}\)需要用到\(\{s_i,y_i\}_{i=0}^k\)，因此若從\(s_0,y_0\)開始連續的存儲m組的話，只能存儲到\(s_{m-1}, y_{m-1}\)，也就是說只能依次計算\(D_1,D_2,...,D_m\)，那么\(D_{m+1}, D_{m+2}\)該如何計算呢？

如果一定要丟掉一些向量，那么肯定優先考慮那些最早生成的限量，具體來說，計算\(D_{m+1}\)時，我們保存\(\{s_i,y_i\}_{i=1}^m\)，丟掉了\(\{s_0,y_0\}\)，計算\(D_{m+2}\)時，我們保存\(\{s_i,y_i\}_{i=2}^{m+1}\)，丟掉\(\{s_1,y_1\}\)...

但是丟掉一些向量后就只能近似計算了，當\(k+1>m\)時，仿照(2.44)可以構造近似公式

\(D_{k+1}= (V_k^TV_{k-1}^T...V_{k-m+2}^TV_{k-m+1}^T)D_0(V_{k-m+1}V_{k-m+2}...V_{k-1}V_k) \\+(V_k^TV_{k-1}^T...V_{k-m+2}^Tv_{k-m+2}^T)(\rho_0s_0s_0^T)(V_{k-m+2}V_{k-m+3}...V_{k-1}V_k) \\+...\\+(V_k^TV_{k-1}^T)(\rho_{k-2}s_{k-2}s_{k-2}^T)(V_{k-1}V_k)\\+V_k^T(\rho_{k-1}s_{k-1}s_{k-1}^T)V_k \\+\rho_ks_ks_k^T\)

若引入\(\hat m = min{k, m-1}\)，則還可以合並寫成（\(k+1\leq m\)也成立）

\(D_{k+1}= (V_k^TV_{k-1}^T...V_{k-\hat m+2}^TV_{k-\hat m+1}^T)D_0(V_{k-\hat m+1}V_{k-\hat m+2}...V_{k-1}V_k) \\+(V_k^TV_{k-1}^T...V_{k-\hat m+2}^Tv_{k-\hat m+2}^T)(\rho_0s_0s_0^T)(V_{k-\hat m+2}V_{k-\hat m+3}...V_{k-1}V_k) \\+...\\+(V_k^TV_{k-1}^T)(\rho_{k-2}s_{k-2}s_{k-2}^T)(V_{k-1}V_k)\\+V_k^T(\rho_{k-1}s_{k-1}s_{k-1}^T)V_k \\+\rho_ks_ks_k^T\)

算法2.4 （\(D_k·g_k\)的快速算法）

Step1 初始化

初始化 \(\delta = \begin{cases} 0 & k\leq m \\ k-m & k>m\end{cases}\)； \(L = \begin{cases} k & k\leq m \\ m & k>m\end{cases}\)； \(q_L=g_k\)

Step2 后向循環

FOR i = L-1, L-2, ...,1,0 DO

{

\[\begin{array}{l} j=i+\delta \\ \alpha_i = \rho_js_j^Tq_{i+1}, \alpha_i需要保存下來，前向循環使用 \\ q_i = q_{i+1} - \alpha_iy_j \end{array} \]

}

Step3 前向循環

\(r_0\) = \(D_0 · q_0\)

FOR i = L-1, L-2, ...,1,0 DO

{

\[\begin{array}{l} j=i+\delta \\ \beta_i = \rho_jy_j^Tt_{i} \\ r_{i+1} = r_i + (\alpha_i - \beta_i)s_j \end{array} \]

}

最后算出的\(r_L\)即為\(H_k · g_k\)的值

二、梯度下降法

在求解節氣學習的模型參數，即無約束優化問題時，梯度下降法是最常用的方法之一，另一種常用的方法是最小二乘法。

2.1 梯度

在微積分里面，對多元函數的參數求偏導，把求得的各個參數的偏導數以向量的形式寫出來，就是梯度。比如函數\(f(x,y)\)分別對\(x,y\)求偏導的梯度向量就是\((\partial f/\partial x, \partial f/\partial y)\)，簡稱\(grad f(x,y)\)或者\(\nabla f(x,y)\)。對於在點\((x_0,y_0)\)的具體梯度向量就是\(\nabla f(x_0,y_0)\)。

這個梯度集合意義上講，就是函數變化最快的地方，具體來說對於函數\(f(x,y)\)在點\((x_0,y_0)\)，沿着梯度向量的方向就是\(\nabla f(x_0,y_0)\)的方向是\(f(x,y)\)增加最快的地方。或者說沿着梯度向量的方向，更容易找到函數的最大值。反過來說，沿着梯度向量相反的方向，也就是\(-\nabla f(x_0,y_0)\)的方向，梯度減少最快，也就是最容易找到函數的最小值。

2.2 梯度下降與梯度上升

在機器學習算法中，在最小化損失函數時，可以通過梯度下降法來一步步的迭代求解，得到最小化的損失函數和模型參數值，反過來，如果我們需要求解損失函數的最大值，就需要使用梯度上升法來迭代了。梯度下降法和梯度上升法是可以互相轉化的。

2.3 梯度下降法算法詳解

首先來看看梯度下降的一個直觀的解釋。比如我們在一座大山上的某處位置，由於我們不知道怎么下山，於是決定走一步算一步，也就是在每走到一個位置的時候，求解當前位置的梯度，沿着梯度的負方向，也就是當前最陡峭的位置向下走一步，然后繼續求解當前位置梯度，向這一步所在位置沿着最陡峭最易下山的位置走一步。這樣一步步的走下去，一直走到覺得我們已經到了山腳。當然這樣走下去，有可能我們不能走到山腳，而是到了某一個局部的山峰低處。

2.3.1 梯度下降法的相關概念

1）步長(Learning rate)：步長決定了在梯度下降的迭代過程中，每一步梯度負方向前進的長度（類似於阻尼牛頓法中的阻尼系數）。用上面下山的例子，步長就是在當前這一步所在位置沿着最陡峭最易下山的位置走的那一步的長度。

2）特征(feature)：指的是樣本中輸入部分，比如兩個單特征的樣本\((x^{(0)},y^{(0)}),(x^{(1)},y^{(1)})\)，牛頓法中就是因為特征值過多可能會導致二階梯度計算復雜，所以引入了擬牛頓法。

3）假設函數(hypothesis function)：在監督學習中，為了擬合輸入樣本，而使用的假設函數，記為\(h_\theta(x)\)，就是假設的函數形式，如\(h_\theta(x)=\theta_0 + \theta_1x\)

4）損失函數(loss function)：為了評估模型擬合的好壞，通常損失函數來度量擬合程度。分類算法中常用的損失函數是交叉熵損失，回歸算法中使用的是平方損失，如\(J(\theta_0, \theta_1)=\sum_{i=1}^m(h_\theta(x_i)-y_i)^2\)

2.3.2 梯度下降法的詳細算法（代數法和矩陣法）

梯度下降法的代數描述

1.先決條件：確認模型的假設函數和損失函數

比如線性回歸，假設函數表示為\(h_\theta(x_1,x_2,...,x_n)=\theta_0 + \theta_1x_1 + ... + \theta_nx_n\)

損失函數為：\(J(\theta_0, theta_1,...,\theta_n) = \frac{1}{2m}\sum_{j=1}^m(h_\theta(x_0^{(j)}, x_1^{(j)},...,x_n^{(j)}), y_j)^2\)

2.算法相關參數初始化

主要是初始化\(\theta_0, \theta_1,...,\theta_n\)，算法終止距離\(\epsilon\)以及步長\(\alpha\)

3.算法過程：

1. 確定當前位置的損失函數的梯度，對於\(\theta_i\)其梯度表達式為：\(\frac{\partial}{\partial\theta_i}J(\theta_0, \theta_1,...,\theta_n)\)

2. 用步長乘以損失函數的梯度，得到當前位置下降的距離：\(\alpha \frac{\partial}{\partial\theta_i}J(\theta_0, \theta_1,...,\theta_n)\)

3. 確定是否所有的\(\theta_i\)梯度下降的距離都小於\(\epsilon\)，如果小於\(\epsilon\)則算法終止，當前所有的\(\theta_i (i=0,1,...,n)\)。否則進入步驟4

4. 更新所有的\(\theta\)，對於\(\theta_i\)其更新表達式為：\(\theta_i = \theta_i - \alpha \frac{\partial}{\partial\theta_i}J(\theta_0, \theta_1,...,\theta_n)\)

以線性回歸的例子來具體描述，假設樣本\((x_1^{(0)},x_2^{(0)},...,x_n^{(0)},y_0),(x_1^{(1)},x_2^{(1)},...,x_n^{(1)},y_1),...,(x_1^{(m)},x_2^{(m)},...,x_n^{(m)},y_m)\)

損失函數為：\(J(\theta_0, theta_1,...,\theta_n) = \frac{1}{2m}\sum_{j=1}^m(h_\theta(x_0^{(j)}, x_1^{(j)},...,x_n^{(j)}), y_j)^2\)

對\(\theta_i\)求偏導：\(\frac{\partial}{\partial\theta_i}J(\theta_0, \theta_1,...,\theta_n)=\frac{1}{m}\sum_{j=0}^m(h_\theta(x_0^{(j)},x_1^{(j)},...,x_n^{(j)})-y_j)x_i^{(j)}\)

由於樣本中沒有\(x_0\)，上式令所有的\(x_0^j\)為1

更新參數：\(\theta_i = \theta_i - \alpha \frac{1}{m}\sum_{j=0}^m(h_\theta(x_0^{(j)},x_1^{(j)},...,x_n^{(j)})-y_j)x_i^{(j)}\)

上文中加上\(\frac{1}{m}\)是為了好理解，由步長為常數，所有\(\alpha \frac{1}{m}\)可以合並到，使用一個常數表示即可(這就是為什么在設置學習率的時候，往往設置的很小，是因為在原學習率的基礎上需要再除以一個樣本數)

梯度下降法的矩陣方式描述

1.先決條件

同上一節，但是可以簡化為矩陣形式\(h_\theta(X) = X\theta\)，其中\(h_\theta(X)\)為\(m \times 1\)的向量，\(\theta\)為\((n+1) \times 1\)的向量，\(X\)為\(m \times (n+1)\)維的矩陣。m代表樣本個數，n+1代表樣本的特征數

損失函數為:\(J(\theta)=\frac{1}{2}(X\theta-Y)^T(X\theta-Y)\)，其中Y是樣本輸出向量，維度為\(m \times 1\)

2.算法相關參數初始化

同上一節

3.算法過程

1. 確定當前損失函數的梯度\(\frac{\partial J(\theta)}{\partial \theta}\)

2. 用步長乘以損失函數的梯度，得到當前位置下降的距離：\(\alpha \frac{\partial}{\partial\theta}J(\theta)\)

3. 確定\(\theta\)中每個值的距離都小於\(\epsilon\)，如果小於\(\epsilon\)則算法終止，當前所有的\(\theta_i (i=0,1,...,n)\)。否則進入步驟4

4. 更新\(\theta\)向量：\(\theta = \theta - \alpha \frac{\partial}{\partial\theta}J(\theta)\)

還是上面的例子，其中對\(\theta\)向量的偏導數計算為\(\frac{\partial}{\partial\theta}J(\theta)=X^T(X\theta-Y)\)

更新：\(\theta = \theta - \alpha \frac{\partial}{\partial\theta}J(\theta) = \theta - \alpha X^T(X\theta - Y)\)

3.3 梯度下降的算法調優

1）算法的步長選擇，可以先大后小。（好像有一個專門的算法，初始時先快速增大，然后再持續減小，查資料）

2）算法參數的初始值選擇。

3）特征歸一化。由於樣本不同特征的取值范圍不一樣，可能導致迭代很慢，為了減少特征取值的影響，可以對特征數據歸一化。（經常看到的一個圖，神經網絡中為什么要對數據做歸一化）

4.梯度下降法和其他無約束優化算法的比較

在機器學習中的無約束優化算法，除了梯度下降以外，還有前面提到的最小二乘法，此外還有牛頓法和擬牛頓法。

　　　　梯度下降法和最小二乘法相比，梯度下降法需要選擇步長，而最小二乘法不需要。梯度下降法是迭代求解，最小二乘法是計算解析解。如果樣本量不算很大，且存在解析解，最小二乘法比起梯度下降法要有優勢，計算速度很快。但是如果樣本量很大，用最小二乘法由於需要求一個超級大的逆矩陣，這時就很難或者很慢才能求解解析解了，使用迭代的梯度下降法比較有優勢。

　　　　梯度下降法和牛頓法/擬牛頓法相比，兩者都是迭代求解，不過梯度下降法是梯度求解，而牛頓法/擬牛頓法是用二階的海森矩陣的逆矩陣或偽逆矩陣求解。相對而言，使用牛頓法/擬牛頓法收斂更快。但是每次迭代的時間比梯度下降法長。

三、梯度下降法的改進

3.1 優化算法通用框架

首先定義：待優化參數\(\theta\)，目標函數\(f(\theta)\)，初始學習率\(\alpha\)，而后開始進行迭代優化，在每個迭代步驟：

1.計算目標函數關於當前參數的梯度 \(g_t = \nabla f(\theta_t)\)

2.根據歷史梯度計算一階動量和二階動量: \(m_t = \phi(g_1, g_2, ...,g_t), V_t=\psi(g_1,g_2,...,g_t)\)

3.計算當前時刻的下降梯度:\(\eta_t = \alpha · m_t / \sqrt{V_t}\) (相當於SGD中介紹的往下走的距離)

4.根據夏季那個梯度進行更新\(\theta_{t+1} = \theta_t - \eta_t\)

注：這里的一階動量和二階動量與牛頓法中的一階導數和二階導數的概念容易混淆；一階導數就是梯度，而一階動量和二階動量是根據歷史梯度計算出來的

根據上述框架，來看看SGD解釋，\(\theta = \theta - \alpha \frac{\partial}{\partial\theta}J(\theta)\)，其中:

下降梯度可以描述為：\(\eta_t=\alpha ·g_t\)

根據通用框架：\(\eta_t = \alpha · m_t / \sqrt{V_t}\)

因此可令\(m_t=g_t, V_t=I^2\)

3.2 SGD with Momentum

為了抑制SGD的震盪，Momentum認為梯度下降過程中可以加入慣性。下坡的時候，如果發現是陡坡，就利用慣性跑的快一些，SGDM（SGD with Momentum）在SGD的基礎上引入了一階動量，

\(m_t = \beta_1m_{t-1} + (1-\beta_1)·g_t\)

一階動量是各個時刻梯度方向的指數移動平均值，約等於\(1/(1-\beta_1)\)個時刻的梯度向量的平均值。也就是說t時刻的下降方向，不僅由當前點的梯度方向決定，而且由此前累積的下降方向決定。\(\beta_1\)的經驗值為0.9，這就意味着下降方向主要是此前積累的下降方向，並略微偏向當前時刻的下降方向。可以想象高速公路上汽車轉彎，在高速向前的同時略微偏向，急轉彎時要出事的。

3.3 SGD with Newterov Acceleration

SGD還有一個問題時困在局部最優的溝壑里面震盪。想象一下你走到一個盆地，四周都是略高的小山，你覺得沒有下坡的方向，那就只能待在這里了。可是如果你爬上高地，就會發現外面的世界還很廣闊。因此，我們不能停留在當前位置去觀察未來的方向，而要向前一步、多看一步、看遠一些。

NAG的全稱Nesterov Accelerated Gradient是在SGDM的基礎上進一步改進，主要是針對梯度的改進。不再計算當前位置的梯度方向，而是按照累積動量走了一步，那個時候下降的方向。即計算梯度的時候，將動量考慮了進來。

\(g_t = \nabla f(\theta_t - \alpha m_{t-1}/\sqrt{V_{t-1}})\)

在此梯度的基礎上，重新計算一階動量和二階動量

3.3 AdaGrad

二階動量的出現，意味着“自適應學習率”優化算法時代的到來。SGD及其變種以同樣的學習率更新每個參數，但神經網絡往往包含大量的參數，這些參數並不是總會用得到，對於經常更新的參數，我們已經積累了大量關於它的知識，不希望被耽擱樣本影響太大，希望學習速率慢一些；對於偶爾更新的參數，我們了解的信息太少，希望能從每個偶然出現的的樣本本身多學一些，即學習速率大一些。利用二階動量，獲取每個維度上迄今為止所有梯度值的平方和：

\(V_t = \sum_{j=1}^t g_j^2\)

根據\(\eta_t = \alpha · m_t / \sqrt{V_t}\)，實際上此時的學習率已經變為\(\alpha / \sqrt{V_t}\)（歷史累積梯度越大，對應維度的學習率就越小）。這一方法在系數數據場景下表現非常好，但也存在一些問題：因為 \(\sqrt{V_t}\)是單調遞增的，會使得學習率單調遞減至0，可能會使訓練過程提前結束，即便后續還有數據也無法學到必要的知識。

3.4 AdaDelta/RMSProp

由於AdaGrad單調遞減的學習率變化過於激進，我們考慮一個改變二階動量計算方法的策略，不累積全部歷史梯度，而只關注過去一段時間穿夠的下降速度。這就是AdaDelta名稱中Delta的來歷。

SGDM中有介紹，指數移動平均大約就是過去一段時間的平均值，因為我們用這一方法來計算二階累積動量：

\(V_t = \beta_2 * V_{t-1} + (1-\beta_2)g_t^2\)

這樣就避免了二階動量的持續累積，導致訓練過程提前結束的問題了。

注意：AdaDelta不需要單獨設置學習率

RMSprop其實是Adadelta的一種特殊情況，其依然靠的是全局的學習率

3.5 Adam

前面介紹，SGDM在SGD的基礎上增加了一階動量，AdaGrad和AdaDelta在SGD基礎上增加了二階動量，把一階和二階動量都用起來，就是Adam了，Adaptive + Nomentum

一階動量：\(m_t = \beta_1 ·m_{t-1} + (1- \beta_1) ·g_t\)

二階動量：\(V_t =\beta_2*V_{t-1}*(1-\beta_2)g_t^2\)

實際使用過程中，取\(\beta_1 = 0.9, \beta_2 = 0.999\)，初始化\(m_0 = 0, V_0 = 0\)

1. \(m_t\)大，\(V_t\)大，梯度大且穩定，遇到一個明顯的大坡，前進方向明確

2)\(m_t\)大，\(V_t\)小，不可能出現

3)\(m_t\)小，\(V_t\)大，梯度不穩定，可能遇到一個峽谷，窯易引起反彈震盪

4)\(m_t\)小，\(V_t\)小，梯度趨於零，可能到這局部最低點，也可能走到一片坡度假韁的平地 ， 此時要避免陷入平原( plateau )

3.6 Nadam

在Nadam是在Adam的基礎上增加了Nestreov，體現在計算梯度的公式上：

\(g_t = \nabla f(\theta_t - \alpha m_{t-1}/\sqrt{V_{t-1}})\)

四、如何選擇優化算法

4.1 Adam存在的幾個問題

4.1.1 可能不收斂

二階動量是固定時間窗口內的累積，隨着時間窗口的變化，遇到的數據可能發生巨變，使得 \(V_t\) 可能會時大時小，不是單調變化。這就可能在訓練后期引起學習率的震盪，導致模型無法收斂。

對二階動量進行修正\(V_t = max(\beta_2 V_{t-1}+(1-\beta_2)g_t^2, V_{t-1})\)，這樣就保證了\(||V_t|| \geq ||V_{t-1}||\)，從而使學習率單調遞減。

4.1.2 可能錯過全局最優解

同樣的一個優化問題，不同的優化算法可能會找到不同的答案，但自適應學習率的算法往往找到非常差的答案。他們通過一個特定的數據例子說明，自適應學習率算法可能會對前期出現的特征過擬合，后期才出現的特征很難糾正前期的擬合效果。

改進Adam的方法：前期用Adam，享受Adam快速收斂的優勢；后期切換到SGD，慢慢尋找最優解。

4.1.3 小結

理解數據對於設計算法的必要性。優化算法的演變歷史，都是基於對數據的某種假設而進行的優化，那么某種算法是否有效，就要看你的數據是否符合該算法的胃口了。

4.2 Adam+SGD

需要注意解決兩個問題：

1. 切換算法以后用什么樣的學習率？——Adam用的是自適應學習率，依賴的是二階動量的累積，SGD接着訓練的話，用什么樣的學習率？
2. 什么時候切換優化算法？——如果切換太晚，Adam可能已經跑到自己的盆地里去了，SGD再怎么好也跑不出來了。

4.2.1 切換算法后用什么樣的學習率

Adam的下降方向是：\(\eta_t^{Adam} = (\alpha /\sqrt{V_t})·m_t\)

SGD的下降方向是：\(\eta_t^{SGD} = \alpha^{SGD}·g_t\)

\(\eta_t^{SGD}\)必定可以分解為\(\eta_t^{Adam}\)所在方向及其正交方向上的兩個方向的和。那么其在\(\eta_t^{Adam}\)上的投影就意味着SGD在Adam算法決定的下降方向上前進的距離，而在\(\eta_t^{Adam}\)的正交方向上的投影是SGD在自己選擇的修正方向上前進的距離。

如果SGD要走完Adam未走完的路，那就首先要結果Adam的大旗，並沿着\(\eta_t^{Adam}\)方向走一步，而后再沿着其正交方向走相應的一步。

這樣我們就知道該如何確定SGD的步長了--SGD在Adam下降方向上的正交投影，應該正好等於Adam的下降方向（含步長），也即：\(proj_{\eta_t^{SGD}} = \eta_t^{Adam}\)

解這個方程可以得到接續進行SGD的學習率：\(\alpha_t^{SGD}=((\eta_t^{Adam})^T\eta_t^{Adam})/((\eta_t^{Adam})^Tg_t)\)

為了減少噪聲的影響，作者使用移動平均值來修正對學習率的預估，即：

\(\lambda_t^{SGD} = \beta_2· \lambda_{t-1}^{SGD} + (1-\beta_2)·\alpha_t^{SGD}\)， \(\tilde \lambda_t^{SGD}=\lambda_t^{SGD}/(1-\beta_2^t)\)

這里直接服用了Adam的\(\beta_2\)參數

4.2.2 什么時候切算法

當SGD的相應學習率的移動平均值基本不變的時候，即\(|\tilde \lambda_t^{SGD} - \alpha_t^{SGD}| < \epsilon\)，每次迭代完都會計算一下SGD接班人的相應學習率；如果發現基本穩定了，那就使用SGD，並以\(\tilde \lambda_t^{SGD}\)為學習率繼續前進

4.3 優化算法常用tricks

1. 各大算法孰優孰劣並無定論，如果剛入門，優先考慮SGD+Nesterov Momentum或者Adam
2. 選擇你熟悉的算法--這樣就可以更加熟練的利用經驗進行調參
3. 充分了解你的數據--如果數據非常系數，那么優先考慮自適應學習率的算法
4. 根據你的需求來選擇--再模型設計實驗過程中，要快速驗證新模型的效果，可以先用Adam進行快速實驗優化，在模型上線或者結果發布前，可以用精調的SGD進行模型的極致優化
5. 先用小數據集進行實驗
6. 考慮不同算法的組合--先用Adam進行快速下降，而后再換到SGD進行充分的調優
7. 數據集一定要充分的打散(shuffle)--這樣再使用自適應學習率算法的時候，可以避免某些特征幾種出現，而導致的有時學習過度、有時學習不足，使得下降方向出現偏差的問題
8. 訓練過程中持續監控訓練數據和驗證數據上的目標函數值以及精度或者AUC等指標變化情況，避免出現過擬合
9. 指定一個合適的學習率衰減策略--可以使用定期衰減，比如每過多少個epoch就衰減一次，或者利用精度或AUC等性能指標來監控，當測試集上的指標不變或者下跌時，就降低學習率。

4.4 為什么深度學習不使用牛頓法或者擬牛頓法

原因1: 牛頓法需要用到梯度和Hessian矩陣，這兩個都難以求解。因為很難寫出深度神經網絡擬合函數的表達式，遑論直接得到其梯度表達式，更不要說得到基於梯度的Hessian矩陣了。

原因2: 即使可以得到梯度和Hessian矩陣，當輸入向量的維度N較大時，Hessian矩陣的大小是N×N，所需要的內存非常大。

原因3: 在高維非凸優化問題中，鞍點相對於局部最小值的數量非常多，而且鞍點處的損失值相對於局部最小值處也比較大。而二階優化算法是尋找梯度為0的點，所以很容易陷入鞍點。