「バッチサイズは小さい方が良い」を勾配ブースティングで検証してみる

名刺の会社のアドベントカレンダー3日目です。

へい社ではR&Dメンバーが論文読み会を定期的に開催しています。私は直近でバッチサイズに関するものを読みましたので、それに関連する数値実験を報告しようと思います。

論文

前者の論文では、Deep Learningにおいて、SGDのバッチサイズを小さくした方が汎化性能の高くなる経験則について、バッチサイズとtrain loss functionの収束先のflatさ/sharpさという観点で分析をしています。バッチサイズが小さいと重みの更新幅の共分散が大きくなり、flatな(局所)解に収束しやすくなることで汎化性能が向上すると主張しています。

後者の論文(の前半部分)では、バッチサイズと重みの更新幅の関係性を学習率で調整した場合には汎化性能の差異がないと反論をしています。

読み会でも、バッチサイズ1にすれば(収束させるのが難しくなりそうですが)汎化性能が最良となるのか議論になりましたので、前者の論文がパラメータチューニングにおいて重要な示唆を与えているものの、そこまで強い主張ができるのか疑わしいと考えています。

今回やること

Deep Learningのパラメータチューニングについては専門家に任せるとして、勾配法でバッチサイズを小さくすれば汎化性能が高くなるというのは勾配ブースティングでも同じなので、XGBoostを用いてこれを検証してみました。 Hyperoptでチューニングして、バッチサイズは小さい方が良いという報告は聞こえてこないので結果は見えていそうですが、とりあえずやってみます。

数値実験

データセットCredit Card Fraud Detectionを用いました。 クレジットカードのトランザクションが不正かどうかをターゲットとしています。正例が0.17%とかなり不均衡なデータとなっています。

今回はバッチサイズ subsample, 学習率 learning_rate 以外はデフォルトのパラメータを用いることとします。 正則化パラメータの効き方がバッチサイズの大小によって異なることが予想されますが、ここでは捨象します。

clf = xgb.XGBClassifier(
    learning_rate=learning_rate,
    subsample=subsample,
    base_score=0.0017,
    n_estimators=num_round, n_jobs=16
)

まずは以下のパラメータで実験してみました。

Model learning_rate subsample
LB (Large Batch) 0.1 1.0
SB (Small Batch) 0.1 0.1

f:id:marugari2:20171204173739p:plain

ここではバッチサイズの小さい方がvalidationの精度が高くなりました。

学習率の調整

XGBoostの更新幅の分散は

 \displaystyle
\mathrm{Var}
\left( \varDelta w_M \right)
= \mathrm{Var}
\left( \eta \frac{G_M}{H_M} \right)
= \frac{\eta ^2}{M} \frac{1}{N} \mathrm{Var}
\left( \varDelta w_N \right)

となります。 NIPS論文に従えば、バッチサイズを1/100にするときには学習率を1/10するのが良いとわかります。 学習率を調整した場合の結果は以下のようになります。

Model learning_rate subsample
LB 0.1 1.0
SB adjusted 0.03 0.1
LB slow 0.03 1.0

f:id:marugari2:20171204173955p:plain

こちらの結果では learning_rate を小さくした効果が本質的でバッチサイズはあまり関係ないように見えます。 バッチサイズの大きい方が収束は速く、より小さい learning_rate とより多い num_round によって精度を上げやすいでしょう。

まとめ

バッチサイズの小さい方が汎化性能に優れるケースも確認できた一方で learning_rate こそが本質的に重要であると思われます。 当然ですが、ICLR論文の模式図のようにバッチサイズを小さくしたからといって劇的に汎化性能の向上は望むべくもないでしょう。

また、Deep LearningではGPUメモリの制約でバッチサイズを小さくしたいことが多い*1一方で、XGBoostでは行方向の圧縮が効いておりバッチサイズを小さくしてもリニアには高速化されず計算効率が悪くなるという事情もあります。

XGBoostのパラメータチューニングにおいては subsample を大きめで固定しておき learning_rate を小さく num_round を大きくする王道が、やはり強いというのを再認識しました。

