从CART到经典GBDT到XGBoost——树模型重温

正如从AdaBoost到GBM到LambdaMART——Boosting方法重温所提到的，GBM中应用最广泛的模型是GBDT，即以CART回归树为基函数的梯度提升树模型。XGBoost作为实现Gradient Boosting的主流框架之一，除了在计算性能上作优化之外，在演算法层面也作了一些改变。这里挑若干值得一提的点稍作对比和总结。

减轻过拟合的手段

机器学习上的优化和一般优化问题最大的区别在于真实数据分布未知，因此无法对真实的优化目标期望泛化误差直接进行优化，只能从分布中抽样出训练集，通过优化训练集上的经验风险函数来间接优化期望泛化误差。这就意味著在经验风险函数上的最优解不一定是期望泛化误差的最优解，通常的情况是在经验风险函数上表现良好，但基于测试集来做估计的泛化误差高，也即所谓的过拟合。从优化的角度来看，可以通过一些约束将经验风险函数上的最优解「往隔壁挪一挪」，让它没那么完美，但在泛化误差上表现更好，从而减轻过拟合。

决策树本身是一个容易高方差的模型，因为它挑选切分特征和切分点是完全基于训练数据的分布来的，也就是多次抽样下训练集之间的差异很容易导致得到的决策树也显著不同，因此方差理论上会比较大。

理论上GBDT作为多个决策树的组合，会得到一个相对单棵决策树更平滑的模型，因此能一定程度上降低方差，也就没那么容易过拟合。但这并非意味著GBDT不存在过拟合，毕竟无论如何它都只能基于训练数据来学习。Friedman在提出GBM时给出了一些降低GBDT过拟合程度的方法，而XGBoost也相应地做了一些补充和改良。

首先将GBDT模型记为 $F(x)=sum_{m=1}^Malpha_mg_m(x)$ ，其中 $g_m(x)=sum_{t=1}^T c_{mt}I(x in R_{mt})$ . 那么 $F(x)=sum_{m=1}^Malpha_mg_m(x)=sum_{m=1}^M sum_{t=1}^T alpha_mc_{mt}I(x in R_{mt})$ . 又记 $w_{mt}=alpha_mc_{mt}$ ，因此 $F(x)=sum_{m=1}^M sum_{t=1}^T w_{mt}I(x in R_{mt})$ .

GBDT模型化简之后可以看作拟合负梯度方向的是决定 $I(x in R_{mt})$ 的一系列规则，而 $w_{mt}$ 是步长。由于 $R_{mt}$ 之间的互斥关系，一个对应于唯一一个前进方向和步长。

1）控制迭代次数和步长

由于Boosting的过程可以看作是在函数空间通过line search方法求最优值解，因此为了减轻过拟合，可以在Boosting层面上可以控制迭代次数和步长，分别对应于和 $w_{mt}$ . 的最优值可以由交叉检验实验得出，而控制 $w_{mt}$ 的常见方法有两种：

shrinkage

即指定收缩因子，使收缩后的模型变为 $F(x)=sum_{m=1}^M etasum_{t=1}^T w_{mt}I(x in R_{mt})$ 。

这种做法可以从多个角度进行解释，比如从Boosting角度可以看成是通过引入学习率对步长进行收缩；从决策树角度看是对叶节点输出值作收缩；若是将 $sum_{t=1}^T w_{mt}I(x in R_{mt})$ 记为一棵新的决策树， $F(x)=sum_{m=1}^M eta f_m(x)$ 又可以看作是降低每个基函数的影响力。

L2正则化

即在损失函数上加上 $w_{mt}$ 对应向量的L2模，使得第次boost的损失函数为 $sum_{i=1}^nL(y_i, F^{m-1}(x_i) +sum_{t=1}^T w_{mt}I(x in R_{mt})) + frac{1}{2}lambda||vec{w}||^2$ 。

