Rest Term

Blog

JavaScriptで機械学習の実装 5 Gradient Boosting

少し間が空いてしまいましたけど、今回はGradient Boosting (勾配ブースティング)と呼ばれる機械学習アルゴリズムを試してみます。機械学習関連のコンペでも大人気の手法ですね。かなり昔(2011年)に決定木をJavaScriptで実装したことはあるので勾配ブースティングも併せて学んでおこうと思います。

Gradient Boosting (勾配ブースティング)

ブースティングについてWikipediaを引用すると、

ブースティング(英: Boosting)とは、教師あり学習を実行するための機械学習メタアルゴリズムの一種。ブースティングは、Michael Kearns の提示した「一連の弱い学習機をまとめることで強い学習機を生成できるか?」という疑問に基づいている[1]。弱い学習機は、真の分類と若干の相関のある分類器と定義される。- ブースティング - Wikipedia

ブースティングは逐次的(弱学習器を1つずつ順番に構築)に弱学習器を構築していく手法で、新しい弱学習器を構築する際にはそれまでに構築された全ての弱学習器の予測結果を利用するという特徴があります。Gradient Boosting (勾配ブースティング)では正解値と弱学習器の予測値との差を負の勾配と見なして、それを最小化するように逐次的に弱学習器を学習させていきます。弱学習器として決定木が用いられるため Gradient Boosting Decision Tree (GBDT) と呼ばれます。ちなみにOpenCVに顔検出機能はAdaBoostと呼ばれるブースティングアルゴリズムを利用していますね。AdaBoostは学生時代に特にお世話になりました。

今回は回帰問題を扱うので弱学習器には回帰木を用います。JavaScript(ES2015)での学習処理部分の実装を抜粋しますが、ここを見るだけでも勾配ブースティングの特徴がよく表れているのではないでしょうか。

検証

データセットはBoston housingを用い、人口 1 人当たりの犯罪発生数やNOxの濃度などの情報からボストン市の住宅価格を予測します。

今回もJavaScript(ES2015)で実装したのですが、いつものようにscikit-learn風なインタフェースにしていて、パラメータのデフォルト値はscikit-learnのGBDT実装(sklearn.ensemble.GradientBoostingRegressor)を参考にしました。もちろん外からパラメータ変更もできます。また、ファイルからのデータ読み込み処理は省略してデータセット自体も1つのモジュールとして埋め込んでいます。ソースコードもいつものようにGitHubに置いてあるので興味があれば。

クライアント側では以下のように使います。

* 出力結果例

モデル評価用にRMSEやR^2スコア(決定係数)計算用のモジュールも併せて作りました。まずは残差プロットでざっくりとモデルの性能を眺めてみます。最近はGoogleスプレッドシートでグラフを描くことが増えてきました。さすがにExcelよりは機能不足ですけど簡単なグラフならすぐ作成できるので便利に使っています。

良い感じに残差が均一に分散しています。それなりに上手く学習できているようです。

定量的な確認も少しだけ。データをシャッフルして3回計測したRMSEとR^2スコア(決定係数)の平均値を載せます。GBDTの学習率は0.1で固定して、回帰木の本数(ブースティングの繰り返し回数)と深さを変更しました。

RMSE R^2 score (決定係数)
回帰木の本数: 100, 深さ: 3 3.503 0.817
回帰木の本数: 200, 深さ: 4 3.156 0.851

ある程度想定した通りの結果が出ているので問題なさそうです。パラメータ調整についてはまだ理論を理解しきれていない部分もあるので適当なのですが、いろいろパラメータを変えて試してみた所、回帰木は深くしすぎてもダメみたいです。3 ~ 5くらいがちょうど良いんでしょうか。

おわりに

勾配ブースティングは以前からずっと作ってみたいと思っていたので、今回実装する機会を得ることができて良かったです。ところでロシアのYandexが7月18日にCatBoostというすばらしい名前の勾配ブースティング実装を公開しました。今はこちらを使っていろいろ勉強中です。まだ日本語のドキュメント類は無いみたいなので情報がまとまったらCatBoostに関するエントリーを書くかもしれません。またよろしくおねがいします。

参考

 

Tags: , ,

GPU対応の類似検索(最近傍探索)ライブラリ Faissの紹介 part2 GPUを利用した検索

