L-BFGS演算法剖析

05-16

在BFGS演算法或者DFP演算法中，需要保存一個 $N imes N$ 的矩陣，該矩陣用於表示Hessian的逆矩陣的近似。其中 $N$ 表示特徵維度。當特徵維度較高時，這個矩陣會超級佔用內存。例如在CTR預估等場景中，特徵動輒百w到千w，存儲這樣一個巨大的矩陣，需要的內存往往一般伺服器都無法滿足。因此研究者們提出了一種近似方法來解決這個問題，這就是L-BFGS，L表示limited-memory或者limited-storage。

其基本思想是，不再存儲Hessian矩陣（其實是逆矩陣），轉而存儲中間的一些信息，當然中間的信息量可能也很大，那就只保存最新的一些，通過他們推導計算出Hessian矩陣。這樣，其與原有的實現邏輯丟失了一些，所以只能是近似演算法，但是怎麼佔用的資源由 $O(N^2)$ 降到了 $O(mN)$ ，其中 $m$ 表示最新的m個中間變數。

前面說過，BFGS的最後一步更新的公式為： $H_{k+1}=(I-frac{{s_k}{{y_k}^T}}{{{y_k}^T}{s_k}})D_k(I -frac{{y_k}{{s_k}^T}}{{{y_k}^T}{s_k}})+frac{{s_k}{{s_k}^T}}{{{y_k}^T}{s_k}}$ 。

可以看到，整個過程只與 $y_k$ 和 $s_k$ 有關，每次保存這倆個變數，就能恢復每一輪的 $H_k$ 。 $y_k$ 和 $s_k$ 就是中間保存的變數，其都是 $n*1$ 的向量。相比巨大的 $N$ ，其一般小很多，一般在實踐中 $m$ 多取3到20之間。

保存的邏輯很簡單，但是實際操作中，由於最後更新的為特徵的權重 $x_k$ ，其一般更新形式為 $x_{k+1}=x_k-{lambda_k}{H_k}{g_k}$ ，其中 ${lambda_k}$ 為更新的一個步長，一般通過一維搜索或信賴域等演算法求得，可參照擬牛頓法之DFP演算法剖析一維搜索部分求得。因此一般，均是採用高效方式求得 ${H_k}{g_k}$ 。

在Numerical Optimization的第7章講述了兩步循環求得 ${H_k}{g_k}$ 的過程，主要是將遞推的方式進行了公式推導。其步驟如下：

這裡的 $Lambda f_k$ 即是我們所的 $g_k$ ，其中 $ho$ 是前面 $H_{k+1}=(I-frac{{s_k}{{y_k}^T}}{{{y_k}^T}{s_k}})D_k(I -frac{{y_k}{{s_k}^T}}{{{y_k}^T}{s_k}})+frac{{s_k}{{s_k}^T}}{{{y_k}^T}{s_k}}$ 中的分母部分，即 $ho_k=frac{1}{{{y_k}^T}{s_k}}$ 。如此更新完後，更新 $x_{k+1}$ ，其它就與BFGS一樣了。

總結

自己在實驗完後，簡單測試幾把，發現其實相較於BFGS，其準確度並沒降低太多。此外，在特徵維度較低時，其實質上每步多了一個恢復 $H_k$ 的操作，因此效率也並不見得比BFGS高效。但是在特徵維度較高時，其優勢體現就明顯很多。最後，代碼在https://github.com/sloth2012/scalaML/blob/master/src/main/scala/com/lx/algos/ml/optim/newton/LBFGS.scala，使用仍在https://github.com/sloth2012/scalaML/blob/master/src/test/scala/com/lx/algos/ml/NewtonTest.scala。熟悉前面BFGS代碼可以發現，僅是在更新 $H_k$ 的部分做了些更改，其它部分完全一致。