本文是对XGBoost的经典论文XGBoost: A Scalable Tree Boosting System的总结。XGBoost 的全称是 eXtreme Gradient Boosting。作者为华盛顿大学研究机器学习的大牛陈天奇，其最大的特点在于，它能够自动利用 CPU 的多线程进行并行，同时在算法上加以改进提高了精度。
本文更多关注算法，主要内容如下：
- XGBoost目标函数推导
- XGBoost防止过拟合策略
- 切分搜索算法
- XGBoost针对稀疏问题的处理
- 系统速度优化
XGboost目标函数推导
预测值
XGboost是基于回归树的集成学习模型。集成学习主要有两种方式：bagging和boosting。XGboost属于boosting, 其预测值通过对各弱分类器的结果求和得到。对于输入\(x_{i}\)： \[ \begin{split}\hat { y } _ { i } = \phi ( x _ { i } ) = \sum _ { k = 1 } ^ { K } f _ { k } ( x _ { i } ) , \quad f _ { k } \in F \qquad (1) \end{split} \] 其中，\(f _ {k}\)即为一个弱分类器，在XGboost中，弱分类器为 CART。F 为 CART 空间。
目标函数
为了学习得到上式中的K个\(f_{k}\), 我们需要通过最小化模型的目标函数来获得。目标函数通常由损失+正则项表示。损失指示离真实值的差距，正则防止过拟合。这里损失为\(\begin{split}\sum _{ i } l ( \hat y _ { i }, y)\end{split}\), 是个可微凸函数。正则项则是K个\(f _ {k}\)中的参数。由于XGBoost的\(f _ {k}\)是回归树，该目标函数如下：
\[ \begin{split}&O b j = \sum _ { i } l ( \hat { y } _ { i } , y _ { i } ) + \sum _ { k } \Omega ( f _ { k } ) \\ &where\quad\Omega ( f ) = \gamma T + \frac { 1 } { 2 } \lambda | w | ^ { 2 } \nonumber \end{split} \qquad (2) \] 其中，$( f ) $ 是惩罚项。\(f_{ k }\)做为回归树，影响其复杂度的参数为叶子节点的个数\(\, T \,\)及叶子节点的权重\(w\).
训练方法
由于\(f _ { k }\)是树，与其他模型\(f\) 可以由数值型向量表示不同，我们不能使用往常的优化算法，如SGD，来寻找\(f\)。这里的解决方案是Additive Training（Boosting）。
\[ \begin{split} &\hat { y } _ { i } ^ { ( 0 ) } = 0 \\ &\hat { y } _ { i } ^ { ( 1 ) } = f _ { 1 } ( x _ { i } ) = \hat { y } _ { i } ^ { ( 0 ) } + f _ { 1 } ( x _ { i } ) \\ &\hat { y } _ { i } ^ { ( 2 ) } = f _ { 1 } ( x _ { i } ) + f _ { 2 } ( x _ { i } ) = \hat { y } _ { i } ^ { ( 1 ) } + f _ { 2 } ( x _ { i } ) \\ & … \\ &\hat { y } _ { i } ^ { ( t ) } = \sum _ { k = 1 } ^ { t } f _ { k } ( x _ { i } ) = \hat { y } _ { i } ^ { ( t - 1 ) } + f _ { t } ( x _ { i } ) \end{split} \] 其中，每一轮添加的\(f\) 是最优化目标函数的结果。
目标函数近似形式
根据Additive Training方式，在第\(t\) 轮，\(\hat { y } _ { i } ^ { ( t ) } = \hat { y } _ { i } ^ { ( t - 1 ) } + f _ { t } ( x _ { i } )\)。 由此在第\(t\)轮的目标函数可以改写为：
\[ \begin{split} O b j ^ { ( t ) } &= \sum _ { i = 1 } ^ { n } l ( y _ { i } , \hat { y } _ { i } ^ { ( t ) } ) + \sum _ { i = 1 } ^ { t } \Omega ( f _ { i } ) \\ &= \sum _ { i = 1 } ^ { n } l ( y _ { i } , \hat { y } _ { i } ^ { ( t - 1 ) } + f _ { t } ( x _ { i } ) ) + \Omega ( f _ { t } ) + constant \end{split} \] 这个目标函数仍然很难求最优。我们用泰勒展开取其近似如下：
\[ O b j ^ { ( t ) } \simeq \sum _ { i = 1 } ^ { n } [ l ( y _ { i } , \hat { y } _ { i } ^ { ( t - 1 ) } ) + g _ { i } f _ { t } ( x _ { i } ) + \frac { 1 } { 2 } h _ { i } f _ { t } ^ { 2 } ( x _ { i } ) ] + \Omega ( f _ { t } ) + constant \] 其中，\(g _ { i } = \partial _ { \hat { y } ^ { ( t - 1 ) } } l ( y _ { i } , \hat { y } ^ { ( t - 1 ) } ) , \quad h _ { i } = \partial _ { \hat { y } ^ { ( t - 1 ) } } ^ { 2 } l ( y _ { i } , \hat { y } ^ { ( t - 1 ) } )\), \(g _ {i}\), \(h _ {i}\)分别为 \(l ( y _ { i } , \hat { y } _ { i } ^ { ( t - 1 ) } )\)对\(\hat y ^ {(t - 1)}\)求偏导的一阶和二阶导数。把上式中的常数项去掉，得：
\[ \tilde { O b j ^ { ( t ) } } = \sum _ { i = 1 } ^ { n } [ g _ { i } f _ { t } ( x _ { i } ) + \frac { 1 } { 2 } h _ { i } f _ { t } ^ { 2 } ( x _ { i } ) ] + \Omega ( f _ { t } ) \qquad (3) \] 这里要注意理解上式是针对每一个样本\(x_{i}\)求损失，然后求所有样本损失和。我们将要在接下来的计算中为了方便计算\(\tilde { O b j ^ { ( t ) } }\)的最优值而换一种对这波操作的形式化表达方式。
目标函数求最优
这里定义 \(I _ { j } = \{ i | q ( x _ { i } ) = j \}\)。\(I_{j}\)为被树分配到叶子节点\(j\) 上的样本集。\(q ( x _ { i } )\)为样本\(x _ {i}\)通过q这个树结构到达叶子节点\(j\)。 \(w_ {j}\)为第\(j\)个叶子节点的权重。则 \(f _ { t } ( x _ { i } ) = w _ { q ( x _ { i } ) }\)。因此式 \((3)\) 可以表示为： \[ \begin{split} O b j ^ { ( t ) } &\simeq \sum _ { i = 1 } ^ { n } [ g _ { i } f _ { t } ( x _ { i } ) + \frac { 1 } { 2 } h _ { i } f _ { t } ^ { 2 } ( x _ { i } ) ] + \Omega ( f _ { t } ) \\ &= \sum _ { i = 1 } ^ { n } [ g _ { i } w _ { q ( x _ { i } ) } + \frac { 1 } { 2 } h _ { i } w _ { q ( x _ { i } ) } ^ { 2 } ] + \gamma T + \lambda \frac { 1 } { 2 } \sum _ { j = 1 } ^ { T } w _ { j } ^ { 2 } \\ &= \sum _ { j = 1 } ^ { T } [ ( \sum _ { i \in I _ { j } } g _ { i } ) w _ { j } + \frac { 1 } { 2 } ( \sum _ { i \in I _ { j } } h _ { i } + \lambda ) w _ { j } ^ { 2 } ] + \gamma T \end{split} \] 中括号里为我们求最优解时最喜欢形式\(ax ^ 2 + bx\)。显然，当
\[ w _ { j } ^ { * } = - \frac { \sum _ { i \in I _ { j } } g _ { i } } { \sum _ { i \in I _ { j } } h _ { i } + \lambda } \] 时，得到最优的\(\tilde {O b j ^ { ( t ) }}\):
\[ \tilde { Obj } ^ { ( t ) } = - \frac { 1 } { 2 } \sum _ { j = 1 } ^ { T } \frac { ( \sum _ { i \in I _ { j } } g _ { i } ) ^ { 2 } } { \sum _ { i \in I _ { j } } h _ { i } + \lambda } + \gamma T \qquad (4) \] 至此，我们完成了求第\(t\)轮时目标函数的最优值。这个值可以用来衡量树结构的好坏。也就是说，判断一棵树的质量，我们只要算出\(g _ {i}\), \(h_{i}\)，然后代入式子\((4)\)即可算出。
防止过拟合策略
- 权重衰减（Shrinkage）: 在每次迭代中对树的每个叶子结点的分数乘上一个缩减权重η，这可以使得每一棵树的影响力不会太大，留下更大的空间给后面生成的树去优化模型。
- 列采样（Column subsampling）: 同对数据集做有放回的采样可以防止拟合一样，对特征做有放回的采样也可以防止过拟合。有两种方式：建树前采样和节点分裂之前采样。不仅如此，此方法还可以加快模型计算速度。
切分搜索算法
当我们把一个节点做了切分后，我们需要一个评估这个切分是好是坏的标准。一个朴素的思想就是切分前节点的目标函数的得分减去切分后左右两节点目标函数之和的得分，如果两者的差距较大，说明切分后的损失较小。计算公式如下： \[ \begin{split} L _ { s p l i t } &= L ^ { ( t ) } - ( L _ { l } ^ { ( t ) } + L _ { r } ^ { ( t ) } ) \\ &= - \frac { 1 } { 2 } \frac { ( \sum _ { i \in I _ { j } } g _ { i } ) ^ { 2 } } { \sum _ { i \in I _ { j } } h _ { i } + \lambda } + \gamma - ( - \frac { 1 } { 2 } \frac { ( \sum _ { i \in I _ { I } } g _ { i } ) ^ { 2 } } { \sum _ { i \in I _ { l } } h _ { i } + \lambda } + - \frac { 1 } { 2 } \frac { ( \sum _ { i \in I _ { r } } g _ { i } ) ^ { 2 } } { \sum _ { i \in I _ { r } } h _ { i } + \lambda } + 2 \gamma )\\ &= \frac { 1 } { 2 } [ \frac { ( \sum _ { i \in I _ { L } } g _ { i } ) ^ { 2 } } { \sum _ { i \in I _ { L } } h _ { i } + \lambda } + \frac { ( \sum _ { i \in I _ { R } } g _ { i } ) ^ { 2 } } { \sum _ { i \in I _ { R } } h _ { i } + \lambda } - \frac { ( \sum _ { i \in I } g _ { i } ) ^ { 2 } } { \sum _ { i \in I } h _ { i } + \lambda } ] - \gamma \end{split} \qquad (5) \]
贪心算法
遍历所有特征的所有可能切分点，计算（5）式，取得分最高的切分点进行切分。
近似算法
近似算法按特征分布的百分位提出候选切分点，相当于将特征的特征值分桶后再遍历。按提出候选切分点的时间不同，有两种切分的方式：全局切分和局部切分。
全局切分就是在建树之初就提出了每个特征的候选切分点，之后对树的所有节点做切分时都从这些节点对应的特征候选切分点中选择。
局部切分是在每次节点切分后优化提出候选切分点。候选切分点提出较全局切分频繁，较全局切分需要较少切分点。适合生成较深的树。
全局切分候选切分点足够多时，可以达到和局部切分差不多的效果。
加权分位点
加权分位点是近似算法的一种。上面提到的近似算法通常是将切分点设在使数据均匀分布的位置。而加权分位点则是将切分点设在使分桶后的数据对loss产生均衡影响的位置。如果分4个桶，加权分位点使得每个桶内的样本集对loss的影响持平为25%。
样本对loss的影响计算方法如下： \[ \begin{split} L ^ { ( t ) } &\simeq \sum _ { i = 1 }^ { n } [ l ( y _ { i } , \hat { y } _ { i } ^ { ( t - 1 ) } ) + g _ { i } f _ { t } ( x _ { i } ) + \frac { 1 } { 2 } h _ { i } f _ { t } ^ { 2 } ( x _ { i } ) ] + \Omega ( f _ { t } ) + constant \\ &= \sum _ { i = 1 } ^ { n } [ \frac { 1 } { 2 } h _ { i } \cdot \frac { 2 \cdot g _ { i } f _ { t } ( x _ { i } ) } { h _ { i } } + \frac { 1 } { 2 } h _ { i } \cdot f _ { t } ^ { 2 } ( x _ { i } ) ] + \Omega ( f _ { t } ) + constant \\ &= \sum _ { i = 1 } ^ { n } \frac { 1 } { 2 } h _ { i } [ 2 \cdot \frac { g _ { i } } { h _ { i } } \cdot f _ { t } ( x _ { i } ) + f _ { t } ^ { 2 } ( x _ { i } ) ] + \Omega ( f _ { t } ) + constant \\ &= \sum _ { i = 1 } ^ { n } \frac { 1 } { 2 } h _ { i } [ ( 2 \cdot \frac { g _ { i } } { h _ { i } } \cdot f _ { t } ( x _ { i } ) + f _ { t } ^ { 2 } ( x _ { i } ) + ( \frac { g _ { i } } { h _ { i } } ) ^ { 2 } ) - ( \frac { g _ { i } } { h _ { i } } ) ^ { 2 } ] + \Omega ( f _ { t } ) + constant \\ &= \sum _ { i = 1 } ^ { n } \frac { 1 } { 2 } h _ { i } [ ( f _ { t } ( x _ { i } ) + \frac { g _ { i } } { h _ { i } } ) ^ { 2 } ] + \Omega ( f _ { t } ) + constant \\ &= \sum _ { i = 1 } ^ { n } \frac { 1 } { 2 } h _ { i } ( f _ { t } ( x _ { i } ) - ( - g _ { i } / h _ { i } ) ) ^ { 2 } + \Omega ( f _ { t } ) \end{split} \] 这里和原文给出的公式差个负号，欢迎想明白的小伙伴留言。
从上式中我们可以看到样本\(x _ {i}\)的预测值\(f _ { t } ( x _ { i })\)对 \(y = - g _ { i } / h _ { i }\)的加权平方损失。权重即为\(h _ {i}\)，表示样本对loss的影响。
节点上特征值小于 \(z\) 的样本集相对于节点上总样本集对loss影响的比重为： \[ r _ { k } ( z ) = \frac { 1 } { \sum _ { ( x , h ) \in D _ { k } } h } \sum _ { ( x , h ) \in D _ { k } , x < z } h \] 我们的目标是要找到这样一组切分点\(\{ s _ { k 1 } , s _ { k 2 } , \cdots s _ { k l } \}\), 使用任意两个切分点之间的样本集对loss影响的比重小于一个数值\(\epsilon\)。比如 \(\epsilon = 10\%\)，意思就是将特征值分为10桶，每桶内对应的样本集对loss的影响小于10%。形式化表示如下： \[ | r _ { k } ( s _ { k , j } ) - r _ { k } ( s _ { k , j + 1 } ) | < \epsilon , \quad s _ { k 1 } = \min _ { i } x _ { i k } , s _ { k l } = \max _ { i } x _ { i k } \]
稀疏问题的处理方法
在实践中经常遇到的数据稀疏，某些样本的某些特征值是缺失的（NULL）或大量特征的特征值为0。对这种情况论文提出统一将这种数据划分到一个默认的分支，或左树或右树，然后分别计算损失，选择较优的那一个。
XGBoost系统优化
针对并行计算的列块
分
block，便于并行化分列存放，便于按列切片，实施 "列采样"
按照特征值事先排好序，便于执行切分
- 在乱序列表上查一个数，只能是线性扫描，复杂度 \(O ( N )\)
- 在有序列表上可以做二分查找， \(O ( \log N )\)
针对 cpu缓存优化
- 贪心算法：缓存预取缓解缓存未命中问题
- 近似算法：减小块的大小来解决缓存未命中问题
针对IO优化
- Block Compression：block按列压缩，在加载到内存的过程中完成解压缩。
- Block Sharding：将数据分片到多个磁盘上。 将预取线程分配给每个磁盘，并将数据取入内存缓冲区中。 然后，训练线程可以从每个缓冲区读取数据。 当有多个磁盘可用时，这有助于提高磁盘读取的吞吐量。
总结
个人觉得XGBoost在面对目标函数很难被计算时的思路很有些意思。在看论文的过程中几次被这种奇思妙想击中，总结一下以便时不时回顾感叹一下：
- 惩罚项是棵树时的处理
- 泰勒展示恰到好处的应用
- \(I_{j}\) 的定义把\(Obj\) 神奇地转化成了二次表达式
此外，在给特征值分桶时，用对loss的影响做为标准也是个不错的想法，并且找到用\(h\)做为权重也是漂亮的一波操作。
以上是个人觉得这篇论文最为出彩的地方，欢迎小伙伴留言分享你的心得。