超初心者が学ぶランダムフォレスト・分類回帰樹木 (CART法)①

[超初心者が学ぶランダムフォレスト・分類回帰樹木 (CART法)⓪](https://qiita.com/LTtoose/items/0a355aa40d2a9cfedcf8)の続き
前回の導入に続いて、やっと本格的になります
初めて数式書きましたが、めっちゃ大変ですね。先人の投稿者の皆様を尊敬します。
#CART法について：簡単に復習
Classification And Regression Tree (CART)
・説明変数の尺度が名義尺度～間隔尺度まで扱える
・応答変数が連続の場合は回帰樹木、分類変数の場合は分類樹木と呼ばれる
・樹木図は上→下へふしが分岐する（条件がTrueの場合は左へが基本）
・終結ふしの値は応答の予測値であり、回帰の場合はふしに属する個体の平均値が、分類の場合には多数決方式が用いられる
#樹木の構築
CART法での樹木構築3ステップ
①前進過程：樹木の成長
・ふしが規則に従って分割していく
・ある停止基準に達するまで分割し続ける
②後進過程：樹木の刈り込み
・①で成長した樹木を刈り込み、様々な樹木列を作る
③最適モデルの選択：最適な樹木の設定
・刈り込んで得られた樹木列の中から、最適な樹木を選択する
##樹木の成長
樹木は、「ふし」から新しい娘ふしを生成することで成長します。その過程で指標となるものは**不均一性測度**です。これは回帰でも分類でも同じ指標ですが、算出方法は異なります。
(以下 n=1,・・・,nとし、ynは応答変数（回帰の場合は連続、分類の場合はクラス）とします。)
###回帰樹木の不均一性測度
回帰樹木では、ふし内の誤差の均一性が基盤になり、ふし内の応答予測値の平均を用いた回帰残差をふし内均一性測度として用います。
つまり、ふしt内の回帰残差 **$\large r\small reg$**は、残差平方和として

$\normalsize r\small reg \normalsize = \frac{1}{N(t)}  \sum_{n∈t}  \\{y\small _n -\normalsize \overline y(t)\\}^2 $

と表せます。ここで$\overline y(t)$はふしt内の応答の標本平均です。ふし内の各個体と平均の残差を表す式ですね。
さて、これを用いてふしtでの**リスクの再代用推定値**と呼ばれる不均一測度

$R(t)=pr(t)r\small reg\normalsize (t)=\frac{1}{N} \sum_{n∈t}  \\{y\small _n -\normalsize \overline y(t)\\}^2$

が定義されます。$pr(t)=N(t)/N$で、データの全個体がふしtに属する確率です。
###分類樹木の不均一測度
分類の場合、応答変数$y\small n$はK個のクラス{1, ・・・, K}のいずれかの値を取るとします。ふしtにおいて、応答変数がクラスkを取る確率は

$pr(k|t)=\sum_{n∈t} 𝟙(y\small n \normalsize =k)/N(t)$

です。左辺はふしtがkの値を取る条件付確率の表記、右辺は応答変数$y\small n$がkの時に1を返す関数です。これをN(t)と、ふし内の個体数で除しているので、確率の推定値になることがわかります。
分類の場合、応答予測値は多数決で決まりますが、$\hat k(t)$が多数決で応答予測値と決まったと仮定しましょう。ここで、ふしtの不均一測度として次の3種類が考えられます

* 誤分類率：$r\small c1\normalsize (t)=1-pr(\hat k(t)|t)$
* ジニ係数：$r\small c2\normalsize (t)=\sum_{1≤k≤K}pr(k|t)(1-pr(k|t))$
* エントロピー基準：$r\small c3\normalsize (t)=-\sum_{1≤k≤K}pr(k|t)\log pr(k|t)$

はい、わけわからないですよね。誤分類率は直感的にわかりますが、ジニ係数とエントロピーについては何してるんだ状態です。（見た感じジニ係数はベルヌーイ分布の分散に似てる？）
ほんとにこいつらで不均一性とやらが測れるのか、具体的に代入してみましょう。

| A | B | C |
|:-----------:|:------------:|:------------:|
| 0.1       | 0.1        | 0.8         |
| 0.3     | 0.3      | 0.4       |
例えばこのようにクラス分けされるふしがあったとしましょう。Cが予測値になりますね。上の3つの式に代入すると

| 誤分類率 | ジニ係数 | エントロピー |
|:-----------:|:------------:|:------------:|
| 0.2       | 0.34        | 0.28         |
| 0.6     | 0.66      | 0.47       |
たしかに、1行目と2行目で値が違います。Cの確率推定値が大きい1行目の方が、指標の値も小さいのでどうやらこれらの式を信用しても良さそうです。

$r\small cm\normalsize (t), m=1, 2, 3$を不均一性測度として、分類樹木でのリスクの再代用推定値は回帰の時と同じく

$R(t)=pr(t)r\small cm\normalsize (t)$　　と定義されます。どの測度を用いるかはモデル次第といったところか。。
###分枝ルール
リスクの再代用推定値R(t)が計算できました。もちろんこの値は小さい方が良いですよね？ということで今度はこの減少量を最大にする分枝規則$S(t)$を適用します。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1031677/f34ba8a4-622a-6493-6ead-ce2ef3e01463.png)
図のように分枝する時、R(t)の減少量は
$ ΔR(s(t), t)=R(t)-(R(t_L)+R(t_R))$
で表されます。根幹ふし{t1}から、各ふしで上式に従って娘ふしが次々と生成されるのですが、停止基準があります。

