Skip to article frontmatterSkip to article content
Site not loading correctly?

This may be due to an incorrect BASE_URL configuration. See the MyST Documentation for reference.

勾配ブースティング決定木

Gradient Boosting Machine

複数の関数に重みをかけて足し合わせた関数

f(x)=f0(x;θ0)+β1f1(x;θ1)+⋯+βMfM(x;θM)\renewcommand{\argmin}{\mathop{\rm arg~min}\limits} f(x) = f_0(x; \theta_0) + \beta_1 f_1(x; \theta_1) + \cdots + \beta_M f_M(x; \theta_M)
(1)

の形で予測モデルを構築することを考える。ここでθ0,…,θM\theta_0, \dots, \theta_Mは関数を形づくるパラメータ（例えば線形回帰の重みや決定木の分岐の閾値）である。

このモデルはβ1,β2,…,βM\beta_1, \beta_2, \dots, \beta_Mとθ0,θ1,…,θM\theta_0, \theta_1, \dots, \theta_Mのパラメータを推定する必要がある。 今回はすべてのパラメータを一度に学習するのではなく、βmfm(x;θm)\beta_m f_m(x; \theta_m)を一つずつ学習していく方法を考える。具体的には次のように行う。

前向き段階的加法モデリング（forward stagewise additive modeling）

  1. f0(x)=0f_0(x) = 0で初期化

  2. m=1m = 1からMMまでについて、

    1. パラメータを推定する：(βm,θm)=arg min⁡β,θ∑i=1NL(yi,fm−1(xi)+βf(xi;θ))(\beta_m, \theta_m) = \argmin_{\beta, \theta} \sum^N_{i=1} L(y_i, f_{m-1}(x_i) + \beta f(x_i; \theta))

    2. 新たなモデルを足す：fm(x)=fm−1(x)+βmf(x;θm)f_m(x) = f_{m-1}(x) + \beta_m f(x; \theta_m)

これを前向き段階的加法モデリング（forward stagewise additive modeling）という。ブースティングはこの方法でアンサンブル学習を行う。

誤差関数が二乗誤差L(y,f(x))=(y−f(x))2L(y, f(x)) = (y - f(x))^2の場合、

L(yi,fm−1(xi)+βf(xi;θ))=(yi−fm−1(xi)−βf(xi;θ))2=(residuali,m−1−βf(xi;θ))2\begin{align} L(y_i, f_{m-1}(x_i) + \beta f(x_i; \theta)) &= (y_i - f_{m-1}(x_i) - \beta f(x_i; \theta))^2 \\ &= (\text{residual}_{i,m-1} - \beta f(x_i; \theta))^2 \end{align}
(2)

となり、m−1m-1回目のモデルの残差residuali,m−1=yi−fm−1(xi)\text{residual}_{i,m-1} = y_i - f_{m-1}(x_i)を近似するようにmm回目のモデルβf(xi;θ)\beta f(x_i; \theta)を学習させていると捉えることができる。残差が大きければそれだけ訓練中に重視されるため「間違えた箇所を重点的に学習する手法」とも捉えることができる。

最適化の観点からの説明

勾配降下法

数理最適化において関数の最小化問題

\renewcommand{\b} when command \b does not yet exist; use \newcommand

