xgboost 正则项_左手论文右手代码深入理解网红算法XGBoost

忘是亡心i 2023-01-04 15:54 145阅读 0赞

这篇文章的起源是上周给同事讲解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:

![Image 1][]

其中

![Image 1][] 表示第k个基学习器 ，特别地，这里的基学习器为树模型。

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

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

![Image 1][]

先来看损失函数

![Image 1][] 。我们不限定损失函数

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

*相应地，xgboost程序包支持自定义损失函数，只需要给 `xgb.train()`的`obj`参数传入一个 signature为 `(ypred, y) -> (grad, hess)` 的函数对象即可（注意不需要计算*

![Image 1][]

*函数本身的值，只需要能计算*

![Image 1][]

*的一阶、二阶导数值即可，后面的推导将会证明这一点）。当然xgboost包也内置了一系列常见的损失函数方便我们直接使用，包括`binary:logistic`, `reg:squarederror`等*

再来看正则项

![Image 1][] , 这里必须设定其为如下具体形式，注意这个具体形式是跟后面的二阶泰勒展开紧密联系的，如果没有这样的而具体形式，后面就无法推导出叶子节点的权重

![Image 1][] 的解析解:

![Image 1][]

加入正则项，是XGBoost相对于传统GBDT的一个改良。后一项容易理解：叶子节点的L2-norm，常见的惩罚项设置方法，不多说了；而前一项

![Image 1][] 则是XGBoost的一个主要Contribution:

**传统GBDT模型中没有正则项，而是通过剪枝等“启发式”(heuristic，其实就是拍脑袋)方法来控制模型复杂度。后面的推导会证明，gamma参数的作用等同于预剪枝，这就将难以解释的启发式方法纳入到了标准的 supervised leanring 框架中，让我们更清楚的知道what you are learning**

接来下的问题是：如何得到

![Image 1][] ？要知道这里的

![Image 1][] 不是数值而是函数，无法使用欧式空间数值优化方法。这里我们采取的训练方式为分步算法 (trained in additive manner) 这是一种贪心算法：

**站在第**

![Image 1][]

**轮 （因此**

![Image 1][]

**轮的训练结果**

![Image 1][]

**可以视为已知的），只考虑本轮的优化目标**

![Image 1][] ：

![Image 1][]

这一步大大简化了优化目标，从一次优化K个函数，变成了每一轮只需要优化一个函数。 但这一步仍然无法让我们解出

![Image 1][] ，因为

![Image 1][] 的具体形式未知。那么一个自然而然的想法就是：对损失函数

![Image 1][] ，在前一轮预测值

![Image 1][] 处进行泰勒展开，将第t轮预测值增量

![Image 1][]作为自变量增量，从而将

![Image 1][] 从

![Image 1][] 中“拿出来”：

![Image 1][]

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

![18bd660a3fe5806e6ec898300544f3b0.png][]

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

![Image 1][]

这里我们套用

![Image 1][] 为“损失函数

![Image 1][] 固定第1项位置参数（真实y label）时的偏函数，即关于

![Image 1][] 的一元函数； 展开点为 上一轮y预测值

![Image 1][] ； 自变量增量为本轮新加入的y预测值增量

![Image 1][]

对于每一个样本i，我们将损失函数在

![Image 1][] 处的一阶、二阶导数值分别记为

![Image 1][]

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

![a8808c069cd8a53cc746889e3caa163f.png][]

现在我们已经成功地将

![Image 1][] 从损失函数中“分离”出来了。这种分离除了下面要看到的数学推导上的作用之外，在工程上还有一个好处，

**那就是 “ isolate loss function from the boosting machine” 即工程实现上的模块化，这也是为什么我们能自定义损失函数的原因。**

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

![Image 1][]

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

![Image 1][]

step 2. 确定树结构: 将

![Image 1][]的解析解代回

![Image 1][] 得到各特征分割点的 loss reduction (gain) 公式，在每一个叶子节点选择gain最大的（特征，分割点）tuple 作为树的下一步生长方向。

将

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

![Image 1][]

也可以这样说：

![Image 1][]

**就是**

![Image 1][]

**， 别忘了f() 是一颗树, 每一个样本 i 最终都要被映射到某一个叶子节点 j上**, 所以我们可以重新组织求和号，从 sum by i 改写成 sum by j, ie. “group by leaf-node”， 同时将正则项中的

![Image 1][] 塞进中括号，与二阶导统计量放在一起。 注意每个叶子节点上所有样本的

