PythonでXGBoostをちゃんと理解する(1)

なんせ、石を投げればxgboostにあたるくらいの人気で、ちょっとググれば解説記事がいくらでも出てくるので、流し読みしただけでなんとなく使えるようになっちゃうので、これまでまとまった時間を取らずに、ノリと勢いだけで使ってきた感があります。が、腹に落とすまで理解して使いこなせるようになるべく、考えながら動かしてみたいと思います。

まずは、下記の記事で大枠を理解。

・勾配ブースティング木の概念について
Gradient Boosting Decision Treeでの特徴選択 in R | 分析のおはなし。

・数式的な説明
パッケージユーザーのための機械学習(12):Xgboost (eXtreme Gradient Boosting) - 銀座で働くデータサイエンティストのブログ

・Python版xgboost本家
Python Package Introduction — xgboost 0.4 documentation


で、次に少し古いですが、このKaggleのフォーラムの順に従って理解を深めて行こうと思います。
とは言いながら、このフォーラムにある内容だけでは、理解したというには程遠い、また記事内ではR版のxgboostを使ってるので、順番には沿ってみたものの、結果的にはやることがかなり異なってきました。

Understanding XGBoost Model on Otto Data

まずはIntroduction

XGBoostについてはブラックボックスでみんなよくわかんないまま使ってるよね−だから理解を助けるためにこの記事書くよってのが主旨らしいです。

Preparatioin of the data

データをロードして、中身の確認。

import numpy as np
import pandas as pd
import xgboost as xgb

df_train = pd.read_csv("…/input/train.csv")
df_test = pd.read_csv("…/input/test.csv")

訓練データとテストデータを確認すると、特徴量は番号が振ってあるだけの味気ないものが93個。

print len(df_train.index)
print len(df_train.columns)
df_train.head()

#61878
#95
# 	id 	feat_1 	feat_2 	feat_3 		... 	 	feat_90 	feat_91 	feat_92 	feat_93 	target
#0 	1 	1 	0 	0 		... 	 	0 	0 	0 	0 	Class_1
#1 	2 	0 	0 	0 	 	... 	 	0 	0 	0 	0 	Class_1
#2 	3 	0 	0 	0 	 	... 	 	0 	0 	0 	0 	Class_1
#3 	4 	1 	0 	0 	 	... 	 	0 	0 	0 	0 	Class_1


print len(df_test.index)
print len(df_test.columns)
df_test.head()

#144368
#94
# 	id 	feat_1 	feat_2 	feat_3 	 	... 	 	feat_90 	feat_91 	feat_92 	feat_93
#0 	1 	0 	0 	0 	 	... 	 	0 	0 	0 	0
#1 	2 	2 	2 	14 	 	... 	 	0 	0 	2 	0
#2 	3 	0 	1 	12 	 	... 	 	0 	0 	0 	1

クラスは9つあり、分布には偏りあり。
xgboostの認識のために、クラスを数値に変換しておきます。

df_train["target"].value_counts()

#出力
#Class_2    16122
#Class_6    14135
#Class_8     8464
#Class_3     8004
#Class_9     4955
#Class_7     2839
#Class_5     2739
#Class_4     2691
#Class_1     1929

Class  = {u'Class_1': 0,
          u'Class_2': 1,
          u'Class_3': 2,
          u'Class_4': 3,
          u'Class_5': 4,
          u'Class_6': 5,
          u'Class_7': 6,
          u'Class_8': 7,
          u'Class_9': 8,}

Model traning

モデルトレーニングですが、本家の説明を見ながらモデリングします。
Scikit-Learnラッパーも用意されていて、馴染みのある文法でも記述できますので、慣れている方はそっちでも良いかもしれません。
が、後に出てくる変数重要度を可視化することなど考えて、通常のモデリングを選択しています。

df_train = df_train.replace(Class)
df_X = df_train.ix[:, 1:94]

df_y = df_train["target"]

#Split data set to train and test data
X_train, X_test, y_train, y_test = train_test_split(
    df_X, df_y, test_size=0.2, random_state=1)

#Get labels
labels = pd.read_csv("…/input/Train_label.csv")
feature_names = list(labels)

df_X_train=pd.DataFrame(X_train)
dtrain = xgb.DMatrix(df_X_train.as_matrix(),label=y_train.tolist(), feature_names = feature_names)

df_X_test=pd.DataFrame(X_test)
dtest=xgb.DMatrix(df_X_test.as_matrix())