参考

追記

DeepLearningの最適化と Hessian Freeの資料によるとHessianを使った最適化だと最適なバッチサイズも大きくなることが解説されていました。XGBoostはHessianを用いているため、論文で言われているような効果はそもそも出ないようです。Hessianを利用すると、flatならflatなりの更新量、sharpならsharpなりの更新量になるのでICLR論文のような状況は生じにくいですね。 Hessianを用いないscikit-learnのGBDTでどうなるのか関心あるところですが、再試する気力はないですね…

*1:たくさんGPUを持っているのでむしろバッチサイズを大きくしたい人達もいるようです。そういう背景があって論文にあるような研究がされてる訳ですね。

Ryzenマシン組んだのでXGBoostのGPU版やってみた

Ryzen 7 1700とGTX 1080 Tiでマシンを組んだので、動作確認がてらXGBoostGPU版を使ってみました。

タイムラインでそういう話題があったのでネタをパクったような形になってしまいましたが、私自身前からやってみたいと思っていたテーマであり、H2O.aiが最近担いでいる話でもあるので許して頂きたい。

構成

  • Ryzen 7 1700 (3.7GHz 8core/16thread)
  • GTX 1080 Ti FOUNDERS EDITION (1582MHz 11GB)

今のところ定格で動かしてます。 記事と前後してしまいますが、Deep Learning以外の機械学習アルゴリズムGPUで気軽に計算を高速化できるような予感があったのでCPUはケチりました。 将来GPUを追加する際の資金の足しにしたい。

インストー

公式のドキュメント

私はUbuntuの環境でmakeを使ってビルドしました。

CUDA

CUDAをインストールして

export PATH="$PATH:/usr/local/cuda-8.0/bin"
export LD_LIBRARY_PATH="$LD_LIBRARY_PATH:/usr/local/cuda-8.0/lib64"

あたりのパスを通しておきます。nvccが動けば問題ないはずなので、cudnnは不要と思われます。

ビルド

make/config.mkファイルの該当行を

CUB_PATH = cub
PLUGIN_UPDATER_GPU = ON

と変更します。

後は通常と同様にmakeするだけです。

インストールについてもpython3 setup.py installとすれば問題ありませんでした。(Anaconda)

使い方

updaterパラメータを変更してGPUでの計算に切り替えることができます。 今までと同じパラメータを使っている限りではGPU版をビルドしても通常のCPUによる計算となります。

param['updater'] = 'grow_gpu'
  • GPU fast histgram
param['updater'] = 'grow_gpu_hist'
param['max_bin'] = 16

max_binLightGBMとかのアレなので、問題に合わせて設定してください。

サンプルをぼんやり眺めて、GPUが1つしかないのに

param['gpu_id'] = 1

とやるとエラーになって死ぬので注意です。

実験

CPU/GPUとexact/fast histgramで計算時間を比較しました。

データは以前の記事と同様にOttoのデータを利用しました。

github.com

f:id:marugari2:20170611104703p:plain

ntread CPU exact CPU hist GPU exact GPU hist
1 143.213904 43.093807 10.793109 7.385136
2 75.984772 24.485343 10.793109 7.385136
3 56.342093 18.662806 10.793109 7.385136
4 43.276191 14.980951 10.793109 7.385136
5 35.087812 12.747035 10.793109 7.385136
6 29.960601 11.476474 10.793109 7.385136
7 26.535089 11.396709 10.793109 7.385136
8 25.082611 9.916659 10.793109 7.385136
9 24.930102 11.761317 10.793109 7.385136
10 23.934886 11.543043 10.793109 7.385136
11 23.353146 11.724661 10.793109 7.385136
12 23.154665 11.528249 10.793109 7.385136
13 22.851917 11.614104 10.793109 7.385136
14 22.556055 11.890477 10.793109 7.385136
15 22.435674 11.309677 10.793109 7.385136
16 22.297691 11.825569 10.793109 7.385136

CPUよりGPU、exactよりもfast histgramが高速ですが、公表されているベンチマーク程は差が出ませんでした。