![Image 1][] 是相同的，而每个样本的一阶、二阶导数值是不同的。

![aafe0275b854904092f216359312bb0d.png][]

到这里我们发现第t步优化目标

![Image 1][] 已经完全化成了关于

![Image 1][] 的一元二次函数，而一元二次函数的最小值是有解析解的！直接套公式即得到上图 (5) 式叶子节点权重解。（回忆初中学过的一元二次函数

![Image 1][] 最小值为

![Image 1][] ）

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

我们将上一步得到的

![Image 1][]

*的解析解代回 (4) 式，得到优化目标*

![Image 1][]

*的值，也称作一个树结构的 score function*:

![1ff7e46a39a9e26c3fa78a36a7202cbe.png][]

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

注意式中最后的那个

![Image 1][] ，进行分割的前提条件是 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轮，将
    
    ![Image 1][]  设为
    
    `params['base_score']` 的值（默认为0.5， 这是因为默认objective为0-1二分类问题）
    
     *  得到第0轮的y预测值 
        
        ![Image 1][]
     *  计算第0轮的损失函数一阶二阶梯度统计量 
        
        ![Image 1][]
 *  在第 
    
    ![Image 1][] 轮

1.  根据上一轮的损失函数一阶二阶梯度统计量
    
    ![Image 1][]  使用贪心算法得到第t轮回归树
    
    ![Image 1][] 
2.  第t轮 y预测值 
    
    ![Image 1][]
3.  计算t轮损失函数梯度统计量 
    
    ![Image 1][]

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

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

![d4c2b9061301b8ff023260f93c702306.png][]

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

![Image 1][] 的生成过程是

*  在每一个叶子节点上，找出最优的（特征，分割点）
    
     *  对每一个特征及其可能的分割，计算 splite gain，注意 split gain 仅需要知道一阶二阶梯度统计量 
        
        ![Image 1][]
     *  选择使得 split gain 最大的那个，作为分割点
 *  若 split gain 小于 `gamma` 或达到树的最大深度，停止分割
 *  注意同一层不同节点、同一个节点的不同特征之间都可以并行，使用的也是同一个
    
    ![Image 1][]  数据

可以用如下图来表示：

![70d6cb42f1b6d8effec9f2ce25185ad0.png][]

--------------------

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

## XGB工程实践与调参 ##

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

### Booster对象 ###

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

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

![55ac323bea85c581ea6cb98e6513f86a.png][]

xgb.Booster

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

![e59d475d328c9ea36914a6fe2e358470.png][]

Booster.tree\_to\_dataframe()

注意叶子节点上的Gain字段值，表示的是叶子节点权重，即解析解

![Image 1][] , 与

`plot_tree()` 中叶子节点的值一致：

![43607520325775438ebc7702d2ce356e.png][]

leaf node score == Wj\*

我们将所有树相加（对于每个样本，在每个树的叶子节点score相加），再加上初始值

![Image 1][] （对应的参数为

`base_score`，默认值0.5），就得到了 `bst.predict()` 返回的值。

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

![Image 1][] 也就是

![Image 1][] , 每一次迭代，我们调整的不再是某个参数向量

![Image 1][] 的值，而是调整一个函数曲线/曲面

Gradient Descent

![Image 1][]

Gradient Boosting

![Image 1][]

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

![Image 1][]

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

*  `weight` **注意这里的权重既不是 leaf-node-score, 也不是二阶导数** 
    
    ![Image 1][]
    
     **而是定义为“该特征用于分裂的次数”。** （这里我也有个疑问请各位指教：XGB论文中多次出现各种"weights"，从定义上看是完全不同的东西，这个只是名称混用的问题，还是说这些weights之间有某种的联系？）
 *  `gain` 前面的公式中的 loss reduction
 *  `cover` 定义为”受其影响的样本数数量“

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

![4327d5d9f2231024ea8a271bc5b1336e.png][]

### XGB调参 ###

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

*  `num_boost_round` and `learning_rate`
 *  直接约束树结构的参数：`max_depth`, `min_child_weight`
 *  正则化参数 `gamma`, `lambda`, `alpha`
 *  随机性参数 `subsample, colsample_by*`.
 *  其他参数：类别不平衡下使用 `scale_pos_weight`

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

![bb1f03f7cd3d5f6381a865ab1050f55b.png][]

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

下图是我做的一个实验：其他参数不变，左图学习率=0.1，右图学习率=1；横轴是迭代次数，利用自定义函数和callback记录下每一轮的 y预测值、损失函数梯度统计量、以及

