论文精读（一）——XGBoost：A Scalable Tree Boosting System

论文精读（一）——XGBoost：A Scalable Tree Boosting System

1. 前导知识：梯度提升树

提升树利用加法模型与前向分步算法实现学习的优化过程；

当损失函数是平方损失和指数损失函数时，每一步优化都是很简单的，但对于一般损失函数而言，往往每一步优化并不那么容易；

针对这一问题，Freidman提出梯度提升算法，利用最速下降法的近似方法，其关键是利用损失函数的负梯度在当前模型的值：
$$-[\frac{\partial L(y,f(x_i))}{\partial f(x_i)}{]}_{f(x)=f_{m-1}(x)}$$作为回归问题提升树算法中的残差的近似值，拟合一个回归树。

注：
为什么说损失函数关于当前模型的负梯度是残差的近似值？
残差是平方损失函数下的特例。对于其他损失函数，其负梯度就不是残差了。至于说负梯度是残差的近似值，可以利用以下过程去理解：
对于模型的损失函数：
$$L(y,f(x))=L(y,y^{t-1}+f_t(x))，f_t(x)为树函数。$$
将其在$L(y,y^{t-1})$邻域做二阶泰勒级数展开，有：
$$L(y,y^{t-1}+f_t(x))=L(y,y^{t-1})+g{f}_t(x)+\frac{1}{2}h{f}_t^2(x)$$
式中g和h是一阶和二阶偏导：
$$g=\frac{\partial L(y,y^{t-1})}{\partial y^{t-1}}$$
$$h=\frac{{\partial }^{2}L(y,y^{t-1})}{\partial y^{t-1}}$$
$f_t(x)$就是泰勒公式里常见的$(x-x_0)$。
   \;
为求极值，对展开后的式子求导。去掉常数$L(y,y^{t-1})$，容易得到:
$$g+h{f}_t(x)$$
让上式等于0， 得到$f_t(x)=-\frac{g}{h}$。
   \;
对于平方损失函数，已知$g=y^{t-1}-y$, 残差等于负梯度，二阶导$h=1$, 且没有更高阶导数，结果就是用负梯度拟合新的树函数。
   \;
对于其他损失函数，仍旧设h=1，舍弃高阶项，相当于拟合残差的近似值。
   \;
如果加上正则化项，变为$f_t(x)=-\frac{g}{h+\lambda }$，就是xgboost的实现啦。
   \;
参考知乎：苏克
链接：https://www.zhihu.com/question/60625492/answer/369366208

泰勒级数：
$$f(x)=f(x_0)+f^{\prime }(x_0)(x-x_0)+\frac{f^{\prime \prime }(x_0)}{2!}(x-x_0{)}^2+\frac{f^{(3)}(x_0)}{3!}(x-x_0{)}^3+\frac{f^{(4)}(x_0)}{4!}(x-x_0{)}^4+...$$
该级数称为$f$关于$x=x_0$的泰勒级数

2. 问题提出

根据原文的介绍：

XGBoost是一种可拓展的树提升系统，XGBoost的可扩展性(scalability)归因于一些重要的系统优化和算法优化。这些优化包括：

• 一种新的tree-learning算法(a novel tree learning algorithm)：用于处理稀疏数据(sparse data)；
• 一种具有合理理论支撑的加权直方图调整框架（a theoretically justified weighted quantile sketch procedure）：用于处理近似的树学习中的实例权重；

由于XGBoost的并行化和分布式计算，使得学习过程比其它模型实现要快。更重要地，XGBoost实现了核外计算(out-of-core computation: 基于外存)，使得数据科学家们可以在pc机上处理上亿的训练实例。最终，会把这些技术结合起来实现一个end-to-end的系统，可以扩展到集群上。

论文列举了文章最大的四点贡献：

• 1.设计和建立了一个高可扩展的端到端的树提升系统（end-to-end tree boosting system）；
• 2.提出了具有合理理论支撑的分布式加权直方图调整框架(weighted quantile sketch procedure)，用于高效地进行预计算；
• 3.介绍了一种新颖的并行自适应稀疏处理树学习算法（sparsity-aware algorithm）；
• 4.提出了一种高效的基于缓存块的结构（cache-aware block structure），用于外存树（out-of-core tree）的学习；

依据以上四点贡献，首先提出下面几个问题，读完论文能够回答这几个问题，则说明对论文有一定理解：

• 1.为什么说XGBoost是端到端的树提升系统，可拓展性体现在哪些方面？
• 2.如何理解分布式加权直方图？该过程为何可以高效预计算？
• 3.如何理解自适应稀疏处理算法？并行化是如何实现的？
• 4.什么是基于缓存块的结构？如何进行外存计算和分布式计算？

下面将结合论文原文做总结，最后回答以上问题，并比较XGBoost与之前的算法。

3. XGBoost算法

XGBoost的方法源自于Friedman的二阶方法，XGBoost在正则化目标函数上做了小改进。

3.1 Regularized Learning Objective

对于一个含n个训练样本，m个features的给定数据集： D = { ( x i , y i ) } ( ∣ D ∣ = n , x i ∈ R m , y i ∈ R ) D = \{(x_i, y_i)\}(|D|=n, x_i \in \R^m, y_i \in \R) D={(xi​,yi​)}(∣D∣=n,xi​∈Rm,yi​∈R)，所使用的树集成模型（tree ensemble model）使用K次求和函数来预测输出（加法模型）：
$$\begin{array}{r}{\hat{y}}_i=\varphi (x_i)=\sum _{k=1}^K{f}_k(x_i),f_k\in F\end{array}\text{\hspace{1em}(1)}$$
其中， F = { f ( x ) = w q ( x ) } ( q : R m → T , w ∈ R T ) F = \{f(x) = w_q(x)\}(q: \R^m \to T, w \in \R^T) F={f(x)=wq​(x)}(q:Rm→T,w∈RT)是回归树（CART）的空间。q表示每棵树的结构，它会将一个训练样本实例映射到相对应的叶子索引上；T是树中的叶子数；每个$f_k$对应于一个独立的树结构$q$和叶子权重$w$。

