左手论文 右手代码 深入理解网红算法XGBoost

本文经授权转载自 知乎@贾茹，原文链接： https://zhuanlan.zhihu.com/p/91817667，也可点击【阅读原文】

以下为正文

这篇文章的起源是上周给同事讲解XGBoost原理+白板手推公式，结果差点翻车，发现有些地方理解还是不够深入，于是又花了一周时间把之前啃过的所有资料和笔记仔细梳理一遍，又想起自己的知乎专栏已经快一年没更新了，所以打算趁着换工作的3天休息间隙完成这篇文章。

考虑到网上关于XGBoost资料笔记已经烂大街了，想要写出新意着实要费一番心思。 一番思考之后终于找到了一个勉强算是创新的点： 将原始论文中的 mathematicl notations 与 XGBoost文档中的类、方法、参数、返回值等一一对应起来。 一方面通过论文中数学推导弄明白XGB代码背后的数学原理，know what you are doing 而不仅仅是当调包侠； 另一方面通过调包调参、研究API设计、甚至是研究xgb源码，反过来加深对XGB底层原理的理解。

这篇文章假设读者已经对树模型和集成学习boosting方法有基本了解。

XGBoost模型与公式推导

mathematics ... is nothing more difficult than notations —— 忘了是谁说的了

首先要说明的一点是，XGBoost并不是一种独立的“算法”， 而是广义GBDT算法的一种高效实现， X 是 extreme 的意思，即梯度提升树的算法优化和高效系统实现。

应该说整个XGBoost模型推导过程中，除了复杂的数学符号表示之外没有任何难点。 首先是模型： 基学习器为树模型、采用boosting集成方法，模型为如下的 additive model:

其中   表示第k个基学习器 ，特别地，这里的基学习器为树模型。

xgboost程序包中，通用参数  booster='gbtree' 是 指定基学习器类型的参数。除了默认的 gbtree 选项之外，还有 dart 和 gblinear 两个选项，前者是改良的树模型；后者则是线性模型。以下推导我们只考虑基学习器为回归树模型的情况，基学习器为线性模型这种情况我自己也暂时没搞懂

接下来设定损失函数和优化目标，我们的优化目标是结构风险最小化：

先来看损失函数   。 我们不限定损失函数   的具体形式，只要其二阶可导。 这样我们不需要为每一个具体的损失函数单独推导一个模型，而是得到一个通用的模型（以及结论）。

相应地，xgboost程序包支持自定义损失函数，只需要给  xgb.train() 的 ob j  参 数传入一个 signature为  (ypred, y) -> (grad, hess)  的函数对象即可（注意不需要计算   函数本身的值，只需要能计算   的一阶、二阶导数值即可，后面的推导将会证明这一点）。当然xgboost包也内置了一系列常见的损失函数方便我们直接使用，包括 binary:logistic, reg:squarederror 等

再来看正则项   , 这里必须设定其为如下具体形式，注意这个具体形式是跟后面的二阶泰勒展开紧密联系的，如果没有这样的而具体形式，后面就无法推导出叶子节点的权重   的解析解:

加入正则项，是XGBoost相对于传统GBDT的一个改良。 后一项容易理解： 叶子节点的L2-norm，常见的惩罚项设置方法，不多说了； 而前一项   则是XGBoost的一个主要Contribution:  传统GBDT模型中没有正则项，而是通过剪枝等“启发式”(heuristic，其实就是拍脑袋)方法来控制模型复杂度。 后面的推导会证明，gamma参数的作用等同于预剪枝，这就将难以解释的启发式方法纳入到了标准的 supervised leanring 框架中，让我们更清楚的知道what you are learning

接来下的问题是： 如何得到   ？ 要知道这里的   不是数值而是函数，无法使用欧式空间数值优化方法。 这里我们采取的训练方式为分步算法 (trained in additive manner) 这是一种贪心算法： 站在第   轮 （因此   轮的训练结果   可以视为已知的），只考虑本轮的优化目标   ：

这一步大大简化了优化目标，从一次优化K个函数，变成了每一轮只需要优化一个函数。 但这一步仍然无法让我们解出   ，因为   的具体形式未知。 那么一个自然而然的想法就是： 对损失函数   ，在前一轮预测值   处进行泰勒展开，将第t轮预测值增量  作为自变量增量，从而将   从   中“拿出来”：