\renewcommand{\b}[1]{ \boldsymbol{#1} }
\min_\b{x} f(\b{x})

を解く方法のひとつに勾配降下法（gradient descent method）あるいは最急降下法（steepest descent method）と呼ばれるものがある。これは目的関数の微分のベクトルである勾配

Undefined control sequence: \b at position 10: \nabla f(\̲b̲{x}) = 
\begin{…

\nabla f(\b{x}) = 
\begin{bmatrix}
    \frac{ \partial f(x_1) }{ \partial x_1 }\\
    \vdots \\
    \frac{ \partial f(x_m) }{ \partial x_m }
\end{bmatrix}

を用いて

Undefined control sequence: \b at position 1: \̲b̲{x}_{m} = \b{x}…

\b{x}_{m} = \b{x}_{m-1} - \alpha_{m-1} \nabla f(\b{x}_{m-1})

という値の更新を何度も繰り返して最適化を行っていく。ここでα\alphaは学習率と呼ばれるパラメータで、値の更新量が多すぎると最適解を通り過ぎてしまうことがあるので小さめの値を乗じて更新幅を抑えるために用いられる。

最終的にMM回反復して得た最適解x∗x^*は

x∗=x0−α1∂f(x1)∂x1−α2∂f(x2)∂x2−⋯−αM∂f(xM)∂xMx^* = x_0 - \alpha_1 \frac{\partial f(x_{1})}{\partial x_{1}} - \alpha_2 \frac{\partial f(x_{2})}{\partial x_{2}} - \cdots - \alpha_M \frac{\partial f(x_{M})}{\partial x_{M}}
(6)

となり、ブースティングにより得られる予測モデル

f(x)=f0(x;θ0)+β1f1(x;θ1)+⋯+βMfM(x;θM)f(x) = f_0(x; \theta_0) + \beta_1 f_1(x; \theta_1) + \cdots + \beta_M f_M(x; \theta_M)
(7)

と同様に重み付き和の形になる。

ブースティング

ブースティングは勾配降下法を機械学習で行っていると捉えることができる。

機械学習においては予測値f(x)f(x)と実測値yyの誤差の最小化問題

min⁡f(x)L(y,f(x))\min_{f(x)} L(y, f(x))
(8)

を解きたいため、勾配は誤差関数の予測モデルによる微分∂L(y,f(x))∂f(x)\frac{ \partial L(y, f(x)) }{ \partial f(x) }によって得られる。

二乗誤差L(y,f(x))=12(y−f(x))2L(y, f(x)) = \frac{1}{2} (y - f(x))^2の場合、負の勾配は残差である

−∂L(y,f(x))∂f(x)=y−f(x)=residual- \frac{ \partial L(y, f(x)) }{ \partial f(x) } = y - f(x) = \text{residual}
(9)

前向き段階的加法モデルの節で「二乗誤差の場合は残差を近似するように学習している」と述べた。

L(yi,fm−1(xi)+βf(xi;θ))=(residuali,m−1−βf(xi;θ))2L(y_i, f_{m-1}(x_i) + \beta f(x_i; \theta)) = (\text{residual}_{i,m-1} - \beta f(x_i; \theta))^2
(10)

これにより、学習されるモデルβf(x;θ)\beta f(x; \theta)は負の勾配を学習するようになり、最終的にそれらの和となるモデルは勾配降下法を解いた状態を近似することになる。

正則化つきGBDT

nn個の観測データがあり、mm次元の特徴量があるとする。

D={(xi,yi)},∣D∣=n,xi∈Rm,y∈R\mathcal{D} = \{ (x_i, y_i) \}, |D| = n, x_i \in \mathbb{R}^m, y \in \mathbb{R}
(11)

勾配ブースティング決定木のモデルは次のように表される

y^i=ϕ(xi)=∑k=1Kfk(xi)\hat{y}_i = \phi(x_i) = \sum^K_{k=1} f_k(x_i)
(12)

ここで

  • fk∈Ff_k \in \mathcal{F}は予測器

  • F={f(x)=wq(x)}\mathcal{F} = \{ f(x) = w_{q(x)} \}は回帰木の空間

    • q:Rm→Tq: \mathbb{R}^m \to Tは入力データxxを木の各葉のインデックスに割り振る写像。TTは各木の葉の数

    • w∈RTw \in \mathbb{R}^Tは葉の重み（weights）と呼ばれ、予測に使われた葉の出力値。予測値はwq(x)w_{q(x)}の和となるので予測値のベースでもある

学習の際は正則化付き誤差関数

L(ϕ)=∑i=1nl(y^,yi)+∑kΩ(fk)\mathcal{L}(\phi) = \sum^n_{i=1} l(\hat{y}, y_i) + \sum_k \Omega(f_k)
(13)

を最小化する。ここでllは微分可能な凸関数である誤差関数で、Ω\Omegaは正則化項

Ω(f)=γT+12λ∣∣w∣∣2\Omega(f) = \gamma T + \frac{1}{2} \lambda ||w||^2
(14)

学習は加法的に行うためtt番目の誤差は次のようになる。

L(t)(ϕ)=∑i=1nl(yi,y^(t−1)+ft(xi))+Ω(ft)\mathcal{L}^{(t)}(\phi) = \sum^n_{i=1} l(y_i, \hat{y}^{(t-1)} + f_t(x_i)) + \Omega(f_t)
(15)

テイラー展開による近似

この誤差を二次近似したものを使うことで計算量を削減することもできることが知られている（Friedman et al., 2000）

L(t)≈∑i=1n[l(yi,y^(t−1))+gift(xi)+12hift2(xi)]+Ω(ft)\mathcal{L}^{(t)} \approx \sum^n_{i=1} [ l(y_i, \hat{y}^{(t-1)}) + g_i f_t(x_i) + \frac{1}{2} h_i f_t^2(x_i) ] + \Omega(f_t)
(16)

ここで

gi=∂l(yi,y^(t−1))∂y^(t−1)hi=∂2l(yi,y^(t−1))∂(y^(t−1))2g_i = \frac{ \partial l(y_i, \hat{y}^{(t-1)}) }{\partial \hat{y}^{(t-1)} }\\ h_i = \frac{ \partial^2 l(y_i, \hat{y}^{(t-1)}) }{\partial (\hat{y}^{(t-1)})^2 }
(17)

定数項を省略すると

L~(t)=∑i=1n[gift(xi)+12hift2(xi)]+Ω(ft)\tilde{\mathcal{L}}^{(t)} = \sum^n_{i=1} [ g_i f_t(x_i) + \frac{1}{2} h_i f_t^2(x_i) ] + \Omega(f_t)
(18)

葉jjにおけるインスタンス（サンプル）の集合をIj={i∣q(xi)=j}I_j = \{i|q(x_i) = j\}と表記すると、次のように書き換えることができる

L~(t)=∑i=1n[gift(xi)+12hift2(xi)]+γT+12λ∑j=1Twj2=∑j=1T[(∑i∈Ijgi)wj+12(∑i∈Ijhi+λ)wj2]+γT\begin{align} \tilde{\mathcal{L}}^{(t)} &= \sum^n_{i=1} [ g_i f_t(x_i) + \frac{1}{2} h_i f_t^2(x_i) ] + \gamma T + \frac{1}{2} \lambda \sum^T_{j=1} w_j^2 \\ &= \sum^T_{j=1} [ (\sum_{i\in I_j} g_i) w_j + \frac{1}{2} (\sum_{i\in I_j} h_i + \lambda ) w^2_j ] + \gamma T \end{align}
(19)

固定した木の構造q(x)q(x)について、葉jjの最適な重みは

wj∗=−∑i∈Ijgi∑i∈Ijhi+λw^*_j = -\frac{ \sum_{i\in I_j} g_i }{ \sum_{i \in I_j} h_i + \lambda }
(20)

となる。

導出

葉jjについての部分だけ取り出して導関数を0とおいて整理する

∂L~j(t)∂wj=∑i∈Ijgi+(∑i∈Ijhi+λ)wj=0  ⟹  (∑i∈Ijhi+λ)wj=−∑i∈Ijgi  ⟹  wj=−∑i∈Ijgi∑i∈Ijhi+λ=wj∗\frac{ \partial \tilde{\mathcal{L}}^{(t)}_j }{\partial w_j } = \sum_{i\in I_j} g_i + (\sum_{i\in I_j} h_i + \lambda ) w_j = 0 \\ \implies (\sum_{i\in I_j} h_i + \lambda ) w_j = -\sum_{i\in I_j} g_i \\ \implies w_j = -\frac{ \sum_{i\in I_j} g_i }{ \sum_{i\in I_j} h_i + \lambda } = w_j^*
(21)

L1の場合

L~(t)=∑i=1n[gift(xi)+12hift2(xi)]+γT+12λ∑j=1Twj2+α∑j=1T∣wj∣=∑j=1T[(∑i∈Ijgi)wj+12(∑i∈Ijhi+λ)wj2]+γT\begin{align} \tilde{\mathcal{L}}^{(t)} &= \sum^n_{i=1} [ g_i f_t(x_i) + \frac{1}{2} h_i f_t^2(x_i) ] + \gamma T + \frac{1}{2} \lambda \sum^T_{j=1} w_j^2 + \alpha \sum^T_{j=1} |w_j| \\ &= \sum^T_{j=1} [ (\sum_{i\in I_j} g_i) w_j + \frac{1}{2} (\sum_{i\in I_j} h_i + \lambda ) w^2_j ] + \gamma T \end{align}
(22)

L1の場合

L~(t)=∑i=1n[gift(xi)+12hift2(xi)]+γT+12λ∑j=1Twj2+α∑j=1T∣wj∣=∑j=1T[(∑i∈Ijgi)wj+α∣wj∣+12(∑i∈Ijhi+λ)wj2]+γT\begin{align} \tilde{\mathcal{L}}^{(t)} &= \sum^n_{i=1} [ g_i f_t(x_i) + \frac{1}{2} h_i f_t^2(x_i) ] + \gamma T + \frac{1}{2} \lambda \sum^T_{j=1} w_j^2 + \alpha \sum^T_{j=1} |w_j| \\ &= \sum^T_{j=1} [ (\sum_{i\in I_j} g_i) w_j + \alpha |w_j| + \frac{1}{2} (\sum_{i\in I_j} h_i + \lambda ) w^2_j ] + \gamma T \end{align}
(23)

∑j=1T[⋅]\sum^T_{j=1} [\cdot]の内側を整理すると

wj∗=−∑i∈Ijgi±α∑i∈Ijhi+λw^*_j = -\frac { \sum_{i\in I_j} g_i \pm \alpha } { \sum_{i \in I_j} h_i + \lambda }
(24)

最適な重みwj∗w_j^*を誤差関数に戻すと

L~(t)(q)=−12∑j=1T(∑i∈Ijgi)2∑i∈Ijhi+λ+γT\tilde{\mathcal{L}}^{(t)}(q) = -\frac{1}{2} \sum^T_{j=1} \frac{ (\sum_{i\in I_j} g_i)^2 } { \sum_{i\in I_j} h_i + \lambda } + \gamma T
(25)

となり、これ木の構造qqの品質をスコアリングする関数として使うことができる。

導出
L~(t)=∑j=1T[(∑i∈Ijgi)(−∑i∈Ijgi∑i∈Ijhi+λ)+12(∑i∈Ijhi+λ)(−∑i∈Ijgi∑i∈Ijhi+λ)2]+γT=∑j=1T[−(∑i∈Ijgi)2∑i∈Ijhi+λ+12(∑i∈Ijgi)2∑i∈Ijhi+λ]+γT=∑j=1T[−12(∑i∈Ijgi)2∑i∈Ijhi+λ]+γT\begin{align} \tilde{\mathcal{L}}^{(t)} &= \sum^T_{j=1} [ (\sum_{i\in I_j} g_i) (-\frac{ \sum_{i\in I_j} g_i }{ \sum_{i \in I_j} h_i + \lambda }) + \frac{1}{2} (\sum_{i\in I_j} h_i + \lambda ) (-\frac{ \sum_{i\in I_j} g_i }{ \sum_{i \in I_j} h_i + \lambda })^2 ] + \gamma T \\ &= \sum^T_{j=1} [ -\frac{ (\sum_{i\in I_j} g_i)^2 }{ \sum_{i \in I_j} h_i + \lambda } + \frac{1}{2} \frac{ (\sum_{i\in I_j} g_i)^2 }{ \sum_{i \in I_j} h_i + \lambda } ] + \gamma T \\ &= \sum^T_{j=1} [ - \frac{1}{2} \frac{ (\sum_{i\in I_j} g_i)^2 }{ \sum_{i \in I_j} h_i + \lambda } ] + \gamma T \end{align}
(26)
機械学習
ブースティング
機械学習
LightGBM