CPU fast histgramでも十分速い。X付きのにしときべきだったかなぁ…

nvidia-smiでモニタリングしたところGPUはまだまだ余裕ある感じだったので、データが巨大になった場合には差が出てくるのかもしれません。

また今更ですが、ntreadをコア数以上に設定しても速くなっていません。

LightGBMは?

LightGBMでもGPU対応されています。 公式のチュートリアル

インストールの仕方の説明でドライバをaptで入れてたので何となく回避しました。

OpenCLを使っているようなので、QuadroAMDのグラボを持ってる人には面白いかもしれません。

最後に

今日、我々はDeep Learningに留まらず勾配ブースティングにおいてもGPUモードをポチっとするだけで計算を高速にできるようになりました。 専門家以外にも扱い易い機械学習アルゴリズムGPUで使われるようになれば、ますますGPU需要が高まることが予想され、AI半導体メーカーの株価動向からは目が離せません。

勾配ブースティング落穂拾い - マルチクラス分類について

勾配ブースティングに特有の問題という訳ではないですが、それはマルチクラス分類をone-vs-allで解いていることに留意しましょうというのが今回の記事です。

2017/04/02追記 LightGBMではマルチクラス分類の際にsoftmaxとone-vs-allのどちらで解くか選択できるようです。 また、XGBoostではsoftmaxのみがサポートされています。 この記事では、マルチクラス分類のためにクラス数と同じ分類器を使っていることを指してone-vs-allとしてしまっています。ご注意ください。

マルチクラス分類について

ここではマルチクラス分類の問題を対象に、アンサンブル学習として代表的なランダムフォレストと勾配ブースティングとを比較します。

ランダムフォレスト

ランダムフォレストによるマルチクラス分類では学習器に分類木を用いています。 scikit-learnではDecisionTreeClassifierです。

コード

分類木はマルチクラスを自然に扱えるので、n_estimatorsと同じ数だけ木を作れば求めるモデルが構築できます。

勾配ブースティング

勾配ブースティングでは弱学習器に回帰木を用いているので事情が変わってきます。 マルチアウトプットの回帰木は一般的ではありませんので、問題をone-vs-allに変換して解くことになります。 そのため勾配ブースティングでは指定したn_estimatorsとクラス数の積の数だけ木を作る必要があります。

コード

勾配ブースティングは速度と精度を高いレベルで両立していますが、マルチクラス分類では

  • クラス数と比例して計算時間が増加する
  • Dark knowledge的なものを利用できない

という欠点が出てくるので、勾配ブースティングのメリットがスポイルされていないか留意する必要がありそうです。

数値例

勾配ブースティング(というかone-vs-all)が苦手そうな例を考えてみました。 セルの中は(クラスAの割合, クラスBの割合)を示しています。

sex/height
(5%, 10%) (15%, 10%) (20%, 20%)
(15%, 15%) (5%, 15%) (20%, 30%)
(20%, 25%) (20%, 25%)

この例にone-vs-allを使うと、クラスBについてsexで分割する以上のものが出てきません。 一方で分類木を使った場合は、クラスBについてsexで分割した後でクラスAの濃い部分を抽出するheightでの分割を発見できます。

マルチクラス分類をXGBoostのみで取組むのは微妙な感じがしますね。 あまり業務でスタッキングしたくないのですが、それも避けられないように思われます。 もしくはDeep Learningするか。

参考

ブースティングやDeep Learningが一般化した時代におけるマルチクラス分類の解き方についてはOtto Group Product Classification ChallengeのForumで盛んな議論があったようです。 実戦的なテクニックも数多いようなので、関心のあるかたはそちらで調査されるのが良いと思われます。

勾配ブースティング落穂拾い - 木の構築について

このシリーズについて

XGBoost芸人を自称してちょこちょこ活動をしてきたのですが、最近になって自分の理解の甘さを痛感するようになりました。

気になった箇所を散発的に復習しているのですが、その成果を備忘録として残しておこうと思います。 今のところ体系的にまとめるつもりはないので、これを読んでも勾配ブースティングの全体像はつかめませんので悪しからず。

今回のテーマ以外にはマルチクラス分類の際の挙動等に関心を持っています。

木の構築について

勾配ブースティングでは  M+1 回目のイテレーションで誤差

 \displaystyle
\sum_{i=1}^N L \left( y_i - \sum_{j=1}^M \eta_j f_j \left( x_i \right) \right)
= \sum_{i=1}^N L \left( y_i - f^M \left( x_i \right) \right)

の勾配を上手く表現した木を構築します。

この部分の処理についてscikit-learnとXGBoostでの違いを確認します。

scikit-learn

カステラ本に準拠した処理になっています。 勾配の計算は

 \displaystyle
\frac{\partial}{\partial f^M \left( x_i \right)} L \left( y_i - f^M \left( x_i \right) \right)

となり、これを各サンプルのラベル扱いにして DecisionTreeRegressor に投げます。

コード

このとき当然 DecisionTreeRegressor は2乗誤差*1を最小にするよう葉に属するサンプルにスコアを付けます。つまり平均値となっています。

ロス関数がRMSEのケースではそれで問題ないのですが、他の関数ではズレが生じます。 そのため、木の分岐は DecisionTreeRegressor の結果を使うのですが、最終的なスコアには補正をかけて対応します。

例えば、ロス関数がMAEの場合は、葉に属するサンプルの中央値をスコアとして採用するのが、木の構造を所与としたときの最適となります。

コード

XGBoost

XGBoostでは正則化項が追加されているので最適化の目的関数が異なります。式の導出は別に譲りますが以下のようになります。

 \displaystyle
\sum_{j=1}^M \left( G_j w_j + \alpha \left| w_j \right| + \frac{1}{2} \left( H_j + \lambda \right) w_j^2 \right) + \gamma T

木の構造を所与としたとき最適な  w_j^*

 \displaystyle
w_j^* = \frac{l \left( G_j , \alpha \right)}{H_j + \lambda} \\
\mathrm{where} \quad l(x, y) =
\begin{array}{l}
x - y \quad \mathrm{if} \ x \geq y \\
x + y \quad \mathrm{else}
\end{array}

となります。これは最適化の1階条件

 \displaystyle
G_j + \alpha \mathrm{sign} \left( w_j^* \right) + \left( H_j + \lambda \right) w_j^* = 0 \\
\displaystyle
w_j^* = -\frac{G_j + \alpha \mathrm{sign} \left( w_j^* \right)}{H_j + \lambda}

によります。ここで  w_j^* \leq 0 となるためには  G_j \geq \alpha でなければなりません。

これを代入すれば目的関数から  w_j が消去できます。 したがってXGBoostでは、ある分割についてはじめから最適なスコアを付ける前提で、最適な分岐を探索していることになります。

ポエム

XGBoostは通常の勾配ブースティングより精度が良いとされており、その原因は主としてL1/L2正則化にあるとされています。 確かに通常の勾配ブースティングでもインサンプルの誤差を0にするのは簡単なので、汎化性能の高さは正則化によるものと考えるのが自然です。

しかしながら、今回見てきたようにXGBoostは通常のものより精度の高い近似計算をしています。 個人的には正則化がそれほど効いていると思えないので、最適化の工夫によって高い汎化性能が実現しているのではないかと疑っています。 XGBoostの収束が速い分だけ、同水準のインサンプル誤差で比較したときに森の複雑性がより小さいと期待できるのではないでしょうか。

参考

*1:オプションで絶対誤差も指定可能

4教科で珍しく良い点取った生徒のカンニングを疑ってはいけない、「竜王を弁護する」について

はじめに

将棋ソフト不正疑惑に関する事件について新たな記事が出ていたので、また計算機を叩いてみました。

三浦九段不正疑惑について、渡辺明竜王を弁護する

2017/02/27追記 当該記事は削除されてしまっています。一時、ハフィントンポストが全文を掲載してくれていたのですが、そちらも削除されています。ご覧なりたい方は自身で魚拓を漁ってください。

三浦九段については晴れて「シロ」判定がされて、これから名誉回復をどうしていくかという段で、このような記事が出てくるのには違和感があります。

ネット上でもあまり評判良くないようなので、今更内容に噛み付く必要もないのですが、データジャーナリズムの手習いとして公開してみます。

前提

ここでは、記事中にある一致率の推定値を真の値として扱います。

一致率については、将棋ソフトや局面の違いによって値に大きな差異が出る可能性が指摘されています。 そのような数字を基に組み立てたロジックは無意味だと一蹴することもできるのでしょうが、一致率を受け入れた上でも反論が可能と私は考えました。

設定

記事中では、ある対局で高い一致率の示される確率について2つの数字が提示されています。

  • 10%: 「当方の計算では、出現率でいえば、良くて10%程度」
  • 25%: 「第三者委員会の報告は、(中略)上位4分の1に入るような出来」

はっきりと「ある対局で高い一致率の示される確率」と書かれている訳ではないですが、文脈からそう判断して問題ないと思われます。

私も以前の記事で計算をしたのですが、ある対局で一致率85%以上となる確率を9.79%と概算しています。 不安定な一致率の推定値から計算したものなので、これが正しいと考えてもいませんが、ここでは三浦九段に不利な10%を採用してみます。

確率計算

弁護記事では「そんな出来の良い将棋が、疑いをもたれた短期間に、4局も」あったと指摘しています。「短期間」が何を指すか明らかではありませんが、私の前回記事では某まとめサイト情報を参考に12局としていました。

12局中4局以上で高い一致率が示される確率

r   = 0.1
num = 12
lim = 3
s   = 1.0 - scipy.stats.binom.cdf(lim, num, r)

2.56%となります。 教科書にあるような有意水準5%だと、通常では起こりにくいことが起こっていて疑わしいとなります。

以前の記事では、12局というのが恣意的に切り取った一部なのでシロともクロとも言えないとコメントしていましたが、よくよく考えれば恣意的な切り取りについても数値的に評価できそうなので、その計算を今回追加します。

一致率の高くなる期間を切り取れる確率

三浦九段の2015年度対局数は42局です。 42局の中に一致率の高い4局を含む12局が出現するか、という問題に変わると、この確率は当然ながら2.56%より高くなります。

解析的に解く方法が思いつかないので、モンテカルロ法を使ってみます。エレガントな解法があれば教えてください。

確率10%のクジを42回引いて、結果のrolling_sumの最大値が3を超えるかで判定しています。

num_ext  = 42
num_iter = 10000
res = np.zeros(num_iter)
for ii in range(num_iter):
    samp = np.random.choice([0, 1], num_ext, replace=True, p=[1.0-r, r])
    seq  = pd.Series(samp).rolling(window=num).sum()
    res[ii] = seq.max()
s = (res > lim).mean()

モンテカルロ法なので結果にブレはありますが、17%程度という結果になりました。

トップ棋士は、ベースの一致率が高くなりがちというだけでなく、対局数つまり試行回数が増えてしまうため、一致率による疑いをかけられやすい構造にあるようです。

計算は以下にまとめています。 github.com

主張

サンプルを恣意的に切り取ることで、17%程度発生するありふれた事象を2%程度に見せることができます。 我々が統計の結果を作る・見るときには、このようなテクニックが使われていないか常に注意する必要があります。

いつもソフトを使ってる訳でないので、一致率の高い4局が問題で、12局だとか42局は関係ないという論を張っている人もいるようですが、これは誤りです。 無実の人が同じ試行をしたときに起こりそうもない結果が出たときのみ疑いをかけるのが正しい姿勢です。 羽生さんの対局から一致率の高い1局を引っ張ってきて、いつもはソフト指ししていないけれど、当該局ではソフト指ししていたという主張は馬鹿げているでしょう。

また、弁護記事では試験のカンニングについて言及されていますが、もし教師の方がいらっしゃたら記事のような手続きで生徒を疑うのはやめてください。個人的に疑いを持って物的証拠を探すのはセーフかもしれませんが、疑いを同僚教師に話したりするのはアウトです。 過去の偏差値からは取れそうもない点数を取る程度の事象であれば、「卒業までに受ける試験の数×カンニングが疑われるような素行の生徒数」の回数の試行があれば、簡単に発生してしまうのが確率というものです。

最期に

最近、話題のBIGについて計算したのですが、そちらでは直感通りのつまらない結果になりました。

MXNetを中心としたCustom Loss functionの話

これはDeepLearning Advent Calendar 2016の12月24日分の記事です。

Deep Learning一般の内容ではないので申し訳ないですが、来年はMXNet Advent Calendarやりたいですねという願いもあり、空き枠に滑り込みました。

TensorFlowとChainerが強くてMXNetは肩身が狭いですが、AWSという錦の御旗も得たのでユーザ増加に期待したいです。 MXNetはドキュメントが弱いとよく言われますが、最近はドキュメントサンプルコードも充実してきています。 9月頃はドキュメントがリンク切れしておりサードパーティーのドキュメントを読んでいたくらいなので随分良くなりましたよ。

さて、今回はCustom Loss functionについて紹介しようと思います。機能自体は前からあったようですが今月になって簡単な資料が公開され、一般ユーザにも存在が知られるようになったものです。

ロス関数について

ここでは、モデル出力値とラベルの差を表現する関数で、勾配法によって最小化する対象をロス関数と呼びます。 どのライブラリでも回帰なら2乗誤差、分類ならクロスエントロピーが使われているはずです。

注意が必要なのは、モデル精度の最終的な評価をするための指標とロス関数は別になるケースもあるということです。 分類問題ではクロスエントロピーよりもAccuracyやAUCで精度評価することが多いと思いますが、そういった場合でもロス関数にはクロスエントロピーが選ばれてクロスエントロピーを最小化するのが通常です。

評価指標をAUCに変えたら、それを最適化してくれると思っているユーザは結構多く、私もそういった方への案内誘導をボランティアでやっています。*1

一方で、ロス関数に拘ったからといって報われるとは限らず、絶対誤差で評価されるタスク等でも2乗誤差を最小化する方が上手くいくケースも多いのが実際と思われます。

ごく限定的な状況ですが、標準のロス関数では対応できないような問題のためにロス関数をカスタムするインターフェースが各ライブラリから提供されています。

MXNetでの事情

MXNetのLinearRegressionOutputSoftmaxOutputといった出力層向けのシンボルは、出力とロス関数が一体化しておりパラメータでロス関数を変更できません。 ロス関数を変えたい場合はMakelossというシンボルを利用する必要があります。

回帰の例

LinearRegressionOutputを使った通常の方法は以下のようになります。これはMSEを最小化します。

data  = mx.sym.Variable('data')
label = mx.sym.Variable('label')
fc    = mx.sym.FullyConnected(data, num_hidden=3)
act   = mx.sym.Activation(fc, act_type='relu')
fco   = mx.sym.FullyConnected(act, num_hidden=1)
net1  = mx.sym.LinearRegressionOutput(fco, label)

同じようにMSEの最小化をCustom Loss functionでやると以下のようになります。

loss = mx.sym.square(
    mx.sym.Reshape(fco, shape=(-1,)) - label
)
net2 = mx.sym.MakeLoss(loss)

mx.sym.squareやelement-wise minusの導関数はMXNetが知っているので自動微分を使って上手く最小化してくれる訳です。 numpy等を使って好き勝手にロス関数を作ることはできません。

注意 MakeLossにmx.sym.squareした結果を渡すとpredictした際にも2乗した結果が出力されるそうです。

custom loss symbol in R/Python · Issue #3368 · dmlc/mxnet · GitHub

確かにどれが出力として求められているかはMXNetで判断しにくいかもしれません。これについては使い易いように修正がされるそうです。

分類の例

MXNetのLogisticRegressionOutputは名前のような挙動をしないことが知られています。

RPubs - MXnet LogisticRegressionOutput test

SoftmaxOutputを使うことが推奨されていますが、これは手前のレイヤーで隠れユニット数をカテゴリ数と同じだけ確保する必要があります。 2値分類のときに隠れユニットを2つ持っておくのは冗長な気がするので、Custom Loss functionでロジスティック回帰を組んでみようというのがこの試みです。

SoftmaxOutputを使った通常の方法は以下のようになります。これはLog Lossを最小化します。

data  = mx.sym.Variable('data')
label = mx.sym.Variable('sm_label')
fc    = mx.sym.FullyConnected(data, num_hidden=3)
act   = mx.sym.Activation(fc, act_type='relu')
fco   = mx.sym.FullyConnected(act, num_hidden=2)
net1  = mx.sym.SoftmaxOutput(fco, name='sm')

ロジスティック回帰らしく書いた場合は以下のようになります。

data  = mx.sym.Variable('data')
label = mx.sym.Variable('label')
fc    = mx.sym.FullyConnected(data, num_hidden=3)
act   = mx.sym.Activation(fc, act_type='relu')
fco   = mx.sym.FullyConnected(act, num_hidden=1)
p     = mx.sym.Activation(fco, act_type='sigmoid')
eps   = 1e-6
p     = mx.sym.minimum(mx.sym.maximum(p, eps), 1.0-eps)
q     = 1.0 - p
lp    = mx.sym.Reshape(mx.sym.log(p), shape=(-1,))
lq    = mx.sym.Reshape(mx.sym.log(q), shape=(-1,))
loss  = label * lp + (1.0 - label) * lq
net2  = mx.sym.MakeLoss(loss)

出力を確率に直すやり方がわからないので、ちゃんと動いているかわからないのですが、 方針自体はこれで良いはず…

中途半端になってしまいましたが、これ以上掘るのは効率が悪いので、本日はここまでとします。 手元のJupyter Notebookはこちらです。

github.com

やっぱり

以前に比べればドキュメントが強化されているとはいえ、最近どういう機能が追加されたか等を個人で調べるのは大変です。 ユーザが増えて互助していければ素敵ですね。

*1:こんなことは働き盛りの人間がやる仕事ではないので、ボット化したいなぁ

LightGBMのPythonパッケージ触ってみた

DMLCひとりアドベントカレンダー0日目の記事です。 強い競合が現れたということで、DMLCとは直接関係ないですがLightGBMについて紹介します。

LightGBMとは

勾配ブースティング木の高速な実装としてXGBoostが有名ですが、Microsoftの開発した更に高速な実装がLightGBMです。 実験によるとXGBoostの数倍高速だということです。

パフォーマンス向上の大きな要因は、決定木の分割をexactに探索するのではなく、ヒストグラムのビン単位にして計算を間引いていることにあるようです。 ビン単位にすると精度が落ちそうな気もしますが、汎化性能を考えたら特徴量が粗いくらいでちょうど良いのだろうと腹落ちしています。

同じアプローチでFastBDTが先に発表されましたが開発が滞ってしまっている一方、LightGBMでは最近Pythonパッケージがベータリリースされました。 今日、我々がPythonで勾配ブースティングをする際にはXGBoostかLightGBMの2択*1となります。

導入

ビルドはこちら、インストールはこちらを参考に。

XGBoostとそれほど変わりません。

使い方

Exampleより

xgboost.DMatrix的なものを作って…

# create dataset for lightgbm
lgb_train = lgb.Dataset(X_train, y_train)
lgb_eval  = lgb.Dataset(X_test, y_test, reference=lgb_train)

パラメータをディクショナリに格納して…

# specify your configurations as a dict
params = {
    'task' : 'train',
    'boosting_type' : 'gbdt',
    'objective' : 'regression',
    'metric' : {'l2', 'auc'},
    'num_leaves' : 31,
    'learning_rate' : 0.05,
    'feature_fraction' : 0.9,
    'bagging_fraction' : 0.8,
    'bagging_freq': 5,
    'verbose' : 0
}

trainするだけ。XGBoostと同じ流れですね。