这里我稍微改了一下原论文中的notation，原论文公式如下：

我相信绝大多数人第一次看到原文中这一步推导时都是一脸懵逼的。 回忆一下高中学过的泰勒展开式 :

这里我们套用   为“损失函数   固定第1项位置参数（真实y label）时的偏函数，即关于   的一元函数； 展开点为 上一轮y预测值   ； 自变量增量为本轮新加入的y预测值增量 

对于每一个样本i，我们将损失函数在   处的一阶、二阶导数值分别记为 

并忽略可视为常数的第一项，即得到了原论文中的：

现在我们已经成功地将   从损失函数中“分离”出来了。 这种分离除了下面要看到的数学推导上的作用之外，在工程上还有一个好处， 那就是 “ isolate loss function from the boosting machine” 即工程实现上的模块化，这也是为什么我们能自定义损失函数的原因。

接下来神奇的两步让我们最终得到

step 1. 假定已知树结构，得到最优叶子节点权重的解析解 

step 2. 确定树结构: 将 的解析解代回   得到各特征分割点的 loss reduction (gain) 公式，在每一个叶子节点选择gain最大的（特征，分割点）tuple 作为树的下一步生长方向。

将   展开，得到下图中（4）式的第一行。 我们说： 对于任意一个样本 i, 一定存在一个叶子节点 j 使得

也可以这样说：    就是   ， 别忘了f() 是一颗树, 每一个样本 i 最终都要被映射到某一个叶子节点 j上, 所以我们可以重新组织求和号，从 sum by i 改写成 sum by j, ie. “group by leaf-node”， 同时将正则项中的   塞进中括号，与二阶导统计量放在一起。 注意每个叶子节点上所有样本的   是相同的，而每个样本的一阶、二阶导数值是不同的。

到这里我们发现第t步优化目标   已经完全化成了关于   的一元二次函数，而一元二次函数的最小值是有解析解的！ 直接套公式即得到上图 (5) 式叶子节点权重解。 （回忆初中学过的一元二次函数   最小值为   ）

现在的问题变成了： 如何得到树结构？ 这是一个NP-hard问题，无法在多项式时间内求解，因此我们采取贪心算法求解： 从根节点开始，遍历所有可能的 (feature, split)，找到使loss reduction最大的那个，作为树的下一步生长方式。

我们将上一步得到的   的解析解代回 (4) 式，得到优化目标   的值，也称作一个树结构的 score function:

对每一个潜在的分割，都可以计算式 (7) 中的 loss reduction, 选择最大的那个作为实际分割点。

注意式中最后的那个   ，进行分割的前提条件是 loss reduction 必须大于0，所以xgb工具包中的 gamma 参数有一个min_split_loss 的alias。

至此我们已经推完了XGBoost数学推导的主干部分。 剩下的诸如缺失值处理、近似分割点选择等等可以去看原始论文自行了解。 其实，理论推导部分到这里就可以结束了，但这里我们还是从头再梳理一遍，完整的走过整个模型训练流程，这样能更好的理解最终的模型是怎么来的。

下面是我自己写的伪代码，去掉诸如 early stopping, customized callback 之类的特性，只保留最主干的部分，boosting过程主循环的逻辑如下：

t = 0

fx[0] = params.base_socre # 第0轮，初始化f0()为param.base_score中设定的值，默认为0.5

yhat[0] = fx[0]

g[0], h[0] = obj(ytrue, yhat[0])

while t < num_boost_round:

t += 1

# tree structure & leaf node weight wj*

fx[t] = generate_tree_with_greedy_algrithom(g[t-1], h[t-1])

# update yhat of round t

yhat[t] = yhat[t-1] + fx[t]

# first and second order gradient statistics of round t

# herer `obj` is customized obj function with signature obj(ytrue, yhat) -> grad, hess

g[t], h[t] = obj(ytrue, yhat[t])

• 第0轮，将    设为 params['base_score'] 的值（默认为0.5， 这是因为默认objective为0-1二分类问题）

  • 得到第0轮的y预测值  

  • 计算第0轮的损失函数一阶二阶梯度统计量  

• 在第    轮

1. 根据上一轮的损失函数一阶二阶梯度统计量    使用贪心算法得到第t轮回归树   

