梯度提升树(GBDT)

上一篇笔记学习了AdaBoost算法，AdaBoost每一轮基学习器训练过后都会更新样本权重，再训练下一个学习器，最后将所有的基学习器加权组合。AdaBoost使用的是指数损失，这个损失函数的缺点是对于异常点非常敏感，因而通常在噪音比较多的数据集上表现不佳。Gradient Boosting在这方面进行了改进，使得可以使用任何损失函数 (只要损失函数是连续可导的)，这样一些比较robust的损失函数就能得以应用，使模型抗噪音能力更强。

Boosting的基本思想是通过某种方式使得每一轮基学习器在训练过程中更加关注上一轮学习错误的样本，区别在于是采用何种方式？AdaBoost采用的是增加上一轮学习错误样本的权重的策略，而在Gradient Boosting中则将负梯度作为上一轮基学习器犯错的衡量指标，在下一轮学习中通过拟合负梯度来纠正上一轮犯的错误。这里的关键问题是：为什么通过拟合负梯度就能纠正上一轮的错误了？Gradient Boosting的发明者给出的答案是：函数空间的梯度下降。

1. 函数空间的梯度下降

首先回顾一下梯度下降(Gradient Descend)相关的内容。说远点，从泰勒公式开始，其基本形式为：

\[f(x) = \sum_{n=0}^\infty \cfrac{f^{(n)}(x_0)}{n!}(x-x_0)^n \]

一阶泰勒展开为