# train
gbm = lgb.train(params,
                lgb_train,
                num_boost_round=20,
                valid_sets=lgb_eval,
                early_stopping_rounds=5)

predictpandas.DataFrameを直接受け入れてくれるのには注意が必要です。

# predict
y_pred = gbm.predict(X_test, num_iteration=gbm.best_iteration)

XGBoostのパラメータとの対応は以下のようになります。

2016/12/18追記 tksさんに指摘を頂いたので記述を修正しました。 github.com

LightGBM XGBoost 候補 備考
objective objective regression, binary, multiclass, ...
boosting_type booster gbdt, dart 線形モデルはなし
max_depth max_depth num_leavesの制約を受けるので注意
learning_rate learning_rate
feature_fraction colsample_bytree
bagging_fraction subsample bagging_freqを1以上にしないと有効にならない
min_gain_to_split gamma
lamnda_l1 alpha
lambda_l2 lambda
metric eval_metrix l1, l2, auc, ...

LightGBMにもmax_depthパラメータが用意されていますが、オフィシャルにはnum_leavesで葉の数を制御します。

また、特徴量をまとめるビンの数をmax_binで設定できます。

詳細は公式のWikiを調べて下さい。

カテゴリ変数のサポート

LightGBMではカテゴリ変数がサポートされています。 One-hot encodingが不要になったので、メモリ的に有利であったり、カラム方向のサンプリングやfeature importanceが直感的になる等のメリットがありそうです。

下の実験からのコード片ですが、traincategorical_featureパラメータにカテゴリ変数として扱うカラムのインデックスをリストで渡してやるだけです。 カテゴリ変数はintのみ対応しているので、文字だったりする場合は0始まりのコードに置き換える必要があります。

obj_lgb = lgb.train(
    prm_lgb, dt_lgb_c, num_boost_round=num_round,
    valid_sets=dv_lgb_c,
    categorical_feature=list(range(len(x.columns.values))))

処理を完全には追っていないのですが、rpartのようなカテゴリ併合をしてくれる訳ではないようで、数値変数であればx <= thresholdで分割されるところが、x is categoryで分割されるようになるだけのようです。

また以下の実験の出力からなのですが、数値変数として扱うと以下のように小数点の閾値で分割されます。

threshold=1.5 3.5 2.5 0.5 0.5 0.5 0.5 0.5 4.5 14.5 0.5 1.5 0.5 0.5 1.5 2.5 0.5 0.5 0.5 0.5 0.5 4.5 0.5 0.5 0.5 0.5 1.5 0.5 2.5 2.5

一方、カテゴリ変数として扱うとデータにある値そのもので分割されます。

threshold=0 0 0 0 1 1 3 2 0 1 0 0 0 1 1 0 0 3 1 1 0 3 0 0 0 0 1 0 1 0

実験

手元でもパフォーマンスを確認してみようと、Ottoのデータで実験してみました。

github.com

モデル MLogLoss 計算時間
XGBoost 0.570366 65.20
LightGBM 0.522378 16.02
LightGBM Cat 0.671286 15.78

この設定ではLightGBMが4倍以上速い結果となりました。精度もLightGBMの方が良好です。

全変数をカテゴリ変数として扱ったLightGBM Catの有り難みがないように見えますが、One-hot encodingによってカラム数が膨らんでいる場合には計算時間の短縮が実感できるはずです。 精度が悪いのは、今回のデータが順序性を持っているためと思われます。

最後に

XGBoost4jのような使い方をすると話は別かもしれませんが、クライアントマシンでせこせこ計算する分にはLightGBMで良さそうです。 開発も活発なので、XGBoostにしかないような機能も続々と取り込まれている上に、コードもスッキリしており改造が比較的容易になっています。

CNTKなんかもパフォーマンス良いようですし、巻き返しに力の入っているMicrosoftさんのLightGBMには将来性も期待できるのではないでしょうか。

次は23日にMXNetのネタを考えています。

*1:scikit-learnのGradientBoostingもMAEで分割ができる等、面白そうではあるんですが…