2. 第t轮 y预测值  

3. 计算t轮损失函数梯度统计量  

• 直到第 num_boost_round 轮，停止。 或者满足early-stopping条件停止

画了一张ppt，看起来更直观：

上图表示boosting过程主循环发生的事情。 而其中在第t轮第1步，回归树    的生成过程是

• 在每一个叶子节点上，找出最优的（特征，分割点）

  • 对每一个特征及其可能的分割，计算 splite gain，注意 split gain 仅需要知道一阶二阶梯度统计量  

  • 选择使得 split gain 最大的那个，作为分割点

• 若 split gain 小于 gamma 或达到树的最大深度，停止分割

• 注意同一层不同节点、同一个节点的不同特征之间都可以并行，使用的也是同一个    数据

可以用如下图来表示：

现在我们做一下调包侠，看一下xbgoost包中的每种对象、函数、方法分别对应的是论文中的哪些部分，加深理解。

XGB工程实践与调参

xgboost包提供两种API: 原生的 Learning API 以及 Sklearn API，这里我们研究原生的Learning API。

Booster对象

xgboost的数据结构（类）只有两个：  DMarix 和 Booster 。  DMatrix是为XGB算法专门设计的、高效存储和计算的数据结构，具体的设计方法可以参考原始论文的第四章 System Design 这里就不展开了（其实是我看不懂:laughing:）

而Booster顾名思义是boosting集成模型对象。 可以通过原生的Learning API xgb.train() 得到。 直观上看，最后得到的模型长这样：

xgb.Booster

booster对象也有一个.trees_to_dataframe() 方法，将集成树模型的信息返回成一个DataFrame，包括树的分裂节点、每个节点的Gain, Cover （cover表示该节点影响的样本数量）

Booster.tree_to_dataframe()

注意叶子节点上的Gain字段值，表示的是叶子节点权重，即解析解   , 与plot_tree() 中叶子节点的值一致：

leaf node score == Wj*

我们将所有树相加（对于每个样本，在每个树的叶子节点score相加），再加上初始值   （对应的参数为base_score，默认值0.5），就得到了 bst.predict() 返回的值。

如果你看过XGBoost其他资料，会发现他们常常提到一个概念叫做 “函数空间的梯度下降” 。 这里我们仔细考察通过 ypred 是如何通过叶子节点score相加得到的，能加深对这个概念的理解。 每一次迭代，我们在上一次的预测基础上，加上一个新的函数/树   也就是   , 每一次迭代，我们调整的不再是某个参数向量   的值，而是调整一个函数曲线/曲面

Gradient Descent

Gradient Boosting

我们上述每一棵树的叶子节点值，都是对上一轮拟合的“修正”，因此将所有树加起来最后得到的就是最终预测函数

对于树模型，我们还可以考察 特征重要性 。 特征重要性有三种不同的计算方式：

• weight 注意这里的权重既不是 leaf-node-score, 也不是二阶导数   而是定义为“该特征用于分裂的次数” 。 （这里我也有个疑问请各位指教：XGB论文中多次出现各种"weights"，从定义上看是完全不同的东西，这个只是名称混用的问题，还是说这些weights之间有某种的联系？）

• gain  前面的公式中的 loss reduction

• cover  定义为”受其影响的样本数数量“

这里我个人认为无论是 gain 和 cover 都好于默认的 weight , 最简单例子，如果一个特征在每个树中都是root node的分裂特征，那显然这个特征是非常有用的，但按照定义它的weight却不一定高，这显然不合理。所以实践中我一般是用更合理、同时解释性也更高的 cover

XGB调参

XGBoost需要调的参数不算多，大概有以下几类：

• num_boost_round  and  learning_rate

• 直接约束树结构的参数： max_depth ,  min_child_weight

• 正则化参数  gamma ,  lambda ,  alpha

• 随机性参数  subsample , colsample_by* .

• 其他参数：类别不平衡下使用  scale_pos_weight

树的个数与学习率： 这两个参数是互相联系的，关于学习率，paper中的解释是这样的：

其实我也不太理解这一段，跑去看了Friedman的论文，也没看懂为啥要乘以一个缩水值，但之一个直观的想法是学习率和树的个数应该是大致呈反比的，之前查资料也见过一种广泛采用的套路： 固定学习率参数，用 xgb.cv 找出当前学习率下的最优树个数。