这种方法相当于在原优化目标中加入 $||vec{w}||^2<C$ 的约束，其中常数与存在某种关联。相比起shrinkage方法的人为指定收缩因子，这种方法对各个分量产生的影响是不均等的，比如说对绝对值大的分量 $w_{mt}$ 上的收缩力度会大一些。

2）控制决策树的结构

以上这些方法能影响的是树的数目、叶节点输出值的平滑程度，若要通过调整树的结构来降低过拟合程度，常用的方法是从随机森林借鉴过来的row subsampling和column subsampling，分别控制每棵树采用的样本集和特征集，从而使得树的结构与使用全集时有所不同。

相对于经典CART树上的后剪枝，以上这些策略更为高效，毕竟GBDT需要拟合多个CART树，假如让每棵CART树都完全生长再剪枝，一是低效，二是让单棵CART树尽可能榨干训练数据的信息对于Boosting方法来说没有必要。

虽然XGBoost也支持shrinkage和subsampling等方法，但使其与经典GBDT区别开来的还是它在经典GBDT的损失函数上加上L2正则项和决策树大小，即损失函数变为 $sum_{i=1}^nL(y_i, F^{m-1}(x_i) +sum_{t=1}^T w_{mt}I(x in R_{mt})) + frac{1}{2}lambda||vec{w}||^2 + u T$ ，其中为当前这棵决策树的叶节点数目。这也使得它的参数估计方法与经典GBDT有所不同。

正则化损失函数下的参数估计

传统GBDT每一步的损失函数为 $sum_{i=1}^nL(y_i, F^{m-1}(x_i) +sum_{t=1}^T w_{mt}I(x in R_{mt}))$

在参数估计时只需要找到最优的 $R_{mt}$ 来拟合负梯度方向 $abla L |_{F^{m-1}(x)}$ ，再求出最优步长 $vec{w_m}^*=argmin_{vec{w_m}} sum_{i=1}^nL(y_i, F^{m-1}(x_i) +sum_{t=1}^T w_{mt}I(x in R_{mt}))$ 。由于 $R_{mt}$ 之间的互斥关系， $w_{mt}^*=argmin_{w_{mt}} sum_{x_i in R_{mt}}L(y_i, F^{m-1}(x_i) + w_{mt})$ ，即各分量的最优解仅与其对应的子区域上的样本有关。

而正如上文所提到的，XGBoost使用正则化的目标函数 $sum_{i=1}^nL(y_i, F^{m-1}(x_i) +sum_{t=1}^T w_{mt}I(x in R_{mt})) + frac{1}{2}lambda||vec{w}||^2 + u T$

这个目标函数对并非连续可导的，因此求最优 $R_{mt}$ 无法直接通过拟合负梯度方向 $abla L |_{F^{m-1}(x)}$ 来完成。但实际上负梯度方向也是通过对损失函数作一阶泰勒近似推导而来的，因此不妨直接通过对损失函数作二阶泰勒近似推导出最优的 $w_{mt}$ 和 $R_{mt}$ 。

1）估计 $w_{mt}$

假设现在处于第次Boosting，当前的树结构是，划分出了个子区域。同样由子区域之间的互斥关系得到在当前树结构下 $w_{mt}$ 的最优解：

$w_{mt}^* = argmin_{w_{mt}}{sum_{x_i in R_{mt}} L(y_i, F^{m-1}(x_i) + w_{mt}) + frac{1}{2}lambda w_{mt}^2 + u}$

对 $L(y_i, F^{m-1}(x_i) + w_{mt})$ 作二阶泰勒近似得到 $w_{mt}^* = argmin_{w_{mt}}{sum_{x_i in R_{mt}} L(y_i, F^{m-1}(x_i)) + w_{mt} * abla L + frac{1}{2}w_{mt}^2* abla^2L + frac{1}{2}lambda w_{mt}^2 + u}\=argmin_{w_{mt}}{sum_{x_i in R_{mt}} w_{mt} * abla L + frac{1}{2}w_{mt}^2* abla^2L + frac{1}{2}lambda w_{mt}^2}\=argmin_{w_{mt}} {sum_{x_i in R_{mt}} w_{mt} * abla L + frac{1}{2}( abla^2L+lambda)w_{mt}^2}$

