pythonでGPUとMCMC(とR) - xiangze's sparse blog

GPUでモンテカルロ法の計算をしたくなったりした場合には普通CUDA,OpenCLを使うことになります。
C++でプログラミングする必要があるのですが、変数の確保、解放などで記述が長くなりがちです。pythonを用いると記述を簡潔にできるところが多いらしいので関連するライブラリを紹介します。

実務的にはPyMCからTheanoが呼べるらしいのでそれを使用するのが良さそうですが、現段階では実装されていない機能もあります。
@beroberoさん
http://heartruptcy.blog.fc2.com/blog-entry-152.html
id:breakbeeさん
http://breakbee.hatenablog.jp/entry/2014/05/20/121722
がインストール、実行に関して解説されているのでここでは実装面を少し説明します。というか備忘録です。

PyCUDA,PyOpenCL

http://documen.tician.de/pyopencl/
http://documen.tician.de/pycuda/
メモリの確保CPU,GPU間のデータの受け渡し、OpenCLのqueueの確保、開放などをpythonで書くことが出来ます。同じ組織が開発しているようです。以下のライブラリの基盤となっていることが多いです。
CUDAとOpenCLではそのGPU側の記述は似ていますが、CPU側の記述はかなり異なります。CUDAは諸々の初期化がCUT_DEVICE_INIT(pyCUDAではそれすら必要なくimport pycuda.autoinit)でほぼ一発なのに対し、OpenCLではcontext, queueの設定が必要です*1。CPU-GPU間の変数の受け渡しは明示する必要があります。
PyCUDA,PyOpenCLともにGPUに行わせる処理(CUDAの__device__,OpenCLの__kernel)はC++で書いたコードを文字列として与える必要があります。

Reikna,CLUDA

http://reikna.publicfields.net/index.html
機械学習のアルゴリズムの実装基盤だそうですが、あまり例は多くないようです(https://github.com/Manticore/reikna/blob/develop/examples/demo_struct_reduce.py)。
Reiknaが基盤としているCLUDAではMakoというテンプレートエンジンをつかってCUDA,OpenCLへの統一的アクセスをしています。reiknaとはアイスランド語で計算という意味らしいです。
少し変わった記述の例

gpustats

pypi
github
PyCUDAをベースにしていてpyMCに確率密度関数のサンプリング、乱数生成の機能を提供しているような形になっています。
使用する場合にはC++の記述を意識する必要はなくなっていますが、現在使用可能な分布関数１変数、多変数の正規分布のみのようです。2011年くらいで開発は停止しているようです。

Theano

Deeplearning(stacked denoising autoencoder)などが実装できるライブラリとして有名です。
Deeplearning.netのサイトも有名です。ニュース、求職情報もあります。
メトロポリス法の実装、そしてなぜかHamilton Monte Carloの実装も公開されています*2。deeplearning関係のpython,matlab実装は大体ここにまとまっているようです。

Theanoでの非常に大まかなプログラミング手順は

shared object(変数)の定義
式の定義 (TはtensorのT)
functionで式のコンパイル
コンパイルされた式の実行

となります(logistic regressionによる手書き文字テストデータMNISTの識別の例)。式の定義とコンパイル、実行が同じpythonで書かれているので分かりやすいと思うか、分かりづらいと思うかは人によるかもしれません。*3
他のpythonライブラリとの変換対応に積極的らしく、数式処理部分ではsympy, GPUでの計算ではPyCUDAとの変換機能が実装されています。
http://deeplearning.net/software/theano/tutorial/gpu_data_convert.html
pure pythonでのコード実行、gccでのコンパイル、nvccでのコンパイルを選択することができて、GPUがない環境でもとりあえず実装を試行することも出来ます。どのアーキテクチャで計算を行うかはpythonの起動時に指定、もしくは.theanorcを記述することで選択できます。(Windows,Linuxともに同じ)
http://deeplearning.net/software/theano/library/config.html
内部の関数のドキュメントもある程度整備されていますが、特にc++,CUDAのコード生成部分のドキュメントが公開されておらず、作成者絶賛募集中と思われます。
http://deeplearning.net/software/theano/library/index.html#libdoc
github *4

Theanoはピタゴラスの妻だそうです。またpylearn2はTheanoに依存しているらしいです。
http://yamitzky.hatenablog.com/entry/2014/04/27/160138

Gnumpy, cudamat

cudamatはcuda関数のラッパーでgnumpyはそのnumpy風のインターフェースを提供しています。
http://www.cs.toronto.edu/~tijmen/gnumpy.html
https://github.com/cudamat/cudamat
http://www.tsc.uc3m.es/~miguel/MLG/adjuntos/slidesCUDA.pdf
cudamatはCUDAのラッパーとしてctypesを使っています。
cudamat、gnumpyともにRestricted Boltzmann Machineの実装を公開しています(確率的勾配降下法)。gnumpyの方が値を返す形になっており、読みやすいです。*5

pygml

https://github.com/slinderman/pyglm
http://nbviewer.ipython.org/github/slinderman/pyglm/blob/master/parallel_test.ipynb
Theanoに依存した一般化線形モデル(GLM)を計算するライブラリでGPUでMAP推定、ベイズ推定もできるそうです。

Numba(pro)

NumbaはpythonをJIT コンパイルすることで高速実行する仕組みだそうです。
ProではCUDA用のAPIが提供されるそうです。
http://docs.continuum.io/numbapro/

Nvidiaによる説明
https://developer.nvidia.com/how-to-cuda-python

R+GPU(gputools)

最後にRの話です。モンテカルロ法というよりは機械学習のGPU実装です。
gputools
http://brainarray.mbni.med.umich.edu/brainarray/rgpgpu/
http://cran.r-project.org/web/packages/gputools/index.html
WindowsバイナリはCRANにはないのでビルドする必要があるらしく、非常にめんどくさそうです。
http://www.slideshare.net/wdkz/rgpu-10228793
こちらではAWSのGPU インスタンスを使ったRpudというライブラリの使用を詳しく解説されています。page 16,17で説明されているコマンドの実行でR, gputoolsのインストールが可能だそうです。Rpudは実装されているアルゴリズムこそ少ないものの引数が基本的なRの関数やe1071パッケージと共通なところが利便性、ベンチマークに優れているようです。
http://www.r-tutor.com/gpu-computing
https://www.youtube.com/watch?v=dov5IUfJlkc