Facebook AI Research (FAIR)が開発したGPU対応の類似検索ライブラリ Faiss に関する2回目の紹介エントリーとなります。前回は、FaissのインストールとC++チュートリアルの説明を行いました。今回はGPUを利用した検索処理を試してみます。

Faissのインストール方法については前回のエントリーを参照してください。

環境

Amazon Linux AMI release 2017.03
Intel Xeon CPU E5-2686 v4 @ 2.30GHz
NVIDIA GK210GL [Tesla K80]
Ubuntu 16.04 LTS (x86_64)
Intel Core i5-2500K CPU @ 3.30GHz
NVIDIA GM204 [GeForce GTX 970]
CUDA Toolkit 8.0
GCC 4.8.3 / 5.4.0
Python 3.5.2
NumPy 1.12

今回からはプログラムのベンチマークも行っていくので自宅のGPUサーバをメインに使って検証します。EC2のp2.xlargeをがっつり使うと請求が怖いし、そもそもレイテンシが大きくてストレスが溜まるので。

前回の導入部分の補足としてCUDAコード部分のコンパイルにはGCCではなくNVCC(CUDAコンパイラ)を利用するので事前にインストール済みであることを確認しておきます。また、GPU製品のCompute Capabilityは3.5以降が必要となっています。

GPUでk-meansクラスタリング

まずは、公式ドキュメントにサンプルとして掲載されているGPUを利用したk-meansクラスタリングプログラムの動作確認をします。日本語での補足も簡単に入れておきます。

faiss::gpuに属するクラスについては後述しますが、faiss::gpu::GpuIndexFlatL2 は前回のエントリーで紹介した faiss::IndexFlatL2 のGPU版となっています。faiss::Clustering::train 実装内部ではイテレーション毎に現在のセントロイドをインデックスしておき、全データをクエリとして検索しています。

* コンパイル

* 実行結果例

GPUを利用したk-meansクラスタリングは過去にもいろいろライブラリ等が出ていますし、サンプルとしてはよく見るものです。ここでの動作確認では5,000,000個のデータを20,000クラスタに分けるのに4分弱で処理できました。いまいちよく分からないのでCPU版と比較してみます。

上記プログラムをCPUでのクラスタリング処理に変更するには、インデックスオブジェクトである faiss::gpu::GpuIndexFlatL2faiss::IndexFlatL2 に変えるだけです。ただし上記プログラムの設定だとCPU版の処理が一向に終わらなかったので以下のようにテストの規模を縮小することで調整しました。

  • データ(128次元): 500,000
  • クラス(分割数): 2,000
  • イテレーション: 20

500,000個のデータを2,000クラスタに分けます。それぞれ5回計測した平均処理時間は以下の通りです。また、参考までに各製品のリリース時期も併記しています。

CPU
Intel Core i5-2500K
(2011/1)
CPU
Intel Xeon E5-2686
(2016/6)
GPU
NVIDIA GM204 [GeForce GTX 970]
(2014/9)
GPU
NVIDIA GK210GL [Tesla K80]
(2014/11)
511.15s 114.10s 5.70s 4.85s

前回のエントリーにも書いた通り、CPU版はOpenMPやSIMDを駆使して実装されており十分に高速です。しかしGPU版はそれと比較しても遙かに高速に動作します。また、CPUからGPUへのデータ転送は必要に応じて都度行われるため(計算に必要なデータだけGPU上に載せる)、VRAM容量を超えるような巨大なデータセットでもクラスタリング可能となっています。

さらに、FP16(半精度浮動小数点値での演算)を有効にしてGTX 970上で比較してみました。こちらも同様に5回計測した平均処理時間を載せます。

  • データ(128次元): 5,000,000
  • クラス(分割数): 20,000
  • イテレーション: 20
FP16無効 FP16有効
228.47s 191.76s

なかなか興味深い結果となりました。確かPascal世代からFP16の性能が格段に上がっていたかと思うので(1SPでFP16二つを同時演算したり)、PascalのGPUでベンチマークしてみたいところです。Facebook ResearchのレポートだとTesla P100でだいたい100秒ほどで完走するらしいです。
* GPU k means example

GPUを利用した検索