Stop1: リスクの再代用推定値が規定値より小さくなる→例：t内の応答変数が完全に均一になる
Stop2: すべての終結ふしの個体数が規定値より小さくなる
Stop3: 終結ふしの数が規定値よりも大きくなる

これで樹木の生成は完了です。ここで得られる樹木を$T_\max$として話を進めましょう。
##樹木の刈り込み
###罰則項
説明変数$x$に対する応答予測値は、$x$が属する終結ふしのある代表値というお話はしました。
樹木$T$のデータへの（見かけの）適合測度は終結ふし$\widetilde{T}$のリスクの再代用推定値の総和が測度になります。
これを $R(T)=\sum_{t\in \widetilde{T}}R(t)$ と表します。
樹木Tは大きくなるほど、この値は小さくなるので良い推定値と見えますが、これは予測性能が良いことを意味するわけではありません。過学習ってやつですね。
これには終結ふしの数$|\widetilde{T}|$が関係しているので、これを罰則項として先ほどの式に加えて新たな測度を定義したものが
$R_α(T)=R(T)+α|\widetilde{T}|$ です。$(α≥0)$
こちらは複雑度コストと呼びます。樹木が大きくなると$R(T)$は大きくなるけど、増えた終結ふしの数の分だけペナルティを与えています。そのペナルティの重みづけしているのがαです。
ここで、Tは$T_{max}$の部分樹木を表します。最大樹木$T_{max}$から刈り込む過程で生成される樹木のことです。（式では$T\preceq T_{max}$ と書く）
この複雑度コスト$R_α(T)$に基づき、
$T(α)= \ argmin_{T\preceq T_{max}} R_α(T)$
で定められる部分樹木$T(α)$が求まります。この過程を最小複雑度コスト刈り込みと呼びます。

さて、ここまでよくわからない言葉と数式ばかりで発狂しかけてるのでアルゴリズムを記述して図に表して見ましょう。

###アルゴリスム
1.$\alpha_0=0$、$k=0$ とする。$R(T_0)＝R (T_{max})$となる部分樹木を$T_0$ とする
2. $T_k$の任意の非終結ふし$t ∈ T_k- \widetilde T_k$として$g_k(t)=\frac{R(t)-R(T_{kt})}{|\widetilde T_{kt}|-1}$が最小となる$t=t^{\*}$を求める
3. $T_k$から枝$T_{kt^{\*}}$を刈り込み、部分樹木$T_{k+1}$、つまり$(T_k-T_{kt^{\*}}) \cup t^{\*}$を得る。また、$\alpha_{k+1}=g_k(t^{*})$とおく
4. $T_{k+1}$が樹木{$t_1$}に到達すれば終了。そうでなければ$k=k+1$として2.に戻る

こちらの論文[樹木構造接近法と最近の発展](https://www.jstage.jst.go.jp/article/jscswabun/18/2/18_KJ00004574188/_article/-char/ja/)の図を拝借して上記アルゴリズムを考えましょう。
1.樹木$T_0$を設定します
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1031677/04a24119-8498-d425-3136-91d8e68f5a8f.png)

2.非終結ふし（〇のふし）の集合は $T_k$ (今はk=0)から非終結ふしを除いたものですね。なのでt=2,3,5が集合となります。今回は例として $t^{*}=t_2$
として進めてみましょう。
 $\alpha_1=g_0(t_2)$として更新されます。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1031677/9a34c9b1-05a6-970a-5080-923c9f32a198.png)


3.次の刈り込まれた樹木$T_1$は$(T_0 -T_{02}*)\cup t_2$と 同等です。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1031677/66949eb1-f815-c726-2423-80e0e5e0eb8a.png)
これを根幹ふし${t_1}$に到達するまで繰り返します。

この流れを一般化し、樹木$T_{max}$から$T_k(k=0,1,・・,K)$が生成されるとすると、適用される式は以下のようになります。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1031677/0832dced-f1da-d2dd-6315-a72987407b95.png)


$α_0<α_1<α_2<・・・<α_K$

