Gradient Boosting 能夠看作是一個整體的算法框架,起始於Friedman 的論文 [Greedy Function Approximation: A Gradient Boosting Machine] 。XGBoost (eXtreme Gradient Boosting) 是於2015年提出的一個新的 Gradient Boosting 實現,由華盛頓大學的 陳天奇 等人開發,在速度和精度上都有顯著提高,於是近年來在 Kaggle 等各大數據科學比賽中都獲得了普遍應用。本文主要對其原理進行闡述,並將其與傳統的 GBDT 進行比較。html
大致來看,XGBoost 在原理方面的改進主要就是在損失函數上做文章。一是在原損失函數的基礎上添加了正則化項產生了新的目標函數,這相似於對每棵樹進行了剪枝並限制了葉結點上的分數來防止過擬合。二是對目標函數進行二階泰勒展開,以相似牛頓法的方式來進行優化(事實上早在 [Friedman, J., Hastie, T. and Tibshirani, R., 1999] 中就已有相似方案,即利用二階導信息來最小化目標函數,陳天奇在論文中也提到了這一點)。python
在上一篇文章中,瞭解到 Gradient Boosting 的思想能夠理解爲經過函數空間的梯度降低來最小化目標函數 \(L(f) = \sum\limits_{i=1}^NL(y_i,f_m(x_i))\),其只使用了一階導信息,而 XGBoost 引入二階導的一大好處是能夠推導出一種新的增益計算方法,事實證實採用新的增益計算方法在優化目標函數上更加有效,精確度上也賽過傳統的 GBDT。因此下面也從目標函數的優化入手進行推導。web
與 GBDT 同樣,XGBoost 一樣採用加法模型,設基學習器爲\(f(x)\),預測值爲\(\hat{y}\):
\[ \hat{y} = \sum\limits_{m=1}^M f_m(x)\;, \quad f_m \in \mathcal{F} \tag{1.1} \]
在第m步,前m-1個基學習器是固定的,於是目標函數爲
\[ \mathcal{L}^{(m)} = \sum\limits_{i=1}^N L(y_i, \;\hat{y_i}^{(m-1)} + f_m(x_i)) + \Omega(f_m) \tag{1.2} \]
在傳統的 GBDT 中,是沒有正則化項 \(\Omega\) 的,而在XGBoost中採用的是
\[ \Omega(f) = \gamma \, J + \frac12 \lambda \sum\limits_{j=1}^J b_j^2 \]
其中 \(J\) 爲葉結點數目,\(b_j\) 爲各個葉結點上的值。該項限制了樹的複雜度,這樣至關於使葉結點的數目變小,同時限制葉結點上的分數,由於一般分數越大學得越快,就越容易過擬合。算法
接下來將 \(\Omega\) 代入\((1.2)\) 式並對目標函數在 \(\hat{y}_i^{(m-1)}\) 處進行二階泰勒展開:
\[ \mathcal{L}^{(m)} \simeq \sum\limits_{i=1}^N \left[ L(y_i, \hat{y}_i^{(m-1)}) + g_if_m(x_i) + \frac12 h_i f_m^2(x_i)\right] + \gamma \, J + \frac12 \lambda \sum\limits_{j=1}^J b_j^2 \tag{1.3} \]
其中 \(g_i = \frac{\partial L(y_i, \,\hat{y}_i^{m-1})}{\partial \, \hat{y}^{(m-1)}}\;\;,\;\; h_i = \frac{\partial^2 L(y_i, \,\hat{y}_i^{(m-1)})}{(\partial\, \hat{y}^{(m-1)})^2}\)多線程
在上一篇 Gradient Boosting 中提到決策樹將特徵空間劃分爲各個獨立區域,每一個樣本只屬於其中一個區域,於是單顆決策樹可表示爲 \(f(x) = \sum\limits_{j=1}^J b_j I(x \in R_j)\),若是將全部樣本加起來,則可表示爲 \(\sum\limits_{i=1}^N f(x_i) = \sum\limits_{j=1}^J \sum\limits _{x \in R_j} b_j\) ,代入 \((1.3)\) 式並將常數項 \(L(y_i, \hat{y}_i^{(m-1)})\) 移去:
\[ \begin{align*} \mathcal{L}^{(m)} &= \sum\limits_{i=1}^N \left[ g_if_m(x_i) + \frac12 h_i f_m^2(x_i)\right] + \gamma \, J + \frac12 \lambda \sum\limits_{j=1}^J b_j^2 \\ &= \sum\limits_{j=1}^J \left[ \sum\limits_{x_i \in R_j}g_ib_j + \frac12 \sum\limits_{x_i \in R_j} h_i b_j^2\right] + \gamma \, J + \frac12 \lambda \sum\limits_{j=1}^J b_j^2 \\ &= \sum\limits_{j=1}^J \left[ \sum\limits_{x_i \in R_j}g_ib_j + \frac12 (\sum\limits_{x_i \in R_j} h_i + \lambda) \,b_j^2 \right] + \gamma \, J \\ & = \sum\limits_{j=1}^J(G_j b_j + \frac12 (H_j + \lambda) \,b_j^2) + \gamma\, J \tag{1.4} \end{align*} \]框架
其中 \(G_j = \sum\limits_{x_i \in R_j} g_i \;\; , \;\; H_j = \sum\limits_{x_i \in R_j} h_i\)機器學習
通過一系列推導後,如今咱們的目標變爲最小化 \((1.4)\) 式,並求得相應的樹結構 \(\{R_j\}^J_1\) 和葉結點上的值 \(\{b_j\}^J_1\) 。分佈式
精確地劃分 \(\{R_j\}^J_1\) 是一個 NP hard 問題,現實中常使用貪心法,遍歷每一個特徵的每一個取值,計算分裂先後的增益,並選擇增益最大的進行分裂。
對於具體增益的衡量標準,在幾種決策樹算法中,ID3 採用了信息增益:
\[ I(D,A) = H(D) - H(D|A) = H(D) - \sum\limits_{v=1}^V\frac{|D^v|}{|D|} H(D^v) \]
其中 V 表示特徵 A 有 V 個可能的取值,\(D^v\) 表示第 v 個取值上的樣本數量。ide
C4.5 採用了信息增益比:
\[ I_R(D,A) = \frac{I(D, A)}{H_A(D)} \]
其中 \(H_A(D) = -\sum\limits_{v=1}^V \frac{|D^v|}{|D|} log_2 \frac{|D^v|}{|D|}\) 。函數
CART 分類樹採用了基尼係數:
\[ Gini(D,A) = \sum\limits_{v=1}^V \frac{|D^v|}{|D|} Gini(D^v) = \sum\limits_{v=1}^V \frac{|D^v|}{|D|} \left(1 - \sum\limits_{k=1}^K \left\lbrace\frac{|C_k|}{|D|}\right\rbrace^2\right) \]
其中 K 爲類別個數,\(|C_k|\) 爲 \(D\) 中屬於第k類的樣本數量。
CART 迴歸樹採用了均方偏差:
\[ \mathop{min}_{A,s}\Bigg[\mathop{min}_{c_1}\sum\limits_{x_i \in R_1(A,s)}(y_i - c_1)^2 + \mathop{min}_{c_2}\sum\limits_{x_i \in R_2(A,s)}(y_i - c_2)^2\Bigg] \]
其中 \(c_1\)爲特徵 A 的切分點 s 劃分出來的 \(R_1\) 區域的樣本輸出均值,\(c_1 = mean(y_i | x_i \in R_1(A,s))\),\(c_2\)爲 \(R_2\) 區域的樣本輸出均值,\(c_2 = mean(y_i | x_i \in R_2(A,s))\) 。
而XGBoost則提出了一種新的增益計算方法。
若是已經肯定了樹的結構 \(\{R_j\}^J_1\) ,則直接對 \((1.4)\) 式求導,最優值爲:
\[ b^*_j = - \frac{G_j}{H_j + \lambda} \]
再代入\((1.4)\)式獲得此時最小的爲目標函數爲:
\[ \mathcal{L^{(m)}} = -\frac12 \sum\limits_{j=1}^J \frac{G_j^2}{H_j+ \lambda} + \gamma\, J \tag{1.5} \]
式 \((1.5)\) 能夠認爲是 XGBoost 的打分函數,該值越小,說明樹的結構越好,下圖示例了該式的計算方法:
該式的優勢是除了能做爲樹分裂的衡量標準外,還能使 XGboost 適用於各類不一樣的損失函數,因此 XGBoost 包中支持自定義損失函數,但前提是一階和二階可導。
從另外一角度看, \((1.5)\) 式就相似於傳統決策樹中的不純度指標,在決策樹中咱們但願一次分裂後兩邊子樹的不純度越低越好,對應到XGBoost中則是但願 \((1.5)\)式 越小越好,即 \(\frac{G_j^2}{H_j+ \lambda}\) 越大越好,這樣分裂先後的增益定義爲:
\[ Gain = \frac{G_L^2}{H_L+ \lambda} + \frac{G_R^2}{H_R+ \lambda} - \frac{(G_L + G_R)^2}{H_L+ H_R + \lambda} - \gamma \tag{1.6} \]
\(\frac{G_L^2}{H_L+ \lambda}\) 爲分裂後左子樹的分數,\(\frac{G_R^2}{H_R+ \lambda}\) 爲分裂後右子樹的分數,$\frac{(G_L + G_R)^2}{H_L+ H_R + \lambda} $ 爲分裂前的分數,\(Gain\) 越大說明越值得分裂。固然 \((1.6)\) 式中 \(Gain\) 可能會變成負的,這個時候分裂後的目標函數不會減小,但這並不意味着不會分裂 。事實上 XGBoost 採用的是後剪枝的策略,建每棵樹的時候會一直分裂到指定的最大深度(max_depth),而後遞歸地從葉結點向上進行剪枝,對以前每一個分裂進行考察,若是該分裂以後的 \(Gain \leqslant 0\),則咔咔掉。 \(\gamma\) 是一個超參數,具備雙重含義,一個是在 \((1.3)\) 式中對葉結點數目進行控制的參數;另外一個是 \((1.6)\) 式中分裂先後 \(Gain\) 增大的閾值,固然兩者的目的是同樣的,即限制樹的規模防止過擬合。
接下來考察決策樹創建的過程。若是是使用貪心法,就是遍歷一個葉結點上的全部候選特徵和取值,分別計算 \(Gain\) ,選擇 \(Gain\) 最大的候選特徵和取值進行分裂,以下樹分裂算法流程 (注意這是單個葉結點的分裂流程):
有了上述單個葉結點上的分裂流程後,咱們能夠總結下整個 XGBoost 的學習算法:
- 初始化 \(f_0(x)\)
- for m=1 to M :
(a) 計算損失函數在每一個訓練樣本點的一階導數 \(g_i = \frac{\partial L(y_i, \,\hat{y}_i^{m-1})}{\partial \, \hat{y}^{(m-1)}}\) 和二階導數 \(h_i = \frac{\partial^2 L(y_i, \,\hat{y}_i^{(m-1)})}{(\partial\, \hat{y}^{(m-1)})^2}\) , $ i = 1,2 \cdots N$
(b) 遞歸地運用樹分裂算法生成一顆決策樹 \(f_m(x)\)
(c) 把新生成的決策樹添加到模型中, $\hat{y_i}^{(m)} = \hat{y_i}^{(m-1)} + f_m(x) $
若是把上述 XGBoost 的學習算法和上一篇中傳統 GBDT 的學習算法做比較,XGBoost 的主要優點是在損失函數中加入正則化項後使得學習出來的樹更加簡單,防止過擬合,但除此之外並不能體現出 XGBoost 的速度優點。XGBoost 之因此快的一大緣由是在工程上實現了 Column Block 方法,使得並行訓練成爲了可能。
對於決策樹來講,對連續值特徵進行劃分一般比較困難,由於連續值特徵每每取值衆多。一般的作法是先按特徵取值對樣本排序,再按順序計算增益選擇劃分點。若每次分裂都要排一次序,那是很是耗時的,因此 XGBoost 的作法是在訓練以前,預先按特徵取值對樣本進行了排序,而後保存爲block結構,採用CSC格式存儲,每一列(一個特徵列)均升序存放,這樣,一次讀入數據並排好序後,之後都可重複使用,大大減少計算量。
因爲已經預先排過序,因此在樹分裂算法中,經過一次線性掃描 (linear scan) 就能找出一個特徵的最優分裂點,以下圖所示,同時能夠看到缺失值並無保存:
陳天奇論文裏有關鍵的一句話:Collecting statistics for each columns can be parallelized 。因爲已經預先保存爲block 結構,因此在對葉結點進行分裂時,每一個特徵的增益計算就能夠開多線程進行,訓練速度也由此提高了不少。並且這種 block 結構也支持列抽樣,只要每次從全部 block 特徵中選擇一個子集做爲候選分裂特徵就能夠了,據個人使用經驗,列抽樣大部分時候都比行抽樣的效果好。
當數據量很是大難以被所有加載進內存時或者在分佈式環境下時,上文的樹分裂算法將再也不適用,於是做者提出了近似分裂算法,不只解決了這個問題,也同時能提高訓練速度,具體流程見下圖:
近似分裂算法其實就是對連續性特徵進行離散化,對於某個特徵 \(k\),算法首先根據特徵分佈找到 \(l\) 個分位點的候選集合 \(S_k = \{s_{k1}, s_{k2}, ... ,s_{kl} \}\) ,而後根據集合 \(S_k\) 中的分位點將相應樣本劃分到桶(bucket)中。在遍歷該特徵時,只需遍歷各個分位點,對每一個桶內的樣本統計值 \(g\)、\(h\) 進行累加統計,尋找最佳分裂點進行分裂。該算法可分爲 global 近似和 local 近似,global 就是在新生成一棵樹以前就對各個特徵計算分位點並劃分樣本,以後在每次分裂過程當中都重複利用這些分位點進行劃分,而 local 就是每一次結點分裂時都要從新計算分位點。
最後總結一下 XGBoost 與傳統 GBDT 的不一樣之處:
傳統 GBDT 在優化時只用到一階導數信息,XGBoost 則對目標函數進行了二階泰勒展開,同時用到了一階和二階導數。另外 XGBoost 工具支持自定義損失函數,只要函數可一階和二階求導。
XGBoost 在損失函數中加入了正則化項,用於控制模型的複雜度,防止過擬合,從而提升模型的泛化能力。
傳統 GBDT 採用的是均方偏差做爲內部分裂的增益計算指標(由於用的都是迴歸樹),而 XGBoost 使用的是通過優化推導後的式子,即式 \((1.6)\) 。
XGBoost 借鑑了隨機森林的作法,支持列抽樣,不只能下降過擬合,還能減小計算量,這也是 XGBoost 異於傳統 GBDT 的一個特性。
XGBoost 添加了對稀疏數據的支持,在計算分裂增益時不會考慮帶有缺失值的樣本,這樣就減小了時間開銷。在分裂點肯定了以後,將帶有缺失值的樣本分別放在左子樹和右子樹,比較二者分裂增益,選擇增益較大的那一邊做爲默認分裂方向。
並行化處理:因爲 Boosting 自己的特性,沒法像隨機森林那樣樹與樹之間的並行化。XGBoost 的並行主要體如今特徵粒度上,在對結點進行分裂時,因爲已預先對特徵排序並保存爲block 結構,每一個特徵的增益計算就能夠開多線程進行,極大提高了訓練速度。
傳統 GBDT 在損失再也不減小時會中止分裂,這是一種預剪枝的貪心策略,容易欠擬合。XGBoost採用的是後剪枝的策略,先分裂到指定的最大深度 (max_depth) 再進行剪枝。並且和通常的後剪枝不一樣, XGBoost 的後剪枝是不須要驗證集的。 不過我並不以爲這是「純粹」的後剪枝,由於通常仍是要預先限制最大深度的呵呵。
說了這麼多 XGBoost 的優勢,其固然也有不完美之處,由於要在訓練以前先對每一個特徵進行預排序並將結果存儲起來,對於空間消耗較大。另外雖然相比傳統的 GBDT 速度是快了不少,但和後來的 LightGBM 比起來仍是慢了很多,不知之後還會不會出現更加快的 Boosting 實現。
XGBoost 推導的關鍵一步是二階泰勒展開,這裏做一下拓展。泰勒公式的簡單解釋就是用多項式函數去逼近原函數,由於用多項式函數每每求解更加容易,而多項式中各項的係數則爲原函數在某一點的n階導數值除以n階乘。這裏力薦 3Blue1Brown 關於泰勒公式的視頻 微積分的本質 - 泰勒級數 ,講得很是形象 。
已知函數 \(f(x)\) 在 \(x=x_0\) 處n階可導,那麼:
\[ \begin{aligned} f(x) &= f(x_0) + f'(x_0)(x-x_0) + \frac{f''(x_0)}{2!}(x-x_0)^2 + \cdots + \frac{f^{(n)}(x_0)}{n!}(x-x_0)^n \\ & = \sum\limits_{n=0}^{\infty}\frac{f^{(n)}x_0}{n!}(x - x_0)^n \end{aligned} \]
例如,\(x_0 = 0, \;\; f(x) = e^x\)時,\(e^x = \sum\limits_{n=0}^{\infty}\frac{x^n}{n!} = 1+x+\frac{x^2}{2!} + \frac{x^3}{3!} + \cdots\)
在機器學習中泰勒公式的一個典型應用就是牛頓法。在牛頓法中,將損失函數 \(L(\theta)\) 在 \(\theta^{k}\) 處進行二階泰勒展開(考慮 \(\theta\) 爲標量):
\[ L(\theta) =L(\theta^k) + L'(\theta^k)(\theta - \theta^k) + \frac{L''(\theta^k)}{2}(\theta - \theta^k)^2 \]
將一階導和二階導分別記爲 \(g_k\) 和 \(h_k\),爲使上式最小化,令
\[ \frac{\partial\left[g(\theta - \theta^k) + \frac{h}{2}(\theta - \theta^k)^2\right]}{\partial \,\theta} = 0, \qquad 求得\; \theta-\theta^k = -\frac{g_k}{h_k}\;, \]
則參數更新公式爲 \(\theta^{k+1} = \theta^k - \frac{g_k}{h_k} \;\),推廣到向量形式則爲 \(\boldsymbol{\theta}^{k+1} = \boldsymbol{\theta}^k - \mathbf{H}^{-1}_k\mathbf{g}_k\;\),\(\mathbf{H}_k\) 爲海森矩陣(Hessian matrix),\(\mathbf{g}_k\) 爲梯度向量:
\[ \nabla_\theta L = \mathbf{g}_k = \begin {bmatrix} \frac{\partial L}{\partial \theta_1}\\ \frac{\partial L}{\partial\theta_2}\\ \vdots \\ \frac{\partial L}{\partial\theta_n} \end {bmatrix} \;\;,\;\; \nabla_{\theta}^2 L = \mathbf{H}_k = \begin {bmatrix} \frac{\partial^2 L}{\partial\, \theta_1^2} & \frac{\partial^2 L}{\partial\,\theta_1 \partial\,\theta_2} & \cdots & \frac{\partial^2 L}{\partial\,\theta_1 \partial\,\theta_n} \\ \frac{\partial^2 L}{\partial\,\theta_2 \partial\,\theta_1} & \frac{\partial^2 L}{\partial\, \theta_2^2} & \cdots & \frac{\partial^2 L}{\partial\,\theta_2 \partial\,\theta_n} \\ \vdots & \vdots & \ddots \\ \frac{\partial^2 L}{\partial\,\theta_n \partial\,\theta_1} & \frac{\partial^2 L}{\partial\,\theta_n \partial\,\theta_2} & \cdots & \frac{\partial^2 L}{\partial\, \theta_n^2} \end {bmatrix} \]
/