次は前回紹介したCPU版チュートリアルにある各種検索プログラムをGPU対応していきます。

GPUリソースの確保・初期化

GPUで処理するために予めGPUリソースの事前確保作業が必要ですが、生のCUDAプログラムと比較すれば十分抽象化されているので特に難しいところは無いかと思います。GPU版のヘッダファイルは faiss/gpu 以下にあります。

  • faiss::gpu::StandardGpuResources
  • faiss::gpu::GpuIndexFlatConfig

前述のk-meansのサンプルコードと同じ設定でだいたい大丈夫です。内部ではcuBLASのハンドル初期化やCUDAのストリーム生成処理等を行っています。デフォルトだとストリーム(cudaStream_t)は2本準備しているようです。

GPU版インデックスクラス

チュートリアルで紹介した各種CPU版インデックス用クラスのGPU対応は簡単なのでまとめて紹介します。

CPU GPU
faiss::IndexFlatL2 faiss::gpu::GpuIndexFlatL2
faiss::IndexIVFFlat faiss::gpu::GpuIndexIVFFlat
faiss::IndexIVFPQ faiss::gpu::GpuIndexIVFPQ

2017/05現在、faiss::gpu::GpuIndexIVFPQ にはGPU版インデックスクラスのポインタを直接受け取るコンストラクタは未実装のようですので別のコンストラクタを利用しています。この場合は引数で渡す faiss::MetricType に応じたGPU版インデックスクラスのインスタンスを内部で生成/保持してくれます。前回紹介していませんでしたが距離尺度は faiss::MetricType 定数で指定するようになっており、チュートリアルではL2ノルムの紹介のみでしたが、もちろん内積(Inner Product)にも対応しています。

また、faiss::gpu::IndicesOptions は以下のように定義されています。これはCPU<->GPU間でのインデックス情報の転送・保存をどのような形式で行うかを指定します。

CPU版とは異なりGPU版の各種インデックスクラスのメンバ変数はなぜかカプセル化されているみたいなので、提供されているgetter/setter経由でデータ操作を行います。GPU対応する際には忘れないようにしましょう。

ベンチマーク

次は最近傍探索処理の測定を行います。チュートリアルのプログラムをベースにパラメータ等を調整して比較しました。論文中ではDeep1BYFCC100Mなどの巨大なデータセットでベンチマークをしていて真似したかったんですけど、自宅のGPU開発環境だと空冷が適当なので kidle_inject が発生し、温度上昇防止のためCPU使用に制限をかけてしまいます。。ここではANN_SIFT1Mと同じ小規模データでベンチマークを行いました。

  • 次元数: 128
  • データ数: 1,000,000 (100万データ)
  • 同時実行クエリ数: 10,000 (1万クエリ)
  • インデックス: faiss::gpu::IndexFlatL2 (全探索) / faiss::gpu::IndexIVFPQ (PQ + 転置インデックス)
  • 転置インデックスのリストサイズ 100、サブベクトル数 8、量子化ビット数 8、FP16有効
  • 検索時のkの値: 1 ~ 1024 の範囲で変更

ベンチマーク用プログラムはここに記載するには少し長くなるのでGistにアップロードしておきます。

僕の個人開発環境だと1万クエリの上位20件を検索するのに、全探索(IndexFlatL2)の場合は約1.1秒、PQ + 転置インデックス(IVFPQ)の場合は約20ミリ秒でした。ただし、Faissの真価はbillion-scaleでの検索処理において発揮されるので、ここでのベンチマークは参考程度ということで。

また、Facebook AI Research公式の詳細なベンチマーク結果は以下を参照してください。

注意点など

GPUでの検索処理にはいくつか制限が設けられています。2017/06現在の仕様なので変更される可能性があります。

GPUテンポラリメモリ

cudaMalloc/cudaFree はコストが高い処理なので事前にある程度のメモリ領域が自動的に確保されるのですが、そのサイズをユーザー側が指定することができます。StandardGpuResources のAPIで指定できます。

テンポラリメモリが足りないと以下のような警告が出力されます。cudaMalloc が新たに必要になってしまうから遅くなるよというヒントであって処理が継続できないわけではありません。

PQ(Product Quantization)のために内部ルックアップテーブルを準備する際にはそこそこ十分な容量を確保しておいた方が良いかと思います。

