[Python]決定木の理論と実装を徹底解説してみた

#はじめに
今回は、決定木の理論についてまとめていきます。

お付き合い頂ければ幸いです。

#理論編
それでは最初に、決定木の理論についてまとめていきます。

##決定木の概要
決定木の可視化すると以下のようになります。

今回はirisデータセットによる分類を`scikit-learn`の`export_graphviz`を用いて可視化しました。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/411258/8caed2cf-dafd-e1a2-60ff-220b97c20bcb.png)

決定木とは、上の画像のようにデータをある条件に従って分割することにより、データの分類または回帰のモデルを作成するアルゴリズムです。分類を行う分類木と回帰を行う回帰木を総称して決定木と呼びます。

かなりシンプルなアルゴリズムのため、他の複雑なアルゴリズムと比較すると精度は出にくい傾向にあります。

しかし、モデルの解釈性が非常に高いです。

モデルの見た目が上の図のように木のように見えることから、決定木と呼ばれています。

代表的なアルゴリズムに`CART`や`C5.0`などがあります。

アルゴリズムについては[こちら](http://halharahal.com/?page_id=1734)の記事を参考にしてください。

ここでは、ニクラス分類を繰り返し行う`CART`のアルゴリズムを取り扱っていきます。

##決定木の基準について
ここまでで決定木のおおよその概要をつかめたと思います。

ここでは、分類木と回帰木の枝分かれの基準について考えていきます。
###分類木の基準
それでは、分類木はどのような判断基準に従って枝分かれさせていくのかを考えていきましょう。

結論から書くと、分類は`不純度`が最小になるような分割をする、という基準に基づいて特徴量や閾値を決めています。

`不純度`とはどれくらい多くのクラスが混じりあっているかを数値にしたもので、`誤分類率`や`ジニ指数`、`交差エントロピー誤差`を用いて表されます。一つのノードに一つのクラスの観測地がある状態のとき、不純度は0になります。

###回帰木の基準
回帰木は、平均二乗誤差で表されるコスト関数を定義し、そのコスト関数の重み付き和が最小となるように特徴量や閾値を選択します。

#数式について
##分類木
不純度を数式を用いて表してみましょう。以下の図について考えます。![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/411258/81eb99e8-4c2a-df21-6e23-068682d965c7.png)

$N_m$個の観測値を持つ領域$R_m$におけるクラスkの観測値の割合を以下のように定義します。

```math
\hat{p}_{mk} = \frac{1}{N_m} \sum_{x_i\in R_m}^{}I(y_i = k)
```

図の上から三番目の段に注目します。mは領域の番号を表しており、m=1で`gini=0.168`の領域を表し、m=2で`gini=0.043`の領域を表しています。kはクラスラベルを表していて、今回は`value`の部分の左からクラス１、クラス２、クラス３と定義しています。

数式にすると難しそうですが、実際に計算すると以下のようになります。

```math
\hat{p}_{11} = \frac{0}{54} \quad \hat{p}_{12} = \frac{49}{54}  \quad \hat{p}_{13} = \frac{5}{54} 
```

なんとなく数式の意味が理解できたでしょうか。

この$\hat{p}$を用いて、以下の三つの関数で不純度を表現します。

####誤分類率
```math
\frac{1}{N_m} \sum_{x_i\in R_m}^{}I(y_i \neq k(m)) = 1-\hat{p}_{mk}
```

####ジニ指数
```math
1 - \sum_{k=1}^{K}\hat{p}_{mk}
```

####交差エントロピー誤差
```math
-\sum_{k=1}^{K}\hat{p}_{mk}log\hat{p}_{mk}
```

`sklearn`で標準で用いる不純度の関数はジニ指数になっているので、実際にジニ指数を計算してみましょう。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/411258/da6f793a-e187-8422-b10b-b1e62a49f4ab.png)

三段目の不純度をジニ指数を用いて計算します。

左の`gini=0.168`のノードのジニ指数は以下の式になりますね。

```math
1 - (\frac{0}{54})^2 - (\frac{49}{54})^2 - (\frac{5}{54})^2 = 0.168
```

当然答えは0.168になります。上の式がsklearnが内部で行っている計算式になります。ついでに右の`gini=0.043`のノードの計算も行いましょう。

```math
1 - (\frac{0}{46})^2 - (\frac{1}{46})^2 - (\frac{45}{46})^2 = 0.043
```