很容易可以求出 $w_{mt}^*=-frac{sum_{x_i in R_{mt}} abla L}{sum_{x_i in R_{mt}} abla^2L+lambda}$

相比起无正则化约束时的经典GBDT的解 $-frac{1}{|R_{mt}|}sum_{x_i in R_{mt}} abla L$

$w_{mt}^*$ 可以直观而不严谨地看作是用Hessian矩阵和正则化系数对负梯度方向作了修正。

2）估计 $R_{mt}$

由于目前是处于第次boost，因此 $L(y_i, F^{m-1}(x_i))$ 是常数，因此在考虑 $sum_{t=1}^T w_{mt}I(x in R_{mt})$ 的生成时可以忽略 $L(y_i, F^{m-1}(x_i))$ ，因此将 $w_{mt}^*$ 代入到 $w_{mt} * abla L + frac{1}{2}w_{mt}^2* abla^2L + frac{1}{2}lambda w_{mt}^2 + u$ 就得到结构下 $R_{mt}$ 这个子区域的损失函数极小值，个子区域的整体损失函数为 $sum_{t=1}^Tw_{mt}^* * abla L + frac{1}{2}w_{mt}^{*2}* abla^2L + frac{1}{2}lambda w_{mt}^{*2} + u T\=sum_{t=1}^T-frac{1}{2}frac{(sum_{x_i in R_{mt}} abla L)^2}{sum_{x_i in R_{mt}} abla^2L + lambda} + u T$

假设现在在结构的基础上要对某一个叶节点 $R_{mt}$ 继续进行切分成两个新的叶节点 $R_{mtL}$ 和 $R_{mtR}$ ，那么切分之后的整体损失函数的下降值为

$frac{1}{2}(frac{(sum_{x_i in R_{mt_L}} abla L)^2}{sum_{x_i in R_{mt_L}} abla^2L + lambda} +frac{(sum_{x_i in R_{mt_R}} abla L)^2}{sum_{x_i in R_{mt_R}} abla^2L + lambda} -frac{(sum_{x_i in R_{mt_R}} abla L)^2}{sum_{x_i in R_{mt_R}} abla^2L + lambda} )- u$

使此下降值最大的即是这一次生长的最优切分特征和最优切分点。使得下降值最大显然也等价于生长后的损失函数最小，这与CART的生成逻辑是一致的。

另一方面，XGBoost利用了损失函数的二阶展开来迭代求解，理论上这会比传统GBDT的收敛速度更快。

对缺失值和无序离散型特征的处理方式

缺失值和离散特征本身也是CART树生成过程中的的重要话题。之所以要强调「无序」是因为有序的离散特征在切分时的做法可以和连续特征一样，而无序离散特征的难点在于它总共有这么多种切分可能，其中为该特征的取值个数。

1）缺失值

先以CART回归树为例，生成过程决定切分特征和切分点的目标函数是 $min_{j,s}{min_{c_1}sum_{x_iin R_1(j,s)} (y_i-c_1)^2 + min_{c_2} sum_{x_i in R_2(j,s)} (y_i - c_2)^2 }$

其中和根据在第维特征上的值与切分点的大小关系来定义。

在训练时假如在第维上的特征缺失，那么目标函数将无法考虑到这个样本的信息，相当于将它丢弃；而在预测时，第维上的特征缺失的样本也无法得到预测结果。

CART处理缺失值的常见方式有两种：

将缺失值看作是这个特征的新取值，填充为「missing」或「NaN」（这是lightGBM采用的方式）
采用代理切分特征。仍然是以回归树为例，每次基于非缺失数据上的目标函数找到最优的之后，再从以外的特征寻找和的「切分行为」很相似的若干个代理特征。在针对维特征缺失的样本做预测时，在这些代理特征里面找一个上不缺失且与相近程度最高的代理特征来完成切分