kおよびnprobe値の上限

kとnprobe値は1024が上限となっていますが、1024もあれば実用上は十分かなと思います。

Product Quantizationのパラメータ

PQの実装上、量子化ビット数(サブベクトル当たりのデータサイズ)は8ビット以下、データの次元数をサブベクトル数で割った余りが0になるように調整する必要があります。今回紹介したベンチマークプログラムではデータは128次元、サブベクトル数は8、量子化ビット数は8で設定しています。

おわりに

今回はFaiss GPU版のインタフェースを紹介しました。CPU版から移行修正しやすいように工夫されているので特に難しいところはなかったかなと思います。GPU版の実装は論文と併せて読み進めていますが、k-selection(選択アルゴリズム)のGPU実装の仕組みが難しくて難航しています。。理論面を調べつつ、そろそろPythonインタフェースを用いたアプリケーション寄りの開発もやってみたいと思っています。

 

Tags: , ,

GPU対応の類似検索(最近傍探索)ライブラリ Faissの紹介 part1 導入/チュートリアル

Facebook AI Research (FAIR)が開発したGPU対応の類似検索ライブラリ Faiss を紹介します。

[06/25追記]
Faiss GPU版の検索についてエントリーを書きました。

論文は以下で公開されています。

論文タイトルの通り、10億スケールの大規模データに対してGPUを駆使した効率的な最近傍探索アルゴリズムの研究事例となっています。分量が多くなりそうなので複数の記事に分けて紹介しようと思います。

環境

Amazon Linux AMI release 2017.03
Intel Xeon CPU E5-2686 v4 @ 2.30GHz
NVIDIA GK210GL [Tesla K80]
Ubuntu 16.04 LTS (x86_64)
Intel Core i5-2500K CPU @ 3.30GHz
NVIDIA GM204 [GeForce GTX 970]
CUDA Toolkit 8.0
GCC 4.8.3 / 5.4.0
Python 3.5.2
NumPy 1.12

最初はTesla系GPUを使いたかったのでEC2のp2.xlargeインスタンス上で動作確認をしてましたけど、やっぱりオレゴンは遠くてターミナルがもっさりするので途中から自宅サーバに切り替えて作業しました。快適です。ちょっと前にGTX 1070を買いましたので、お下がりのGTX 970をLinuxサーバに挿しています。GTX 760は使い途がないのでメルカリかヤフオク行きですかね。

導入

※ GPU環境がある場合はNVIDIA DriverとCUDA Toolkit 8.0以降の事前導入が必要です

ここでは通常のC++ライブラリ(GPU非対応)、GPU対応のC++ライブラリ、Python用のライブラリをそれぞれソースコードからコンパイル/インストールします。プラットフォーム毎にMakefileから読み込む設定(makefile.inc.xx)が用意されているのでコピー・編集して使います。2017/05現在、CMake対応も進められているようですけれどまだ作りかけなので普通にmakeでコンパイルした方がいいかと思います。

CPU版C++/Pythonライブラリはビルド時にBLASライブラリをリンクしますが、ここではBLAS実装としてOpenBLASを利用しました。Intel MKLやAtlasでもOK、macOSならAccelerateフレームワークにBLAS実装が含まれているのでそれも利用できます。対応ライブラリの中ではIntel MKLが最速ですので、もしプロダクション環境等で可能な限り高速化したい場合はMKLをリンクさせましょう。GPU版はcuBLASが利用されます。こちらはCUDA Toolkitに含まれているので別途インストールする必要はありません。最近のCUDAだとBLAS drop-inライブラリとしてNVBLASが付属しているのでそちらをリンクするプロダクトも増えてきたように思いますが、FaissではcuBLASでこつこつと実装されているようです。

GTX 970のCompute Capability(CC)は5.2なのでそれ以外のコード生成オプションは削っておきます。デフォルトのままでも問題ないですけど、無駄を省いてコンパイル時間を短くしたいならハードウェア環境に合わせて設定すると良いです。

PythonインタフェースはSWIG(3.x以上)でビルドされます。依存ライブラリとしてNumPyも必要なので事前にインストールしておきます。

