GitHubコードから読み解く！DecisionTreeClassifier@scikit-learnのfeature_importances_計算方法
はじめに
昨今の人工知能ブームにて、機械学習によるプロダクトやサービスが、珍しくなくなってきました。
従来、人間が下してきた判断を、機械が代替することで、自動化を促す形です。
自動運転の歩行者識別や、スマホの顔認証などが、それですね。
人間の判断は、次々と機械に置き換えられていくことでしょう。
尚、機械学習製品のリリースには2つの壁があります。
それは、「精度の高さ」と「説明性の高さ」です。
また、それらには、一般にトレードオフの関係性があります。
共にクリアできれば問題ないですが、難しい場合は、そのトレードオフと向き合わなくてはなりません。
製品の目的を意識して、調整する必要があるのです。
以下は、よく知られた記事で、「もしも、火災保険の引受査定をAIが行ったら」というものです。
絵の女性は、保険加入をDENIED（拒否）されてしまいました。
理由を尋ねますが、AIは「私の判断アルゴリズムは、説明性の低いものですので何とも…。」と答えます。
果たしてユーザーは、それで納得するでしょうか？
答えは当然NOですよね。
最も直感的で、説明性の高い機械学習手法「決定木」
ここで、紹介したいのは「決定木」という機械学習手法と、その説明性についてです。
決定木は、幾つかの条件分岐を経ることにより、判断を下すアルゴリズムです。
先ほどの火災保険引受を判断する決定木を作るとすると、例えば、以下のようになります。
図から、女性が査定を通らなかったのは…
① 海岸までの距離が100m未満であったから
② 海抜が3m未満であったから
という理由からだと分かります。
判断基準が正当性はさておき、説明性は高いですね。
決定木は、このように判断基準がハッキリした、説明性の高いアルゴリズムとなります。
決定木は、どの項目が効いて判断に至ったかを「feature_importance_（重要度）」で示す
先ほどの決定木の例ですと、パッと見で判断基準が把握できますが、例えば次のパターンはどうでしょう？
scikit-learnにてロードが可能なbreast-cancerデータに対して、scikit-learnのDecisionTreeClassifierで学習した決定木の、分岐条件を可視化したものです。
各Boxには、分岐条件や、その条件が学習データにおいてもたらした結果などが記載されています。
尚、最下段以外の各Boxに記載されている、「≦」記号で表現されるものが、分割条件です。
その条件を満たすデータは左に、満たさないデータは右にと流れていきます。
この位の複雑さで、割と判断条件の把握が難しくなってきますね。
追うのが一苦労です。
決定木は、解こうとする問題によっては、10段以上もの条件を積み重ねる必要があったりもしますので、そうなると更に大変です。
試しに、breast-cancerデータに対して、判断の精度を高めるために、もう少し条件を複雑にしてみましょう。
うーん、これを解釈しようとすると、気が遠くなってきますね…。
さて、どう整理すると、判断基準が上手く解釈できるでしょうか…。
ここで役立つのが、scikit-learnのDecisionTreeClassifierクラスに組み込まれている「feature_importances_（重要度）」の計算機能です。
判断の精度向上に対して、どの説明変数が、どの位の貢献をしたであろうかを、定量的に示してくれるものです。
その定量は、決定木学習時の収集情報から、計算します。
例えば、先程紹介した、条件が最大で7段存在するものですと、下記のような「feature_importances_」の図となります。
breast cancerのデータにおいて存在する説明変数は、30種存在しますが、その各種毎に、貢献度合いが割り振られます。
貢献度合いは、貢献割合として求められ、合計値が1.0になるようスケーリングされています。
こうして、割合で示してくれると、とても分かりやすいですね。
判断根拠が分かると、予測器への信頼も高まるかと思います。
「feature_importances_」の計算方法を抽象的に理解する
さて、上記の便利な指標ですが、これはどうやって求められているか。
先ずは、抽象的に説明します。
「feature_importances_」は、breast-cancerか否かの判断に、寄与している説明変数ほど高い値が示されます。
その計算は…
① 上段にある説明変数ほど高くなる
② 登場回数が多い条件項目が高い
という根拠の下に計算されています。
①については、条件が多段になっている場合の、より上段での条件に使用される説明変数ほど、重要度が高くなるという意味です。
決定木は、より効率の良い分割条件を探索しながら、上段から下段に向けて、分割条件を重ねていくアルゴリズムです。
よって、上段から下段に条件推移していくに連れて、条件適用のデータ数が段々と少なくなっていく形となります。
つまり、上段の方が、より多くのデータを精度高く分割する、マクロ向きの条件が設定されています。
逆に、下段の条件は、一部のデータに効果的な、ミクロ向きな条件であるとも言えます。
そうした場合、全体の分割に寄与しているのは、上段の条件と言えるでしょう。
その為、より上段にある条件項目の方が重要度が高くなるのです。
②については、複数の条件がある内、より多くの条件に使用される説明変数ほど、重要度が高くなるという意味です。
重要度の計算としては、先ず、各条件について、重要度（分割貢献度）を計算します。
次に、条件への使用が重複する項目については、重要度の合算を行います。
尚、各条件毎の重要度は全て、ZERO以上の値が算出されます。
よって、登場回数が多いほど、多く合算されて、重要度が大きくなりやすいのです。
その為、仮に、決定木上の全ての条件が、ある1種類の説明変数だけでまかなわれてしまった場合、その説明変数の「feature_importanes_」値が1.0となり、他の説明変数の重要度がZEROとなります。
以上が、「feature_importanes_」の抽象的な説明となります。
尚、この「feature_importanes_」の計算アルゴリズムは、決して絶対的に正しい指標という訳ではありません。
あくまでも、ヒューリスティックに考え出された形でしょう。
例えば、お客さんに対して、この指標を報告するような場合は、そういった前提を忘れずに、慎重に向き合うと良いかと思います。
「feature_importances_」の計算方法を具体的に理解する
さて、次に、「feature_importances_」の計算方法を具体的に説明させていただこうと思います。
尚、計算方法は色々な流儀があるようで、今回はあくまで、scikit-learnのソースコードに実装されているものを紹介します。
scikit-learnでの計算方法は、以下となります。
【feature_importances_の計算方法】
① 条件分割前のデータ集合について、「不純度」を計算する
② ①に条件分割前のデータ集合のデータ数を掛ける
③ 「項目X ≦ 値XX.XX」条件を満たすデータ集合について、「不純度」を計算する
④ ③に「項目X ≦ 値XX.XX」条件を満たすデータ集合のデータ数を掛ける
⑤ 「項目X > 値XX.XX」条件を満たすデータ集合について、「不純度」を計算する
⑥ ⑤に「項目X > 値XX.XX」条件を満たすデータ集合のデータ数を掛ける
⑦ 「②－④－⑥」を計算し、その値を条件の対象となっている説明変数に加算していく
（※加算前は、各説明変数について、ZEROポイントから開始する）
⑧ ①〜⑦までの計算工程を、全条件に対して実施しながら、加算を繰り返す
⑨ ⑧の各説明変数毎の加算結果を、全説明変数の加算ポイント合計で割る
（※そうして求められた各説明変数毎のポイントは、合計が1となる）
以上。
「不純度」については、条件適用対象データにおける目的変数の値（今回の例だと、breast cancerか否か）の混ざり具合を、定量的に示す指標のことです。
数値が低いほど混ざりが少なく、数値が高いほど混ざりが多い、という風な感じです。
最も「不純度」が低くなるのは、グループデータ内に1種類の目的変数値のみが存在する状態です。
つまり、決定木の分岐条件は、条件を経て、極力「不純度」が低い分割を行っていく形です。
尚、scikit-learnのDecisionTreeClassifierクラスにおいては、「不純度」が「Gini係数」や「Cross Entropy」などから選ぶことができます。
一般的には（Default設定では）、「Gini係数」が用いられます。
尚、Gini係数の計算方法は下記となります。
【Gini係数の計算方法】
ｋは、目的変数の持つ値の種類数、つまり、識別のクラス数を示します。
例えば、2クラスの属性を持つデータが、100件あったとします。
内訳は、Aクラスのデータが70件、Bクラスのデータが30件だったとします。
すると、Gini係数で説明される不純度は、「1-(70/100)²-(30/100)²=0.42」となります。
もし、Aクラスのデータ数が85件、Bクラスのデータ15件だったとすると、Gini係数の値は「1-(85/100)²-(15/100)²=0.255」です。
Aクラスのデータ数が100件、Bクラスのデータ0件ならば、「1-(100/100)²-(0/100)²=0」です。
Aクラスのデータ数が50件、Bクラスのデータ50件ならば、「1-(50/100)²-(50/100)²=0.5」です。
属性別のデータ数が均衡している程、不純度が高くなり。
一方の属性に片寄っている程、不純度が低くなる。
そんな計算式となっています。
「feature_importances_」の計算方法を厳密に再現してみる
さて、上記の具体的な計算方法を、実装してみようと思います。
プログラムは、幾つかのパート別に説明を実施していきます。
ご自分のPC上で実行してみようという方は、恐縮ですが、パート別のプログラムを、最終的に連結する形で、コピー踏襲していただければと思います。
【パート1：データ準備】
先ずは、ライブラリのインポートから、データ読み込みまで実施するプログラムです。
# import lib
import numpy as np
import matplotlib.pyplot as plt
import graphviz
import pydotplus as pdp
from sklearn import tree
from sklearn.datasets import load_breast_cancer# data load
data = load_breast_cancer()
X = data.data
y = data.target[:, np.newaxis]
X_title = data.feature_names# check
print('np.shape(X) = [%d, %d]' % np.shape(X))
print('np.shape(y) = [%d, %d]' % np.shape(y))
print('len(X_title) = [%d]' % len(X_title))
print('np.unique(y) = %s' % np.unique(y))
上記コード実施によって、以下ログが出力されます。
np.shape(X) = [569, 30]
np.shape(y) = [569, 1]
len(X_title) = [30]
np.unique(y) = [0 1]scikit-learnに保持されている「breast cancer」か否かを識別するための2クラスのデータで、件数が569件、説明変数の種類数が30列の構成となっています。
【パート2：決定木学習と「feature_importances_」確認】
次に、ロードしたデータを決定木で学習し、「feature_importance_」を可視化します。
model = tree.DecisionTreeClassifier(random_state=0, max_depth=10)
model = model.fit(X, y)feat_imp = model.feature_importances_
sort_idx = np.argsort(-feat_imp)plt.figure(figsize=(12,5),dpi=100)
plt.rcParams["font.size"] = 20plt.bar(np.arange(len(sort_idx)), feat_imp[sort_idx])
plt.ylim(-0.05, 1.05)
plt.xlim(-0.9, (len(sort_idx) -0.1))
plt.xticks(np.arange(len(sort_idx)))
pld.gca().set_xticklabels(X_title[sort_idx], rotation=90)
plt.ylabel('feature_importances_')
plt.grid()plt.show()
「feature_importances_」が、バーグラフで可視化されました。
先に説明したように、30項目についての「feature_importances_」値を合計すると1.0になるようにスケールされている為、グラフは必ず0から1の間の値域に収まります。
【パート3：「不純度」等の計算の再現】
次に、「feature_importances_」を計算する前に、その計算元となる値を計算してみたいと思います。
ここからが本題です。
ロードしたデータと、モデルが保持している分割条件から、学習時の集計情報を再現してみます。
以下が整理です。
（再現計算の元となるINPUT変数）
① model.tree_.node_count：条件の数（「2^(深さ+1)-1」で求まる値、条件の無いノードも含む）
② model.tree_.feature：条件対象項目を条件毎に管理
③ model.tree_.threshold：条件値を条件毎に管理
④ model.tree_.children_left：条件を満たした場合の、次の条件インデックスを管理
⑤ model.tree_.children_right：条件を満たさなかった場合の、次の条件インデックスを管理
↓
（再現計算の対象となるOUTPUT変数）
⑥ model.tree_.impurity：学習時に、学習データより算出された「不純度」値を管理
⑦ model.tree_.n_node_samples：学習時に、学習データより算出された条件分岐前後のデータ数を管理
⑧ model.tree_.value：学習時に、学習データより算出された条件毎の満たしたデータ数／満たさなかったデータ数を管理
より理解を深めるための補足計算になります。
尚、「不純度」値はGini係数により求めます。
# prepare variable to stock result of reproduction
class_num = len(np.unique(y))
n_node_samples_rep = np.zeros(model.tree_.node_count)
impurity_rep = np.zeros(model.tree_.node_count)
value_rep = np.zeros([model.tree_.node_count, class_num])
y_work = -np.ones([len(y), model.tree_.node_count])# reproduction at top node
n_node_samples_rep[0] = len(y)
for class_i in range(class_num):
value_rep[0, class_i] = np.sum(y == class_i)
impurity_rep[0] = 1
for class_i in range(class_num):
impurity_rep[0] -= ((value_rep[0, class_i] / len(y)) ** 2)
y_work[:, 0] = y.ravel()# reproduction at all node
for node_i in range(model.tree_.node_count):
y_tmp = y_work[:, node_i].copy()
if (model.tree_.children_left[node_i] > 0):
feature_tmp = model.tree_.feature[node_i]
threshold_tmp = model.tree_.threshold[node_i]
X_tmp = X[:, feature_tmp]
y_left_idx = (X_tmp <= threshold_tmp)
y_right_idx = (X_tmp > threshold_tmp)
y_left = y_tmp.copy()
y_left[y_right_idx] = -1
y_right = y_tmp.copy()
y_right[y_left_idx] = -1
node_j = model.tree_.children_left[node_i]
node_k = model.tree_.children_right[node_i]
n_node_samples_rep[node_j] = np.sum(y_left >= 0)
n_node_samples_rep[node_k] = np.sum(y_right >= 0)
for class_i in range(class_num):
value_rep[node_j, class_i] = np.sum(y_left == class_i)
value_rep[node_k, class_i] = np.sum(y_right == class_i)
impurity_rep[node_j] = 1
impurity_rep[node_k] = 1
for class_i in range(class_num):
impurity_rep[node_j] -= ((value_rep[node_j, class_i] / n_node_samples_rep[node_j]) ** 2)
impurity_rep[node_k] -= ((value_rep[node_k, class_i] / n_node_samples_rep[node_k]) ** 2)
y_work[:, node_j] = y_left.ravel()
y_work[:, node_k] = y_right.ravel()# check value
print('model.tree_.impurity = \n%s' % model.tree_.impurity)
print('')
print('impurity_rep = \n%s' % impurity_rep)
print('')
print('model.tree_.n_node_samples = \n%s' % model.tree_.n_node_samples)
print('')
print('n_node_samples_rep = \n%s' % n_node_samples_rep.astype('int'))
print('')
print('model.tree_.value = \n%s' % model.tree_.value.reshape([np.shape(model.tree_.value)[0], np.shape(model.tree_.value)[2]]))
print('')
print('value_rep = \n%s' % value_rep)# check diff
print()
print('diff of impurity = %.8f' % np.sum(np.abs(model.tree_.impurity - impurity_rep)))
print('diff of n_node_samples = %.8f' % np.sum(np.abs(model.tree_.n_node_samples - n_node_samples_rep)))
print('diff of value = %.8f' % np.sum(np.abs(model.tree_.value - value_rep[:, np.newaxis, :])))
最後のprint文によって、以下のLOGが確認できると思います。
diff of impurity = 0.00000000
diff of n_node_samples = 0.00000000
diff of value = 0.00000000scikit-learnの算出値と、再現値に差分が無いことが確認できました。
【パート4：「feature_importances_」の計算の再現】
最後に、「feature_importances_」を計算してみたいと思います。
パート3で行った計算よりも単純です。
feat_imp = np.zeros([np.shape(X)[1]])for node_i in range(model.tree_.node_count):
if (model.tree_.children_left[node_i] > 0):
node_j = model.tree_.children_left[node_i]
node_k = model.tree_.children_right[node_i]
root_n_node_samples = model.tree_.n_node_samples[node_i]
left_n_node_samples = model.tree_.n_node_samples[node_j]
right_n_node_samples = model.tree_.n_node_samples[node_k]
root_impurity = model.tree_.impurity[node_i]
left_impurity = model.tree_.impurity[node_j]
right_impurity = model.tree_.impurity[node_k]
root_feature = model.tree_.feature[node_i]
feature_importances_rep[root_feature] += ((root_n_node_samples * root_impurity ) -
(left_n_node_samples * left_impurity ) -
(right_n_node_samples * right_impurity))
normalizer = np.sum(feature_importances_rep)
if normalizer > 0.0:
# Avoid dividing by zero (e.g., when root is pure)
feature_importances_rep /= normalizer# check value
print()
print('model.feature_importances_ = \n%s' % model.feature_importances_)
print()
print('feature_importances_rep = \n%s' % feature_importances_rep)# check diff
print()
print('diff of feature_importances_ = %.8f' % np.sum(np.abs(model.feature_importances_ - feature_importances_rep)))
最後のprint文にて、以下のログが出力されると思います。
diff of feature_importances_ = 0.00000000これにて、「feature_importances_」の計算が再現できました。
（お付き合いいただいた方、ご苦労さまでございました…。🙇♂💦）
おまけ
ちなみに、上記コードの後に、下記コードを追加すると、決定木の木構造を可視化することができます。
（Jupyter上で見たい場合）
dot_data = tree.export_graphviz(model,
out_file=None,
feature_names=X_title,
class_names='target',
filled=True,
rounded=True,
special_characters=True)
graph = graphviz.Source(dot_data)
graph（画像ファイルに出力したい場合）
dot_data = tree.export_graphviz(model,
out_file=None,
feature_names=X_title,
class_names='target',
filled=True,
rounded=True,
special_characters=True)graph = pdp.graph_from_dot_data(dot_data)
graph.write_png('./graph.png')おわりに
以上で、DecisionTreeClassifier@scikit-learnにおける、「feature_importances_」値の計算方法を説明となります。
ここまで読んで下さった方、本当にありがとうございます。
個人的に、人が書いたコードをリバースして、キッチリ確認するのは気持ちが良いです…。
今後、他アルゴリズムでもやっていこうと思います。
P.S
記事の修正を指摘してくれた同僚の土屋さんに感謝します。🙇♂