![Image 1][] 。从图中可以看出，当学习率较小时，模型是”逐渐“学习到最终的结果，而学习率=1在前几轮就已经充分学习了，剩下的时间都在以较大幅度”震荡“。

![70f77a4026e86af928084496ea5163a4.png][]

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

`min_child_weight`这个参数虽然直接调的不多，但它的推导和数学意义却很有趣。它的含义是节点上所有样本的在当前模型输出处的二阶导数值之和

![Image 1][] ，这个值必须大于一定的阈值才继续分裂。但为什么是

![Image 1][] 呢？

回忆文章最开始部分的泰勒展开式是一个关于

![Image 1][] 的一元二次函数，运用初中的简单代数知识将其整理成一元二次函数的平方和形式：

![ca77502bc319f01af62e626bbff44382.png][]

会发现，（lambda=0情况下）我们在做的事情是 “以二阶导

![Image 1][] 为权重，以

![Image 1][] 为真实label，以RMSE最小化为目标，拟合

![Image 1][] "。而

![Image 1][] 的解，跟我们之前推导出来的

![Image 1][] 解析解也是一致的。另外再多说一句，如果你之前熟悉牛顿法，会发现这个形式

![Image 1][]

跟牛顿法参数更新公式

![Image 1][]

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

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

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

![524ef90a1788db516cdc2252b80c75b9.png][]

## 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`是怎么来的、特征排序、块计算等。
 *  陈天奇会议pre https://www.youtube.com/watch?v=Vly8xGnNiWs&list=PLSz28ynAy1RohdgsPfC4l4t3lHu863PGx&index=3&t=1786s 演讲中的有些细节是ppt上没有的
 *  A Gentle Introduction to GBDT https://machinelearningmastery.com/gentle-introduction-gradient-boosting-algorithm-machine-learning/ 这篇文章从历史沿革的角度，讲解boosting算法是如何一步一步发展起来的，用一句话概括出每篇论文的contribution, 可以当文献综述看。
 *  李航统计学习方法 Chapter 7 AdaBoost, Gradient Tree boosting

[Image 1]: 
[18bd660a3fe5806e6ec898300544f3b0.png]: https://img-blog.csdnimg.cn/img_convert/18bd660a3fe5806e6ec898300544f3b0.png
[a8808c069cd8a53cc746889e3caa163f.png]: https://img-blog.csdnimg.cn/img_convert/a8808c069cd8a53cc746889e3caa163f.png
[aafe0275b854904092f216359312bb0d.png]: https://img-blog.csdnimg.cn/img_convert/aafe0275b854904092f216359312bb0d.png
[1ff7e46a39a9e26c3fa78a36a7202cbe.png]: https://img-blog.csdnimg.cn/img_convert/1ff7e46a39a9e26c3fa78a36a7202cbe.png
[d4c2b9061301b8ff023260f93c702306.png]: https://img-blog.csdnimg.cn/img_convert/d4c2b9061301b8ff023260f93c702306.png
[70d6cb42f1b6d8effec9f2ce25185ad0.png]: https://img-blog.csdnimg.cn/img_convert/70d6cb42f1b6d8effec9f2ce25185ad0.png
[55ac323bea85c581ea6cb98e6513f86a.png]: https://img-blog.csdnimg.cn/img_convert/55ac323bea85c581ea6cb98e6513f86a.png
[e59d475d328c9ea36914a6fe2e358470.png]: https://img-blog.csdnimg.cn/img_convert/e59d475d328c9ea36914a6fe2e358470.png
[43607520325775438ebc7702d2ce356e.png]: https://img-blog.csdnimg.cn/img_convert/43607520325775438ebc7702d2ce356e.png
[4327d5d9f2231024ea8a271bc5b1336e.png]: https://img-blog.csdnimg.cn/img_convert/4327d5d9f2231024ea8a271bc5b1336e.png
[bb1f03f7cd3d5f6381a865ab1050f55b.png]: https://img-blog.csdnimg.cn/img_convert/bb1f03f7cd3d5f6381a865ab1050f55b.png
[70f77a4026e86af928084496ea5163a4.png]: https://img-blog.csdnimg.cn/img_convert/70f77a4026e86af928084496ea5163a4.png
[ca77502bc319f01af62e626bbff44382.png]: https://img-blog.csdnimg.cn/img_convert/ca77502bc319f01af62e626bbff44382.png
[524ef90a1788db516cdc2252b80c75b9.png]: https://img-blog.csdnimg.cn/img_convert/524ef90a1788db516cdc2252b80c75b9.png