Makefileの準備と依存ライブラリの導入が終わったら後はコンパイルするだけですが、最低限 SSE4.2/POPCOUNT/OpenMP に対応しているCPU/C++コンパイラが必要となっています。AVXは使ってないようなのでSandy Bridge以前のかなり古いモデルでも大丈夫なはずです。

GPU版 C++ライブラリはコンパイルにかなり時間がかかると思います。最近だとこういう手順でプロダクトをインストールするのは逆に新鮮に感じます。リリース初期の頃のCaffeみたいだ。

ちなみに付属のMakefileはinstallターゲットがないので、手動でライブラリとヘッダファイル群を任意の場所に配置します。ここでは既にパスの通ったディレクトリに置きました。

動作確認

CPU版 C++ライブラリのビルド時にテストプログラムもコンパイルされて実行可能バイナリが生成されているのでそのまま実行できます。GPU版も用意されているのですが、これは2017/05現在、バグがあって特定のGPUだとSEGVしてしまうので完走できません。どうやらGeForce系GPU(GTX 970や1080)でのクラッシュ報告が多いようです。EC2 p2-xlargeインスタンスのTesla-K80なら大丈夫でした。

テストプログラム内部で行われていることはチュートリアルを進めていけばわかるようになるのでここでは省略します。公式のチュートリアルよりもコードは綺麗な気がしますけど。。

チュートリアル

公式のチュートリアルに補足を加えながら順番に動作確認していきます。

Faissが提供するクラスや各機能が扱うデータの型を意識しながら学んだ方がより理解しやすいと思うので、最初はC++インタフェースで試します。Pythonインタフェースは前述の通りSWIGで生成されているため、C++インタフェースが使えるようになればPythonインタフェースを使うための学習コストは少なくなります。

データの準備

Faissが扱うベクトルデータは基本的に単精度(32bits)浮動小数点値を要素とした配列を用います。C++ならfloat型配列 、Python(NumPy)だと numpy.ndarray オブジェクトとなり、1つのベクトルを1行としたrow-majorな行列データとして扱うことになります。

64次元ベクトルデータをデータベース登録用に10万個、検索クエリ用に1万個それぞれ作成しています。各ベクトルの1次元目の値だけ修正してますが、検索時にその意図がわかるのでここでは気にしなくても良いです。このチュートリアルでは通常のfloat型配列で入力データを管理していますが、Faissの実装を読んでみたところ、インデックス作成時に std::vector のコンテナにデータを詰め替えています。もしOpenCVと併用するなら cv::Mat に変換するアダプタを作っておきたいですね。

インデックスの作成

データベースのインデックスを作成します。ここからはFaissのAPIを使っていきますが、インデックスの種類はたくさんあるのでここでは基本的なものをピックアップして紹介していこうと思います。

faiss::IndexFlatL2

まずは基本となる距離尺度にL2ノルム(L2距離)を利用したインデックスの作り方を確認します。

IndexFlatL2 はL2ノルムでの全探索用(Brute-force search)インデックスとなっており、入力ベクトルデータは未加工のまま全て登録します。publicなメンバとして is_trainedntotal を参照していますが、それぞれインデックスの学習済フラグ(データ分布等を事前解析、IndexFlatL2では利用されない)、登録データ数を示しています。

検索

データの登録が終わったら検索してみます。前手順で作成した検索クエリ用データを利用してk-近傍探索を行います。

* 出力結果例

データベースに登録したデータ、検索クエリ用のデータも全てランダムに生成したにもかかわらず、類似検索結果のインデックスに偏りがあることがわかります。ここで、データ生成時に各ベクトルの1次元目のデータだけ少し加工したことを思い出しましょう。データ配列のインデックス値に依存する数値が加算されており、検索結果にもその影響が表れていることが確認できます。

もし、データを加工しない場合は以下のような検索結果になります。

2017/05現在のCPU版の実装では、次元数が4の倍数かつ同時検索クエリ数が20件以下の場合はSSEを使って距離計算をしてくれるようです。つまり、上記コード中の最初のsanity checkではSSEを利用、検索クエリデータ xq での検索ではインストール時に採用したBLAS実装(ここではOpenBLAS)に計算処理が委譲されます。また、マルチスレッドで処理される箇所はOpenMPを使って実装されているので大変読みやすくて助かります。