与决策树不同的是，每棵回归树包含了在每个叶子上的一个连续分值（梯度提升算法中的平均值c），我们使用$w_i$来表示第$i$个叶子上的分值。对于一个给定样本实例，我们会使用树上的决策规则(由q给定)来将它分类到叶子上，并通过将相应叶子上的分值(由$w$给定)做求和，计算最终的预测值。为了在该模型中学到这些函数集合，我们最小化下面的正则化目标函数：
$$\begin{array}{r}L(\varphi )=\sum _{i}l({\hat{y}}_i,y_i)+\sum _i\Omega (f_k)\end{array}\text{\hspace{1em}(2)}$$
其中：$\Omega (f)=\gamma T+\frac{1}{2}\lambda ||\omega |{|}^2$

其中，$l$是一个可微凸函数（differentiable convex loss function），可以计算预测值$\widehat{y_i}$与目标值$y_i$间的微分。第二项$\Omega$会惩罚模型的复杂度。正则项可以对最终学到的权重进行平滑，避免overfitting。相类似的正则化技术也用在RGF模型(正则贪婪树)上。XGBoost的目标函数与相应的学习算法比RGF简单，更容易并行化。当正则参数设置为0时，目标函数就相当于传统的gradient tree boosting方法。

注：
正则化项$\Omega (f)$为什么写成$\gamma T+\frac{1}{2}\lambda ||\omega |{|}^2$ ？
   \;

• 因为损失函数的第一项是预测值与实际值之间的偏差，在训练时，为了使模型拟合的更好，会尽可能减小这一部分，从而造成模型的过拟合。
• 为了避免过拟合，加入惩罚项$\Omega (f)$，惩罚项的第二部分就是L2正则化，直观上理解就是为使${\sum }_{i}l({\hat{y}}_i,y_i)$这一项更小，会将数据集中出现更多的样本的权重$w$加大，为了避免这些权重过大导致过拟合，损失函数中加入$||w|{|}^2$，对过大的$w$做惩罚，减小这些$w$，从而得到合适的权重。
• 惩罚项的第一项加入$T$，也就是树中的叶子数。直观上理解就是为使模型拟合的更好，会生成更复杂的决策树，因而树的叶子节点也就更多，这里对树的叶子节点做惩罚，限制叶子结点的个数，也就是限制树的复杂程度，当然可以一定程度上防止过拟合。

3.2 Gradient Tree Boosting

式（2）中的tree ensemble模型将函数作为参数，通常这是一个复杂的优化问题，不能使用在欧拉空间中的传统优化方法进行优化。模型以一种叠加的方式进行训练，也就是前向分步算法。

前向分步算法：
加法模型的损失函数极小化问题：
$$mi{n}_{{\gamma }_k}\sum _{i=1}^{N}L(y_i,\sum _{k=1}^{K}f(x_i,{\gamma }_k))$$
因为学习的是加法模型，如果能够从前到后，每一步只学习一个基函数及其参数，逐步逼近优化目标函数，那么就可以简化优化的复杂度。具体的，每步只需优化如下损失函数：
$$mi{n}_{\gamma }\sum _{i=1}^{N}L(y_i,f(x_i,\gamma ))$$

由（2）式，第一部分是训练误差，第二部分是每棵树的复杂度的和。

我们采用前向分步算法进行模型学习，每一次保留原来的模型不变，加入一个新的函数$f(x)$到我们的模型中，如下图。

那么，我们如何选择每一轮加入什么$f(x)$呢？
可以选取一个$f(x)$来使得我们的目标函数尽可能大地降低(加入$f(x)$后的预测值与实际值误差减少)。

这样我们引入$f_t(x_i)$，${\hat{y}}_i^t={\hat{y}}_i^{t-1}+f_t(x_i)$，其中${\widehat{y_i}}^{(t)}$为第i个实例在第t次迭代时的预测值，然后最小化下面的目标函数：
$$L^{(t)}=\sum _{i=1}^{n}l(y_i,{\hat{y}}_i^{t-1}+f_t(x_i))+\Omega (f_t)$$
对于$l$是平方误差时：

对于$l$不是平方误差的情况：
采用如下的泰勒展开定义一个近似的目标函数：

注：

• 其中，$f_t(x_i)$就是泰勒展开中的$\Delta x$，${\hat{y}}_i^{(t-1)}$是泰勒展开中的$x_0$。
• $g_i={\partial }_{{\hat{y}}^{(t-1)}}l(y_i,{\hat{y}}^{(t-1)}),h_i={\partial }_{{\hat{y}}^{(t-1)}}^{2}l(y_i,{\hat{y}}^{(t-1)})$分别是损失函数上的一阶梯度和二阶梯度。

移除常数项(真实值与上一轮的预测值之差)，目标函数只依赖于每个数据点的在误差函数上的一阶导数和二阶导数：
$$\begin{array}{r}{L}^{(t)}=\sum _{i=1}^n[g_i{f}_t(x_i)+\frac{1}{2}{h}_i{f}_t^2(x_i)]+\Omega (f_t)\end{array}\text{\hspace{1em}(3)}$$
定义$q$函数将输入$x$映射到某个叶节点上，则$f_t(x)=w_{q(x)}$（每一个样本所在叶子索引的分数 ），此外，定义每个叶子节点$j$上的样本集合为 I j = { i ∣ q ( x i ) = j } I_j= \{i|q(x_i)=j\} Ij​={i∣q(xi​)=j}，将(3)式进行重写，并展开$\Omega$项：