下图是我做的一个实验： 其他参数不变，左图学习率=0.1，右图学习率=1； 横轴是迭代次数，利用自定义函数和callback记录下每一轮的 y预测值、损失函数梯度统计量、以及   。 从图中可以看出，当学习率较小时，模型是”逐渐“学习到最终的结果，而学习率=1在前几轮就已经充分学习了，剩下的时间都在以较大幅度”震荡“。

直接约束树结构的参数，最常用的就是树的深度max_depth 。 因为XGB所有基学习器都是二叉树，限制了树深度，也就相应的限制了叶子节点的个数，从而限制了最终预测值ypred 的uniqueN() 个数 （eg. max_depth=3, 每棵树叶子节点最多 2^3 = 8个，若num_boster=2则最终集成模型最多有 8^2 = 64个不同的预测值）

min_child_weight这个参数虽然直接调的不多，但它的推导和数学意义却很有趣。 它的含义是节点上所有样本的在当前模型输出处的二阶导数值之和   ，这个值必须大于一定的阈值才继续分裂。 但为什么是   呢？

回忆文章最开始部分的泰勒展开式是一个关于   的一元二次函数，运用初中的简单代数知识将其整理成一元二次函数的平方和形式：

会发现，（lambda=0情况下）我们在做的事情是 “以二阶导   为权重，以   为真实label，以RMSE最小化为目标，拟合   "。 而   的解，跟我们之前推导出来的   解析解也是一致的。 另外再多说一句，如果你之前熟悉牛顿法，会发现这个形式

跟牛顿法参数更新公式

有神奇的相似之处。 因此也有人把XBG这种方法称之为 Newton Boosting

随机性参数包括两种，一是关于样本的随机采样率 subsample 另一种是关于特征的随机采样率 colsample_by* ，思想来源于随机森林。 按论文中的说法， 经验上colsample比rowsample更能防止过拟合, 并且最好是用比较激进的colsample, ie < 0.5。 这种对样本或特征进行随机抽取的方法，思想上类似于数值空间最优化中的 stochastic gradient descent，于是有些资料将其称之为 ”stochastic gradient boosting“。 除了控制过拟合之外，在工程上，对样本或特征进行随机选取也能大大降低计算耗时。

总结一下，论文中的Notations与文档中params的对应关系：

Fancy features of XGB

XGB包有许多其他有趣的特性，这里限于篇幅（实在是写不动了o(╥﹏╥)o）只简单列举，有兴趣的可以自行探索：

1. XGB程序包中，每一轮迭代的 for loop 是显式地暴露在 native language (Python, R) 中，并且留下来用户自定义callback的接口，其中最有用的当属 early stopping。用户可以自定义callback 回调函数接收一个 CallbackEvn 作为参数，而 CallbackEnv 是一个NamedTuple()，包含每一轮的booster对象、迭代轮数、CVfold、evaluation_result 等

2. XGB可以用来训练随机森林。随机森林模型与boosting tree 形式上完全相同，只需要设置三个参数 num_parallel_trees=100 ，  num_boosting_round=1 , learning_rate=1 ， colsample_by_* < 1 即可

3. 特征交叉。 可以设置哪些特征之间可以做交叉，哪些不可以

参考文献

• XGB原始论文。 如果只看model前三个SECTION基本就够了。 section4是系统设计方面的东西，包括数据如何优化存储 ie DMatrix是怎么来的、特征排序、块计算等。

• 陈天奇会议演讲中的有些细节是ppt上没有的

  https://www.youtube.com/watch?v=Vly8xGnNiWs&list=PLSz28ynAy1RohdgsPfC4l4t3lHu863PGx&index=3&t=1786s

• A Gentle Introduction to GBDT 这篇文章从历史沿革的角度，讲解boosting算法是如何一步一步发展起来的，用一句话概括出每篇论文的contribution, 可以当文献综述看。

  https://machinelearningmastery.com/gentle-introduction-gradient-boosting-algorithm-machine-learning/

• 李航统计学习方法 Chapter 7 AdaBoost, Gradient Tree boosting

封面图来源： Photo by Radek Grzybowski on Unsplash