faiss::IndexIVFFlat

IndexFlatL2 をそのまま利用すると、いくらSSEやOpemMPを駆使して計算をしたとしても検索速度が遅くなってしまうため、Billion-scaleの大規模データに対応することは困難です。ここでは、入力ベクトルデータをクラスタリングしておき、検索時にはクエリデータに対応するセルおよび近隣のセルをいくつか選び出して、それらのセルに属するデータのみを最近傍点の候補として距離計算を行うことで検索時間を削減します。

* 出力結果例

IndexIVFFlat は各セルの代表ベクトル(母点)を保存しておくための転置ファイル(転置インデックス)用のクラスとなっており、処理手順としてはベクトル量子化処理の学習ステップが一つ増えています(IndexIVFFlat::train)。ここではデータベース登録対象データをそのまま学習に利用していますが、別途用意したデータで学習しても良いです。publicなメンバである nprobe で検索対象のセル数を管理しており、この値を増やすと検索精度の向上と引き替えに検索時間が増えます。上記チュートリアルでは全探索と比較してnprobe=1で15 ~ 20倍、nprobe=10で3 ~ 5倍くらい検索速度が向上しました。

faiss::IndexIVFPQ

IndexFlatL2 および IndexIVFFlat を併用した場合でも、入力ベクトルデータは全て未加工でインデックス登録していました。これでは数億オーダーのデータを全てメモリ上に載せるのは現実問題として困難です(NVMe等ストレージハードウェアも急激に進歩しているのですがここでは考えないこととします)。次は、入力ベクトルデータをProduct Quantization(PQ: 直積量子化)により圧縮することでインデックスするデータ量を削減します。PQはハッシュベースの手法に替わってここ数年の論文等でよく見かけますね。Product Quantizationに関する日本語の概要説明については以下の映像奮闘記さんのエントリーが読みやすくて参考になるかと思います。

  1. クエリベクトルをM分割し、N/M次元のサブベクトルに分解する
  2. クエリのサブベクトルと、それぞれのコードブックに記された代表ベクトル間の距離を予め計算しておく
  3. それぞれのエントリがどのインデックスに位置しているのか確認し、2で計算した距離を参照して足し算を行う
  4. 最小距離をもつエントリを候補として提出する

* 出力結果例

IndexIVFPQ はPQを伴う転置ファイル用クラスです。上記チュートリアルではサブベクトル数を8、各サブベクトルは8bitsにベクトル量子化されます。つまり、元データは64(次元数)*32(float)=2048bits、PQすると8*8=64bitsとなり32分の1のデータサイズとなります。また、各サブベクトル毎の距離計算は独立しているため並列処理とも相性が良さそうです。IndexFlatL2 による全探索(線形探索)結果を比較すると量子化の影響は見られますが、高い精度を維持できていることが確認できます。

注意点として、C++版チュートリアルのコードはGPUではなくCPUで動作します。GPU版実装の詳細は次回パート以降に紹介する予定です。

おわりに

今回はチュートリアルレベルの内容でしたので特に目新しい部分は無かったかと思います。Faissの論文は一通り読みましたが、主にGPU上で効率的に計算するアルゴリズムの記述が主でしたので、次パート以降はその理論や実装面についても深掘りしていきたいと思います。

SSEやOpenMP、CUDAを活用した処理は学生時代にけっこう勉強したんですけど、すっかり知識が抜けていたので焦りました。普段のWebアプリケーション開発の現場だと触れる機会がないですから、この分野の技術はプライベートの時間を使ってじっくり学んでいきたいです。

参考

 

Tags: , ,

NumPy入門記事をアップデートしました

4年以上前の記事ですけど、今もたくさんのアクセスがあるようなので内容を少し更新しました。さすがに古すぎる記事をずっと見続けてもらうのは気が引けますから。

更新内容

  • CentOS 7.0, Python 3.6, NumPy 1.12で動作確認
  • Record Arrays, Automatic Reshaping, 乱数生成, Linear Algebra(線形代数)について新規説明を追加

NumPyは多くの機械学習系フレームワークや画像処理ライブラリ内部で利用されているのでインタフェースの互換性には十分気を遣っているようです。今後も必要に応じて内容更新していければと思います。

 

Tags: