DecisionTreeClassifierの中を見てみる
概要
機械学習をするためのメソッドは数多くありますが、kaggleなどで人気が高いのはXGBMやLightGBMなどの決定木を利用するものです。今回はそんな決定木が学習をするアルゴリズムを調べてみました。
詳細の前に
決定木とはなんぞや、についてはこの記事などを読んでみてください。
また、以下の前提を設けます。
・scikit-learnのDecisionTreeClassifierのfit()メソッドのソースコード を調べてみる。
・予測する変数は1種類
・サンプル毎の重みはつけない
・不純度の指標はエントロピー
全体像
1 stackに学習データを入れる
2 stackの先頭のデータ（X）を取り出す。
3 Xの特徴量（fi）を1つ選ぶ
4 Xを区切る位置（pos）を選ぶ
5 そのように区切った場合の不純度の減少量を求める
6 区切り位置をずらす
7 4に戻る
8 特徴量を変える
9 3に戻る
10 不純度の減少量が最大となる特徴量、区切り位置でXを分割
11 これ以上分割できない or 不純度が減少しなくなったら終了。
そうでない場合は分割した各々をstackに入れる
12 2に戻る
13 stackが空になったら終了
不純度というのはデータがどれだけ綺麗に分かれているかの指標で、小さい方が綺麗に分かれているということになります。不純度の計算方法はいくつかありますが、今回はエントロピーを例にしたいと思います。（エントロピーについては上記記事を参照）
全体像をコードで
BaseDecisionTreeのfit()メソッドを見てみます。
最初の方はパラメータのチェックで、実質TreeBuilderのbuild()メソッドを呼んでいるだけです。
sklearn/tree/_classes.py
class BaseDecisionTree(MultiOutputMixin, BaseEstimator, metaclass=ABCMeta):def fit(self, X, y, sample_weight=None, check_input=True,X_idx_sorted="deprecated"):
//略
builder.build(self.tree_, X, y, sample_weight)TreeBuilderのbuild()メソッドで上記の処理が行われます。
ここでは、DepthFirstTreeBuilderクラスのbuild()メソッドを見てみます。
sklearn/tree/_tree.pyx
cdef class DepthFirstTreeBuilder(TreeBuilder):
"""Build a decision tree in depth-first fashion."""cpdef build(self, Tree tree, object X, np.ndarray y, np.ndarray sample_weight=None):
"""Build a decision tree from the training set (X, y)."""
//略 変数の定義
with nogil:
# push root node onto stack
rc = stack.push(0, n_node_samples, 0, _TREE_UNDEFINED, 0, INFINITY, 0)
//1 stackに学習データを入れる
if rc ==-1:
# got return code -1 - out-of-memory
with gil:
raiseMemoryError()
while not stack.is_empty():
stack.pop(&stack_record) //2 stackの先頭データを取り出す
//略 データを取り出す処理先頭でstackにデータ全体をpushし、stackにデータがある間はそのデータを取り出し、処理するというループが走ります。
//略 データを取り出す処理
is_leaf = (is_leaf or
(impurity <= min_impurity_split))
// 3-10 データを分割
if not is_leaf:
splitter.node_split(impurity, &split, &n_constant_features)
# If EPSILON=0 in the below comparison, float precision
# issues stop splitting, producing trees that are
# dissimilar to v0.18
is_leaf = (is_leaf or split.pos >= end or
(split.improvement + EPSILON <
min_impurity_decrease))
//これ以上分割できるかチェック
node_id = tree._add_node(parent, is_left, is_leaf,
split.feature, split.threshold,
impurity, n_node_samples,
weighted_n_node_samples)splitter.node_splitが上記3-10のデータを分割する特徴量と位置を探して最適なところで分割する処理です（詳細は後述します）。
tree.add_nodeで決定木を更新します。
#略
if not is_leaf:
# Push right child on stack
rc = stack.push(split.pos, end, depth +1, node_id, 0,
split.impurity_right, n_constant_features)
if rc ==-1:
break
# Push left child on stack
rc = stack.push(start, split.pos, depth +1, node_id, 1,
split.impurity_left, n_constant_features)
if rc ==-1:
break
1個上のコードでis_leafの確認を行います。より良い分割位置が見つからなかった（split.pos >= end）もしくは不純度の減少量が小さい
（split.improvement + EPSILON < min_impurity_decrease）場合はstackには何も追加しません。そうでない場合は分割した2つのデータをそれぞれstackに追加します（11）。
これでwhileループは終了です。is_leafでない場合は分割したデータそれぞれについて同じ手順で分割を繰り返していきます。
データの分割
それではnode_splitメソッドの処理を見てみようと思います。分割を行うsplitterクラスも複数あるのですが、今回は最適な分割を探すBestSplitterのnode_splitメソッドを見てみたいと思います。
sklearn/tree/_splitter.pyx
cdef class BestSplitter(BaseDenseSplitter):
"""Splitter for finding the best split."""本メソッドでは「全体像」のところで説明したように、データを分割するのに最適な特徴量と分割位置を探します。
特徴量のループ
まずは特徴量についてのループを回しますが、ここでは計算量節約のために1つずつ順番に調べるのではありません。与えられたデータサンプルで、データの値がほとんど変わらない場合はその特徴量については分割位置の探索はスキップして次の特徴量に移ります。
（結構長いです。ここはあくまで計算量を減らすテクニックなので次の章に移ってもいいです。）
while (f_i > n_total_constants and
# Stop early if remaining features
# are constant
(n_visited_features < max_features or
# At least one drawn features must be non constant
n_visited_features <= n_found_constants +
n_drawn_constants)):
n_visited_features +=1
f_j = rand_int(n_drawn_constants, f_i - n_found_constants,
random_state)
if f_j < n_known_constants:
# f_j in the interval [n_drawn_constants,
# n_known_constants[
features[n_drawn_constants], features[f_j] =
features[f_j], features[n_drawn_constants]
n_drawn_constants +=1 //14
else:
# f_j in the interval [n_known_constants, f_i -
# n_found_constants[
f_j += n_found_constants
# f_j in the interval [n_total_constants, f_i[
current.feature = features[f_j]
# Sort samples along that feature; by
# copying the values into an array and
# sorting the array in a manner which utilizes the
cache more effectively.
for i in range(start, end):
Xf[i] =self.X[samples[i], current.feature]
//15 データを取り出し小さい順にソート
sort(Xf + start, samples + start, end - start)
if Xf[end -1] <= Xf[start] + FEATURE_THRESHOLD:
features[f_j], features[n_total_constants] =
features[n_total_constants], features[f_j]
// 16 Xが定数
n_found_constants +=1
n_total_constants +=1
else:
f_i -=1
// 17 Xが定数でない
features[f_i], features[f_j] =
features[f_j], features[f_i]探索中、特徴量は以下のように並んでいます。
: n_drawn_constants ・・・ 前回の調査で定数だと分かっている、かつ今回すでに当たったもの
n_drwan_constants: n_known_constants ・・・ 前回の調査ですでに定数だと分かっているもの
n_known_constants: n_found_constants ・・・ 今回の調査で定数だと分かったもの
n_found_constants: fi ・・・前回の調査では定数ではなかったもので、今回まだ当たっていないもの
fi: ・・・ 今回すでに当り、定数でないと分かったもの
各特徴量を当たったら、この規則が成り立つように並び替えます。
まず（省略してしまいましたが）fiは最後の特徴量からスタートします。n_drawn_constantsからfi-n_found_constantsの範囲で（未調査の範囲）特徴量fjを選びます。
・fjがn_known_constantsより小さい場合は、その特徴量は定数なのでfjをn_drawn_constantsの次に移動して次のループへ行きます。（14）
・それ以外の場合はデータXのfjの列を取り出し、値の小さい順に並べます。(15)
⭐︎ Xの最大が最小とほぼ同じならば定数とみなしてn_found_constantsに追加し、次のル
ープへ行きます。(16)
⭐︎ Xが定数でなければfiを1つ前に進めてfjと入れ替えます。その後、分割位置の探索に入
ります。(17)
分割位置のループ
self.criterion.reset()
p = start
while p < end:
while (p +1< end and
Xf[p +1] <= Xf[p] + FEATURE_THRESHOLD):
//18 Xの値が変わらない場合はスキップ
p +=1
p +=1
if p < end:
current.pos = p
# Reject if min_samples_leaf is not guaranteed
if (((current.pos - start) < min_samples_leaf) or
((end - current.pos) < min_samples_leaf)):
continue
self.criterion.update(current.pos) //19 ターゲットの分割
# Reject if min_weight_leaf is not satisfied
if ((self.criterion.weighted_n_left < min_weight_leaf) or
(self.criterion.weighted_n_right < min_weight_leaf)):
continue
current_proxy_improvement
=self.criterion.proxy_impurity_improvement()
//20 不純度の増減分を計算
if current_proxy_improvement > best_proxy_improvement:
best_proxy_improvement = current_proxy_improvement
# sum of halves is used to avoid infinite value
current.threshold = Xf[p -1] /2.0+ Xf[p] /2.0
if ((current.threshold == Xf[p]) or
(current.threshold == INFINITY) or
(current.threshold ==-INFINITY)):
current.threshold = Xf[p -1]
best = current # copy //21 最適な分割やっと分割位置（ここでは変数p）の探索が始まります。
まず、前後でXの値がほとんど変わらない場合は分割位置とはならないので飛ばします（18）。
19で選んだ分割位置でターゲット（y）を分割します。（詳細は後述）
20ではその分割によって不純度がどれだけ変わるかを計算しています。これまでの最適な分割位置よりも変化量が大きい場合はその分割位置を最適な分割位置として更新します（21）
if best.pos < end:
partition_end = end
p = start
//22 最適な分割位置で分割
while p < partition_end:
if self.X[samples[p], best.feature] <= best.threshold:
p +=1
else:
partition_end -=1
samples[p], samples[partition_end] =
samples[partition_end], samples[p]
self.criterion.reset()
self.criterion.update(best.pos)
self.criterion.children_impurity(&best.impurity_left,
&best.impurity_right)
// 23 木全体で不純度を計算
best.improvement =self.criterion.impurity_improvement(
impurity, best.impurity_left, best.impurity_right)
# Respect invariant for constant features: the original
order of
# element in features[:n_known_constants] must be preserved
for sibling
# and child nodes
memcpy(features, constant_features, sizeof(SIZE_t) *
n_known_constants)
# Copy newly found constant features
memcpy(constant_features + n_known_constants,
features + n_known_constants,
sizeof(SIZE_t) * n_found_constants)
# Return values
split[0] = best
n_constant_features[0] = n_total_constants
return0ループが終わったら不純度の変化量が最大になる分割位置でデータを分割します（22）。
そして、その分割位置で木全体の不純度を計算します（23）。
これでデータの分割処理は終了です。以下ではターゲットを分割し、不純度を計算する過程を見ていきます。
ターゲットの分割
19のself.criterion.updateの中を見ていきます。
criterionはCriterion型のインスタンスです。Criterionクラスからは複数のクラスが派生していますが、今回はEntropyクラスを例にしたいと思います。（EntropyはClassificationCriterionクラスから派生します）
sklearn/tree/_criterion.pyx
cdef class ClassificationCriterion(Criterion):
"""Abstract criterion for classification."""cdef class Entropy(ClassificationCriterion):
r"""Cross Entropy impurity criterion.
This handles cases where the target is a classification taking values 0, 1, ... K-2, K-1. If node m represents a region Rm with Nm observations, then let
count_k = 1 / Nm \sum_{x_i in Rm} I(yi = k)
be the proportion of class k observations in node m.
The cross-entropy is then defined as
cross-entropy = -\sum_{k=0}^{K-1} count_k log(count_k)
"""updateメソッドはClassificationCriterionクラスのものをそのまま使います。
cdef int update(self, SIZE_t new_pos) nogil except-1:
//略 変数の定義など
if (new_pos - pos) <= (end - new_pos):
for p in range(pos, new_pos):
i = samples[p]
if sample_weight !=NULL:
w = sample_weight[i]
for k in range(self.n_outputs):
label_index = k *self.sum_stride +
<SIZE_t>self.y[i, k]
sum_left[label_index] += w
self.weighted_n_left += w
else:
self.reverse_reset()
for p in range(end -1, new_pos -1, -1):
i = samples[p]
if sample_weight !=NULL:
w = sample_weight[i]
for k inrange(self.n_outputs):
label_index = k *self.sum_stride +
<SIZE_t>self.y[i, k]
sum_left[label_index] -= w
self.weighted_n_left -= w
# Update right part statistics
self.weighted_n_right =self.weighted_n_node_samples -
self.weighted_n_left
for k in range(self.n_outputs):
for c in range(n_classes[k]):
sum_right[c] = sum_total[c] - sum_left[c]
sum_right +=self.sum_stride
sum_left +=self.sum_stride
sum_total +=self.sum_stride
self.pos = new_pos
return0n_outputsはyの次元数（つまりターゲットデータの種類数）です。ここではn_output=1として考えます。（多次元の場合も同じループが回るだけです）
sum_stridesはカテゴリ数です（多次元の場合はその最大値）。
sum_leftはカテゴリ数のサイズを持つ配列で、各要素は分割より左（つまりインデックスが小さい）のデータのうち、yがその値になる個数を表しています。sum_rightは右側、sum_totalはデータ全体について同様です。
weighted_n_left、weighted_n_right、weighted_n_node_samplesは左側、右側、全体のデータ数です。
ここでは一から計算するのではなくて、前回からの差分を計算しています。if文内は左が増えた場合、else文内は右が増えた場合です。
不純度の増減の計算
ターゲットを分割したのでそれに基づき不純度の増減を計算します。Entropyクラスのchildren_impurityメソッドです。
cdef void children_impurity(self, double* impurity_left,double* impurity_right) nogil:
//略 変数の定義など
for k in range(self.n_outputs):
for c in range(n_classes[k]):
count_k = sum_left[c]
if count_k >0.0:
count_k /=self.weighted_n_left
entropy_left -= count_k * log(count_k)
count_k = sum_right[c]
if count_k >0.0:
count_k /=self.weighted_n_right
entropy_right -= count_k * log(count_k)
sum_left +=self.sum_stride
sum_right +=self.sum_stride
impurity_left[0] = entropy_left /self.n_outputs
impurity_right[0] = entropy_right /self.n_outputsこれは定義通りエントロピーを計算しています。
参考文献
最後に
機械学習の本質に関わる部分なのか計算量を減らすためのテクニックなのか判別するのが若干大変でした。
コメントとか変数の定義とか一部の処理を省略してしまったので気になる方は是非ご自身で読んでみてください。
あと、pythonだと思っているとポインタが出てきたり、C++だと思っていると定義せずに変数が使われたりとなかなか難しかったです。
最新記事
すべて表示概要 pythonでデータ解析を行っている。解析自体はpandasを用いて行い、最終結果はpandas.DataFrameの形式で保持されている。 この結果を他のアプリケーションで利用するため、json形式でファイル出力したい。 やり方 1...
現象 raspberry piでfirestoreをimportしようとするとタイトルのエラーが発生。 from from firebase_admin import firestore ImportError: Failed to import the Cloud...
概要 フィッティングを行いたい場合、pythonならばscipy.optimize.leastsqなどでできます。 しかし、フィッティングを行う場合、フィッティングパラメータに条件を付けたい場合も多々あります。 例えば、下記のようにパラメータa、bは共に正の範囲で最適な値を...
コメント