$T_{max} \succeq T_0 \succ  ... \succ T_K = $$ \{ t_1 \}$
kが大きくなるほどにパラメーターαは大きくなり、終結ふしの重みは大きくなる、というパラメーターαの変化をしっかりと理解しましょう。
樹木の大きさとパラメーターαはトレードオフの関係にあって、バランスを取ることが大切なんですね！
##最適樹木の選択
###交差検証法
ここまでで部分樹木列{$T_j$} ($j=0,1,,,J$)が得られました。最良の樹木の決定方法は、機械学習ではお馴染みの交差検証法（クロスバリデーション）です。
CART法では通常、交差検証法に基づいて、
$R^{*}(T_j)$= $E[L(y, d_{ \widetilde T_j }(x))]$のより良い推定値を構成し、予測性能の高い樹木を選択します。$d_{ \widetilde T_j }(x)$は、$x_n$が落ちる終結ふしにおける$y_n$の予測子です。
ここで$L(y_n, a)$について
回帰問題では$L(y_n, a)=(y_n-a)^2$, 分類問題では$L(y_n, a)= 𝟙 (y_n\neq a)$, クラス確率問題では$L(y_n, a)= \sum^K_{k=1} (𝟙 (y_n=a)-a_k)^2$で定義されます。
$y_n$と予測子の差が大きいほど、これらの値も大きくなりそうですね。

ここで用いるのはV重交差確認法というもので、以下の流れで進みます

1.全体標本$\psi=${${(y_n, x_n), n=1,2,,,N}$}を部分標本$\psi^1_{sub}...\psi^V_{sub}$にランダム分割する
2.$\psi_{cv}^{(v)}$ = $\psi$$-$$\psi_{sub}^{(v)}$をv番目の学習標本として、$\psi^1_{sub}...\psi^V_{sub}$を検証用のテスト標本にする
3.$\psi_{cv}^{(v)}$を用いて樹木$T_{max}^{(v)}$を構築、刈り込みし、部分樹木列{$T_j^{(v)}$}を得る
4.複雑度$\alpha_j'$で指定される$T_{max}^{(v)}(\alpha_j')$に対応する予測子を$d_j^{(v)}=d_{ \widetilde T_{max(\alpha_j')} }(x)$とする
5.予測子$d_j^{(v)}$に、テスト標本の$d_{sub}^{(v)}$を当てはめ、$R^{*}(T_j)$の交差確認推定値、
$R^{cv}(T_j)=\frac{1}{N}\sum^V_{v=1}\sum_{(y_n, x_n)\in d_{sub}^{v}}L(y_n, d_{j}^{(v)}(x_n)), j=0,,,,J$を得る
ここで、予測子$d_j^{(v)}(x_n)$は、個体nの属する終結ふし$t \in \widetilde T_{max(\alpha_j')}$において、回帰では$\overline y(t)$、分類では多数決で決められたクラス$\hat k(t)$のことを指します。

さて、この交差検証法は良い予測性能を持つ最適樹木を選ぶことが目的でした。これにおいて部分樹木列{$T_j$}$(j=0,,,J)$の中から$R^{cv}(T_j)$の最小値
$R^{cv}(T_{j0})=\min_jR^{CV}(T_j)$を持つ樹木$T_{j0}$が選択されます。前項とやっていることは同じですね。

※交差確認法では、標本サイズが**200~300**くらいあればおおよそうまくいくとされているので、サンプル選びの時は気を付けましょう

###1SEルール
$R^{cv}(T_{j0})$により選択される樹木が不安定な時、1SEルールというものがしばしば推奨されます。SEは標準誤差です。
1SEルールを適用した樹木を$T_{j1}$とすると、樹木番号$j1$は
$R^{cv}(T_{j1})\leq R^{cv}(T_{j0}) + \widehat {SE}(R^{cv}(T_{j0}))$
を満たす最大の整数になり、1SE分足したことで保守的な樹木が出来上がります。他にも2SEルールというのもあり、樹木が煩雑な場合にこれらのルールを用いて、樹木の簡明化をしているようです。

##まとめ
CART法での樹木構築3ステップ
①前進過程：樹木の成長
・ふしが規則に従って分割していく→リスクの再代用推定値が小さくなる方向へ進む（誤分類率、ジニ係数、エントロピー基準）
・ある停止基準に達するまで分割し続ける→3つの停止基準（リスクの再代用推定値、ふし内の個体数、ふし数）
②後進過程：樹木の刈り込み
・①で成長した樹木を刈り込み、様々な樹木列を作る→罰則パラメーターαと終結ふしの関係に注意
③最適モデルの選択：最適な樹木の設定
・刈り込んで得られた樹木列の中から、最適な樹木を選択する→交差検証法を用いる。標本サイズは200~300以上が理想


体力が尽きたのでここまでです。数式、説明ばかりで複雑ですが、ポイントを押さえて自分で樹木図を書いてみたり、数字を入れて計算してみると理解が深まると思います。

###参考
[樹木構造接近法と最近の発展](https://www.jstage.jst.go.jp/article/jscswabun/18/2/18_KJ00004574188/_article/-char/ja/)