こちらも一致しましたね。それでは、それぞれのにデータの重みをかけることで全体の不純度を計算しましょう。以下の式になります。

```math
\frac{54}{100} ×0.168 + \frac{46}{100} ×0.043 = 0.111
```

これで全体の不純度が導出できました。決定木はこの不純度を小さくするような特徴量や閾値を選択することにより、モデルを構築しています。

##回帰木
回帰木において、コスト関数を以下のように定義します。

```math
\hat{c}_m = \frac{1}{N_m}\sum_{x_i \in R_m}^{}y_i\\
Q_m(T) = \frac{1}{N_m} \sum_{x_i \in R_m}(y_i - \hat{c}_m)^2
```

$\hat{c}_m$がそのノードに含まれる観測値の平均を表しているので、コスト関数は平均二乗誤差となります。

このコスト関数はそれぞれのノードについて計算されるため、そのコスト関数の重み付き和が最小となるように特徴量や閾値を設定します。

# 実装編
ここまでで、分類木と回帰木の理論が理解できたと思います。

それでは、これから決定木の実装を行っていきましょう。

##分類木の実装


それでは分類木を実装しています。

まず、分類するデータを作成します。

```.py
from sklearn.datasets import make_moons
from sklearn.model_selection import train_test_split
from sklearn.tree import DecisionTreeClassifier
from sklearn.tree import export_graphviz
from matplotlib.colors import ListedColormap
import graphviz

moons = make_moons(n_samples=300, noise=0.2, random_state=0)
X = moons[0]
Y = moons[1]
X_train, X_test, Y_train, Y_test = train_test_split(X, Y, random_state=0, stratify=Y)

plt.figure(figsize=(12, 8))
mglearn.discrete_scatter(X[:, 0], X[:, 1], Y)
plt.show()
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/411258/0b5dcdf9-35a9-698d-dd6f-923b4b8213dc.png)

このデータを分類するモデルを作成していきます。

以下のコードです。

```.py
clf_model = DecisionTreeClassifier(max_depth=3)
clf_model.fit(X_train, Y_train)
print(clf_model.score(X_test, Y_test))
```

>0.8933333333333333

そこそこの精度がでましたね。

この分類木のモデルを可視化してみましょう。以下のコードです。

```.py
dot_data = export_graphviz(clf_model)
graph = graphviz.Source(dot_data)
graph.render('moon-tree', format='png')
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/411258/5e41eb76-82d4-97ac-182e-8f561a7b516c.png)

`graphviz`を用いた可視化については、[こちら](https://analytics-note.xyz/machine-learning/graphviz-decision-tree/)の記事を参考にしてください。

このように、決定木はモデルの解釈性が非常に高いです。今回のモデルでは`max_depth=3`としたため、深さが3のモデルになっています。

分類後のモデルを可視化してみましょう。以下のコードです。

```.py
plt.figure(figsize=(12, 8))
_x1 = np.linspace(X[:, 0].min() - 0.5, X[:, 0].max() + 0.5, 100)
_x2 = np.linspace(X[:, 1].min() - 0.5, X[:, 1].max() + 0.5, 100)
x1, x2 = np.meshgrid(_x1, _x2)
X_stack = np.hstack((x1.ravel().reshape(-1, 1), x2.ravel().reshape(-1, 1)))
y_pred = clf_model.predict(X_stack).reshape(x1.shape)
custom_cmap = ListedColormap(['mediumblue', 'orangered'])
plt.contourf(x1, x2, y_pred, alpha=0.3, cmap=custom_cmap)
mglearn.discrete_scatter(X[:, 0], X[:, 1], Y)
plt.show()
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/411258/0b611ef9-c244-37c2-7b09-2e6ba6812681.png)


格子点を用いて色を変化させる方法については[こちら](https://qiita.com/kazukiii/items/8f0e9559524c0bb0f5ba)の記事を参考にしてください。

`_x1`と`_x2`でx軸方向とy軸方向の格子点の領域を指定し、`x1, x2 = np.meshgrid(_x1, _x2)`で格子点を作成しています。

`X_stack = np.hstack((x1.ravel().reshape(-1, 1), x2.ravel().reshape(-1, 1)))`の部分で100×100の格子点を一次元配列に変化させた後、10000×1の二次元配列に変換し、水平方向に結合して10000×2のデータに変換しています。

`y_pred = clf_model.predict(X_stack).reshape(x1.shape)`の部分で10000×2のデータを0と1のデータに変換し、それを100×100のデータに変換しています。データを分離する線の片側が0でもう一方が1になります。

`plt.contourf(x1, x2, y_pred, alpha=0.3, cmap=custom_cmap)`の部分で、等高線を描画します。色を`cmap=custom_cmap`で指定しています。

ここまでで分類木の実装は終了です。


##回帰木の実装
それでは回帰木の実装を行いましょう。

データを準備して描画してみましょう。

```.py
import mglearn
from sklearn.tree import DecisionTreeRegressor
import numpy as np
import matplotlib.pyplot as plt
from sklearn.tree import export_graphviz
import graphviz

X, Y = mglearn.datasets.make_wave(n_samples=200)

plt.figure(figsize=(12, 8))
plt.plot(X, Y, 'bo', ms=15)
plt.show()
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/411258/03a3f65c-75c5-a8c4-957b-290d1a42ae1a.png)

このようなデータになっています。それではモデルを作成しましょう。

```.py
tree_reg_model = DecisionTreeRegressor(max_depth=3)
tree_reg_model.fit(X, Y)
print(tree_reg_model.score(X, Y))
```


>0.7755211625482443

`score`で$R^2$を表示させることができます。

>score(self, X, y[, sample_weight])	Returns the coefficient of determination R^2 of the prediction.

あまり精度はよくありませんね。

次のコードでモデルを可視化してみましょう。

```.py
dot_data = export_graphviz(tree_reg_model)
graph = graphviz.Source(dot_data)
graph.render('wave-tree', format='png')
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/411258/82cb830a-bc65-c20e-0733-043be48f0439.png)

このように、回帰木は(決定木もそうですが)モデルの解釈性がとても良いですね。

それでは、次のコードで回帰直線を図示してみましょう。

```.py
X1 = np.linspace(X.min() - 1, X.max() + 1, 1000).reshape(-1, 1)
y_pred = tree_reg_model.predict(X1)
plt.xlabel('x', fontsize=10)
plt.ylabel('y', fontsize=10, rotation=-0)
plt.plot(X, Y, 'bo', ms=5)
plt.plot(X1, y_pred, 'r-', linewidth=3)
plt.show()
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/411258/1bc5338d-0101-3ecf-88a1-1d49b290866a.png)


図示してみてわかる通り、あまり正しくはありませんね。

###深さを指定せずに実装
それでは、次は深さを実装せずに実装してみましょう。

深さを指定しない以外は同じなので、同じ手順は関数にまとめましょう。

```.py
import mglearn
from sklearn.tree import DecisionTreeRegressor
import numpy as np
import matplotlib.pyplot as plt
from sklearn.tree import export_graphviz
import graphviz

X, Y = mglearn.datasets.make_wave(n_samples=200)

tree_reg_model_2 = DecisionTreeRegressor()

tree_reg_model_2.fit(X, Y)

print(tree_reg_model_2.score(X, Y))
```

>1.0

$R^2$が1になりました。いったいどのようなモデルになっているのでしょうか。以下で分岐を図示しましょう。

```.py
def graph_export(model):
    dot_data = export_graphviz(model)
    graph = graphviz.Source(dot_data)
    graph.render('test', format='png')

graph_export(tree_reg_model_2)
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/411258/ee4b010f-92eb-7dc0-dceb-d8abad71b7fd.png)

恐ろしいほどの分岐になりました。もう、どのように分岐しているか読めませんね。

次のコードでモデルを図示しましょう。

```.py
def plot_regression_predictions(tree_reg, x, y):
    x1 = np.linspace(x.min() - 1, x.max() + 1, 500).reshape(-1, 1)
    y_pred = tree_reg.predict(x1)
    plt.xlabel('x', fontsize=10)
    plt.ylabel('y', fontsize=10, rotation=-0)
    plt.plot(x, y, 'bo', ms=5)
    plt.plot(x1, y_pred, 'r-', linewidth=1)
    plt.show()


plot_regression_predictions(tree_reg_model_2, X, Y)
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/411258/78ccc094-c72c-4705-2d7d-f6ab6b0c1cce.png)


上の図を見て頂ければ分かりますが、これは明らかに過学習になっていますね。

これでは未知のデータを予測することができないので、適切に木の深さを設定する大切さが理解できると思います。

#終わりに
今回はここまでになります。

ここまでお付き合い頂きありがとうございました。