パラメーターに関しては今回はデフォルト値でやってみてます。

# Set parameteres
params={"objective": "multi:softprob", 
        "eval_metric": "mlogloss",
        "eta": 0.3, 
        "max_depth": 6, 
        "min_child_weight": 1, 
        "subsample": 1, 
        "colsample_bytree": 1,
        "num_class": 9
        }

記事では、ここでクロスバリデーションの説明をしてますが、その結果をどう使うか的な記述は一切なく、この記事だけでは消化不良。
なので、別のフォーラムから、「XGBOOST and cross validation」なんてのを拾ってきました。

www.kaggle.com

大したこと書いてあるわけじゃないのですが、num_boost_roundは少なくとも200でやって、テストデータのerror rateの一番少ない回数選びなさいよとのことです。
このnum_boost_roundについては、納得出来てない部分もあるので、後ほど少し考察してみます。
とりあえずで、やってみました。133回が良さそうです。

cv=xgb.cv(params,dtrain,num_boost_round=200,nfold=10)
print(cv)

#出力
#[0]cv-test-mlogloss:1.552386+0.008590     cv-train-mlogloss:1.538197+0.001742
#[1]cv-test-mlogloss:1.301161+0.010400     cv-train-mlogloss:1.279041+0.001476
#[2]cv-test-mlogloss:1.137894+0.011714     cv-train-mlogloss:1.108999+0.001610
#[3]cv-test-mlogloss:1.023316+0.012670     cv-train-mlogloss:0.988626+0.001223
#   ・
#      ・
#[131]cv-test-mlogloss:0.489264+0.008465     cv-train-mlogloss:0.219077+0.001969
#[132]cv-test-mlogloss:0.489191+0.008445     cv-train-mlogloss:0.218099+0.001892
#[133]cv-test-mlogloss:0.489147+0.008295     cv-train-mlogloss:0.216975+0.002032
#[134]cv-test-mlogloss:0.489262+0.008136     cv-train-mlogloss:0.215958+0.002184
#   ・
#   ・
#[197]cv-test-mlogloss:0.493435+0.006741     cv-train-mlogloss:0.162068+0.001473
#[198]cv-test-mlogloss:0.493468+0.006863     cv-train-mlogloss:0.161414+0.001497
#[199]cv-test-mlogloss:0.493638+0.006781     cv-train-mlogloss:0.160707+0.001603

#Training
bst=xgb.train(params,dtrain,num_boost_round=133)

Model understanding

最近の記事に、Python版で特徴量重要度をプロットできるようプルリク出して改良して下さった方がいらっしゃったので、使ってみることにしました。
こういう方々本当に頭が下がります、ありがとうございいます。

sinhrks.hatenablog.com

import matplotlib.pyplot as plt

# 変数重要度を求める
imp=bst.get_fscore()
xgb.plot_importance(bst)

Ottoのデータだと特徴量が多すぎて、上手くプロットが見えません。泣
図のサイズを変えようと、matplotlibの方でも色々頑張ってみましたが、結果上手くいかず。
時間があれば改良のお手伝いをするかもです。

id:sinhrksさんご本人にご教示頂き、サイズを大きくすることができました。
本当にありがとうございます。

fig, ax = plt.subplots(1, 1, figsize=(7, 25))
xgb.plot_importance(bst, ax=ax)

f:id:tkzs:20151003030504p:plain


で、気を取り直して、pandasのデータフレームで、重要な特徴量をソートして出すようにしてみました。
xgboostのpython版は特徴量のラベルを引き継がないので、自分で再度作りなおして貼り付けてます。
上記図でも同様にDMatrixへのラベル引き継ぎも出来てます。コード修正済。

#Get labels created
labels = pd.read_csv("…/input/Train_label.csv")

#Sort ascending by importance
imp = pd.Series(imp)
imp = pd.DataFrame(imp)
imp.sort(columns = 0, ascending=False).head(5)

#出力
# 	0
#f67 	1426
#f24 	1347
#f86 	1192
#f25 	1139
#f40 	1125

で、記事はこの情報をfeature selectionやモデル解釈に役立てなさい、みたいに書いてあるけど、具体的に何をすべきかについては情報薄でした。。。


とここまで、ざっとフォーラムの投稿にそってやってみましたが、不完全燃焼感が半端無いので、次の3点について次回もう少し詳しく考察してみようと思います。