\[f(x) \approx f(x_0) + f^{'}(x_0)(x-x_0) \]

假设\(x^t = x^{t-1} + \Delta x\)，将\(f(x^t)\)在\(x^{t-1}\)处进行泰勒展开，则

\[\begin{aligned} f(x^t) &= f(x^{t-1} + \Delta x) \\ &\approx f(x^{t-1})+f^{'}(x^{t-1})\Delta x \end{aligned}\]

在机器学习任务中，需要最小化损失函数\(L(\theta)\)，其中\(\theta\)是要求解的模型参数。梯度下降法常用来求解这种无约束最优化问题，它是一种迭代方法：选取初值\(\theta_0\)，不断迭代，更新\(\theta\)的值，进行损失函数的极小化。其迭代公式可以表示为

\[\theta^t = \theta^{t-1} + \Delta \theta \]

根据泰勒公式，将\(L(\theta^t)\)在\(\theta^{t-1}\)处进行一阶泰勒展开，即

\[\begin{aligned} L(\theta^t) &= L(\theta^{t-1} + \Delta \theta) \\ &\approx L(\theta^{t-1})+L^{'}(\theta^{t-1})\Delta \theta \end{aligned}\]

要使得\(L(\theta^t) \lt L(\theta^{t-1})\)。可取：\(\Delta\theta = -\alpha L^{'}(\theta^{t-1})\)，则

\[\theta^t = \theta^{t-1} - \alpha L^{'}(\theta^{t-1}) \]

这里的\(\alpha\)是步长，可以通过line search确定，但是一般直接赋一个小的数。如果将\(\theta\)推广到向量的话，\(L^{'}(\theta^{t-1})\)使用偏导数替换，即

\[\theta^t = \theta^{t-1} - \alpha \cdot \frac{\partial}{\partial \theta}L(\theta) \]

上式中\(\theta^t\)表示第\(t\)次迭代后的参数，\(\theta^{t-1}\)表示第\(t-1\)次迭代后的参数，参数更新方向是负梯度方向，我们用\(\theta_t\)表示\(- \alpha \cdot \frac{\partial}{\partial \theta}L(\theta)\)，并记作\(-\alpha_tg_t\)，则\(\theta_t\)表示第\(t\)次迭代参数的增量，我们发现最终参数等于每次迭代的增量的累加和

\[\theta = \sum_{t=0}^T\theta_t \]

其中，\(\theta_0\)为初始值。这是梯度下降在参数空间（\(\theta\)是模型参数）的描述。

那么，如果我们将\(\theta^t\)换成函数\(f^t(x)\)，即把\(f^t(x)\)看作参数，则

\[f^t(x) = f^{t-1}(x) + f_t(x) \]

其中，\(f^t(x)\)表示第\(t\)次迭代后的函数，\(f^{t-1}(x)\)表示第\(t-1\)次迭代后的函数，同样的拟合负梯度，\(f_t(x) = - \alpha \cdot \frac{\partial}{\partial f^{t-1}(x)}L(f(x))\)，记作\(f_t(x) = -\alpha_tg_t(x)\)，则\(\theta_t\)表示第\(t\)次迭代参数的增量，最终函数等于每次迭代的增量的累加和

\[F(x) = \sum_{t=0}^Tf_t(x) \]

其中，\(f_0(x)\)为模型初始值，通常为常数。这就是在函数空间（\(f(x)\)实际是一个函数）的梯度下降。

具体的，Gradient Boosting采用和AdaBoost同样的加法模型，在第\(t\)次迭代中，前\(t-1\)个基学习器都是固定的，即

\[f_t(x) = f_{t-1}(x) + \alpha_t h_t(x) \]

因而在第\(t\)步我们的目标是最小化损失函数\(L(f) = \sum\limits_{i=1}^mL(y_i,f_t(x_i))\)，进而求得相应的基学习器。如前所述，若将\(f(x)\)当成参数，则同样可以使用梯度下降法：

\[f_t(x) = f_{t-1}(x) - \alpha_t \cdot \frac{\partial}{\partial f_{t-1}(x)}L(y,f_{t-1}(x)) \]

对比可以发现，若将\(h_t(x) \approx -\frac{\partial L(y,f_{t-1}(x))}{\partial f_{t-1}(x)}\)，即用基学习器\(h_t(x)\)拟合前一轮模型损失函数的负梯度，就是通过梯度下降法最小化\(L(f)\)。

负梯度也被称为“响应 (response)”或“伪残差 (pseudo residual)”，从名字可以看出是一个与残差接近的概念。直觉上来看，残差\(r=y−f(x)\)越大，表明前一轮学习器\(f(x)\)的结果与真实值\(y\)相差较大，那么下一轮学习器通过拟合残差或负梯度，就能纠正之前的学习器犯错较大的地方。

2. GBDT的负梯度拟合

GBDT基于函数空间的梯度下降解决损失函数的拟合问题，用损失函数的负梯度来拟合本轮损失的近似值，进而拟合一个CART回归树。

第\(t\)轮的第\(i\)个样本的损失函数的负梯度表示为

\[r_{ti} = -\bigg[\frac{\partial L(y_i, f(x_i)))}{\partial f(x_i)}\bigg]_{f(x) = f_{t-1}\;\; (x)} \]

利用\((x_i,r_{ti}), i=1,2,\cdots,m\)，我们可以拟合一棵CART回归树，得到了第\(t\)棵回归树，其对应的叶结点区域\(R_{tj}, j=1,2,\cdots,J\)。其中，\(J\)表示叶子结点的个数。

针对每一个叶子结点里的样本，我们求出使损失函数最小，也就是拟合叶结点最好的输出值\(c_{tj}\)，如下

\[c_{tj} = \smash{\arg\min_{c}}\sum\limits_{x_i \in R_{tj}} L(y_i,f_{t-1}(x_i) +c) \]

这样我们就得到了本轮的决策树拟合函数如下：

\[h_t(x) = \sum\limits_{j=1}^{J}c_{tj}I(x \in R_{tj}) \]

从而本轮最终得到的强学习器的表达式如下：

\[f_{t}(x) = f_{t-1}(x) + \sum\limits_{j=1}^{J}c_{tj}I(x \in R_{tj}) \]

通过损失函数的负梯度来拟合，我们找到了一种通用的拟合损失误差的办法，这样无轮是分类问题还是回归问题，我们通过其损失函数的负梯度的拟合，就可以用GBDT来解决我们的分类回归问题。区别仅仅在于损失函数不同导致的负梯度不同而已。

3. GBDT回归算法

好了，有了上面的思路，下面我们总结下GBDT的回归算法。为什么没有加上分类算法一起？那是因为分类算法的输出是不连续的类别值，需要一些处理才能使用负梯度，我们在下一节讲。

输入：是训练集样本\(T=\{(x_1,y_1),(x_2,y_2), \cdots, (x_m,y_m)\}\)，最大迭代次数\(T\), 损失函数\(L\)；

输出：是强学习器\(f(x)\)；

1）初始化弱学习器

\[f_0(x) = \smash{\arg\min_{c}}\sum\limits_{i=1}^{m}L(y_i, c) \]

2）对迭代轮数\(t=1,2,\cdots,T\)，有：

对样本\(i=1,2,\cdots,m\)，计算负梯度

\[r_{ti} = -\bigg[\frac{\partial L(y_i, f(x_i)))}{\partial f(x_i)}\bigg]_{f(x) = f_{t-1}(x)} \]

利用\((x_i,r_{ti}), i=1,2,\cdots,m\)，拟合一棵CART回归树，得到第\(t\)棵回归树，其对应的叶结点区域为\(R_{tj}, j=1,2,\cdots,J\)。其中，\(J\)表示叶子结点的个数。
对叶子区域\(j=1,2,\cdots,J\)，计算最佳拟合值

\[c_{tj} = \smash{\arg\min_{c}}\sum\limits_{x_i \in R_{tj}} L(y_i,f_{t-1}(x_i) +c) \]

更新强学习器

\[f_{t}(x) = f_{t-1}(x) + \sum\limits_{j=1}^{J}c_{tj}I(x \in R_{tj}) \]

3）得到强学习器\(f(x)\)的表达式

\[f(x) = f_T(x) =f_0(x) + \sum\limits_{t=1}^{T}\sum\limits_{j=1}^{J}c_{tj}I(x \in R_{tj}) \]

第1步初始化，估计使损失函数极小化的常数值，它是只有一个根结点的树。第2.1步计算损失函数的负梯度在当前模型的值，将它作为残差的估计。对于平方损失函数，它就是通常所说的残差；对于一般损失函数，它就是残差的近似值。第2.2步估计回归树叶结点区域，以拟合残差的近似值。第2.3步利用线性搜索估计叶结点区域的值，使损失函数极小化。第2.4步更新回归树。第3步得到输出的最终模型\(f(x)\)。

4. GBDT分类算法

我们再看看GBDT分类算法，GBDT的分类算法从思想上和GBDT的回归算法没有区别，但是由于样本输出不是连续的值，而是离散的类别，导致我们无法直接从输出类别去拟合类别输出的误差。

为了解决这个问题，主要有两个方法，一个是用指数损失函数，此时GBDT退化为Adaboost算法。另一种方法是用类似于逻辑回归的对数似然损失函数的方法。也就是说，我们用的是类别的预测概率值和真实概率值的差来拟合损失。本文仅讨论用对数似然损失函数的GBDT分类。而对于对数似然损失函数，我们又有二元分类和多元分类的区别。

4.1 二元GBDT分类算法

对于二元GBDT，如果用类似于逻辑回归的对数似然损失函数，则损失函数为：

\[L(y, f(x)) = log(1+ e^{-yf(x)}) \]

其中\(y\in \{-1,+1\}\)。此时的负梯度误差为

\[r_{ti} = -\bigg[\frac{\partial L(y, f(x_i)))}{\partial f(x_i)}\bigg]_{f(x) = f_{t-1}\;\; (x)} = \cfrac{y_i}{1+e^{y_if(x_i)}} \]

对于生成的决策树，我们各个叶子节点的最佳负梯度拟合值为

\[c_{tj} = \smash{\arg\min_{c}}\sum\limits_{x_i \in R_{tj}} log(1+e^{-y_i(f_{t-1}(x_i) +c))} \]

由于上式比较难优化，我们一般使用近似值代替

\[c_{tj} = \cfrac{\sum\limits_{x_i \in R_{tj}}r_{ti}}{\sum\limits_{x_i \in R_{tj}}|r_{ti}|(1-|r_{ti}|)} \]

除了负梯度计算和叶子节点的最佳负梯度拟合的线性搜索，二元GBDT分类和GBDT回归算法过程相同。

4.2 多元GBDT分类算法

多元GBDT要比二元GBDT复杂一些，对应的是多元逻辑回归和二元逻辑回归的复杂度差别。假设类别数为\(K\)，则此时我们的对数似然损失函数为：

\[L(y, f(x)) = - \sum\limits_{k=1}^{K}y_klog\;p_k(x) \]

其中，如果样本输出类别为\(k\)，则\(y_k=1\)。第\(k\)类的概率\(p_k(x)\)的表达式为：

\[p_k(x) =\cfrac{e^{f_k(x)}}{\sum\limits_{l=1}^{K} e^{f_l(x)}} \]

结合上两式，我们可以计算出第\(t\)轮的第\(i\)个样本对应类别\(l\)的负梯度误差为

\[r_{til} = -\bigg[\frac{\partial L(y_i, f(x_i)))}{\partial f(x_i)}\bigg]_{f_k(x) = f_{l, t-1}(x)} = y_{il} - p_{l, t-1}(x_i) \]

观察上式可以看出，其实这里的误差就是样本i对应类别\(l\)的真实概率和\(t−1\)轮预测概率的差值。

对于生成的决策树，我们各个叶子节点的最佳负梯度拟合值为

\[c_{tjl} = \smash{\arg\min_{c_{jl}}}\sum\limits_{i=0}^{m}\sum\limits_{k=1}^{K} L(y_k, f_{t-1, l}(x) + \sum\limits_{j=0}^{J}c_{jl} I(x_i \in R_{tjl})) \]

由于上式比较难优化，我们一般使用近似值代替

\[c_{tjl} = \frac{K-1}{K} \; \frac{\sum\limits_{x_i \in R_{tjl}}r_{til}}{\sum\limits_{x_i \in R_{til}}|r_{til}|(1-|r_{til}|)} \]

除了负梯度计算和叶子节点的最佳负梯度拟合的线性搜索，多元GBDT分类和二元GBDT分类以及GBDT回归算法过程相同。

这里补充一个更直观的GBDT多分类算法流程的说明：

我们在训练的时候，是针对样本\(X\)每个可能的类都训练一个分类回归树。举例说明，目前样本有三类，也就是\(K = 3\)。样本\(x\)属于第二类。那么针对该样本\(x\)的分类结果，其实我们可以用一个三维向量\([0,1,0]\)来表示。\(0\)表示样本不属于该类，\(1\)表示样本属于该类。由于样本已经属于第二类了，所以第二类对应的向量维度为1，其他位置为\(0\)。

针对样本有三类的情况，我们实质上是在每轮的训练的时候是同时训练三棵树。第一棵树针对样本\(x\)的第一类，输入为\((x,0)\)。第二棵树输入针对样本\(x\)的第二类，输入为\(x,1)\)。第三棵树针对样本\(x\)的第三类，输入为\((x，0)\)。

在这里每棵树的训练过程其实就是就是我们之前已经提到过的CATR决策树的生成过程。在此处我们参照之前的生成树的程序即可解出三棵树，以及三棵树对\(x\)类别的预测值\(f_1(x),f_2(x),f_3(x)\)。那么在此类训练中，我们仿照多分类的逻辑回归，使用softmax来产生概率，则属于类别\(1\)的概率

\[p_{1}=\cfrac{e^{f_{1}{(x)}}}{\sum_{k= 1}^{3}e^{f_k(x)}} \]

并且我们可以针对类别\(1\)求出残差\(y_{11}(x) = 0-p_1(x)\)，类别\(2\)求出残差\(y_{22}(x)= 1-p_2(x)\)，类别\(3\)求出残差\(y_{33}(x)= 0-p_{3}(x)\)。

然后开始第二轮训练，针对第一类，输入为\((x,y_{11}(x))\), 针对第二类输入为\((x,y_{22}(x))\), 针对第三类输入为\((x,y_{33}(x))\)。继续训练出三棵树。一直迭代\(M\)轮。每轮构建3棵树。

所以当\(K=3\)。我们其实应该有三个式子

\[F_{1M}{(x)}=\sum_{m=1}^{M}{\hat{C_{1m}}I(x\epsilon R_{1m})} \]

\[F_{2M}{(x)}=\sum_{m=1}^{M}{\hat{C_{2m}}I(x\epsilon R_{2m})} \]

\[F_{3M}{(x)}=\sum_{m=1}^{M}{\hat{C_{3m}}I(x\epsilon R_{3m})} \]

当训练完毕以后，新来一个样本\(x_1\)，我们需要预测该样本的类别的时候，便可以有这三个式子产生三个值，\(f_1(x),f_2(x),f_3(x)\)。样本属于某个类别\(c\)的概率为

\[p_{c}=\cfrac{e^{f_{c}{(x)}}}{\sum_{k= 1}^{3}e^{f_k(x)}} \]

5. GBDT的正则化

和Adaboost一样，我们也需要对GBDT进行正则化，防止过拟合。GBDT的正则化主要有三种方式。

1）和Adaboost类似的正则化项，即步长(learning rate)。定义为\(ν\),对于前面的弱学习器的迭代

\[f_k(x)=f_{k−1}(x)+h_k(x) \]

如果我们加上了正则化项，则有

\[f_{k}(x) = f_{k-1}(x) + \nu h_k(x) \]

\(ν\)的取值范围为\(0\lt ν \le 1\)。对于同样的训练集学习效果，较小的\(ν\)意味着我们需要更多的弱学习器的迭代次数。通常我们用步长和迭代最大次数一起来决定算法的拟合效果。

2）正则化的方式是通过子采样比例（subsample）。取值为(0,1]。注意这里的子采样和随机森林不一样，随机森林使用的是放回抽样，而这里是不放回抽样。如果取值为1，则全部样本都使用，等于没有使用子采样。如果取值小于1，则只有一部分样本会去做GBDT的决策树拟合。选择小于1的比例可以减少方差，即防止过拟合，但是会增加样本拟合的偏差，因此取值不能太低。推荐在[0.5, 0.8]之间。

使用了子采样的GBDT有时也称作随机梯度提升树(Stochastic Gradient Boosting Tree, SGBT)。由于使用了子采样，程序可以通过采样分发到不同的任务去做boosting的迭代过程，最后形成新树，从而减少弱学习器难以并行学习的弱点。

3）对于弱学习器即CART回归树进行正则化剪枝。在决策树原理篇里我们已经讲过，这里就不重复了。

6. 特征重要性

单棵决策树的可解释性很强，GBDT则继承了这个优点。

对于单棵树\(T\)，用下式来衡量每个特征\(X_l\)的重要性：

\[\mathcal{I}_{l}^2(T) = \sum\limits_{t=1}^{J-1}\hat{i}_t^2I(v(t) = l) \]

其中\(J\)表示叶结点 (leaf node) 数量，\(J−1\)表示内部结点 (internal node) 数量，\(X_{v(t)}\)是与内部结点\(t\)相关联的分裂特征。对于每个内部结点\(t\)，用特征\(X_{v(t)}\)来模拟划分特征空间，得到一个分裂后的平方误差减少量，即\(\hat{i}_t^2\)，最后将所有内部节点上的误差减少量加起来，就是特征\(X_l\)的重要性。总误差减少地越多，该特征就越重要。

对于\(M\)棵树的集成而言，特征重要性就是各棵树相应值的平均：

\[\mathcal{I}_{l}^2 = \frac1M\sum\limits_{m=1}^M\mathcal{I}_l^2(T_m) \]