k_oomoriであります!これまで数値計算といえばPythonのNumPy一択かなと思っていたのですが、Scalaでも似たようなライブラリはないものかと調べてみたところ、ありました、breezeが。ここではベクトルや行列の演算について性能を調べてみたいと思います。なお、使用したサーバはAWSのc3.large EC2インスタンス、OSはUbuntu 14.04、ソフトウェアバージョンはPython 2.7.6, NumPy 1.8.2, Scala 2.10.4, breeze 0.10となります。
ベクトルの内積
まずはベクトルの内積 から。NumPyでは例えばランダムな値を持つ3次元ベクトルは
import numpy as np v = np.random.random(size=(3)) v array([ 0.8716417 , 0.01981258, 0.25848169])
で生成できます。自分自身との内積は
np.dot(v,v) # NumPyの関数 v.dot(v) # ndarrayのメソッド (v * v).sum() # 成分ごとに積を取ってから和を取る (v ** 2).sum() # 成分ごとに2乗してから和を取る
などといった方法で計算できるので、いろいろな次元のベクトルに対して速度を計測してみます。なお、IPythonの%timeitを使うと計測が捗ります。
%timeit np.dot(v,v) 1000000 loops, best of 3: 672 ns per loop
これは1回あたり数百ナノ秒(1ナノ秒は10億分の1秒)で処理が完了するため、100万回試行してみたら上位3回の平均が672ナノ秒でした!という意味のようです。
対してbreezeでは
import breeze.linalg._ val v = DenseVector.rand[Double](3) v dot v //DenseVectorのメソッド v.t * v //縦ベクトルvを転置して得られた横ベクトルとvとの行列としての積 sum(v:*v) //成分ごとに積を取ってから和を取る
のようになります。では早速結果を。
ベクトルの次元 | 10 | 100 | 1000 | 10000 | 100000 | |
NumPy | np.dot(v,v) | 557 | 666 | 1,560 | 10,300 | 98,000 |
v.dot(v) | 918 | 1,020 | 1,920 | 10,700 | 98,700 | |
(v*v).sum() | 5,910 | 6,220 | 8,140 | 23,000 | 171,000 | |
(v**2).sum() | 6,450 | 6,910 | 8,790 | 23,800 | 172,000 | |
breeze | v dot v | 106 | 188 | 1,102 | 10,297 | 102,300 |
v.t * v | 116 | 198 | 1,104 | 10,210 | 101,400 | |
sum(v:*v) | 68 | 416 | 4,373 | 42,110 | 401,500 |
時間の単位はナノ秒です。これから読み取れるのは
- NumPyでは成分ごとの積を取るよりdotに任せた方が速い
- メソッド版よりnp.dotの方が若干速い
- 2乗を求めるにはv**2よりv*vの方が若干速い
- breezeではベクトルの次元が小さい時は成分ごとに掛けて和を取った方が速いが、数十次元以上になる場合はdotや*の方が速い
- v dot vとv.t*vはあまり差がなく、次元によって優位性が異なるらしい
- NumPy最速 対 breeze最速の比較では、次元が10000くらいまでではbreezeが速いが、それ以上になるとNumPyが逆転する
行列の積
行列の生成とその積を求める方法は、NumPyでは
m = np.mat(np.random.random(size=(5, 3))) m.T*m
breezeでは
val m = DenseMatrix.rand(5, 3) m.t*m
とします。(転置を取っているのは、そうしないと非正方行列の積が定義できないからです)
行列サイズ | 10x10 | 100x100 | 1000x1000 | 60000x500 |
NumPy | 5 | 893 | 752,000 | 19,500,000 |
breeze | 2 | 901 | 743,600 | 15,800,000 |
時間単位はマイクロ秒(s, 100万分の1秒)です。サイズが大きくなったときにbreezeが少し速いかなといった程度で、それほど大きな差は見られません。60000x500という行列のサイズは、たまたま手元にあったデータがそのくらいの大きさでこの検証を始めたことに由来します。
行列の要素ごとの積
2つの行列
に対し要素ごとの積
を考えます。こんな積を使うことがあるのか?と思うかもしれませんが、非負値行列因子分解で使うのです。2つの行列
を受け取り、要素ごとの積を取った結果を行列として返すためには、NumPyでは
np.mat(np.array(m) * np.array(n))
breezeでは
m:*n
とします。
行列サイズ | 10x10 | 100x100 | 1000x1000 | 60000x500 |
NumPy | 7.44 | 29 | 35,300 | 988,000 |
breeze | 0.43 | 26 | 2,800 | 114,000 |
単位はマイクロ秒です。np.arrayへの変換等がない分有利なのか、こちらはbreezeの方がかなり速いですね。
行列ノルム(の2乗)
最後に行列の各成分の2乗和であるフロベニウスノルムを考えます。
これを2つの行列の差に対して適用すれば、それらがどれだけ離れているかの基準となります。breezeではシンプルにsum(m:*m)として計算できますが、内積計算では成分ごとに計算するよりdotで計算した方が早かったことを思い出すと、行列を
次元の縦ベクトルが横に
個並んだものとして計算してみたくなります。
これを実現するScalaコードは
(for (i <- 0 until m.cols) yield {val v = m(::, i); v dot v}).sum
と書けます。逆に次元横ベクトルが縦に
個並んだものとみなすこともでき、その時のコードは
(for (i <- 0 until m.rows) yield {val v = m.t(::, i); v dot v}).sum
となります。
NumPyの成分版はma = np.array(m); (ma*ma).sum()、ベクトルのdot積版はScalaのように一息で書く方法が思いつかなかったので
def column_vectors(m): norm = 0 for i in xrange(m.shape[1]): v = np.array(m[:, i]) norm += np.dot(v.T,v)[0][0] return norm def row_vectors(m): norm = 0 for i in xrange(m.shape[0]): v = np.array(m[i, :]) norm += np.dot(v,v.T)[0][0] return norm
としました。結果は
行列サイズ | 3x2 | 100x10 | 1000x100 | 60000x500 | |
NumPy | ma = np.array(m); (ma*ma).sum() | 8 | 10 | 378 | 631,000 |
column_vectors(m) | 27 | 128 | 1,730 | 360,000 | |
row_vectors(m) | 40 | 1,240 | 12,700 | 838,000 | |
breeze | sum(m:*m) | 1 | 4 | 369 | 143,593 |
(for (i <- 0 until m.cols) yield {val v = m(::, i); v dot v}).sum |
1 | 3 | 126 | 36,900 | |
(for (i <- 0 until m.rows) yield {val v = m.t(::, i); v dot v}).sum |
1 | 16 | 284 | 177,360 |
単位はマイクロ秒です。これは圧倒的にbreezeの方が速いですね。また、あえて縦に長い行列を用いて実験したのですが、その場合は縦ベクトル、つまり長い方を一つのベクトルとみなしてベクトル演算(dot)した方が効率的であることがわかりました。
結論
私はコレで、Pythonやめました。(昭和を生きた人向けのネタです)
おまけ:インストール方法
breeze→ こちらを参照してください。
NumPy→ python-numpyをyumやapt-getでinstallしてください。RPMパッケージを作成したい方はこちら。