GBDT,全称为梯度提升决策树,即Gradient Boosting Decision Tree,它与Friedman等人的《Additive logistic regression: a statistical view of boosting》这篇文章有极大的联系,基本上可以说是由这篇文章发展而来的,但是由于GBDT自己有着很好的实践效果,因此对于学习集成学习来说读一读原本的论文也是很必要的。因此这篇博客主要是对GBDT的原论文进行精读和理解。
《Greedy function approximation: A gradient boosting machine》是Jerome H. Friedman 于2001年在《Annals of Statistics》(统计四大金刚之一)发表的。下面就是论文的概览以及全文翻译。
文章结构:
- 函数估计:提出了一般的函数估计的概念以及经常使用的损失函数,并给出了问题的解决思路,即把函数优化问题变为参数优化问题,然后再用数值优化算法例如最速下降求解;
- 函数空间中的数值优化:将最速下降法在函数空间中运用,得到增量函数的方向 g m ( x ) g_m(x) gm(x)以及线性搜索的步长 ρ m \rho_m ρm;
- 有限数据:在有限数据样本中,借助参数化的求解方法,用贪婪分步算法求得的 h ( x ; a m ) h(x;a_m) h(x;am)来替代最速下降中的 g m ( x ) g_m(x) gm(x),并求得相应的线性搜索的步长 ρ m \rho_m ρm,由此得到梯度提升Gradient Boost的一个通用算法框架(算法1);
- 应用:加性模型:将Gradient Boost应用于几种常用的损失准则:运用到最小二乘(LS)得到算法2的 L S LS LS_Boost;后面分析了回归树的特性,将它作为基学习器,再运用到最小绝对偏差(LAD),这样得到了算法3的 L A D LAD LAD_TreeBoost;再运用到Huber(M)损失函数就得到了 M M M_TreeBoost;运用到logistic二项对数似然(L)得到了 L K L_K LK_TreeBoost(二分类和多分类)。然后再讨论了权重修剪,即删除那些 w i w_i wi值小于阈值的所有样本,由此达到减少计算量但不损失精度的目的;
- 正则化:讨论了两种防止过拟合的参数,即学习速率 ν \nu ν和迭代次数 M M M,通过模拟实验来阐述这样的 ν − M \nu-M ν−M之间的权衡,说明降低学习速率可以显著地提高性能,但迭代次数可能会上升,但实际中迭代次数应该在计算方便或者可行的时候尽可能大;
- 模拟研究:通过蒙特卡罗模拟评估算法在许多不同情况下的表现,即产生100个随机生成的目标函数,并设定误差分布。这样对比了 L S LS LS_TreeBoost、 M M M_TreeBoost和 L A D LAD LAD_TreeBoost三种算法,在正常误差和误差极端情况中 M M M_TreeBoost都表现较好;对比了 L S LS LS_TreeBoost与MARS, M M M_TreeBoost对于正态误差几乎与 L S LS LS_TreeBoost一样准确,而且对输出 y y y的离群值具有很强的抵抗力;对比了 L k L_k Lk_TreeBoost、K分类LogitBoost和Adaboost.MH,当对这三种方法中的每一种都仔细调整收缩参数时,它们之间的性能差别可能很小;
- 提升树:树的大小 J J J是基学习器的主要的超参数,它表示输入变量中占主导的交互作用的最高阶数,在实践中,这可以通过使用独立的测试集来调整;
- 解释:利用输入变量的相对重要性以及部分依赖图(可视化的解释)定量地翻译了目标函数 F ∗ ( x ) F^*(x) F∗(x)的性质,提供了关于输入变量 x x x和输出变量 y y y的潜在关系的信息;
- 真实数据:将上述几种提升树回归算法将在两个实际数据集上进行说明,并分析了变量的相对重要性以及部分依赖图,结果和之前模拟时类似;
- 数据挖掘:总结了提升树的一些特性,继承了决策树的良好特性,例如不用考虑变量转换、对长尾分布和异常值不敏感、树的内部特征选择导致稳健性好、处理缺失值很方便等,也改善了单一树不准确不稳定的缺点,解释性也有一定的工具,而且在计算方面也有一定优势,可以作为预处理工具,且可以在先前停止的地方重新启动,不会造成计算损失。
1. 函数估计
在函数估计或预测学习的问题中,我们有一个由随机的“输出”或“响应”变量 y y y和一组随机的“输入”或“解释”变量 x = { x 1 , … , x n } x=\{x_1,\dots,x_n\} x={x1,…,xn}组成的系统。使用已知的训练样本 { y i , x i } i N \{y_i,x_i\}_i^N {yi,xi}iN,目标是得到从 x x x到 y y y的映射函数 F ∗ ( x ) F^*(x) F∗(x)的一个估计或者近似 F ^ ( x ) \hat{F}(x) F^(x),这个估计要在所有的样本的联合分布上最小化某些特定的损失函数 L ( y , F ( X ) ) L(y,F(X)) L(y,F(X))的期望 F ∗ = arg min F E y , x L ( y , F ( x ) ) = arg min F E x [ E y ( L ( y , F ( x ) ) ) ∣ x ] F^*=\arg\min_FE_{y,x}L(y,F(x))=\arg\min_FE_{x}\bigg[E_y\bigg(L\big(y,F(x)\big)\bigg)\bigg|x\bigg] F∗=argFminEy,xL(y,F(x))=argFminEx[Ey(L(y,F(x)))∣∣∣∣x]经常使用损失函数 L ( y , F ( x ) ) L(y,F(x)) L(y,F(x))包括平方误差 ( y − F ( x ) ) 2 (y-F(x))^2 (y−F(x))2、对 y ∈ R 1 y\in R^1 y∈R1(回归)时的绝对误差 ∣ y − F ∣ |y-F| ∣y−F∣以及当 y ∈ { − 1 , 1 } y\in \{-1,1\} y∈{−1,1}(分类)时的负伯努利对数似然 l o g ( 1 + e − 2 y F ) log(1+e^{-2yF}) log(1+e−2yF)。
常见的是将 F ( x ) F(x) F(x)限制为参数类函数 F ( x ; P ) F(x;P) F(x;P)中的一个,其中 P = { P 1 , P 2 , … } P=\{P_1,P_2,\dots\} P={P1,P2,…}是一个有限的参数集,其联合值标识各个类成员。本文中我们关注以下加性展式 F ( x ; { β m , a m } 1 M ) = ∑ m = 1 M β m h ( x ; a m ) (2) F(x;\{\beta_m,a_m\}_1^M)=\sum_{m=1}^M\beta_mh(x;a_m)\tag{2} F(x;{βm,am}1M)=m=1∑Mβmh(x;am)(2)
上式中生成函数 h ( x ; a ) h(x;a) h(x;a)通常是一个输入变量 x x x的简单的参数函数,特征参数为 a = { a 1 , a 2 , … } a=\{a_1,a_2,\dots\} a={a1,a2,…}。每一项在为这些参数选择的联合值 a m a_m am上有所不同。诸如此类的扩展是对诸如神经网络[Rumelhart, Hinton, and Williams(1986)]、径向基函数[Powell(1987)]、MARS [Friedman(1991)]、小波[Donoho(1993)]和支持向量机[Vapnik(1995)]等许多函数逼近方法的验证。这里有趣的是每个函数 h ( x ; a m ) h(x;a_m) h(x;am)都是一个小的回归树,例如CART决策树[Breiman, Friedman, Olshen and Stone (1983)]。对于回归树来说,参数 a m a_m am是单个树的分裂变量、分裂位置和终端节点均值。
1.1 数值优化
一般来说,选择一个参数模型 F ( x ; P ) F(x;P) F(x;P)把函数优化问题变为参数优化问题 P ∗ = arg min P Φ ( P ) P^*=\arg\min_P\Phi(P) P∗=argPminΦ(P)其中 Φ ( P ) = E y , x L ( y , F ( x ; P ) ) \Phi(P)=E_{y,x}L(y,F(x;P)) Φ(P)=Ey,xL(y,F(x;P))则 F ∗ ( x ) = F ( x ; P ∗ ) F^*(x)=F(x;P^*) F∗(x)=F(x;P∗)
对大部分的 F ( x ; P ) F(x;P) F(x;P)和 L L L,求解上式一定会用到数值优化方法。这通常涉及用下式表示参数的解 P ∗ = ∑ m = 0 M p m P^*=\sum_{m=0}^Mp_m P∗=m=0∑Mpm其中 p 0 p_0 p0是初始化猜测, { p m } 1 M \{p_m\}_1^M {pm}1M是连续的增量(“步骤”或“提升”),每个步骤基于前面步骤的顺序。利用优化方法确定了计算各步 p m p_m pm的方法。
1.2 最速下降
最速下降是常用的数值最小化方法中最简单的一种。它如下定义了增量 { p m } 1 M \{p_m\}_1^M {pm}1M。首先计算当前梯度 g m g_m gm: g m = { g j m } = { [ ∂ Φ ( P ) ∂ P j ] P = P m − 1 } g_m=\{g_{jm}\}=\bigg\{\bigg[{\partial\Phi(P)\over\partial P_j}\bigg]_{P=P_{m-1}}\bigg\} gm={gjm}={[∂Pj∂Φ(P)]P=Pm−1}这里 P m − 1 = ∑ i = 1 m − 1 p i P_{m-1}=\sum_{i=1}^{m-1}p_i Pm−1=i=1∑m−1pi然后每步取为 p m = − ρ m g m p_m=-\rho_mg_m pm=−ρmgm其中 ρ m = arg min ρ Φ ( P m − 1 − ρ g m ) \rho_m=\arg\min_\rho\Phi(P_{m-1}-\rho g_m) ρm=argρminΦ(Pm−1−ρgm)负梯度 − g m -g_m −gm称为最速下降方向,且上式的 ρ \rho ρ称为沿着这个方向的线性搜索。
2. 函数空间中的数值优化
在此,我们采取“非参”的方法,并将数值优化应用于函数空间。即我们将在每个点 x x x估计的 F ( x ) F(x) F(x)作为一个“参数”,并最小化 Φ ( F ) = E y , x L ( y , F ( x ) ) = E x [ E y ( L ( y , F ( x ) ) ) ∣ x ] \Phi(F)=E_{y,x}L(y,F(x))=E_x\bigg[E_y\bigg(L\big(y,F(x)\big)\bigg)\bigg|x\bigg] Φ(F)=Ey,xL(y,F(x))=Ex[Ey(L(y,F(x)))∣∣∣∣x]或等价地,在每个点 x x x直接关于 F ( x ) F(x) F(x)最小化 ϕ ( F ( x ) ) = E y [ L ( y , F ( x ) ) ∣ x ] \phi(F(x))=E_y\bigg[L\bigg(y,F(x)\bigg)\bigg|x\bigg] ϕ(F(x))=Ey[L(y,F(x))∣∣∣∣x]在函数空间中有无数这样的参数,但在(下面讨论的)数据集中只有有限个 { F ( x i ) } 1 N \{F(x_i)\}_1^N {F(xi)}1N。按照数值优化的范例,我们把解取为 F ∗ ( x ) = ∑ m = 1 M f m ( x ) F^*(x)=\sum_{m=1}^Mf_m(x) F∗(x)=m=1∑Mfm(x)其中 f 0 ( x ) f_0(x) f0(x)是初始化猜测, { f m ( x ) } 1 M \{f_m(x)\}_1^M {fm(x)}1M是由优化方法定义的增量函数(“步骤”或“提升”)。
对于最速下降 f m ( x ) = − ρ m g m ( x ) (6) f_m(x)=-\rho_mg_m(x)\tag{6} fm(x)=−ρmgm(x)(6)其中 g m ( x ) = [ ∂ ϕ ( F ( x ) ) ∂ F ( x ) ] F ( x ) = F m − 1 ( x ) = [ ∂ E y [ L ( y , F ( x ) ) ∣ x ] ∂ F ( x ) ] F ( x ) = F m − 1 ( x ) g_m(x)=\bigg[{\partial\phi(F(x))\over\partial F(x)}\bigg]_{F(x)=F_{m-1}(x)}=\bigg[{\partial E_y\big[L\big(y,F(x)\big)\big|x\big]\over\partial F(x)}\bigg]_{F(x)=F_{m-1}(x)} gm(x)=[∂F(x)∂ϕ(F(x))]F(x)=Fm−1(x)=[∂F(x)∂Ey[L(y,F(x))∣∣x]]F(x)=Fm−1(x)且 F m − 1 ( x ) = ∑ i = 0 m − 1 f i ( x ) F_{m-1}(x)=\sum_{i=0}^{m-1}f_i(x) Fm−1(x)=i=0∑m−1fi(x)
假设有充分的正则性可以交换微分和积分,它就变成了 g m ( x ) = E y [ ∂ L ( y , F ( x ) ) ∂ F ( x ) ∣ x ] F ( x ) = F m − 1 ( x ) (7) g_m(x)=E_y\bigg[{\partial L\big(y,F(x)\big)\over\partial F(x)}\bigg|x\bigg]_{F(x)=F_{m-1}(x)}\tag{7} gm(x)=Ey[∂F(x)∂L(y,F(x))∣∣∣∣x]F(x)=Fm−1(x)(7)由线性搜索给出乘数 ρ m \rho_m ρm ρ m = arg min ρ E y , x L ( y , F m − 1 ( x ) − ρ g m ( x ) ) (8) \rho_m=\arg\min_\rho E_{y,x}L\big(y,F_{m-1}(x)-\rho g_m(x)\big)\tag{8} ρm=argρminEy,xL(y,Fm−1(x)−ρgm(x))(8)
3. 有限数据
这种非参数方法在用有限数据样本 { y i , x i } 1 N \{y_i,x_i\}_1^N {yi,xi}1N来估计 ( y , x ) (y,x) (y,x)的联合分布时失效,在这种情况下 E y [ ⋅ ∣ x ] E_y[\cdot|x] Ey[⋅∣x]无法通过每个样本点 x x x的数据值来准确估计,即使可能,我们也希望在 x x x值处估计 F ∗ ( x ) F^*(x) F∗(x)而不是在训练样本上。必须从附近的数据点得到一些帮助,通过在解上施加平滑度。一种方法是采用参数化形式,如(2),并进行1.1节中讨论的参数优化,以最小化相应的数据的期望损失估计, { β m , a m } 1 M = arg min { β m ′ , a m ′ } 1 M ∑ i = 1 N L ( y i , ∑ m = 1 M β m ′ h ( x i ; a m ′ ) ) \{\beta_m,a_m\}_1^M=\arg\min_{\{\beta'_m,a'_m\}_1^M}\sum_{i=1}^NL\bigg(y_i,\sum_{m=1}^M\beta'_mh(x_i;a'_m)\bigg) {βm,am}1M=arg{βm′,am′}1Mmini=1∑NL(yi,m=1∑Mβm′h(xi;am′))在这种不可行的情况下我们可以试着用“贪婪分步”方法,对 m = 1 , 2 , … , M m=1,2,\dots,M m=1,2,…,M ( β m , a m ) = arg min β , a ∑ i = 1 N L ( y i , F m − 1 ( x i ) + β h ( x i ; a ) ) (9) (\beta_m,a_m)=\arg\min_{\beta,a}\sum_{i=1}^NL\big(y_i,F_{m-1}(x_i)+\beta h(x_i;a)\big)\tag{9} (βm,am)=argβ,amini=1∑NL(yi,Fm−1(xi)+βh(xi;a))(9)则 F m ( x ) = F m − 1 ( x ) + β m h ( x ; a m ) (10) F_m(x)=F_{m-1}(x)+\beta_mh(x;a_m)\tag{10} Fm(x)=Fm−1(x)+βmh(x;am)(10)
注意,这个分步策略与逐步方法不同,后者在添加新项时重新调整先前输入的项。
在信号处理中这样的分步策略被称为“匹配的追求”(Mallat和张(1993)),其中 L ( y , F ) L(y,F) L(y,F)是均方误损失,且 { h ( x , a m ) } 1 M \{h(x,a_m)\}_1^M {h(x,am)}1M被称为基函数,通常是从一个过完备的类小波字典取得。在机器学习中,(9)、(10)被称为“boosting”,其中 y ∈ { − 1 , 1 } y\in\{-1,1\} y∈{−1,1}且 L ( y , F ) L(y,F) L(y,F)要么是一个指数损失准则 e − y F e^{−yF} e−yF[Freund and Schapire (1996), Schapire and Singer (1998)]要么是负二项对数似然[Friedman, Hastie and Tibshirani (2000)]。函数 h ( x ; a ) h(x;a) h(x;a)被称为“弱学习器”或“基学习器”,通常是一个分类树。
假设对于一个特定的损失 L ( y , F ) L(y,F) L(y,F)和/或基学习器 h ( x ; a ) h(x;a) h(x;a),(9)式的解很难求。给定任意逼近 F m − 1 ( x ) F_{m-1}(x) Fm−1(x),(9)式 β m h ( x , a m ) \beta_mh(x,a_m) βmh(x,am)和(10)式可以看做是在步方向 h ( x ; a m ) h(x;a_m) h(x;am)为参数化的函数类 h ( x ; a ) h(x;a) h(x;a)中的一个这样的限制下,对(1)式中基于数的估计 F ∗ ( x ) F^*(x) F∗(x)的最佳贪婪的步骤。因此它可以看做是在这个限制下(6)式中的最速下降步。通过限制,基于数据的无约束负梯度(7)模拟 − g m ( x i ) = − [ ∂ L ( y i , F ( x i ) ) ∂ F ( x i ) ] F ( x ) = F m − 1 ( x ) -g_m(x_i)=-\bigg[{\partial L\big(y_i,F(x_i)\big)\over\partial F(x_i)}\bigg]_{F(x)=F_{m-1}(x)} −gm(xi)=−[∂F(xi)∂L(yi,F(xi))]F(x)=Fm−1(x)给出了在 F m − 1 ( x ) F_{m-1}(x) Fm−1(x)的N维数据空间中最佳的最速下降步方向 − g m = { − g m ( x i ) } 1 N -g_m=\{-g_m(x_i)\}_1^N −gm={−gm(xi)}1N。然而,这种梯度只有在数据点 { x i } 1 N \{x_i\}^N_1 {xi}1N上定义,不能推广到其他的 x x x值。推广的可能性是选择参数化的类的成员 h ( x , a m ) h(x,a_m) h(x,am),它产生最平行于 − g m ∈ R N -g_m\in R^N −gm∈RN的 h m = { h ( x , a m ) } 1 N h_m=\{h(x,a_m)\}_1^N hm={h(x,am)}1N。这是在数据分布上和 − g m ( x ) -g_m(x) −gm(x)最高度相关的 h ( x , a ) h(x,a) h(x,a)。它可以由解下式得到 a m = arg min a , β ∑ i = 1 N [ − g m ( x i ) − β h ( x i ; a ) ] 2 (11) a_m=\arg\min_{a,\beta}\sum_{i=1}^N[-g_m(x_i)-\beta h(x_i;a)]^2\tag{11} am=arga,βmini=1∑N[−gm(xi)−βh(xi;a)]2(11)
这个约束的负梯度 h ( x ; a m ) h(x;a_m) h(x;am)用来替代最速下降策略中无约束的(7) − g m ( x ) -g_m(x) −gm(x)。特别地,线性搜索(8)为 ρ m = arg min ρ ∑ i = 1 N L ( y i , F m − 1 ( x i ) + ρ h ( x i ; a m ) ) (12) \rho_m=\arg\min_\rho\sum_{i=1}^NL(y_i,F_{m-1}(x_i)+\rho h(x_i;a_m))\tag{12} ρm=argρmini=1∑NL(yi,Fm−1(xi)+ρh(xi;am))(12)且更新的近似为 F m ( x ) = F m − 1 ( x ) + ρ m h ( x ; a m ) F_m(x)=F_{m-1}(x)+\rho_mh(x;a_m) Fm(x)=Fm−1(x)+ρmh(x;am)
基本上,都是通过用 h ( x ; a ) h(x;a) h(x;a)拟合“伪响应” { y ~ i = − g m ( x i ) } 1 N \{\tilde{y}_i=-g_m(x_i)\}_1^N {y~i=−gm(xi)}1N将约束应用于无约束(粗糙)的解,而不是平滑约束下获得解(9)。通过最小二乘函数的最小化(11),这允许替换比较难的最小化问题(9)的函数,这样就只有一个基于最初的标准(12)的参数优化。因此,对于任意求解(11)存在可行的最小二乘算法的 h ( x ; a ) h(x;a) h(x;a),可以使用这种方法结合前向分步加性模型来最小化任何可微的损失 L ( y , F ) L(y,F) L(y,F)。这就引出了下面使用最速下降的(通用)算法。
算法1(梯度提升Gradient Boost)
- F 0 ( x ) = arg min ρ ∑ i = 1 N L ( y i , ρ ) F_0(x)=\arg\min_\rho\sum_{i=1}^NL(y_i,\rho) F0(x)=argminρ∑i=1NL(yi,ρ)
- 对 m = 1 m=1 m=1到 M M M执行:
- y ~ i = − [ ∂ L ( y i , F ( x i ) ) ∂ F ( x i ) ] F ( x ) = F m − 1 ( x ) , i = 1 , N \tilde{y}_i=-\bigg[{\partial L\big(y_i,F(x_i)\big)\over\partial F(x_i)}\bigg]_{F(x)=F_{m-1}(x)},\ \ i=1,N y~i=−[∂F(xi)∂L(yi,F(xi))]F(x)=Fm−1(x), i=1,N
- a m = arg min a , β ∑ i = 1 N [ y ~ i − β h ( x i ; a ) ] 2 a_m=\arg\min_{a,\beta}\sum_{i=1}^N[\tilde{y}_i-\beta h(x_i;a)]^2 am=argmina,β∑i=1N[y~i−βh(xi;a)]2
- ρ m = arg min ρ ∑ i = 1 N L ( y i , F m − 1 ( x i ) + ρ h ( x i ; a m ) ) \rho_m=\arg\min_\rho\sum_{i=1}^NL(y_i,F_{m-1}(x_i)+\rho h(x_i;a_m)) ρm=argminρ∑i=1NL(yi,Fm−1(xi)+ρh(xi;am))
- F m ( x ) = F m − 1 ( x ) + ρ m h ( x ; a m ) F_m(x)=F_{m-1}(x)+\rho_mh(x;a_m) Fm(x)=Fm−1(x)+ρmh(x;am)
- 结束循环,结束算法
注意,任何估计条件期望(给定 x x x)的拟合准则原则上都可以用来估计算法1中第4步的(平滑的)负梯度(7)。由于许多最小二乘算法具有优越的计算性能,因此最小二乘(11)是一种自然选择。
在 y ∈ { − 1 , 1 } y\in\{-1,1\} y∈{−1,1}以及只通过乘积 L ( y , F ) = L ( y F ) L(y,F)=L(yF) L(y,F)=L(yF)取决于 y y y和 F F F的损失函数 L ( y , F ) L(y,F) L(y,F)的特殊情况下,(9)、(10)类似的用最速下降最小化的boosting已经在机器学习文献中指出[Ratsch, Onoda and Muller (1998), Breiman (1999)]。Duffy和Helmbold(1999)巧妙地利用这个类似来激发他们的GeoLev和GeoArc程序。 y F yF yF称为“边界”(margin),并在边界值的空间而不是在函数值 F F F的空间内进行最速下降。后一种方法允许应用于更一般的损失函数,在这种损失函数中不出现边界。Drucker(1997)在AdaBoost算法的背景下采用了一种不同的将回归转换为到分类框架的策略[FreundandSchapire(1996)]。
4. 应用:加性模型
在本节中,梯度增强方法应用于几种常用的损失准则:最小二乘(LS)、最小绝对偏差(LAD)、Huber(M)和logistic二项对数似然(L)。第一个是作为“实现验证”,而其他得到了新的boosting算法。
4.1 最小二乘回归
这里
L
(
y
,
F
)
=
(
y
−
F
)
2
/
2
L(y,F)=(y-F)^2/2
L(y,F)=(y−F)2/2,算法1中第3行的伪响应为
y
~
i
=
y
i
−
F
m
−
1
(
x
i
)
\tilde{y}_i=y_i-F_{m-1}(x_i)
y~i=yi−Fm−1(xi)。所以,第4行只是拟合了当前的残差,并且第5行中的线性搜索得到
ρ
m
=
β
m
\rho_m=\beta_m
ρm=βm,这里
β
m
\beta_m
βm是第4行中最小化的
β
\beta
β。因此,在均方误差损失上的梯度提升产生了一个通常的分步算法,它迭代地拟合了当前的残差。
4.2 最小绝对偏差(LAD)回归
对于损失函数 L ( y , F ) = ∣ y − F ∣ L(y,F)=|y-F| L(y,F)=∣y−F∣,有 y ~ i = − [ ∂ L ( y i , F ( x i ) ) ∂ F ( x i ) ] F ( x ) = F m − 1 ( x ) = s i g n ( y i − F m − 1 ( x i ) ) (13) \tilde{y}_i=-\bigg[{\partial L\big(y_i,F(x_i)\big)\over\partial F(x_i)}\bigg]_{F(x)=F_{m-1}(x)}=sign(y_i-F_{m-1}(x_i))\tag{13} y~i=−[∂F(xi)∂L(yi,F(xi))]F(x)=Fm−1(x)=sign(yi−Fm−1(xi))(13)这意味着 h ( x ; a ) h(x;a) h(x;a)通过最小二乘来拟合算法1中的第4行的当前误差的sign函数。线性搜索变为 ρ m = arg min ρ ∑ i = 1 N ∣ y i − F m − 1 ( x i ) − ρ h ( x i ; a m ) ∣ = arg min ρ ∑ i = 1 N ∣ h ( x i ; a m ) ∣ ⋅ ∣ y i − F m − 1 ( x i ) h ( x i ; a m ) − ρ ∣ = m e d i a n w { y i − F m − 1 ( x i ) h ( x i ; a m ) } 1 N , w i = ∣ h ( x i ; a m ) ∣ (14) \begin{aligned} \rho_m&=\arg\min_\rho\sum_{i=1}^N|y_i-F_{m-1}(x_i)-\rho h(x_i;a_m)|\\ &=\arg\min_\rho\sum_{i=1}^N|h(x_i;a_m)|\cdot |{y_i-F_{m-1}(x_i)\over h(x_i;a_m)}-\rho |\\ &=median_w\bigg\{{y_i-F_{m-1}(x_i)\over h(x_i;a_m)}\bigg\}_1^N\ ,\ \ \ w_i=|h(x_i;a_m)|\\ \end{aligned}\tag{14} ρm=argρmini=1∑N∣yi−Fm−1(xi)−ρh(xi;am)∣=argρmini=1∑N∣h(xi;am)∣⋅∣h(xi;am)yi−Fm−1(xi)−ρ∣=medianw{h(xi;am)yi−Fm−1(xi)}1N , wi=∣h(xi;am)∣(14)
这里 m e d i a n w { ⋅ } median_w\{\cdot\} medianw{⋅}表示关于权重 w i w_i wi的加权中位数。将(13)、(14)的结果带入算法1,再使用基学习器 h ( x ; a ) h(x;a) h(x;a)就得到了对于最小绝对偏差的boosting算法。
4.3 回归树
这里我们考虑一个特殊的情况,即每个基学习器是一个J个叶结点的回归树[Breiman, Friedman, Olshen andStone(1983)]。每一个回归树模型都有自己的加性形式 h ( x ; { b j , R j } 1 J ) = ∑ j = 1 J b j 1 ( x ∈ R j ) (15) h(x;\{b_j,R_j\}_1^J)=\sum_{j=1}^Jb_j1(x\in R_j)\tag{15} h(x;{bj,Rj}1J)=j=1∑Jbj1(x∈Rj)(15)
这里 { R j } 1 J \{R_j\}_1^J {Rj}1J不相交的区域,共同覆盖了预测变量 x x x的所有联合取值空间。这些区域由对应的决策树的叶结点来表示。示性函数 1 ( ⋅ ) 1(\cdot) 1(⋅)表示为真时取1,否则取0。(15)中这些基学习器的参数为系数 { b j } 1 J \{b_j\}_1^J {bj}1J以及定义了区域边界的量 { R j } 1 J \{R_j\}_1^J {Rj}1J。这些是分裂变量和代表树的非终端节点上的分裂的那些变量的值。因为这些区域是不相交的,(15)等价于这样的预测规则:若 x ∈ R j x\in R_j x∈Rj,则 h ( x ) = b j h(x)=b_j h(x)=bj
对于一个回归树,算法1中的第6行的更新就变成了 F m ( x ) = F m − 1 ( x ) + ρ m ∑ j = 1 J b j m 1 ( x ∈ R j m ) (16) F_m(x)=F_{m-1}(x)+\rho_m\sum_{j=1}^Jb_{jm}1(x\in R_{jm})\tag{16} Fm(x)=Fm−1(x)+ρmj=1∑Jbjm1(x∈Rjm)(16)这里的 { R j m } 1 J \{R_{jm}\}_1^J {Rjm}1J是在 m m m次迭代中由树的叶结点定义的区域。通过第4行中的最小二乘来构建这些区域去得到预测的伪响应 { y ~ i } 1 N \{\tilde{y}_{i}\}_1^N {y~i}1N, { b j m } \{b_{jm}\} {bjm}为对应的最小二乘的系数 b j m = a v e x i ∈ R j m y ~ i b_{jm}=ave_{x_i\in R_{jm}}\tilde{y}_{i} bjm=avexi∈Rjmy~i比例因子 ρ m \rho_m ρm为第5行中的线性搜索的解。
(16)式的更新可以等价的表示为 F m ( x ) = F m − 1 ( x ) + ∑ j = 1 J γ j m 1 ( x ∈ R j m ) (17) F_m(x)=F_{m-1}(x)+\sum_{j=1}^J\gamma_{jm}1(x\in R_{jm})\tag{17} Fm(x)=Fm−1(x)+j=1∑Jγjm1(x∈Rjm)(17)其中 γ j m = ρ j m b j m \gamma_{jm}=\rho_{jm}b_{jm} γjm=ρjmbjm。我们可以将(17)式看做是在每步 { 1 ( x ∈ R j m ) } 1 J \{1(x\in R_{jm})\}_1^J {1(x∈Rjm)}1J增加了J个独立的基函数,而不是像(16)式中单独的一个加性函数。因此,在这种情况下,可以通过使用这些独立基函数(17)的最优系数进一步提高拟合的质量。这些最优系数是 { γ j m } 1 J = arg min { γ j } 1 J ∑ i = 1 N L ( y i , F m − 1 ( x i ) + ∑ j = 1 J γ j 1 ( x ∈ R j m ) ) \{\gamma_{jm}\}_1^J=\arg\min_{\{\gamma_{j}\}_1^J}\sum_{i=1}^NL\bigg(y_i,F_{m-1}(x_i)+\sum_{j=1}^J\gamma_j 1(x\in R_{jm})\bigg) {γjm}1J=arg{γj}1Jmini=1∑NL(yi,Fm−1(xi)+j=1∑Jγj1(x∈Rjm))由于回归树产生的区域具有不相交的性质,这就简化为 γ j m = arg min γ ∑ x i ∈ R j m L ( y i , F m − 1 ( x i ) + γ ) (18) \gamma_{jm}=\arg\min_{\gamma}\sum_{x_i\in R_{jm}}L\bigg(y_i,F_{m-1}(x_i)+\gamma \bigg)\tag{18} γjm=argγminxi∈Rjm∑L(yi,Fm−1(xi)+γ)(18)这只是在给定当前近似 F m − 1 ( x ) F_{m-1}(x) Fm−1(x)的情况下,根据损失函数 L L L,每个终端节点区域进行的最优常数更新。
对于LAD回归,(18)式变为
γ
j
m
=
m
e
d
i
a
n
x
i
∈
R
j
m
{
y
i
−
F
m
−
1
(
x
i
)
}
\gamma_{jm}=median_{x_i\in R_{jm}}\{y_i-F_{m-1}(x_i)\}
γjm=medianxi∈Rjm{yi−Fm−1(xi)}即第
m
m
m次迭代时第
j
j
j个终端节点的当前残差的中位数。在每次迭代,回归树基于最小二乘准则来构建当前残差
y
i
−
F
m
−
1
(
x
i
)
y_i-F_{m-1}(x_i)
yi−Fm−1(xi)的sign函数的最佳预测。然后通过在每个派生的终端节点中添加残差的中位数来更新这个近似。
这个算法高度地稳健,决策树只使用了每一个输入变量 x j x_j xj的阶数信息,伪响应 y ~ i \tilde{y}_i y~i只有两个值 y ~ i ∈ { − 1 , 1 } \tilde{y}_i\in\{-1,1\} y~i∈{−1,1},终端结点的更新是基于中位数。另一种方法是建立一个树来直接最小化损失准则 t r e e m ( x ) = arg min J − n o d e t r e e ∑ i = 1 N ∣ y i − F m − 1 ( x i ) − t r e e ( x i ) ∣ tree_m(x)=\arg\min_{J-node\ tree}\sum_{i=1}^N|y_i-F_{m-1}(x_i)-tree(x_i)| treem(x)=argJ−node treemini=1∑N∣yi−Fm−1(xi)−tree(xi)∣ F m ( x ) = F m − 1 ( x ) + t r e e m ( x ) F_m(x)=F_{m-1}(x)+tree_m(x) Fm(x)=Fm−1(x)+treem(x)然而,算法3要快得多,因为它使用最小二乘来诱导树。当构建树的过程中搜索分割时,平方误差损失比平均绝对偏差更新得快得多。
4.4 M回归
M回归技术试图抵抗长尾误差分布和异常值,同时保持对正态分布误差的高效率。考虑Huber损失函数[Huber (1964)] L ( y , F ) = { 1 2 ( y − F ) 2 , ∣ y − F ∣ ⩽ δ δ ( ∣ y − F ∣ − δ / 2 ) , ∣ y − F ∣ > δ (19) L(y,F)=\begin{cases} {1\over2}(y-F)^2, & |y-F|\leqslant\delta\\ \delta(|y-F|-\delta/2), & |y-F|>\delta \end{cases}\tag{19} L(y,F)={21(y−F)2,δ(∣y−F∣−δ/2),∣y−F∣⩽δ∣y−F∣>δ(19)这里的伪响应为 y ~ i = − [ ∂ L ( y i , F ( x i ) ) ∂ F ( x i ) ] F ( x ) = F m − 1 ( x ) = { y i − F m − 1 ( x i ) , ∣ y i − F m − 1 ( x i ) ∣ ⩽ δ δ ⋅ s i g n ( y i − F m − 1 ( x i ) ) , ∣ y i − F m − 1 ( x i ) ∣ > δ \begin{aligned} \tilde{y}_i&=-\bigg[{\partial L\big(y_i,F(x_i)\big)\over\partial F(x_i)}\bigg]_{F(x)=F_{m-1}(x)}\\ &=\begin{cases} y_i-F_{m-1}(x_i), & |y_i-F_{m-1}(x_i)|\leqslant\delta\\ \delta\cdot sign(y_i-F_{m-1}(x_i)), & |y_i-F_{m-1}(x_i)|>\delta \end{cases} \end{aligned} y~i=−[∂F(xi)∂L(yi,F(xi))]F(x)=Fm−1(x)={yi−Fm−1(xi),δ⋅sign(yi−Fm−1(xi)),∣yi−Fm−1(xi)∣⩽δ∣yi−Fm−1(xi)∣>δ线性搜索变为 ρ m = arg min ρ ∑ i = 1 N L ( y i , F m − 1 ( x i ) + ρ h ( x i ; a m ) ) (20) \rho_m=\arg\min_\rho\sum_{i=1}^NL(y_i,F_{m-1}(x_i)+\rho h(x_i;a_m))\tag{20} ρm=argρmini=1∑NL(yi,Fm−1(xi)+ρh(xi;am))(20)这里的 L L L由(19)式给出。(19)、(20)式得解可以由标准的迭代方法得到[see Huber (1964)]。
边界点 δ \delta δ的值受绝对损失函数而不是平方误差损失控制,定义了这些被认为是“离群点”的残差值。最优值将取决于 y − F ∗ ( x ) y-F^*(x) y−F∗(x)的分布,其中 F ∗ F^* F∗是真正的目标函数(1)。选择 δ \delta δ值的常见的做法是选为分布 ∣ y − F ∗ ( x ) ∣ |y-F^*(x)| ∣y−F∗(x)∣的 α \alpha α分位数,这里 ( 1 − α ) (1-\alpha) (1−α)控制了程序的断点。这个“断点”把样本分为一部分,这部分可以任意修改而不会降低结果的质量。因为 F ∗ ( x ) F^*(x) F∗(x)是未知的,我们在第 m m m次迭代中使用当前估计 F m − 1 ( x ) F_{m-1}(x) Fm−1(x)作为近似。 ∣ y − F m − 1 ( x ) ∣ |y-F_{m-1}(x)| ∣y−Fm−1(x)∣的分布由当前残差来估计,这样 δ m = q u a n t i l e α { ∣ y i − F m − 1 ( x i ) ∣ } 1 N \delta_m=quantile_\alpha\{|y_i-F_{m-1}(x_i)|\}_1^N δm=quantileα{∣yi−Fm−1(xi)∣}1N
使用回归树作为基学习器,我们使用4.3节的策略,即在每个终端节点
R
j
m
R_{jm}
Rjm中单独更新(18)。对于Huber损失(19),(18)的解可以近似为从中位数开始的标准迭代过程[Huber(1964)]单个步骤
r
~
j
m
=
m
e
d
i
a
n
x
i
∈
R
j
m
{
r
m
−
1
(
x
i
)
}
\tilde{r}_{jm}=median_{x_i\in R_{jm}}\{r_{m-1}(x_i)\}
r~jm=medianxi∈Rjm{rm−1(xi)}其中
{
r
m
−
1
(
x
i
)
}
1
N
\{r_{m-1}(x_i)\}_1^N
{rm−1(xi)}1N为当前残差
r
m
−
1
(
x
i
)
=
y
i
−
F
m
−
1
(
x
i
)
r_{m-1}(x_i)=y_i-F_{m-1}(x_i)
rm−1(xi)=yi−Fm−1(xi)近似为
γ
j
m
=
r
~
j
m
+
1
N
j
m
∑
x
i
∈
R
j
m
s
i
g
n
(
r
m
−
1
(
x
i
)
−
r
~
j
m
)
⋅
min
(
δ
m
,
a
b
s
(
r
m
−
1
(
x
i
)
−
r
~
j
m
)
)
\gamma_{jm}=\tilde{r}_{jm}+{1\over N_{jm}}\sum_{x_i\in R_{jm}}sign(r_{m-1}(x_i)-\tilde{r}_{jm})\cdot \min(\delta_m,abs(r_{m-1}(x_i)-\tilde{r}_{jm}))
γjm=r~jm+Njm1xi∈Rjm∑sign(rm−1(xi)−r~jm)⋅min(δm,abs(rm−1(xi)−r~jm))这里的
N
j
m
N_{jm}
Njm是在第
j
j
j个终端结点中的样本数。这样就给出了以下基于Huber损失的boosting回归树算法。
根据稳健回归背后的动机,该算法应该具有与正态分布误差的最小二乘boosting算法(Algorithm2)相似的性质,与极长尾分布的最小绝对偏差回归(算法3)相似的性质。对于只有中等长度尾部的误差分布,它的性能可能优于两者(参见第6.2节)。
4.5 二分类逻辑回归和分类
这里的损失函数为负二项对数似然(FHT00)
L
(
y
,
F
)
=
l
o
g
(
1
+
e
x
p
(
−
2
y
F
)
)
,
y
∈
{
−
1
,
1
}
L(y,F)=log(1+exp(-2yF)),\ \ y\in\{-1,1\}
L(y,F)=log(1+exp(−2yF)), y∈{−1,1}其中
F
(
x
)
=
1
2
l
o
g
[
P
r
(
y
=
1
∣
x
)
P
r
(
y
=
−
1
∣
x
)
]
(21)
F(x)={1\over2}log\bigg[{Pr(y=1|x)\over Pr(y=-1|x)}\bigg]\tag{21}
F(x)=21log[Pr(y=−1∣x)Pr(y=1∣x)](21)伪响应为
y
~
i
=
−
[
∂
L
(
y
i
,
F
(
x
i
)
)
∂
F
(
x
i
)
]
F
(
x
)
=
F
m
−
1
(
x
)
=
2
y
i
/
(
1
+
e
x
p
(
2
y
i
F
m
−
1
(
x
)
)
)
(22)
\tilde{y}_i=-\bigg[{\partial L(y_i,F(x_i))\over \partial F(x_i)}\bigg]_{F(x)=F_{m-1}(x)}=2y_i\big/\bigg(1+exp\big(2y_iF_{m-1}(x)\big)\bigg)\tag{22}
y~i=−[∂F(xi)∂L(yi,F(xi))]F(x)=Fm−1(x)=2yi/(1+exp(2yiFm−1(x)))(22)线性搜索变为
ρ
m
=
arg
min
ρ
∑
i
=
1
N
l
o
g
(
1
+
e
x
p
(
−
2
y
i
(
F
m
−
1
(
x
i
)
+
ρ
h
(
x
i
;
a
m
)
)
)
)
\rho_m=\arg\min_\rho\sum_{i=1}^Nlog\bigg(1+exp\big(-2y_i(F_{m-1}(x_i)+\rho h(x_i;a_m))\big)\bigg)
ρm=argρmini=1∑Nlog(1+exp(−2yi(Fm−1(xi)+ρh(xi;am))))用回归树作为基学习器,我们再次使用4.3节中的办法,单独更新每个终端结点
R
j
m
R_{jm}
Rjm:
γ
j
m
=
arg
min
γ
∑
x
i
∈
R
j
m
l
o
g
(
1
+
e
x
p
(
−
2
y
i
(
F
m
−
1
(
x
i
)
+
γ
)
)
)
(23)
\gamma_{jm}=\arg\min_\gamma\sum_{x_i\in R_{jm}}log\bigg(1+exp\big(-2y_i(F_{m-1}(x_i)+\gamma)\big)\bigg)\tag{23}
γjm=argγminxi∈Rjm∑log(1+exp(−2yi(Fm−1(xi)+γ)))(23)(23)式没有闭式解。在负二项对数似然的损失函数下,我们用一个Newton–Raphson迭代步来替代,即
γ
j
m
=
∑
x
i
∈
R
j
m
y
~
i
/
∑
x
i
∈
R
j
m
∣
y
~
i
∣
(
2
−
∣
y
~
i
∣
)
\gamma_{jm}=\sum_{x_i\in R_{jm}}\tilde{y}_i\bigg/\sum_{x_i\in R_{jm}}|\tilde{y}_i|(2-|\tilde{y}_i|)
γjm=xi∈Rjm∑y~i/xi∈Rjm∑∣y~i∣(2−∣y~i∣)其中
y
~
i
\tilde{y}_i
y~i由(22)式给出。这就给出了下面用回归树的似然梯度boosting算法。
最后的近似
F
M
(
x
)
F_M(x)
FM(x)通过(21)式和对数几率联系起来,这样可以转化得到概率估计
p
+
(
x
)
=
P
r
^
(
y
=
1
∣
x
)
=
1
/
(
1
+
e
−
2
F
M
(
x
)
)
p
−
(
x
)
=
P
r
^
(
y
=
−
1
∣
x
)
=
1
/
(
1
+
e
2
F
M
(
x
)
)
p_+(x)=\hat{Pr}(y=1|x)=1/(1+e^{-2F_M(x)})\\ p_-(x)=\hat{Pr}(y=-1|x)=1/(1+e^{2F_M(x)})
p+(x)=Pr^(y=1∣x)=1/(1+e−2FM(x))p−(x)=Pr^(y=−1∣x)=1/(1+e2FM(x))这些反过来又可以用于分类
y
^
(
x
)
=
2
⋅
1
[
c
(
−
1
,
1
)
p
+
(
x
)
>
c
(
1
,
−
1
)
p
−
(
x
)
]
−
1
\hat{y}(x)=2\cdot1[c(-1,1)p_+(x)>c(1,-1)p_-(x)]-1
y^(x)=2⋅1[c(−1,1)p+(x)>c(1,−1)p−(x)]−1其中
c
(
y
^
,
y
)
c(\hat y,y)
c(y^,y)是当真实值为
y
y
y预测值为
y
^
\hat y
y^时有关的损失。
4.5.1 影响修剪(权重修剪)
对于二分类的逻辑回归问题在
m
m
m次迭代时的经验损失函数为
ϕ
m
(
ρ
,
a
)
=
∑
i
=
1
N
l
o
g
[
1
+
e
x
p
(
−
2
y
i
F
m
−
1
(
x
i
)
)
⋅
e
x
p
(
−
2
y
i
ρ
h
(
x
i
;
a
)
)
]
(24)
\phi_m(\rho,a)=\sum_{i=1}^Nlog[1+exp(-2y_iF_{m-1}(x_i))\cdot exp(-2y_i\rho h(x_i;a))]\tag{24}
ϕm(ρ,a)=i=1∑Nlog[1+exp(−2yiFm−1(xi))⋅exp(−2yiρh(xi;a))](24)
若
y
i
F
m
−
1
(
x
i
)
y_iF_{m-1}(x_i)
yiFm−1(xi)非常大,则对于接近于零的很小到一般的数(24)式都几乎不相关于
ρ
h
(
x
i
;
a
)
\rho h(x_i;a)
ρh(xi;a)。这意味着第
i
i
i个样本
(
y
i
,
x
i
)
(y_i,x_i)
(yi,xi)对损失函数几乎没有影响,因此在其解也没有影响
(
ρ
m
,
a
m
)
=
arg
min
ρ
,
a
ϕ
m
(
ρ
,
a
)
(\rho_m,a_m)=\arg\min_{\rho,a}\phi_m(\rho,a)
(ρm,am)=argρ,aminϕm(ρ,a)这表明可以从第
m
m
m次迭代的所有计算中删除所有的那些
y
i
F
m
−
1
(
x
i
)
y_iF_{m-1}(x_i)
yiFm−1(xi)相对很大的样本
(
y
i
,
x
i
)
(y_i,x_i)
(yi,xi),而不会对结果产生实质性的影响。因此
w
i
=
e
x
p
(
−
2
y
i
F
m
−
1
(
x
i
)
)
(25)
w_i=exp(-2y_iF_{m-1}(x_i))\tag{25}
wi=exp(−2yiFm−1(xi))(25)可以看作是第
i
i
i个样本在估计
ρ
h
(
x
;
a
m
)
\rho h(x;a_m)
ρh(x;am)上的影响或者权重的一种度量。
更一般地,从第2节的非参数函数空间角度来看,参数为样本函数值 { F ( x i ) } 1 N \{F(x_i)\}_1^N {F(xi)}1N。对于“参数”值 F ( x i ) F(x_i) F(xi)(保持所有其他参数固定)的变化,估计的影响可以用损失函数关于这个参数求二阶导数来衡量。这里在第 m m m次迭代中的二阶导数为 ∣ y ~ i ∣ ( 2 − ∣ y ~ i ∣ ) |\tilde{y}_i|(2-|\tilde{y}_i|) ∣y~i∣(2−∣y~i∣),其中 ∣ y ~ i ∣ |\tilde{y}_i| ∣y~i∣由(22)给出。因此,另一种衡量第 m m m次迭代中第 i i i个样本在估计 ρ h ( x ; a m ) \rho h(x;a_m) ρh(x;am)上的影响或者权重的方法为 w i = ∣ y ~ i ∣ ( 2 − ∣ y ~ i ∣ ) w_i=|\tilde{y}_i|(2-|\tilde{y}_i|) wi=∣y~i∣(2−∣y~i∣)
影响修剪就是删除那些 w i w_i wi值小于 w l ( α ) w_{l(\alpha)} wl(α)的所有样本,这里 l ( α ) l(\alpha) l(α)是下式的解 ∑ i = 1 l ( α ) w i = α ∑ I = 1 N w i (27) \sum_{i=1}^{l(\alpha)}w_i=\alpha\sum_{I=1}^Nw_i\tag{27} i=1∑l(α)wi=αI=1∑Nwi(27)这里 { w ( i ) } 1 N \{w_{(i)}\}_1^N {w(i)}1N是将权重 { w i } 1 N \{w_{i}\}_1^N {wi}1N以升序排列,特殊值是 α ∈ [ 0.05 , 0.2 ] \alpha\in[0.05,0.2] α∈[0.05,0.2]。注意,基于(25)、(27)的影响修剪与Real AdaBoost使用的“权重修剪”策略相同,而(26)、(27)与FHT00中使用的LogitBoost相同。在那里可以看到,使用任何一种影响度量,90%到95%的样本常被删除,而不牺牲估计的准确性。这使得计算量相应减少了10到20倍。
4.6 多分类逻辑回归和分类
这里我们得到了一个 K K K分类问题的梯度下降boosting算法。损失函数为 L ( { y k , F k ( x ) } 1 K ) = − ∑ k = 1 K y k log p k ( x ) (28) L(\{y_k,F_k(x)\}_1^K)=-\sum_{k=1}^Ky_k\log p_k(x)\tag{28} L({yk,Fk(x)}1K)=−k=1∑Kyklogpk(x)(28)其中 y k = 1 ( c l a s s = k ) ∈ { 0 , 1 } y_k=1(class=k)\in\{0,1\} yk=1(class=k)∈{0,1}且 p k ( x ) = P r ( y k = 1 ∣ x ) p_k(x)=Pr(y_k=1|x) pk(x)=Pr(yk=1∣x)。在FHT00之后,我们使用对称的多重逻辑变换 F k ( x ) = log p k ( x ) − 1 K ∑ l = 1 K log p l ( x ) (29) F_k(x)=\log p_k(x)-{1\over K}\sum_{l=1}^K\log p_l(x)\tag{29} Fk(x)=logpk(x)−K1l=1∑Klogpl(x)(29)或者等价地 p k ( x ) = exp ( F k ( x ) ) / ∑ l = 1 K exp ( F l ( x ) ) (30) p_k(x)=\exp(F_k(x))\bigg/\sum_{l=1}^K\exp(F_l(x))\tag{30} pk(x)=exp(Fk(x))/l=1∑Kexp(Fl(x))(30)将(30)代入(28),求一阶导数得到 y ~ i k = − [ ∂ L ( { y i l , F l ( x i ) } l = 1 K ) ∂ F k ( x i ) ] { F l ( x ) = F l , m − 1 ( x ) } 1 K = y i k − p k , m − 1 ( x i ) (31) \tilde{y}_{ik}=-\bigg[{\partial L(\{y_{il},F_l(x_i)\}_{l=1}^K)\over \partial F_k(x_i)}\bigg]_{\{F_l(x)=F_{l,m-1}(x)\}_1^K}=y_{ik}-p_{k,m-1}(x_i)\tag{31} y~ik=−[∂Fk(xi)∂L({yil,Fl(xi)}l=1K)]{Fl(x)=Fl,m−1(x)}1K=yik−pk,m−1(xi)(31)这里 p k , m − 1 ( x ) p_{k,m-1}(x) pk,m−1(x)是通过(30)式从 F k , m − 1 ( x ) F_{k,m-1}(x) Fk,m−1(x)中推得。因此,在每个迭代 m m m中引入 K K K树来预测每个类在概率尺度上的相应的当前残差。这些树每个都有 J J J个终端结点,对应区域为 { R j k m } j = 1 J \{R_{jkm}\}_{j=1}^J {Rjkm}j=1J。模型更新 γ j k m \gamma_{jkm} γjkm相当于这些区域是下式的解 { γ j k m } = arg min { γ j k } ∑ i = 1 N ∑ k = 1 K ϕ ( y i k , F k , m − 1 ( x i ) + ∑ j = 1 J γ j k 1 ( x i ∈ R j m ) ) \{\gamma_{jkm}\}=\arg\min_{\{\gamma_{jk}\}}\sum_{i=1}^N\sum_{k=1}^K\phi\bigg(y_{ik},F_{k,m-1}(x_i)+\sum_{j=1}^J\gamma_{jk}1(x_i\in R_{jm})\bigg) {γjkm}=arg{γjk}mini=1∑Nk=1∑Kϕ(yik,Fk,m−1(xi)+j=1∑Jγjk1(xi∈Rjm))其中由(28)式 ϕ ( y k , F k ) = − y k log p k \phi(y_k,F_k)=-y_k\log p_k ϕ(yk,Fk)=−yklogpk, F k F_k Fk由(30)式和 p k p_k pk有关。它没有闭式解,此外,不同类别的树所对应的区域是重叠的,这样解就不会像(18)那样,在每棵树的每个区域中进行单独的计算。在FHT00之后,我们用单个Newton–Raphson步骤来近似解,使用对海森的对角近似。这将问题分解为每个树的每个终端节点的单独计算。结果为 γ j k m = K − 1 K ∑ x i ∈ R j k m y ~ i k ∑ x i ∈ R j k m ∣ y ~ i k ∣ ( 1 − ∣ y ~ i k ∣ ) (32) \gamma_{jkm}={K-1\over K}{\sum_{x_i\in R_{jkm}}\tilde{y}_{ik}\over\sum_{x_i\in R_{jkm}}|\tilde{y}_{ik}|(1-|\tilde{y}_{ik}|)}\tag{32} γjkm=KK−1∑xi∈Rjkm∣y~ik∣(1−∣y~ik∣)∑xi∈Rjkmy~ik(32)
这样就得到了下面的
K
K
K分类问题的逻辑梯度下降boosting算法。
最后的估计
{
F
k
M
(
x
)
}
1
K
\{F_{kM}(x)\}_1^K
{FkM(x)}1K可以通过(30)式用来得到对应的概率估计
{
p
k
M
(
x
)
}
1
K
\{p_{kM}(x)\}_1^K
{pkM(x)}1K。这些反过来又可以用于分类
k
^
(
x
)
=
arg
min
1
⩽
k
⩽
K
∑
k
′
=
1
K
c
(
k
,
k
′
)
p
k
′
M
(
x
)
\hat k(x)=\arg\min_{1\leqslant k\leqslant K}\sum_{k'=1}^Kc(k,k')p_{k'M}(x)
k^(x)=arg1⩽k⩽Kmink′=1∑Kc(k,k′)pk′M(x)这里
c
(
k
,
k
′
)
c(k,k')
c(k,k′)表示当真实是
k
′
k'
k′预测为
k
k
k类时的损失。注意对
K
=
2
K=2
K=2,算法6等价于算法5。
算法6与FHT00的 K K K分类的LogitBoost过程相似,该过程是基于Newton–Raphson而不是函数空间中的梯度下降。在该算法中, K K K棵树被诱导,每棵树都对每个样本 ( y ~ i k , x i ) (\tilde{y}_{ik},x_i) (y~ik,xi)使用相应的伪响应 y ~ i k = K − 1 K y i k − p k ( x i ) p k ( x i ) ( 1 − p k ( x i ) ) (33) \tilde{y}_{ik}={K-1\over K}{y_{ik}-p_k(x_i)\over p_k(x_i)(1-p_k(x_i))}\tag{33} y~ik=KK−1pk(xi)(1−pk(xi))yik−pk(xi)(33)以及权重 w k ( x i ) = p k ( x i ) ( 1 − p k ( x i ) ) (34) w_k(x_i)=p_k(x_i)(1-p_k(x_i))\tag{34} wk(xi)=pk(xi)(1−pk(xi))(34)终端结点的更新为 γ j k m = ∑ x i ∈ R j k m w k ( x i ) y ~ i k ∑ x i ∈ R j k m w k ( x i ) \gamma_{jkm}={\sum_{x_i\in R_{jkm}}w_k(x_i)\tilde{y}_{ik}\over \sum_{x_i\in R_{jkm}}w_k(x_i)} γjkm=∑xi∈Rjkmwk(xi)∑xi∈Rjkmwk(xi)y~ik这和(32)式等价。这两种算法的区别在于用于诱导树和最终区域的分裂准则 { R j k m } 1 J \{R_{jkm}\}_1^J {Rjkm}1J。
用于评估一个当前终端区域 R R R可能分裂为两个子区域 ( R l , R r ) (R_l,R_r) (Rl,Rr)的最小二乘改进准则为 i 2 ( R l , R r ) = w l w r w l + w r ( y ˉ l − y ˉ r ) 2 (35) i^2(R_l,R_r)={w_lw_r\over w_l+w_r}(\bar y_l-\bar y_r)^2\tag{35} i2(Rl,Rr)=wl+wrwlwr(yˉl−yˉr)2(35)其中 y ˉ l , y ˉ r \bar y_l,\bar y_r yˉl,yˉr分别为左、右两边的子响应均值, w l , w r w_l,w_r wl,wr为对应的权重之和。对于一个给定的分裂,利用带单位权重的(31)式或者带权重(34)式的(33)式,得到 y ˉ l , y ˉ r \bar y_l,\bar y_r yˉl,yˉr有相同的值。然而,权重之和 w l , w r w_l,w_r wl,wr是不同的。单位权重( L K L_K LK_TreeBoost)偏向在每个子节点中关于样本数对称的分裂,然而(34)(LogitBoost)偏向使得当前估计的响应方差 v a r ( y i k ) = p k ( x i ) ( 1 − p k ( x i ) ) var (y_{ik})=p_k(x_i)(1-p_k(x_i)) var(yik)=pk(xi)(1−pk(xi))之和更相等的分裂。
在数值稳定性上 L K L_K LK_TreeBoost具有实际优势,当对于任何样本 x i x_i xi,(34)的值接近于0时,LogitBoost在数值上变得不稳定,这种现象经常发生。这是Newton–Raphson所面临的二阶导数消失的困难的结果。它的性能受到处理这个问题的方式的强烈影响(参见FHT00)。 L K L_K LK_TreeBoost只有当(34)对一个终端节点的所有样本值接近于0时才会出现这种困难。这种情况发生的频率要低得多,发生时也更容易处理。
多分类过程的影响修剪与第4.5.1节中概述的二分类情况的影响修剪的实现方法相同。与每个“样本” ( y i k , x i ) (y_{ik},x_i) (yik,xi)有关的是一个影响 w i k = ∣ y ~ i k ∣ ( 1 − ∣ y ~ i k ∣ ) w_{ik}=|\tilde{y}_{ik}|(1-|\tilde{y}_{ik}|) wik=∣y~ik∣(1−∣y~ik∣),当诱导第 k k k颗树时它用于在当前迭代 m m m中删除样本(27)。
5. 正则化
在预测问题中,过于紧密地拟合训练数据可能会适得其反。将除了一些点的训练数据上的期望损失减少,将导致分布的期望损失停止下降,并常常开始增加。正则化方法试图通过约束拟合过程来防止这种“过拟合”。对于加性展式(2),一个自然的正则化参数是项的数量 M M M,这类似于分步回归,其中 { h ( x ; a m ) } 1 M \{h(x;a_m)\}_1^M {h(x;am)}1M被认为是依次加入的解释变量。控制 M M M的值可以调节训练数据上期望损失最小化的程度。 M M M的最佳值可以通过一些模型选择方法来估计,比如使用一个独立的“测试”集,或者交叉验证。
通过控制展式中的项的数量来进行正则化,这就隐含了一种先验信念,即“稀疏”近似的项数越少,预测效果越好。然而,人们经常发现,通过收缩实现的正则化比通过限制项的数量实现的正则化效果更好[Copas(1983)]。在以前向分步规则(9)(10)构造的可加性模型(2)中,一个简单的收缩策略是将通用算法(算法1)的第6行替换为 F m ( x ) = F m − 1 ( x ) + ν ⋅ ρ m h ( x ; a m ) , 0 < ν ⩽ 1 (36) F_m(x)=F_{m-1}(x)+\nu\cdot\rho_mh(x;a_m),\ \ \ 0<\nu\leqslant 1\tag{36} Fm(x)=Fm−1(x)+ν⋅ρmh(x;am), 0<ν⩽1(36)且在所有特定的算法(算法2-6)中作相应的等价变化。每次更新都是简单地根据“学习率”参数 ν \nu ν的值进行缩放。
以这种方式引入收缩到梯度提升(36),提供了两种正则化参数,学习速率 ν \nu ν和项的数量 M M M。任何一个都可以控制拟合的程度,从而影响到另一个的最佳值。减少 ν \nu ν的值会增加 M M M的最佳值。理想地,我们应该同时最小化关于这两个参数的模型选择准则来估计这两个的最优值。这里还有计算方面的考虑,增加 M M M的大小会使计算量按比例增加。
我们通过模拟实验来阐述这样的 ν − M \nu-M ν−M之间的权衡。训练集有5000个样本 { y i , x i } \{y_i,x_i\} {yi,xi}组成,其中 y i = F ∗ ( x i ) + ϵ i y_i=F^*(x_i)+\epsilon_i yi=F∗(xi)+ϵi目标函数 F ∗ ( x i ) , x ∈ R 10 F^*(x_i),\ x\in R^{10} F∗(xi), x∈R10是由6.1节中描述那样随机生成的,噪音 ϵ \epsilon ϵ是由零均值的正态分布生成,方差调整以满足在给定2/1的信噪比下 E ∣ ϵ ∣ = 1 2 E x ∣ F ∗ ( x ) − m e d i a n x F ∗ ( x ) ∣ E|\epsilon|={1\over2}E_x|F^*(x)-median_xF^*(x)| E∣ϵ∣=21Ex∣F∗(x)−medianxF∗(x)∣对于这个阐述,基学习器 h ( x ; a ) h(x;a) h(x;a)取为11个终端节点的以最佳优先的方式诱导(FHT00)的回归树。关于树大小选择的一般性讨论在第7节。
图1显示了对于不同的收缩参数 ν ∈ { 1.0 , 0.25 , 0.125 , 0.06 } \nu\in\{1.0,0.25,0.125,0.06\} ν∈{1.0,0.25,0.125,0.06}, L S LS LS_TreeBoost、 L A D LAD LAD_TreeBoost和 L 2 L_2 L2_TreeBoost的拟合损失(LOF)作为项的数量(迭代次数) M M M的函数。对前两种方法,LOF是由估计 F ^ M ( x ) \hat F_M(x) F^M(x)相对于最优常数解的平均绝对误差来度量的 A ( F ^ M ( x ) ) = E x ∣ F ∗ ( x ) − F ^ M ( x ) ∣ E x ∣ F ∗ ( x ) − m e d i a n x F ∗ ( x ) ∣ A(\hat F_M(x))={E_x|F^*(x)-\hat F_M(x)|\over E_x|F^*(x)-median_xF^*(x)|} A(F^M(x))=Ex∣F∗(x)−medianxF∗(x)∣Ex∣F∗(x)−F^M(x)∣对于逻辑回归的 y y y值,是在 x x x值的分布上把 F ∗ ( x ) F^*(x) F∗(x)的中位数作为阈值来得到的; F ∗ ( x i ) F^*(x_i) F∗(xi)的值比中位数大则 y i = 1 y_i=1 yi=1,比中位数小则 y i = − 1 y_i=-1 yi=−1。因此贝叶斯错误率为零,但是决策边界相当的复杂。对 L 2 L_2 L2_TreeBoost有两种LOF计算方法;负二倍对数似然(偏差)和分类错误率 E x [ 1 ( y ≠ s i g n ( F ^ M ( x ) ) ) ] E_x[1(y\neq sign(\hat F_M(x)))] Ex[1(y=sign(F^M(x)))]。所有LOF测量值都是通过使用一个包含10,000个观测值的独立验证集来计算的。
如图1所示,小的收缩参数 ν \nu ν(收缩)被认为导致更好的性能,尽管存在收益递减的最小值。对于较大的值,可以观察到过拟合的行为特征;性能达到 M M M的某个值时达到最优,之后随着 M M M的增加而降低。这种效果在 L A D LAD LAD_TreeBoost和 L 2 L_2 L2_TreeBoost的错误率标准下不太明显。对较小的 ν \nu ν,可以减少过度拟合,正如所预期的那样。
尽管除了 ν = 1 \nu=1 ν=1很难看到分类错误率(右下)在逻辑似然可能已经达到最佳(左下)继续减少。因此,通过过拟合来降低似然实际上改善了分类错误率。尽管这可能违反直觉,但这并不矛盾;似然和错误率度量了拟合质量的不同方面。错误率只取决于 F ^ M ( x ) \hat F_M(x) F^M(x)的符号而似然受其符号和量级的影响。显然,过拟合降低了量级估计的质量而不影响(有时改进)符号。因此,分类错误率对过度拟合的敏感性要小得多。
表1总结了包括如图1所示的几个 ν \nu ν值的模拟结果。显示了每个 ν \nu ν值(行)到达最低LOF的的迭代次数和相应的最小值(双列)。
ν − M \nu-M ν−M之间的权衡一览无遗; ν \nu ν的较小产生较大的最佳 M M M值,他们还提供更高的精度, ν < 0.125 \nu<0.125 ν<0.125时收益递减。分类错误率当 M ≳ 200 M\gtrsim200 M≳200非常平缓,因此最佳 M M M值不稳定。
虽然这里只阐述了一个目标函数和基学习器(11个终端节点树),但这些结果的定性性质是相当普遍的,其他的目标函数和树的大小(未显示)会导致相同的行为。这表明 ν \nu ν的最佳值取决于迭代次数 M M M。后者应该在计算方便或者可行的时候尽可能大。 ν \nu ν值应该调整使得LOF达到其最小值,接近选择的 M M M值。如果LOF在最后迭代时仍在减小, ν \nu ν值或迭代次数 M M M应该增加,最好是后者。考虑到算法的顺序性,它可以很容易地在以前完成的地方重新启动,因此不需要重复计算。LOF作为迭代次数的函数,最方便的估计方法是使用遗漏的测试样本。
如这里所示,降低学习速率可以明显地提高性能,通常是显著地提高性能。其原因尚不清楚。在每次迭代中收缩模型更新(36)产生的效果比整个模型的成比例收缩更复杂 F ^ ν ( x ) = y ˉ + ν ⋅ ( F ^ M ( x ) − y ˉ ) \hat F_\nu(x)=\bar y+\nu\cdot(\hat F_M(x)-\bar y) F^ν(x)=yˉ+ν⋅(F^M(x)−yˉ)其中 F ^ M ( x ) \hat F_M(x) F^M(x)是没有收缩的模型得到的。每个迭代中的更新 ρ m h ( x ; a m ) \rho_mh(x;a_m) ρmh(x;am)取决于前一个迭代的特定更新序列。增量收缩(36)产生的模型与全量收缩(38)非常不同。经验证据(未显示)表明,全量收缩(38)在没有收缩的情况下最多只能提供边际改善,远远达不到增量收缩的显著效果。增量收缩成功背后的神秘目前正在研究中。
6. 模拟研究
任何函数估计方法的性能取决于它所应用的特定问题。重要的影响性能的问题的特征包括训练样本大小 N N N,真正的潜在的“目标”函数 F ∗ ( x ) F^*(x) F∗(x)(1)以及来自 F ∗ ( x ) F^*(x) F∗(x)的 y ∣ x y|x y∣x的 ϵ \epsilon ϵ的分布。对于任何给定的问题, N N N总是知道,有时 ϵ \epsilon ϵ的分布也是知道的,例如当 y y y是二元(伯努利)。当 y y y是一个广义的实值变量 ϵ \epsilon ϵ的分布很少可知。在几乎所有情况下, F ∗ ( x ) F^*(x) F∗(x)的本质是未知的。
为了评估任何一种估计方法,有必要准确地评估它在许多不同情况下的表现。这是通过蒙特卡罗模拟最方便的实现,在蒙特卡罗模拟中,根据惯例的种类和精确计算的结果性能来生成数据。在本节中,描述了这样的几种研究,以试图理解前面几节中介绍的梯度增强程序的特性。尽管这样的研究远比仅仅通过几个选定的例子来评价这些方法要深入,但真实的或模拟的结果都具有一定的启发性。
6.1 随机函数生成器
任何问题影响性能的最重要的特征之一是真正的潜在目标函数 F ∗ ( x ) F^*(x) F∗(x)(1),每种方法都有特定的目标,它是最合适的而其他就不是。由于目标函数的性质在不同的问题上可能有很大的不同,而且很少可知,所以我们比较了回归树梯度增强算法在各种随机生成目标函数上的优点。每个都有自己的形式 F ∗ ( x ) = ∑ l = 1 20 a l g l ( z l ) (39) F^*(x)=\sum_{l=1}^{20}a_lg_l(z_l)\tag{39} F∗(x)=l=1∑20algl(zl)(39)系数 { a l } 1 20 \{a_l\}_1^{20} {al}120是从均匀分布中随机生成的,即 a l ∼ U [ − 1 , 1 ] a_l\sim U[-1,1] al∼U[−1,1]。每个 g l ( z l ) g_l(z_l) gl(zl)都是一个大小为 n l n_l nl,有 n n n个输入变量 x x x的随机选择的子集的函数。特别地, z l = { x P l ( j ) } j = 1 n l z_l=\{x_{P_l(j)}\}_{j=1}^{n_l} zl={xPl(j)}j=1nl这里 P l P_l Pl是一个独立的随机排列的整数 { 1 , 2 , … , n } \{1,2,\dots,n\} {1,2,…,n},子集的大小 n l {n_l} nl也是一个随机数, n l = ⌊ 1.5 + r ⌋ {n_l}=\lfloor1.5+r\rfloor nl=⌊1.5+r⌋, r r r是从均值 λ = 2 \lambda=2 λ=2的指数分布中抽取的。因此,对每个 g l ( z l ) g_l(z_l) gl(zl)的输入变量的预期数目在3到4之间。然而,大多数情况下会更少。这反映了对强的非常高阶交互作用影响的偏见。然而,对于已实现的 F ∗ ( x ) F^*(x) F∗(x),20个函数 g l ( z l ) g_l(z_l) gl(zl)中很大可能至少有几个将涉及高阶交互作用。在任何情况下 F ∗ ( x ) F^*(x) F∗(x)都会是所有,或者几乎所有输入变量的函数。
每个 g l ( z l ) g_l(z_l) gl(zl)都是一个 n l {n_l} nl维的高斯函数 g l ( z l ) = exp ( − 1 2 ( z l − μ l ) T V l ( z l − μ l ) ) g_l(z_l)=\exp\big(-{1\over2}(z_l-\mu_l)^TV_l(z_l-\mu_l)\big) gl(zl)=exp(−21(zl−μl)TVl(zl−μl))其中每个均值向量 { μ l } 1 20 \{\mu_l\}_1^{20} {μl}120是和输入变量 x x x相同的分布中随机生成的。 n l × n l n_l\times n_l nl×nl的协方差矩阵 V l V_l Vl也是随机生成的。特别地, V l = U l D l U l T V_l=U_lD_lU_l^T Vl=UlDlUlT这里 U l U_l Ul是随机正交矩阵(均匀或者Haar测度下)且 D l = d i a g { d 1 l … d n l l } D_l=diag \{d_{1l}\dots d_{n_ll}\} Dl=diag{d1l…dnll},特征根的平方根是从均匀分布中随机生成的,即 d j l ∼ U [ a , b ] \sqrt{d_{jl}}\sim U[a,b] djl∼U[a,b],其中 a , b a,b a,b取决于输入变量 x x x的分布。
对于这里提出的所有研究,输入变量的数量取为 n = 10 n = 10 n=10,且他们的联合分布取为标准正态分布 x ∼ N ( 0 , I ) x\sim N(0,I) x∼N(0,I)。特征值的限制为 a = 0.1 , b = 2.0 a=0.1,b= 2.0 a=0.1,b=2.0。虽然正态分布的尾端往往比实际遇到的数据尾端要短,但仍然比模拟研究中常用的均匀分布输入更具有真实感。此外,回归树不受长尾输入变量分布的影响,因此短尾在比较中具有相对优势。
在下面的模拟研究中,根据上述规则(39),(40)随机生成100个目标函数 F ∗ ( x ) F^*(x) F∗(x)。性能是根据这些不同目标的近似误差分布[相对近似误差(37)或误分类风险]来评估的。这种方法允许在十维输入空间中根据轮廓的形状生成各种不同的目标函数。虽然低阶交互作用是有利的,但这些函数并不特别适合加性回归树。决策树生成的是张量积基函数,但目标函数 F ∗ ( x ) F^*(x) F∗(x)的组成 g l ( z l ) g_l(z_l) gl(zl)不是张量积函数。使用第8节中描述的方法,第一个随机生成的函数对一些更重要的参数的相关关系的可视化显示在第8.3节中。
虽然只有10个输入变量,但每个目标都是它们的函数。在许多数据挖掘应用程序中,有许多超过10个输入。然而,相关的维度是输入空间的内在维度以及实际影响输出响应变量 y y y的输入的数量。在许多输入变量的问题中,通常存在很多变量之间的共线性度高的问题,并且独立变量的大致数量(近似内在维度)要小得多。而且,目标函数通常只依赖于所有输入的一小部分。
6.2 误差分布
在本节中, L S LS LS_TreeBoost、 L A D LAD LAD_TreeBoost和 M M M_TreeBoost在100个目标函数上通过两种不同的误差分布对它们的性能进行了比较。所有的算法都使用了拥有11个终端节点的最佳优先回归树。 M M M_TreeBoost的参数分解设置为默认值 α = 0.9 \alpha=0.9 α=0.9,学习速率参数(36)在这里的模拟中都设置为 ν = 0.1 \nu=0.1 ν=0.1。
根据下式生成一百个数据集 { y i , x i } 1 N \{y_i,x_i\}_1^N {yi,xi}1N y i = F ∗ ( x i ) + ϵ i y_i=F^*(x_i)+\epsilon_i yi=F∗(xi)+ϵi其中 F ∗ ( x ) F^*(x) F∗(x)表示100个6.1节中描述的随机生成的每一个目标函数。对第一个实验,误差有一个零均值的正态分布生成,方差调整以满足给定1/1的信噪比下 E ∣ ϵ ∣ = E x ∣ F ∗ ( x ) − m e d i a n x F ∗ ( x ) ∣ E|\epsilon|=E_x|F^*(x)-median_xF^*(x)| E∣ϵ∣=Ex∣F∗(x)−medianxF∗(x)∣对第二个实验,误差由一个slash分布生成, ϵ i = s ⋅ ( u / v ) \epsilon_i=s\cdot(u/v) ϵi=s⋅(u/v),这里 u ∼ N ( 0 , 1 ) , v ∼ U [ 0 , 1 ] u\sim N(0,1),\ v\sim U[0,1] u∼N(0,1), v∼U[0,1]。尺度因子 s s s调整得到1/1的信噪比(41)。slash分布有非常厚的尾常被用来测试稳健性。训练集大小取为 N = 7500 N=7500 N=7500,其中5000个作为训练,2500个剩下来作为测试样本估计最优的项目数量 M M M。对100次试验的每次试验都生成额外5000个的验证样本(没有错误的)用来评估该试验的近似误差(37)。
图2中的左边显示了三种方法的两种误差分布在100个目标上近似误差(37)的箱形图。每个箱形图的阴影部分显示了四分位范围内的分布,其中白色条是中位数。外部链条表示距离(上/下)四分位数最近的点(±)1.5个四分位范围单位的点。孤立的条代表了这个范围之外的个别点(离群值)。
这些图允许对总体分布进行比较,但是没有提供关于单个目标函数的相对性能的信息。图2的右侧两个图像尝试提供这样的总结。它们显示了错误率的分布,而不是错误本身。对于每个目标函数和方法,该目标上方法的误差除以该目标上获得的最小误差,所有方法(这里有三种)都进行比较。因此,对于100个试验中的每一个,最佳方法的值为1.0,其他方法的值更大。如果一个特定的方法对于所有100个目标函数来说是最好的(误差最小),那么它的结果分布(箱线图)将是1.0的点质量。请注意,此比率的对数在右下角的图像中绘制。
从图2的左侧图像可以看出,100个目标对于所有三种方法来说都代表了相当广泛的难度范围;近似值误差相差两倍以上。对于正态分布误差, L S LS LS_TreeBoost的性能更好,这是可以预料到的。在73次试验中,它的误差最小,而在另外27次试验中, M M M_TreeBoost最好。平均 L S LS LS_TreeBoost比最佳组差0.2%, M M M_TreeBoost差0.9%, L A D LAD LAD_TreeBoost比最佳组差7.4%。
对于slash分布的误差,情况正好相反。 L S LS LS_TreeBoost的平均近似误差为0.95,因此只能解释5%的目标变化。然而在单个试验中,情况可能更好,也可能更糟。无论是 L A D LAD LAD_TreeBoost还是 M M M_TreeBoost,两者的性能都要好得多,具有可比性。 L A D LAD LAD_TreeBoost最好32次, M M M_TreeBoost最好68次。平均来说, L A D LAD LAD_TreeBoost比最佳的差4.1%, M M M_TreeBoost差1.0%, L S LS LS_TreeBoost比best差364.6%,超过100个目标。
在这三种情况中,建议选择 M M M_TreeBoost方法。在表现很好的(正常)和表现很差的(slash)误差的极端情况中,它的表现都非常接近最好的情况。相比之下, L A D LAD LAD_TreeBoost在正常误差情况下略有损失, L S LS LS_TreeBoost在slash误差情况下损失惨重。
6.3 L S LS LS_TreeBoost与MARS
所有的梯度提升树算法都会产生分段常数逼近。尽管这种分段的数量通常比由单独的一颗树生成的更大,这方面的近似函数 F ^ M ( x ) \hat F_M(x) F^M(x)可能关于提供连续近似的方法存在缺点,特别是当真正的潜在目标函数 F ∗ ( x ) F^*(x) F∗(x)(1)是连续和相当光滑时。所有随机生成的目标函数(39)、(40)都是连续且非常平滑的。在这一节中,我们在100个以上的目标函数中通过对比梯度树提升和MARS的精确度[Friedman (1991)] ,以此来研究分段常数的劣势。和TreeBoost一样,MARS也推出了一个基于张量积的近似。然而,它使用连续函数作为乘积因子,从而产生一个连续的近似。它还使用了一个更复杂的(逐步的)策略来得到张量积。
由于MARS基于最小二乘拟合,我们使用正态分布误差将其与 L S LS LS_TreeBoost进行比较,同样采用1/1信噪比(41)。实验设置与6.2节相同。值得注意的是,在这里,MARS的性能通过使用2500个观察测试集进行模型选择而大大提升,而不是使用默认的广义交叉验证(GCV)准则[Friedman(1991)]。
图3的左上角图像在100个随机生成的目标函数(39)和(40)中比较了MARS的平均绝对近似误差的分布,与图2中的 L S LS LS_TreeBoost的平均绝对近似误差分布。MARS的分布范围要大得多,变化幅度几乎是三分之一。在许多目标上,MARS做得比 L S LS LS_TreeBoost好得多,而在许多目标上,MARS做得比 L S LS LS_TreeBoost差得多。这进一步说明了目标函数的性质对不同方法的相对性能有很大的影响。图3的右上角图像显示了误差的分布,相对于每个目标的最佳情况。基于平均绝对误差,两种方法具有相似的性能。在许多目标函数中,每个目标的表现都明显优于另一个目标。
图3的底部两个图像显示了基于均方根误差的结果。在评估性能缺陷时,这给较大的错误增加了相应的权重。对于 L S LS LS_TreeBoost来说,这两个误差的度量值对于100个目标来说都接近于相同的值。然而对于MARS,均方根误差通常比平均绝对误差高30%。这表明MARS的预测要么非常接近目标,要么离目标很远。 L S LS LS_TreeBoost进行的误差分布更加均匀。它有更少的非常大的错误或非常小的错误。后者可能是由于分段逼近的结果,使得用有限大小的逼近法很难任意地接近非常平滑变化的目标函数。如图3所示,相对性能对用于度量它的标准非常敏感。
这些结果表明,提升树近似的分段常数方面不是一个严重的缺点。在正态误差和正态输入变量分布相当原始的环境中,它与MARS是竞争的。分段常数法的优点是鲁棒性;具体来说,它对宽尾和输入变量 x x x中存在异常值这些负面影响免疫。产生连续近似的方法,比如MARS,对这类问题非常敏感。同样,如6.2节所示, M M M_TreeBoost(算法4)对于正态误差几乎与 L S LS LS_TreeBoost一样准确,而且对输出 y y y的离群值具有很强的抵抗力。因此,在数据挖掘应用中,如果数据的清洁度不能保证且 x x x和/或 y y y存在离群值,那么 M M M_TreeBoost的相对较高的精度、一致的性能和鲁棒性可能是一个很大的优势。
6.4 L k L_k Lk_TreeBoost、K分类LogitBoost和Adaboost.MH
在本节中,将 L k L_k Lk_TreeBoost的性能与K分类LogitBoost (FHT00)和Adaboost.MH [Schapire和Singer(1998)]在100个随机生成的目标(第6.1节)上进行了比较。这里 K = 5 K = 5 K=5个类是通过对输入 x x x值分布上每个目标的0.2、0.4、0.6和0.8分位数进行阈值化而生成的。这里有 N = 7500 N=7500 N=7500个训练样本,其中每次试验(每个类1500)分5000个作为训练样本和2500样本用来模型选择(迭代次数,M)。使用独立生成的5000个验证样本来估计每个目标函数的错误率。
对于所有的目标函数,贝叶斯错误率为零,但诱导的策边界可以变得相当复杂,取决于每个目标函数 F ∗ ( x ) F^*(x) F∗(x)的性质。每个方法使用11个终端节点的回归树。
图4显示了这三种方法在100个目标函数上的错误率分布(左)和除以最小值的比率(右)。三种方法的误差率在这些目标上有实质性的不同。 L k L_k Lk_TreeBoost被认为是普遍表现最好的。在78个试验中,它的误差最小,平均误差比每次实验中最好的高0.6%。LogitBoost在21个目标中表现最好,有一个并列。它的错误率比最好的平均高出3.5%。Adaboost.MH从来都不是最好的,平均比最好的还要差15%。
图5显示了与LogitBoost和Adaboost.MH对应的比较。Adaboost.MH程序修改为增量收缩(36),与 L k L_k Lk_TreeBoost一样收缩参数设置为相同的(默认)值 ν = 0.1 \nu= 0.1 ν=0.1。在这里,我们看到的是有些区别的画面。LogitBoost和AdaBoost.MH都大大得益于收缩。这三种程序的性能现在几乎相同,LogitBoost可能有一点优势,平均来说,它的错误率比最好的差0.5%; L K L_K LK_TreeBoost和AdaBoost.MH的对应值分别为2.3%和3.9%。这些结果表明,这些方法的相对性能更依赖于它们的aggressiveness,例如以学习率作为参数,而不是它们的结构差异。当分母接近0时,LogitBoost有一个额外的内部收缩,以稳定其伪响应(33)(FHT00),这也许可以解释为什么它在这种比较中略胜一筹。事实上,和图5所示的LogitBoost一样,当收缩增加应用于 L K L_K LK_TreeBoost( ν = 0.05 \nu= 0.05 ν=0.05)它的性能改善。当对这三种方法中的每一种都仔细调整收缩参数时,它们之间的性能差别可能很小。
7. 提升树
梯度提升算法(算法1)有两个主要的超参数,迭代次数 M M M和学习速率 ν \nu ν(36)。这些都在第5节中有讨论。除了这些以外,还有一些和估计基学习器 h ( x ; a ) h(x;a) h(x;a)有关的超参数。本文的主要研究的重点是使用 J J J个终端节点数量固定的最佳优先回归树,因此, J J J是基学习器的主要的超参数,它的最佳值十分依赖于目标函数的性质,换句话说,就是在这些变量中占主导的交互作用的最高阶数。
考虑函数的方差分析展式 F ( x ) = ∑ j f j ( x j ) + ∑ j , k f j k ( x j , x k ) + ∑ j , k , l f j k l ( x j , x k , x l ) + ⋯ (42) F(x)=\sum_jf_j(x_j)+\sum_{j,k}f_{jk}(x_j,x_k)+\sum_{j,k,l}f_{jkl}(x_j,x_k,x_l)+\cdots\tag{42} F(x)=j∑fj(xj)+j,k∑fjk(xj,xk)+j,k,l∑fjkl(xj,xk,xl)+⋯(42)第一个和称为 F ( x ) F(x) F(x)的“主效应”,它是由一些函数的和组成,这些函数每个都只依赖于一个输入变量。特殊的函数 { f j ( x j ) } 1 N \{f_j(x_j)\}_1^N {fj(xj)}1N是那些在这种加性限制下提供了最近的 F ( x ) F(x) F(x)的逼近。这有时候被称为加性模型,因为每个 x j , f j ( x j ) x_j,f_j(x_j) xj,fj(xj)的贡献都会其他的加在一起。这是一个与(2)不同的、限制性更强的“加性”定义。第二个加和是由一对输入变量的函数组成,它们被称为双变量的“交互作用效应”。这样选择是在主效应下,不超过两个变量的限制下提供最接近 F ( x ) F(x) F(x)的近似。第三个加和表示三变量交互作用效果,等等。
最高的交互阶数可能由输入变量的个数 n n n限制。然而,特别是对很大的 n n n,许多在实际中遇到的目标函数 F ∗ ( x ) F^*(x) F∗(x)可以通过很低阶数的方差分析模型来很好地近似。(42)中只有前几项需要捕捉到 F ∗ ( x ) F^*(x) F∗(x)的主要变化。事实上,相当大的成功往往是通过单独的加性项实现的[Hastie and Tibshirani(1990)]。纯粹的加性近似也可由“朴素”贝叶斯方法产生[Warner, Toronto, Veasey and Stephenson(1961)],这种方法通常在分类上非常成功。这些考虑激发了用于模拟研究的随机生成目标函数(第6.1节)中对低阶交互的偏向。
函数估计的目的是产生一个非常匹配目标函数 F ∗ ( x ) F^*(x) F∗(x)的近似 F ^ ( x ) \hat F(x) F^(x)。这通常要求 F ^ ( x ) \hat F(x) F^(x)的主要交互阶数要和 F ∗ ( x ) F^*(x) F∗(x)的相似。在boosting回归树中,交互阶数可以通过限制每次迭代生成的单颗树的大小来控制。一个有着 J J J个终端结点的树产生交互阶数最多为 min ( J − 1 , n ) \min(J-1,n) min(J−1,n)的函数,boosting程序是加性的,因此整个近似的交互阶数应该不超过单颗树的最大交互阶数。因此,对于任何一种提升树程序,树的最佳大小 J J J都是被目标函数 F ∗ ( x ) F^*(x) F∗(x)的有效交互阶数控制。这通常是未知的,因此 J J J就变成了需要用模型选择准则例如交叉验证或者留出法来调整的超参数。然而,正如上面所讨论的,我们不太需要或者想要很大的树。
图6显示了对于100个(6.1节)用来模拟研究的随机生成的函数,树的大小对于逼近精确性的影响。实验的初始设置和6.2节的一样,对于 J ∈ { 2 , 3 , 6 , 11 , 21 } J\in\{2,3,6,11,21\} J∈{2,3,6,11,21},绝对误差(37)的分布(左图)以及对于每个目标误差相对于最小误差(右图)给出。第一个 J = 2 J=2 J=2只生成了加性的主效应部分; J = 3 J=3 J=3生成了加性的和双变量的交互项,等等。一个 J J J个终端结点的树可以产生的交互阶数最大为 min ( J − 1 , n ) \min(J-1,n) min(J−1,n),一般来说都小于这个值,特别是当 J − 1 ≲ n J-1\lesssim n J−1≲n。
图6可以看出,最小的树 J ∈ { 2 , 3 } J\in\{2,3\} J∈{2,3}平均产生了较低的精度,但是他们的分布比其他的宽很多,这意味着它们会产生更精确,甚至更不准确的近似。较小的树被限制在低阶交互中,因此能够更好地利用低交互阶数低的目标函数。然而,当他们尝试近似高阶交互作用的目标时表现就十分不好。 J ∈ { 6 , 11 , 21 } J\in\{6,11,21\} J∈{6,11,21}的较大的树会更一致,他们牺牲了在低阶交互目标上的精度,但是在高阶函数中做得更好。在这些较大的树中也有一些差异,对 J = 21 J=21 J=21来说可能有轻微的退化。 J = 2 J=2 J=2的树产生了8次最精确的近似; J ∈ { 3 , 6 , 11 , 21 } J\in\{3,6,11,21\} J∈{3,6,11,21}分别对应产生了2、30、31、29次。平均上, J = 2 J=2 J=2的树比每个目标函数的最低误差要高23.2%,而其他的分别对应16.4%、2.4%、2.2%和3.7%。当对每个目标分别估计最佳树大小 J J J时,应获得更高的精度。在实践中,这可以通过使用独立的测试集来评估不同树大小的使用情况来实现,如第9节所示。
8. 解释
在许多应用中,能够解释推导出来的近似 F ^ ( x ) \hat F(x) F^(x)是非常有用的,这包括增加对于那些最影响变化的特定的输入变量的理解,以及在这些有影响的输入上 F ^ ( x ) \hat F(x) F^(x)的独立性。某种程度上 F ^ ( x ) \hat F(x) F^(x)至少定量地翻译了目标函数 F ∗ ( x ) F^*(x) F∗(x)(1)的性质,这样的工具能够提供关于输入变量 x x x和输出变量 y y y的潜在关系的信息。在本节中,将会用一些工具来解释提升树逼近。尽管他们被用来解释单颗决策树,但是他们似乎在提升树(特别是很小的树)中更有效。这些解释工具在第9节的实际数据的例子中阐述。
8.1 输入变量的相对重要性
在一个逼近 F ^ ( x ) \hat F(x) F^(x)中最重要的描述就是单个输入变量 x j x_j xj在联合输入变量的分布上,对 F ^ ( x ) \hat F(x) F^(x)变化的相对影响 I j I_j Ij。我们这样定义 I j = ( E x [ ∂ F ^ ( x ) ∂ x j ] 2 ⋅ v a r x [ x j ] ) 1 / 2 (43) I_j=\bigg(E_x\bigg[{\partial \hat F(x)\over \partial x_j}\bigg]^2\cdot var_x[x_j]\bigg)^{1/2}\tag{43} Ij=(Ex[∂xj∂F^(x)]2⋅varx[xj])1/2(43)对于由决策树产生的分段常数近似,(43)并不严格存在,它必须由反映其性质的替代度量来近似。Breiman、Friedman、Olshen和Stone(1983)提出 I ^ j 2 ( T ) = ∑ t = 1 J − 1 l ^ t 2 1 ( v t = j ) (44) \hat I_j^2(T)=\sum_{t=1}^{J-1}\hat l_t^21(v_t=j)\tag{44} I^j2(T)=t=1∑J−1l^t21(vt=j)(44)其中这个和是对 J J J个终端结点的决策树 T T T上的非终端节点 t t t求和, v t v_t vt是和结点 t t t有关的分裂变量, l ^ t 2 \hat l_t^2 l^t2是作为分裂结果对应的平方误差的经验改进(35)。(44)式右端是和平方误差影响有关,因此它的单位对应了(43)的部分。Breiman、Friedman、Olshen和Stone(1983)直接使用(44)作为影响的度量,而不是平方影响。对于一个由boosting得到的决策树集 { T m } 1 M \{T_m\}_1^M {Tm}1M,(44)可以用序列中所有树的平均值来推广,即 I ^ j 2 = 1 M ∑ m = 1 M I ^ j 2 ( T m ) (45) \hat I_j^2={1\over M}\sum_{m=1}^{M}\hat I_j^2(T_m)\tag{45} I^j2=M1m=1∑MI^j2(Tm)(45)
(44)、(45)的动机纯粹是基于启发式的论证。作为部分的证明,我们证明了在最简单的环境中使用它会产生预期的结果。考虑一个线性目标函数 F ∗ ( x ) = a 0 + ∑ j = 1 n a j x j (46) F^*(x)=a_0+\sum_{j=1}^na_jx_j\tag{46} F∗(x)=a0+j=1∑najxj(46)其中输入变量的协方差矩阵是单位矩阵的倍数 E x [ ( x − x ˉ ) ( x − x ˉ ) T ] = c I n E_x\big[(x-\bar x)(x-\bar x)^T\big]=cI_n Ex[(x−xˉ)(x−xˉ)T]=cIn在这种情况下,影响度量(43)得到 I j = ∣ a j ∣ (47) I_j=|a_j|\tag{47} Ij=∣aj∣(47)表2显示了与第六节类似的一个小的模拟实验的结果,但由于系数的取值 a j = ( − 1 ) j j (48) a_j=(-1)^jj\tag{48} aj=(−1)jj(48)以及1/1的信噪比, F ∗ ( x ) F^*(x) F∗(x)(46)是线性的。在10个随机样本上给出了(44)和(45)的均值和标准差, F ∗ ( x ) F^*(x) F∗(x)由(46)(48)给出。估计出的最优影响的变量的影响 x j ∗ x_{j^*} xj∗人为地赋予值 I j ∗ = 100 I_{j^*}=100 Ij∗=100,其他变量的估计值响应地进行放缩。输入变量的估计的重要性排序在10次试验中每一次都是正确的。正如表2所示,估计的相对重要性的值和(47)与(48)式给出的一致。
在Breiman, Friedman, Olshen和Stone1983年的研究中,影响度量(44)被一种涉及替代分割的策略所扩张,该策略旨在发现与影响变量高度相关的其他隐藏的重要变量。这种策略对于单个决策树最有帮助,因为由于(44)中树的大小 J J J的限制,单个决策树参与分裂的变量的机会有限。然而,在boosting中,分裂的机会大大增加了(45),而相应地,用替代方法去揭露就不那么重要了。
在 K K K分类的逻辑回归和分类(4.6节)中,存在 K K K个逻辑回归函数 { F k M ( x ) } k = 1 K \{F_{kM}(x)\}_{k=1}^K {FkM(x)}k=1K,每一个都是由一列 M M M决策树描述。这样(45)就推广为 I ^ j k 2 = 1 M ∑ m = 1 M I ^ j 2 ( T k m ) (49) \hat I_{jk}^2={1\over M}\sum_{m=1}^{M}\hat I_j^2(T_{km})\tag{49} I^jk2=M1m=1∑MI^j2(Tkm)(49)其中 T k m T_{km} Tkm是在第 m m m次迭代中的第 k k k类引导出来的树。 I ^ j k \hat I_{jk} I^jk的值可以解释为预测变量 x j x_j xj将 k k k类与其他类分隔开来的相关性。 x j x_j xj的所以的相关性可以由所有类求平均得到 I ^ j = 1 K ∑ k = 1 K I ^ j k \hat I_j={1\over K}\sum_{k=1}^K\hat I_{jk} I^j=K1k=1∑KI^jk然而,单个 I ^ j k \hat I_{jk} I^jk也非常有用。通常情况下,不同变量的子集与不同类的子集高度相关。这种更详细的知识可能导致仅通过检查整体相关性无法获得。
8.2 部分依赖图
可视化是最重要的解释工具之一。作为其参数的函数 F ^ ( x ) \hat F(x) F^(x)的图像描述提供了在输入变量的联合分布上的依赖性的综合总结。但不幸的是,这样的可视化只限制在低维参数中,单个实值变量 x x x的函数 F ^ ( x ) \hat F(x) F^(x)可以用 F ^ ( x ) \hat F(x) F^(x)的值与 x x x的每个对应值的关系图来表示。单个分类变量的函数可以由条形图来表示,每条表示其中一个值,且条形的高表示函数的值。两个实值变量的函数可以用轮廓或透视网格图来画出。一个分类变量和另外一个变量(实值或者分类)最好由一系列图来总结,每一个显示了限制在第一个变量的各自的值上,在第二个变量上 F ^ ( x ) \hat F(x) F^(x)的依赖性[Becker and Cleveland (1996)]。
观察高维参数的函数更加困难,因此观察逼近 F ^ ( x ) \hat F(x) F^(x)在选择的输入变量的小的子集上的部分依赖性是很有用的。尽管这样的图像的集合不能提供一个全面的逼近的描述,但是它常常可以提供有用的线索,特别是当 F ^ ( x ) \hat F(x) F^(x)由一些低阶交互作用主导时(第7节)。
令 z l z_l zl表示输入变量 x x x的被选中的目标子集,大小为 l l l z l = { z 1 , … , z l } ⊂ { x 1 , … , x n } z_l=\{z_1,\dots,z_l\}\subset \{x_1,\dots,x_n\} zl={z1,…,zl}⊂{x1,…,xn}且 z \ l z_{\backslash l} z\l表示补集 z \ l ∪ z l = x z_{\backslash l}\cup z_l=x z\l∪zl=x逼近 F ^ ( x ) \hat F(x) F^(x)主要取决于这两个子集的变量 F ^ ( x ) = F ^ ( z l , z \ l ) \hat F(x)=\hat F(z_l,z_{\backslash l}) F^(x)=F^(zl,z\l)如果我们对 z \ l z_{\backslash l} z\l里的变量规定一些特殊的值,那么 F ^ ( x ) \hat F(x) F^(x)可以看做是仅作为所选子集 z l z_l zl中变量的函数 F ^ z \ l ( z l ) = F ^ ( z l ∣ z \ l ) (50) \hat F_{z_{\backslash l}}(z_l)=\hat F(z_l|z_{\backslash l})\tag{50} F^z\l(zl)=F^(zl∣z\l)(50)
一般地, F ^ z \ l ( z l ) \hat F_{z_{\backslash l}}(z_l) F^z\l(zl)的函数形式将取决于为 z \ l z_{\backslash l} z\l所选的特定的值。然而,如果这种依赖性太强,那么平均函数 F ˉ l ( z l ) = E z \ l [ F ^ ( x ) ] = ∫ F ^ ( z l , z \ l ) p \ l ( z l ) d z \ l (51) \bar F_l(z_l)=E_{z_{\backslash l}}[\hat F(x)]=\int\hat F(z_l,{z_{\backslash l}}){p_{\backslash l}}(z_l)d{z_{\backslash l}}\tag{51} Fˉl(zl)=Ez\l[F^(x)]=∫F^(zl,z\l)p\l(zl)dz\l(51)就可以表示在选择的变量子集 z l z_l zl上 F ^ ( x ) \hat F(x) F^(x)的部分依赖的有用的总结。这里 p \ l ( z l ) {p_{\backslash l}}(z_l) p\l(zl)为 z l z_l zl的边际概率密度 p \ l ( z l ) = ∫ p ( x ) d z l (52) {p_{\backslash l}}(z_l)=\int p(x)dz_l\tag{52} p\l(zl)=∫p(x)dzl(52)其中 p ( x ) p(x) p(x)是所有输入 x x x的联合密度。根据训练数据估计补集的边际密度(52),则(51)为 F ˉ l ( z l ) = E z \ l [ F ^ ( x ) ] = 1 N ∑ i = 1 N F ^ ( z l , z i , \ l ) (53) \bar F_l(z_l)=E_{z_{\backslash l}}[\hat F(x)]={1\over N}\sum_{i=1}^N\hat F(z_l,{z_{i,\backslash l}})\tag{53} Fˉl(zl)=Ez\l[F^(x)]=N1i=1∑NF^(zl,zi,\l)(53)
在特殊情况下, F ^ ( x ) \hat F(x) F^(x)在 z l z_l zl上的依赖是加性的 F ^ ( x ) = F ^ l ( z l ) + F ^ \ l ( z \ l ) (54) \hat F(x)=\hat F_l(z_l)+\hat F_{\backslash l}(z_{\backslash l})\tag{54} F^(x)=F^l(zl)+F^\l(z\l)(54)或者乘性的 F ^ ( x ) = F ^ l ( z l ) ⋅ F ^ \ l ( z \ l ) (55) \hat F(x)=\hat F_l(z_l)\cdot\hat F_{\backslash l}(z_{\backslash l})\tag{55} F^(x)=F^l(zl)⋅F^\l(z\l)(55)(50)式 F ^ z \ l ( z l ) \hat F_{z_{\backslash l}}(z_l) F^z\l(zl)的形式不依赖于补变量 z \ l z_{\backslash l} z\l的联合值。则(51)式 F ˉ l ( z l ) \bar F_l(z_l) Fˉl(zl)给出了在选择的输入变量子集 z l z_l zl上 F ^ ( x ) \hat F(x) F^(x)变化的性质的完整描述。
在子集 z l z_l zl上 F ^ ( x ) \hat F(x) F^(x)的依赖性的另外一种总结方式是直接将 F ^ ( x ) \hat F(x) F^(x)作为训练数据上 z l z_l zl的函数模型 F ~ l ( z l ) = E x [ F ^ ( x ) ∣ z l ] = ∫ F ^ ( x ) p ( z \ l ∣ z l ) d z \ l (56) \tilde F_l(z_l)=E_x[\hat F(x)|z_l]=\int\hat F(x){p}(z_{\backslash l}|z_l)d{z_{\backslash l}}\tag{56} F~l(zl)=Ex[F^(x)∣zl]=∫F^(x)p(z\l∣zl)dz\l(56)然而,在条件密度(56)而不是边际密度(51)上取平均,导致 F ~ l ( z l ) \tilde F_l(z_l) F~l(zl)不仅反映在选择的输入变量子集 z l z_l zl上 F ^ ( x ) \hat F(x) F^(x)的依赖性,此外,还有仅仅通过它们和补变量 z \ l z_{\backslash l} z\l之间的联系诱导的明显的依赖性。例如,如果 z l z_l zl恰好是加性的(54)或乘性的(55),(56) F ~ l ( z l ) \tilde F_l(z_l) F~l(zl)不会估计对应的项或因子 F ^ l ( z l ) \hat F_l(z_l) F^l(zl),除非联合密度 p ( x ) p(x) p(x)碰巧是乘积 p ( x ) = p l ( z l ) ⋅ p \ l ( z \ l ) (57) p(x)=p_l(z_l)\cdot p_{\backslash l}(z_{\backslash l})\tag{57} p(x)=pl(zl)⋅p\l(z\l)(57)
部分依赖函数(51)可用于解释任何“黑盒”预测方法生成的模型,如神经网络、支持向量机、最近邻、径向基函数等。当有大量的预测变量时,有一个依赖性度量是非常有用的(第8.1节),以减少潜在的大量变量和要考虑的变量组合。同样,需要对数据进行一次传递(53),以便为参数的每一组联合值 z l z_l zl对每个 F ˉ l ( z l ) \bar F_l(z_l) Fˉl(zl)求值。对于大型数据集来说,这可能非常耗时,不过子集抽样可能会有所帮助。
然而对于基于单变量分裂的回归树, F ^ ( x ) \hat F(x) F^(x)在特定目标变量子集 z l z_l zl(51)上的部分依赖性直接计算给定的树,而不参考数据本身(53)。对变量 z l z_l zl的一组特定值,将执行树的加权遍历。在树的根节点,赋值为1。对于所访问的每个非终端节点,如果其分割变量位于目标子集 z l z_l zl中,则访问适当的左子节点或右子节点,并且不修改权重。如果节点的分类变量在补子集 z \ l z_{\backslash l} z\l中,则访问两个子节点,并将当前权重乘以在该节点分别向左或向右移动的训练样本的比例。
遍历过程中访问的每个终端节点都被赋予当前的权值。当树遍历完成, F ˉ l ( z l ) \bar F_l(z_l) Fˉl(zl)的值是在树遍历访问的终端结点上 F ^ ( x ) \hat F(x) F^(x)相应值的加权平均。对于通过boosting得到的 M M M个回归树的集合,对每个树的结果进行简单平均。
通过图的展示来达到解释的目的,如果输入变量子集是低势的( l ⩽ 2 l\leqslant 2 l⩽2),大部分情况下是有用的。这样的包含信息最多的子集可能由那些被认为是导致 F ^ ( x ) \hat F(x) F^(x)变化的最具影响力的(44),(45)的输入变量组成,第8.3节和第9节中提供了依赖说明。
在子集 z l z_l zl上 F ^ ( x ) \hat F(x) F^(x)的依赖性越是加性(54)或者乘性(55)的,那么部分依赖函数 F ˉ l ( z l ) \bar F_l(z_l) Fˉl(zl)(51)在推导出来的逼近 F ^ ( x ) \hat F(x) F^(x)上就越是能捕捉到 z l z_l zl中变量的影响的性质。因此,将有影响的输入组成起来的子集 z l z_l zl提供了最真实的部分依赖图像。作为诊断, F ˉ l ( z l ) \bar F_l(z_l) Fˉl(zl)和 F l ( z \ l ) F_l(z_{\backslash l}) Fl(z\l)都可以分别计算候选子集。在训练集上 F ^ ( x ) \hat F(x) F^(x)和 { F ˉ l ( z l ) , F l ( z \ l ) } \{\bar F_l(z_l),F_l(z_{\backslash l})\} {Fˉl(zl),Fl(z\l)}之间,或者和 F ˉ l ( z l ) ⋅ F l ( z \ l ) \bar F_l(z_l)\cdot F_l(z_{\backslash l}) Fˉl(zl)⋅Fl(z\l)之间的多种依赖系数的值都可以用来衡量 F ^ ( x ) \hat F(x) F^(x)关于一个选择的子集 z l z_l zl的加性或者因子可分解性的程度。作为加性诊断,(50)式 F ^ z \ l ( z l ) \hat F_{z_{\backslash l}}(z_l) F^z\l(zl)可以用来计算从训练集中随机选择的比较少的 z \ l z_{\backslash l} z\l值。 z l z_l zl的结果函数可以比作 F ˉ l ( z l ) \bar F_l(z_l) Fˉl(zl)来判断在 z l z_l zl上关于 z \ l z_{\backslash l} z\l的改变值的部分依赖性。
在 K K K分类逻辑回归和分类(4.6节)中,存在 K K K个逻辑回归函数 { F k ( x ) } k = 1 K \{F_k(x)\}_{k=1}^K {Fk(x)}k=1K,由(29)式每个都是和 p k ( x ) = P r ( y = k ∣ x ) p_k(x)=Pr(y=k|x) pk(x)=Pr(y=k∣x)对数依赖。 F k ( x ) F_k(x) Fk(x)的较大值意味着更高的在 x x x处观察到 k k k类的概率。每个 F k ( x ) F_k(x) Fk(x)在变量子集 z l z_l zl上对最依赖的那个类的部分依赖图(49)提供了输入变量如何影响各自的类概率的信息。
8.3 随机生成的函数
在本节中,前面两节中描述的解释工具适用于第6节蒙特卡罗研究中使用的第一个(100个)随机生成的函数(第6.1节)。
图7显示了10个输入预测变量的估计相对重要性(44)、(45)。有些变量被认为比其他的更有影响力,但似乎没有一个小的子集占主导地位。这与用来生成这些函数的机制是一致的。
图8显示了在6个最具影响力的变量上的单变量( l = 1 l = 1 l=1)部分依赖图(53)。每个图底部的散列标记表示相应预测变量分布的十分位数。近似的分段常数性质是明显的。与大多数近似方法不同,TreeBoost模型没有明确的平滑度约束。任意尖锐的不连续点可以容纳。这些图中显示的一般平滑趋势表明,平滑近似最能描述这个目标。这与生成这些函数的方式是一致的。
图9显示了两个变量( l = 2 l = 2 l=2)对一些更有影响力的变量的部分依赖图。这些变量对之间存在不同程度的交互作用。这与这些目标函数实际生成的方式一致(39),(40)。
考虑到这些生成的目标作为其参数的函数的一般复杂性,不太可能通过一系列这样的部分依赖图来揭示它们完整的详细的函数形式。我们的目标是获得对函数关系的一些重要方面的可理解的描述。在本例中,目标函数是由已知的设定生成的,因此至少我们可以定性地验证这里的情况。
9. 真实数据
在本节中,提升树回归算法将在两个中等大小的数据集上进行说明。6.4节的结果表明,分类算法 L k L_k Lk_TreeBoost的性质和LogitBoost的很像,都是把FHT00广泛应用于数据。第一个(科学)数据集包含对岩石样本的化学浓度测量,第二个(人口统计学)数据集是样本调查问卷数据。两个数据集都被分割为一个占三分之二的学习样本,剩下的作为测试样本来选择模型大小(迭代次数 M M M)。收缩参数(36)设置为 ν = 0.1 \nu=0.1 ν=0.1。
9.1石榴石数据
这个数据集包含
N
=
13317
N=13317
N=13317个从世界各地收集的石榴石样本 [Griffin, Fisher, Friedman, Ryan and O’ Reilly (1997)]。石榴石是一种复杂的钙-镁-铁-铬硅酸盐复合物,通常作为组成地幔的岩石中的一个小相位出现。与每一个石榴石相关的变量是各种化学物质的浓度和岩石收集的构造板块环境:
(
T
i
O
2
,
C
r
2
O
3
,
F
e
O
,
M
n
O
,
M
g
O
,
C
a
O
,
Z
n
,
G
a
,
S
r
,
Y
,
Z
r
,
t
e
c
)
(TiO_2, Cr_2O_3, FeO, MnO, MgO, CaO, Zn, Ga, Sr, Y, Zr, tec)
(TiO2,Cr2O3,FeO,MnO,MgO,CaO,Zn,Ga,Sr,Y,Zr,tec)
代表浓度的前11个变量是实值。最后一个变量(tec)有三个分类值:“古代稳定盾构”、“元古代盾构区”和“年轻造山带”。这些数据中没有缺失的值,但是许多变量的分布倾向于更大的值,有许多异常值。
本实验的目的是估计钛( T i O 2 TiO_2 TiO2)的浓度作为其他化学物质的联合浓度和构造板块指数的函数。
表3是基于测试样本,对几个不同大小值(终端结点数) J J J个构成树的 L S LS LS_TreeBoost、 L A D LAD LAD_TreeBoost以及 M M M_TreeBoost的输出变量 y y y相对于最优常数预测的平均绝对误差 A ( y , F ^ ( x ) ) = E y , x ∣ y − F ^ ( x ) ∣ E y ∣ y − m e d i a n ( y ) ∣ (58) A(y,\hat F(x))={E_{y,x}\big|y-\hat F(x)\big|\over E_y\big|y-median(y)\big|}\tag{58} A(y,F^(x))=Ey∣∣y−median(y)∣∣Ey,x∣∣y−F^(x)∣∣(58)注意,这个预测误差(58)包括关于潜在(未知的)目标函数 F ∗ ( x ) F^*(x) F∗(x)加性的不可约的错误。这不可约错误在表3所有条目中增加相同的数量。因此,这些项的差异反映了对目标函数本身的近似误差(37)按比例更大的改进。
对于这三种方法,加法(J = 2)近似明显不如使用更大的树,这表明输入变量之间存在交互作用(第7节)。六个终端节点树被认为是足够的,而仅使用三个终端节点树被认为提供了最好的10%的准确率。 L A D LAD LAD_TreeBoost和 M M M_TreeBoost的误差比 L S LS LS_TreeBoost的小,且彼此相似,可能 M M M_TreeBoost有微弱的优势。这些结果与模拟研究结果一致,如图2和图6所示。
图10显示了用六个终端结点树基于 M M M_TreeBoost逼近来预测 T i O 2 TiO_2 TiO2浓度的11个输入变量的相对重要性(44)和(45)。对于表3中的其他模型,结果非常相似,并且有类似的错误。Ga和Zr被认为是最有影响力的,而MnO则不那么重要。图11的上面三个小图展示了在这三个最优影响的变量上,近似 F ^ ( x ) \hat F(x) F^(x)的部分相关性(51) 。下面三个底部的图像显示了这些变量三个配对上的 F ^ ( x ) \hat F(x) F^(x)的部分相关性(51) 。Ga与Zr之间存在明显的强交互作用。当其他变量取值很小时, F ^ ( x ) \hat F(x) F^(x)几乎和这些变量没有相关性。随着他们其中之一的值增加, F ^ ( x ) \hat F(x) F^(x)在其他变量上的相关性也相应地增大。MnO与Zr之间的交互效应较小。
9.2 人口统计数据
该数据集由旧金山湾区的购物中心的顾客填写的 N = 9409 N = 9409 N=9409份问卷组成[Impact Resources, Inc, Columbus, Ohio(1987)]。在这里,我们使用前14个问题的答案来举例说明,这些问题与人口统计学有关。这些问题列在表4中。这些数据由若干个连续和分类的变量组成,每个变量有一个小数目的不同的值。有许多缺失的值。
我们通过将收入建模为其他13个变量的函数来说明提升树。表5显示了三种回归提升树算法在预测收益时相对于最佳常数预测器(58)的平均绝对误差。
这三种方法的性能差别不大。由于这些数据的高度离散性,实值输入或输出 y y y之间不存在异常值或长尾分布。随着组成树大小 J J J的增加,误差几乎没有减少,表明输入变量之间缺乏交互作用;在单独的输入变量( J = 2 J = 2 J=2)中添加一个加法逼近似乎是足够的。
图12显示了基于( J = 2 J = 2 J=2) L S LS LS_TreeBoost逼近的输入变量在预测收入中的相对重要性。没有一个小的子集是占主导地位的。图13显示了6个最具影响力的变量的部分依赖图。分类变量用条形图表示,所有的条形图都以数据均值为零为中心。因为逼近只由主要影响组成[(42)中的第一个和],这些图像完全描述了每个输入的相应的贡献 f j ( x j ) f_j(x_j) fj(xj)。
在图13中似乎没有任何令人惊讶的结果。相关性在很大程度上证实了先前的猜测,并表明这种近似在直觉上是合理的。
10. 数据挖掘
作为预测数据挖掘的“现成”工具,提升树程序具有一些吸引人的特性。他们继承了决策树良好的性质,而减轻了一些不好的。所有提升树过程在所有单个输入变量的(严格)单调变换下都是不变的。例如,用 x j , log x j , e x j x_j,\log x_j,e^{x_j} xj,logxj,exj或 x j a x^a_j xja作为第 j j j个输入变量都会得到相同的结果。因此,考虑输入变量转换的必要性被消除了。作为这种不变性的结果,对长尾分布和异常值的敏感性也被消除了。此外, L A D LAD LAD_TreeBoost对于输出变量 y y y中的异常值也是完全健壮的。 M M M_TreeBoost对输出异常值的健壮性也有一定的度量。
决策树诱导的另一个优点是内部特征选择。对于添加不相关的输入变量,树往往非常健壮。此外,基于树的模型以统一而优雅的方式处理缺失的值[Breiman, Friedman, Olshen and Stone(1983)],没有必要考虑额外的处理方案。提升树显然也继承了这些特性。
单一树模型的主要缺点是不准确。这是由于他们的分段常数近似粗糙性(特别是对于较小的树)和不稳定性(特别是对于较大的树),以及他们所涉及的主要是高阶相互作用的事实的结果。所有这些都通过boosting得到了缓解。提升树程序产生分段常数近似,但是间隔尺寸要细得多。提升树通过使用小树和在许多小树上取平均值的效果来提高稳定性。通过限制各组成树的大小,有效地控制了提升树近似的交互水平。
单一树模型的最大优点之一是可解释性,而提升树则被认为缺乏这一特性。小树容易解释,但由于不稳定,这种解释应该谨慎对待。较大的树的可解释性是有问题的[Ripley(1996)]。提升树近似可以与输入变量相对重要性度量一起,使用部分依赖图进行解释,如第8.3节和第9节所示。虽然没有提供一个完整的描述,但至少提供了对输入-输出关系本质的一些了解。虽然这些工具可以用于任何近似方法,但基于树的模型的特殊特性使它们能够快速计算。部分依赖图也可以用于单个回归树,但是如上所述,由于更大的不稳定性,需要更加谨慎。
在对输入变量进行排序之后,回归提升树程序( L S LS LS_, L A D LAD LAD_和 M M M_TreeBoost)的计算量与样本数量 N N N、输入变量的数量 n n n和迭代次数 M M M成线性关系。分类算法 L K L_K LK_TreeBoost的规模约为构成树的大小 J J J的对数。此外,分类算法 L K L_K LK_TreeBoost的规模与类的数量 K K K成线性关系;但是,如果使用影响修整(章节4.5.1),则与迭代次 M M M呈高度的亚线性关系。作为一个参考,对第9.1节的石榴石数据( N = 13317 , n = 11 , J = 6 , M = 500 N = 13317,n = 11,J = 6,M = 500 N=13317,n=11,J=6,M=500)应用 M M M_TreeBoost,在933Mh奔腾III的电脑上需要20秒。
正如第5节所述,基于小的收缩值参数 ν \nu ν(36),许多boosting迭代( M ≃ 500 M\simeq 500 M≃500)可以获得最佳提升树近似。这在一定程度上由于在每次迭代中引入的树非常小而得到缓解。然而,正如图1所示,改进最初非常快,然后逐渐变慢。因此,几乎可以在相当早期( M ≃ 100 M\simeq 100 M≃100)以及相应的更少的计算来实现最优近似。这些接近最优的近似可用于最初的探索,并表明最终的近似是否具有足够的准确性以保证继续下去。如果在最初的几次迭代中(比如100次),拟合不匹配的改善非常小,那么之后就不太可能有显著的改善。如果判断延续是合理的,则可以在先前停止的地方重新启动过程,这样就不会损失计算的投资。此外,我们还可以使用更大的收缩参数来加速初始化。如图1所示,使用 ν = 0.25 \nu= 0.25 ν=0.25在迭代20次后就提供了精度在10%的最优( ν = 0.1 \nu= 0.1 ν=0.1)解。然而,在这种情况下,如果随后使用较小的收缩参数值,则必须重新启动boosting。
提升树程序能够提供潜在可预测性的快速指示,加上其强大的鲁棒性,使其成为一个有用的预处理工具,可以应用于不完美的数据。如果有足够的预指示能力,可以进行进一步的数据清理,使其适合于更复杂、更不健壮的建模过程。
如果在建模完成后有更多的数据可用,则可以从以前的解开始对新数据进行增强。如果使用独立的测试集来监视改进,以防止对新数据的过度拟合,则通常可以提高准确性。虽然精度的提高一般小于对合并数据重新进行整个分析所得到的精度,但是节省了大量的计算。
当随机存取内存不足以存储整个数据集时,可以在连续子集的数据上运用boosting。boosting可以应用于“arcbites”[Breiman(1997)]数据顺序读入内存,每次从当前的解开始,在时间允许的情况下在先前的子集上循环。同样,重要的是使用一个独立的测试集来停止对每个单独子集的训练,此时合并近似的估计精度开始下降。