定义$G_j$(每个叶子节点里面一阶梯度的和)和$H_j$(每个叶子节点里面二阶梯度的和)：
$$G_j=\sum _{i\in I_j}{g}_i{H}_j=\sum _{i\in I_j}{h}_i$$
目标函数改写为：
$$\begin{array}{rl}{L}^{(t)}& =\sum _{i=1}^n[(\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\\ & =\sum _{i=1}^n[G_j{w}_j+\frac{1}{2}(H_j+\lambda )w_j^2]+\gamma T\end{array}\text{\hspace{1em}(4)}$$

因此，现在要做的是两件事：

• 确定树的结构, 这样，这一轮的目标函数就确定了下来。
• 求使得当前这一轮(第$t$轮)的目标函数最小的叶结点分数$w$。(Obj代表了当我们指定一个树的结构的时候，我们在目标上面最多减少多少，也称为结构分数，structure score）。

假设已经知道了树的结构，那么第2件事情是十分简单的，对于一个确定的结构$q(x)$，我们可以计算最优的权重$w_j^{\ast }$，对目标函数求偏导得到：
$$\begin{array}{r}{w}_j^{\ast }=-\frac{G_j}{H_j+\lambda }\end{array}\text{\hspace{1em}(5)}$$
将式（5）带入目标函数（3）中，得到：
$$\begin{array}{r}Obj=L^{(t)}(q)=-\frac{1}{2}\sum _{j=1}^T\frac{G_j^2}{H_j+\lambda }+\gamma T\end{array}\text{\hspace{1em}(6)}$$
Obj代表了当我们指定一个树的结构的时候，我们在目标上面最多减少多少，可叫做结构分数(structure score)，Obj计算示例：

我们只需要在每个叶子上对梯度和二阶梯度统计求和，然后应用得分公式（scoring formula）来获得质量分（quality score）。

确定树的结构

接下来要解决的就是上面提到的第一个问题，即如何确定树的结构。

暴力枚举所有的树结构，然后选择结构分数最小的。 树的结构太多了，这样枚举一般不可行。

通常采用贪心法，每次尝试分裂一个叶节点，计算分裂后的增益，选增益最大的。这个方法在之前的决策树算法中大量被使用。而增益的计算方式比如ID3的信息增益，C4.5的信息增益率，CART的Gini系数等。

回想式子（6），衡量了每个叶子节点对总体损失的贡献，我们希望目标函数越小越好，因此 $\frac{({\sum }_{i\in I_j}{g}_i{)}^2}{{\sum }_{i\in I_j}{h}_i+\lambda }$越大越好。

从而得到如下打分公式：

通常，不可能枚举所有可能的树结构$q$。而贪婪算法会从单个叶子出发，迭代添加分枝到树中。
假设$I_L$和$I_R$是一次划分(split)后的左节点和右节点所对应的实例集合。$I=I_L\bigcup I_R$，比较分割后的增益。利用这个打分函数来寻找出一个最优结构的树，加入到我们的模型中，再重复这样的操作。

引入新叶子的惩罚项

优化这个目标对应了树的剪枝， 当引入的分割带来的增益小于一个阀值的时候，我们可以剪掉这个分割。
这样根据推导引入了选取分裂节点的计算分数和叶子的惩罚项，替代了回归树的基尼系数与剪枝操作。

当我们正式地推导目标的时候，像计算分数和剪枝这样的策略都会自然地出现，而不再是一种因为heuristic(启发式)而进行的操作了。

3.3 Shrinkage and Column Subsampling

除了采用前面提到的正则化目标函数外，还会使用两种额外的技术来进一步阻止过拟合：

• 收缩(Shrinkage)：Shrinkage会在每一步树提升时，将新加入的weights通过一个因子$\eta$进行缩放。与随机优化中的学习速率相类似，对于用于提升模型的新增树(future trees)，shrinkage可以减少每棵单独的树和叶子空间（leaves space）的影响。
• 列特征子抽样(column feature subsampling)：该技术借鉴自随机森林。根据用户的反馈，比起传统的行子抽样（row sub-sampling：同样也支持），使用列子抽样可以阻止过拟合。列子抽样的使用可以加速并行算法的计算(后面会描述)。

4. 划分点查找算法

4.1 精确的贪心算法（exact greedy algorithm）

树学习算法中的一个关键问题就是找到最好的分裂点(best split)。为了达到这个目标，分裂算法会在所有特征上，枚举所有可能的划分，我们称它为“精确贪婪算法(exact greedy algorithm)”。该算法如下图所示。它会对连续型特征枚举所有可能的split，为了更高效，该算法必须首先根据特征值对数据进行排序，以有序的方式访问数据来枚举打分公式中的结构得分（structure score）的梯度统计(gradient statistics)。

算法原理：
对于每棵树m，针对样本的每个特征j，按照样本进行排序，然后计算该特征下的Score，选取最好的划分点，然后对每个节点循环执行贪婪算法的过程，得到决策树。

4.2 近似算法（approximate algorithm）

完全贪婪算法使得XGBoost的每一步都按照分裂后增益最大的分裂点进行分裂，而分裂点的选取是枚举所有分割点。当数据量十分庞大，以致于不能全部放入内存时，贪婪算法就会很慢，因此XGBoost引入了近似的算法，如下。

该算法首先会根据特征分布的百分位数(percentiles of feature distribution)，提出候选划分点(candidate splitting points)。接着，该算法将连续型特征映射到由这些候选点划分的分桶(buckets)中，聚合统计信息，基于该聚合统计找到proposal中的最优解。

简单的说，就是根据特征$k$的分布来确定$l$个候选切分点，然后根据这些候选切分点把相应的样本放入对应的桶中，对每个桶的$G,H$（一阶梯度，二阶梯度）进行累加。最后在候选切分点集合上采用贪心算法。

算法原理：
第一个for循环：对特征K按照特征分布的分位数找到切割点的候选集合 $S_k=s_{k1},s_{k2},...,s_{kl}$ ，这样做的目的是提取出部分的切分点不用遍历所有的切分点。其中获取某个特征K的候选切割点的方式叫proposal。主要有两种proposal方式：global proposal和local proposal。
第二个for循环：将每个特征的取值映射到由这些特征对应的候选点集划分的分桶(buckets)区间 { s k , v ≥ x j k > s k , v − 1 } \{s_{k,v} ≥ x_{jk} \gt s_{k,v−1}\} {sk,v​≥xjk​>sk,v−1​}中，对每个桶（区间）内的样本统计值 $G,H$进行累加统计，最后在这些累计的统计量上寻找最佳分裂点。这样做的主要目的是获取每个特征的候选分割点的$G,H$量。

给定了候选切分点后，计算打分函数的例子如下：

分桶当然是越多越好，但越多意味计算量也就更大，极端情况下就是一个数据一个桶，但这违背了分桶的初衷，所以分桶数量需要折中。

那么，现在有两个问题：

• 如何选取候选切分点 S k = { s k 1 , s k 2 , … , s k l } S_k=\{s_{k1}, s_{k2}, \dots,s_{kl}\} Sk​={sk1​,sk2​,…,skl​} ？
• 什么时候进行候选切分点的选取？

对于候选切分点的选取时机，原文中给了两种方式：global proposal和local proposal。

• Global： 学习每棵树前， 提出候选切分点
• Local： 每次分裂前， 重新提出候选切分点

原文对两种方法做了比较，最终得出的结论是：

• 由于global方式后续不再改变，因此往往需要更多的候选点；local方式不用像global方式那样设定分桶分的很细，但是每次分裂前都要重新提出候选切分点，因此带来的复杂度也更大。通常情况下，都会优先选用global的形式。
• 全局切分点的个数足够多时，和Exact greedy算法性能相当。
• 局部切分点个数不需要那么多，因为每一次分裂都重新进行选择。

对于如何选取候选切分点，可以采用分位数，也可以直接构造梯度统计的近似直方图等。文中介绍了加权直方图(Weighted Quantile Sketch)方法。

4.3 加权分位法 / 直方图(Weighted Quantile Sketch)

直方图 (Quantile Sketch)：

原因：
采用直方图是因为，当一个序列无法全部加载到内存时，常常采用分位数的直方图近似计算分位点，以近似获取特定的查询。

方式：
使用随机映射 (Random projections) 将数据流投射在一个小的存储空间内作为整个数据流的概要，这个小空间存储的概要数据称为Sketch，可用于近似回答特定的查询。需要保留原序列中的最小值和最大值。

近似算法中的加权直方图 (Weighted Quantile Sketch)：

通常一个特征的百分位数可以被用来让候选在数据上进行均匀地分布。对特征$k$构造一个 multi-set： D k = { ( x 1 k , h 1 ) , ( x 2 k , h 2 ) , . . . , ( x n k , h n ) } D_k=\left\{ \left( x_{1k},h_1 \right),\left( x_{2k},h_2 \right),...,\left( x_{nk},h_n \right) \right\} Dk​={(x1k​,h1​),(x2k​,h2​),...,(xnk​,hn​)}， 其中$x_{ik}$表示样本$i$的特征$k$的取值，$n$是样本个数。

$D_k$表示第$k$个特征$featur{e}_k$与二阶偏导$H$之间的集合【每个训练实例的第$k$个特征值以及它的二阶梯度值统计。其中$h_i$表示第$i$个实例的第$k$个特征值对应的二阶梯度值统计，可看作第$i$个实例的第$k$个特征值的权重。

一般来说对于加权直方图的取值是根据 Rank 值进行的，Rank 计算公式如下：
$$r_k(z)=\frac{{\sum }_{(x,h)\in D_k,x<z}h}{{\sum }_{(x,h)\in D_k}h}$$

该Rank函数，输入为某个特征值z，计算的是该特征所有可取值中小于z的特征值的总权重占总的所有可取值的总权重和的比例，输出为一个比例值。

于是就使用下面这个不等式寻找候选分割点 { s k 1 , s k 2 , . . . , s k l } \{ s_{k1}, s_{k2}, ..., s_{kl} \} {sk1​,sk2​,...,skl​}：
$$|r_k(s_{k,j})-r_k(s_{k,j+1})|<ϵ,s_{k1}=\underset{i}{\min }{x}_{ik},s_{kl}=ma{x}_i{x}_{ik}$$

其中：$s_{k1}$是特征$k$的取值中最小的值$x_{ik},s_{kl}$是特征$k$的取值中最大的值$x_{ik}$，这是分位数略图的要求，需要保留原序列中的最小值和最大值。$ϵ$是一个近似比例，或者说是扫描步幅。可以理解为在特征$k$的取值范围上，按照步幅，挑选出特征 $ϵ$的取值候选点，组成候选点集。起初是从$s_{k1}$起，每次增加$ϵ\times (s_{kl}-s_{k1})$作为候选点，加入到候选集中。如此计算的话，这意味着大约是$\frac{1}{ϵ}$个候选点。

此时特征$k$的取值中${\min }_i{x}_{ik},{\max }_i{x}_{ik}$来自 multi-set $D_k$，对于$D_k$的数据集有两种定义：

• (1) 一开始选好，然后每次树切分都不变，也就是说是在总体样本里选 ${\min }_i{x}_{ik},{\max }_i{x}_{ik}$，这就是我们之前定义的global proposal；
• (2) 是树每次确定好切分点的分割后样本也需要进行分割，${\min }_i{x}_{ik},{\max }_i{x}_{ik}$来自子树的样本集$D_k$，这就是local proposal。

一个例子如下：

上图中，要切分为3个，总和为1.8，因此第1个在0.6处，第2个在1.2处。

$ϵ$是近似因子（approximation factor），直觉上，这意味着大约有$\frac{1}{ϵ}$个候选点。这里，每个数据点通过$h_i$加权。

为什么$h_i$可以表示权重呢？将前面泰勒二阶展开后的目标函数进行配方：
$$\begin{array}{rl}& \sum _{i=1}^N\left(g_i{f}_t({\mathbf{x}}_{\mathbf{i}})+\frac{1}{2}{h}_i{f}_t^2({\mathbf{x}}_{\mathbf{i}})\right)+\Omega (f_t)\\ & =\sum _{i=1}^N\frac{1}{2}{h}_i\left(2\frac{g_i}{h_i}{f}_t({\mathbf{x}}_{\mathbf{i}})+f_t^2({\mathbf{x}}_{\mathbf{i}})\right)+\Omega (f_t)\\ & =\sum _{i=1}^N\frac{1}{2}{h}_i\left(\frac{g_i^2}{h_i^2}+2\frac{g_i}{h_i}{f}_t({\mathbf{x}}_{\mathbf{i}})+f_t^2({\mathbf{x}}_{\mathbf{i}})\right)+\Omega (f_t)\\ & =\sum _{i=1}^N\frac{1}{2}{h}_i{\left(f_t({\mathbf{x}}_{\mathbf{i}})–(-\frac{g_i}{h_i})\right)}^2+\Omega (f_t)\end{array}\text{\hspace{1em}(7)}$$

推导第三行可以加入$\frac{g_i^2}{h_i^2}$是因为$g_i$和$h_i$是上一轮的损失函数求导，是常量，对优化没有影响。

从式（7）可以看出，就像是标签为$-g_i/h_i$，权重为$h_i$的平方损失，因此用$h_i$加权。

4.4 稀疏自适应分裂算法（Sparsity-aware Split Finding）

在许多现实问题中，输入x是稀疏的。有多种可能的情况造成稀疏：

• 数据中的missing values
• 统计中常见的0条目
• 特征工程：比如one-hot encoding

让算法意识到数据中的稀疏模式很重要。为了这么做，我们提出了在每个树节点上增加一个缺省的方向（default direction），如下图所示。当稀疏矩阵x中的值缺失时，样本实例被归类到缺省方向上。在每个分枝上，缺省方向有两种选择，最优的缺省方向可以从数据中学到。

算法如下图所示，关键的改进点是：只访问非缺失的条目$I_k$。上述算法会将未出现值（non-presence）当成是一个missing value，学到最好的方向来处理 missing values。

XGBoost能对缺失值自动进行处理，其思想是对于缺失值自动学习出它该被划分的方向（左子树or右子树）。

其实算法的实现过程很简单：

• 让特征k的所有缺失值的都到右子树，然后和之前的一样，枚举划分点，计算最大的gain
• 让特征k的所有缺失值的都到左子树，然后和之前的一样，枚举划分点，计算最大的gain

这样最后求出最大增益的同时，也知道了缺失值的样本应该往左边还是往右边。使用了该方法，相当于比传统方法多遍历了一次，但是它只在非缺失值的样本上进行迭代，因此其复杂度与非缺失值的样本成线性关系。论文在Allstate-10k数据集上实验，比传统方法快了50倍。

5. 系统设计（SYSTEM DESIGN）

5.1 适应于并行学习的列块（Column Block for Parallel Learning）

在建树的过程中，最耗时是找最优的切分点，而这个过程中，最耗时的部分是将数据排序。为了减少排序的时间，提出分片结构来存储数据。

• Block中的数据以稀疏格式CSC进行存储
• Block中的特征进行排序（不对缺失值排序）
• Block 中特征还需存储指向样本的索引，这样才能根据特征的值来取梯度。
• 一个Block中存储一个或多个特征的值
  
  可以看出，只需在建树前排序一次，后面节点分裂时可以直接根据索引得到梯度信息。

注：开始时候老是搞不明白：对于每个节点进行分裂时，都需要计算对应增益，而本轮树在不停“生长”，如何能够实现一次排序，以后就直接索引了呢？梯度左右梯度难道不会变化的吗？
事实上，先从下面 XGBoost 增益计算公式说起：
$$Gain=\frac{1}{2}[\underset{左子树分数}{\underbrace{\frac{G_L^2}{H_L+\lambda }}}+\underset{右子树分数}{\underbrace{\frac{G_R^2}{H_R+\lambda }}}–\underset{分裂前分数}{\underbrace{\frac{(G_L+G_R{)}^2}{H_L+H_R+\lambda }}}]–\underset{新叶节点复杂度}{\underbrace{\gamma }}$$
XGBoost 提出 block 结构，如上图，计算出各特征对应值下的 Gradient statistics of each example，然后依照特征值进行排序，这样，每次进行节点分裂时候，对于 exact greedy alrogithm 时，直接对每个特征从小到大的取值，选取左右增益进行统计，按照上式计算分割增益，由于由 black 结构，可针对 block 结构进行并行计算选定特征及特征值下对应的最大增益再比较选出最大增益对应特征及特征值完成节点分裂。对于 XGBoost 提出的 approximate algorithm 也是类似，依旧是针对各 block 进行各加权分桶完成梯度统计计算，选择出最大 gain 对应的分桶特征及分桶值。

事实上，需要明白的是，XGBoost 作为 GBDT 的升级版，虽然不是采用负梯度，而是采用了二阶泰勒展开，但是其依旧是 Boosting 架构，也就是说这里所求的一阶、二阶梯度均是针对之前已经建好的树模型即$f_1+f_2+\cdots +f_{t-1}$，与本轮树无关，故而本轮树进行节点分类时候各节点针对所有特征对应下的特征值进行求解梯度计算分裂增益的时候，其梯度值已经固定，故而 block 将各特征按照特征值进行排序，然后对该特征下特征值与每个样本的梯度统计建立索引。这样，当进行切分计算增益时， 只需将切分点前后梯度进行累计，带入计算公式即可求得对应特征、特征值下的切分增益，选出最大的增益，该增益对应的切分特征及特征值即为我们所求，至此，完成切分。

其本质是解耦了特征排序与梯度累积，将排序过程提前完成，这样计算增益的时候只需按照切分点进行统计即可；同时，收益于该 block 结构，其寻找切分点的过程可并行处理，这也是为什么一直说 XGBoost 是并行算法的原因！

时间复杂度分析：
d为树的最大深度，K为树的总树目。对于exact greedy algorithm，原始的稀疏感知算法的时间复杂度：
$$O(Kd|x{|}_{0}logn)$$

这里，我们使用$|x{|}_0$来表示在训练数据中未缺失条目（non-missing entries）的数目。其中$|x{|}_{0}logn$是特征排序时间复杂度，而$Kd$则是 树的总数目∗树的深度=counts of split finding，当不做任何优化时，对每棵树，每层，每个节点进行split finding时，对每个节点的所有特征对应的样本进行该特征下的大小排序，故而总的时间复杂度为上述步骤连乘。
而 $Kd|x{|}_0$则是进行split finding 的时间复杂度开销，最多进行$Kd|x{|}_0$次查找（对应每个样本为一个节点）。

另一方面，块结构上的tree boosting的开销为：
$$O(Kd{|x|}_0+{|x|}_{0}logn)$$
这里， $O({\Vert x\Vert }_{0}logn)$是一次预处理开销(one time preprocessing cost)(即排序开销)，可以分期(be amortized)。

因为采用了block structure之后，对于exact greedy algorithm 就是采用的是先对所有特征进行列排序，然后每次进行split finding时候，直接在block spliting 中进行查找线性遍历，得到各节点所有特征所有划分的对应增益，选择最大增益对应的特征及划分即为所求。
而$Kd|x{|}_0$则是进行split finding 的时间复杂度开销，最多进行$Kd|x{|}_0$次查找（对应每个样本为一个节点）。

该分析展示了块结构可以帮助节省一个额外的$logn$因子，其中当n非常大时就很大。

对于近似算法，使用二分查找的原始算法时间复杂度为：
$$O(Kd{|x|}_{0}logq)$$
这里的q是在数据集中建议候选的数目。其中，q通常为32~100之间，log因子仍会引入间接开销。使用块结构，我们可以将时间减小到：
$$(Kd{|x|}_0+{|x|}_{0}logB)$$
其中B是在每个块中的行的最大数。同样的，我们可以在计算中节约额外的$logq$因子。

解释同上，只是采用了近似算法后，采用了本论文提出的 Weighted Quantile Sketch 进行样本点候选筛选，这样即将原始样本点 n 降到候选点 q。

• 在Exact greedy算法中，将整个数据集存放在一个Block中。这样，复杂度从原来的$O(Hd||x|{|}_0\log n)$降为$O(Hd||x|{|}_0+||x|{|}_0\log n)$，其中$||x|{|}_0$为训练集中非缺失值的个数。这样，Exact greedy算法就省去了每一步中的排序开销。
• 在近似算法中，使用多个Block，每个Block对应原来数据的子集。不同的Block可以在不同的机器上计算。该方法对Local策略尤其有效，因为Local策略每次分支都重新生成候选切分点。

Block结构还有其它好处，数据按列存储，可以同时访问所有的列，很容易实现并行的寻找分裂点算法。此外也可以方便实现之后要讲的out-of score计算。

缺点是空间消耗大了一倍。

利用列块进行并行计算：在我们训练过程中我们主要是做分支处理，分支处理就要对每一列（特征）找出适合的分裂点。通常来说，我们更青睐使用csc存储，这样我们就方便取出来。再者我们在分支的时候都会预先对数据按照其特征值进行排序。所以我们将数据按照列存储成一个数据块方便我们在分支的时候并行处理。所以我们要知道XGB的并行计算的粒度不在树上，而是在特征上，尤其是不同分支节点上（leaf-wise）。当然这也成为XGB的一个问题所在，需要额外的空间存储pre-sort的数据。而且每次分支后，我们都要找处落在下一个子节点上的样本，并组织好它。 后来就有了LightGBM。

5.2 自适应缓存的访问（Cache-aware Access）

使用Block结构的一个缺点是取梯度的时候，是通过索引来获取的，而这些梯度的获取顺序是按照特征的大小顺序的。这将导致非连续的内存访问，可能使得CPU cache缓存命中率低，从而影响算法效率。

因此，对于exact greedy算法中, 使用缓存预取。具体来说，对每个线程分配一个连续的buffer，读取梯度信息并存入Buffer中（这样就实现了非连续到连续的转化），然后再统计梯度信息。该方式在训练样本数大的时候特别有用。

在approximate 算法中，对Block的大小进行了合理的设置。定义Block的大小为Block中最多的样本数。选择一个过小的block size会导致每个thread会小负载（small workload）运行，并引起低效的并行化(inefficient parallelization)。在另一方面，过大的block size会导致cache miss，梯度统计将不能装载到CPU cache中。block size的好的选择会平衡两者。设置合适的大小是很重要的，设置过大则容易导致命中率低，过小则容易导致并行化效率不高。经过实验，发现2^16比较好。

缓存处理能力：对于有大量数据或者说分布式系统来说，我们不可能将所有的数据都放进内存里面。因此我们都需要将其放在外存上或者分布式存储。但是这有一个问题，这样做每次都要从外存上读取数据到内存，这将会是十分耗时的操作。因此我们使用预读取（prefetching）将下一块将要读取的数据预先放进内存里面。其实就是多开一个线程，该线程与训练的线程独立并负责数据读取。此外，我还要考虑block的大小问题。如果我们设置最大的block来存储所有样本在k特征上的值和梯度的话，cache未必能一次性处理如此多的梯度做统计。如果我们设置过少block size，这样不能充分利用的多线程的优势，也就是训练线程已经训练完数据，但是prefetching thread还没把数据放入内存或者cache中。经过测试，作者发现block size设置为2^16个examples最好：

5.3 核外计算的块（Blocks for Out-of-core Computation）

当数据量太大不能全部放入主内存的时候，为了使得out-of-core计算成为可能，将数据划分为多个Block并存放在磁盘上。

• 计算的时候，使用独立的线程预先将Block放入主内存，因此可以在计算的同时读取磁盘
• Block压缩，貌似采用的是近些年性能出色的LZ4 压缩算法，按列进行压缩，读取的时候用另外的线程解压。对于行索引，只保存第一个索引值，然后用16位的整数保存与该block第一个索引的差值。
• Block Sharding， 将数据划分到不同硬盘上，提高磁盘吞吐率

块压缩（Block Compression）：块通过列(column)进行压缩，当加载进主存时可以由一个独立的线程即时解压(decompressed on the fly)。它会使用磁盘读开销来获得一些解压时的计算。我们使用一个通用目的的压缩算法来计算特征值。对于行索引（row index），我们从块的起始索引处开始抽取行索引，使用一个16bit的整数来存储每个偏移(offset)。这需要每个块有216216个训练样本，这证明是一个好的设置。在我们测试的大多数数据集中，我们达到大约26% ~ 29%的压缩率。

块分片（Block Sharding）：在多个磁盘上以一种可选的方式共享数据。一个pre-fetcher thread被分配到每个磁盘上，取到数据，并装载进一个in-memory buffer中。训练线程（training thread）接着从每个bufer中选择性读取数据。当提供多个磁盘时，这可以帮助增加磁盘读(disk reading)的吞吐量。

数据块以外的计算力提高：对于超大型的数据，我们不可能都放入放入内存，因此大部分都放入外存上。假如我们将数据存于外存上将给我们带来读写速度受限的问题。文中有两种方法，一种是对数据进行压缩存于外存中，到内存中需要训练时再解压，这样来增加系统的吞吐率，尽管消耗了一些时间来做编码和解码但还是值得的。**另一种就是多外存存储，其实本质上就是分布式存储。**这样说有多个线程对分布式结构管理，吞吐率自然高啦。

6. 总结

6.1 问题分析

本篇的开头提出了如下几个问题，读完整篇论文，下面来分析这几个问题：

（1）为什么说XGBoost是端到端（end-to-end）的树提升系统，可拓展性（scalable）体现在哪些方面？

end-to-end的好处：通过缩减人工预处理和后续处理，尽可能使模型从原始输入到最终输出，给模型更多可以根据数据自动调节的空间，增加模型的整体契合度。

从这个方面来分析，XGBoost是端到端的树提升系统，应该是体现在它的自适应上的，比如：稀疏自适应（Sparsity-aware）、自适应缓存的访问（Cache-aware Access）。Sparsity-aware可以处理稀疏数据，自动通过数据学习缺失值的默认划分方向；Cache-aware Access通过选择合适的块大小将梯度信息存入buffer，解决索引导致的非连续内存访问。通过自适应的方式减少人工处理，使模型根据数据自动调节，因此可以达到端到端的学习。

由于XGBoost只接受数值特征输入，不直接支持类别特征，因此对于原始数据，只需要将类别特征经过one-hot，然后就能将数据给XGBoost进行训练，因此基本可以实现端到端的学习。

可拓展性的体现：

• 模型的scalability：弱分类器可以支持cart也可以支持lr和linear， 但其实这是Boosting算法做的事情，XGBoost只是实现了而已。
• 目标函数的scalability： 支持不同的loss function, 支持自定义loss function，只要一、二阶可导。有这个特性是因为泰勒二阶展开，得到通用的目标函数形式。
• 学习方法的scalability：Block结构支持并行化，支持 Out-of-core计算。

（2）如何理解分布式加权直方图？该过程为何可以高效预计算？

加权直方图简单来说就是将二阶梯度作为每个数据点的权重，对数据点尽量均匀划分，找到候选点，将连续型特征映射到由这些候选点划分的分桶(buckets)中，聚合统计信息，这样按照特征对数据排序的时候，就可以将样本复杂度从原始样本点 n 降到候选点q，当然也就降低了时间复杂度。

（3）如何理解自适应稀疏处理算法？并行化是如何实现的？

自适应稀疏处理算法简单来说就是按照某个特征将缺失值同时划分到左子树和右子树，比较哪个方向的增益更大，就将该方向设为默认划分方向。

XGBoost并行计算的理解：

初始学习XGBoost，一个很容易困惑的点就是XGBoost的并行计算，毕竟其采用的加法模型架构，怎么看都是一种串行计算；事实上，XGBoost的并行不是在计算各轮添加的回归树上，而是体现在上述Split Finding中的排序上，即并行对各个样本，依据其所有m个特征值进行m轮从小到大排序，然后计算各个特征下各样本从左到右做划分对应的score。这才是并行计算的实质所在，由于XGBoost主要的计算量都击中在找Split Finding上，故而对这一步的并行计算大大加速了XGBoost的模型运算速度。

（4）什么是基于缓存块的结构？如何进行外存计算和分布式计算？

首先，缓存块是用于减少数据排序的时间，采用以下方式进行存储：缓存块中的数据以稀疏格式进行存储，一个Block中存储一个或多个特征的值，按照特征对数据进行排序，块中的特征需要存储指向样本的索引，这样可以根据特征的值来取梯度。

外存计算使用了两种方法：块压缩（Block Compression） 和块分片（Block Sharding） 。块压缩是采用不同的线程进行压缩和解压，减少读取数据的时间；块分片就是采用分布式存储，使用多个线程对分布式结构进行管理，提高数据吞吐率。

6.2 算法比较

6.2.1 XGBoost与传统GBDT

• 传统GBDT以CART作为基分类器，XGBoost还支持线性分类器，这个时候XGBoost相当于带L1和L2正则化项的Logistic回归（分类问题）或者线性回归（回归问题）。

注：
Boosting 模型本质还是加性模型，而XGBoost做为boosting算法的一种并行实现，除了对传统的GBDT进行了很多改进之外，又增加了其扩展性，其中一个就是支持不同的基分类器，不再局限于CART回归树；
&ThickSpace; \;
因为XGBoost可以自己设定损失函数，只需要满足该设定损失函数一阶、二阶可导，而损失函数旨在计算真实值与分类器预测值之间的差异。故而，当采用线型模型（包括线型回归和Logistic Regression）采用线性加和的方式进行预测${\hat{y}}_i={\sum }_j{w}_j{x}_{ij}$，这里的预测值 y yy 可以由不同的解释，比如我们可以把它作为回归目标的输出，或者进行 sigmoid 变换得到概率（即用$\frac{1}{1+e^{-{\hat{y}}_i}}$来预测正例的概率），或者作为排序的指标等。此时，基模型即采用的线型模型。
&ThickSpace; \;
此时，XGBoost的正则项$\Omega (f_t)=\gamma T+\frac{1}{2}\gamma {\sum }_{j=1}^T{w}_j^2$，此时$T$用于限制线型分类器个数，而$w_j^2$则是用于限制线型分类器预测分数（相较于传统的L2正则化，多了一项基分类器的个数正则项）。

• 传统的GBDT只用了一阶导数信息（使用牛顿法的除外），而XGBoost对损失函数做了二阶泰勒展开。并且XGBoost支持自定义损失函数，只要损失函数一阶、二阶可导。

• XGBoost的目标函数多了正则项， 相当于预剪枝，使得学习出来的模型更加不容易过拟合。

• XGBoost还有列抽样，进一步防止过拟合。

• 对缺失值的处理。对于特征的值有缺失的样本，XGBoost可以自动学习出它的分裂方向。

• XGBoost工具支持并行。当然这个并行是在特征的粒度上，而非tree粒度，因为本质还是boosting算法。

xgboost工具支持并行。boosting不是一种串行的结构吗?怎么并行的？注意xgboost的并行不是tree粒度的并行，xgboost也是一次迭代完才能进行下一次迭代的（第t次迭代的代价函数里包含了前面t-1次迭代的预测值）。xgboost的并行是在特征粒度上的。我们知道，决策树的学习最耗时的一个步骤就是对特征的值进行排序（因为要确定最佳分割点），xgboost在训练之前，预先对数据进行了排序，然后保存为block结构，后面的迭代中重复地使用这个结构，大大减小计算量。这个block结构也使得并行成为了可能，在进行节点的分裂时，需要计算每个特征的增益，最终选增益最大的那个特征去做分裂，那么各个特征的增益计算就可以开多线程进行。

6.2.2 XGBoost 与 线性回归

• 线性回归模型的解释性是决策树、随机森林、xgboost无法比拟的，也无法取代。

• 线性回归可以建立线性模型，而xgboost是不可以的。举个例子，即使是简单的$y=x+1$的线性关系，xgboost也无法做到。

• 线性模型计算简单，适用于快速部署。

• 决策树/森林回归是把数据空间非线性地分成$N$份（N为叶节点数），相当于做 N 聚类。然后再每个叶节点上求平均，其实是一种条件期望（conditional mean），$y^{=}E(y|c)$，$c$是聚类标签，可以看做是内插（interpolation）。所以预测点周围需要被训练点包围，这样才会有好的效果。

• 线性模型回归是用线 / 平面去拟合所有训练数据。当引入$x^2$等高维变量后，原有 x 空间是非线性曲线 / 曲面。线性模型回归可以作内插（interpolation，预测点曾经见过），也可以做外插（extrapolation，预测点从没见过，比如$y=x+1$）。

• 决策树优点：决策树用垂直线段去拟合曲线，不用显示的写出曲线的表达式。而线性模型会尝试写出曲线的表达式。所以线性模型更难去拟合复杂的曲线。

• 决策树缺点：决策树是基于数据（统计），而线性模型是基于模型（学习）。决策树需要大量数据去尽量覆盖所有可能的输入，适合于大数据。所谓的“大”是指训练样本的覆盖面要大。线性模型可以处理训练样本小的情况。

• 决策树和线性模型的组合。在每个叶节点上，单纯的决策树回归是用一个点代表所有值。此时可以再用线性模型，用一条线去拟合子空间的数据。比如：GBDT-LR，先采用GBDT构建高维组合特征，再采用LR训练模型。

• 最重要的一点是两者的训练方式不同，LR 模型是一个基模型，其可采用梯度下降法进行参数拟合，而 XGBoost 是一个组合模型，其是有多个基学习器组合（boosting）得到，因此，不能采用梯度下降法进行求解；

6.3 XGBoost文档

XGBoost Python Package

参考文献：

[1] XGBoost：A Scalable Tree Boosting System
[2] 《统计学习方法》 李航
[3] CTR预估 论文精读(一)–XGBoost
[4] xgboost公式推导
[5] 机器之心：XGBoost