乱択アルゴリズムな打ち切り特異値分解のgemを公開した

※Truncated SVD (Singular Value Decomposition) って和訳ありますか？勝手に打ち切り特異値分解としました。

はじめに

機械学習アルゴリズムで特異値分解（もしくは固有値分解）をするときは、だいたい、大きい（もしくは小さい）方からk個だけ特異値/特異ベクトルを取り出すことが多い。NMatrixに実装されている特異値分解は、LAPACKのラッパーで、全ての特異値を求めるスタイル（Full SVD）なので、分解したい行列が大きい場合、実行時間が大きなものとなる。そこで、Randomized AlgorithmなTruncated SVDを実装した。

randsvd | RubyGems.org | your community gem host

NMatrixがARPACKとかに対応すると、存在意義を失うかも。

インストール

行列の扱いにはSciRubyのNMatrixを使用する。また、アルゴリズム中で特異値分解とQR分解を行うため、NMatirx-Lapackeに依存する。

$ gem install randsvd

使い方

内部の特異値分解で、NMatrix::LAPACK.gesvdを使うものと、NMatrix::LAPACK.gesddを使うもので二種類のメソッドを用意した。わかりやすさのため、メソッド名はNMatrixのものと同じにした。gesddの方が、分割統治法で速いらしいっすよ。

require 'randsvd'

# 適当に、分解する行列として、m x nの大きさの行列を用意する。
m = 1000
n = 100
input_matrix = NMatrix.rand([m, n])

# 欲しい特異値の数を定義する。
k = 10

# 特異値分解を実行する。分解したい行列と、欲しい特異値の数を指定する。
# 特異値が大きいほうからk個の特異値と特異ベクトルが得られる。
# ここが u, s, vt = RandSVD.gesdd(input_matrix, k) でも良い。
u, s, vt = RandSVD.gesvd(input_matrix, k)

# 上の例では、
#   u は m x kの大きさの左特異ベクトルが並ぶ行列、
#   s は k個の特異値が並ぶベクトル、
#   vt は k x nの大きさの右特異ベクトルが並ぶ行列となる。
# これら特異値と特異ベクトルで元の行列を復元する場合、以下の様になる。
reconstructed_matrix = u.dot(NMatrix.diag(s).dot(vt))

ベンチマーク

Core m3 1.1GHzで8GBメモリのMacBookで、ベンチマークとってみた。 10,000 x 1,000の大きさでランクが10な行列の特異値分解を行う。欲しいのは特異値の大きい方から10個となる。 NMatrix::LAPACKの特異値分解は、実際には、必要な特異値を取り出す処理（例えば u=u[0..m - 1, 0..k - 1] みたいなの）が必要だけど、速度に大した影響ないだろうと思って省略した。

require 'benchmark'
require 'randsvd'

m = 10000
n = 1000
k = 10

mat = NMatrix.rand([m,k]).dot(NMatrix.rand([k,n]))

Benchmark.bm 10 do |r|
  r.report 'NMatrix' do
    u, s, vt = mat.gesvd
  end

  r.report 'RandSVD' do
    u, s, vt = RandSVD.gesvd(mat, k)
  end
end

実行結果は以下のとおり。totalを見ると、RandSVDの方が約10倍くらい速い。

$ ruby bench_svd.rb
                 user     system      total        real
NMatrix     64.910000   1.300000  66.210000 ( 31.076849)
RandSVD      6.790000   0.860000   7.650000 (  5.485305)

アルゴリズム

乱択アルゴリズムな特異値分解は、色んなバリエーションが提案されている。 NIPSのワークショップの論文に、良い感じに要約されたアルゴリズムが掲載されていたので、主にこれを参考にして実装した。

P.-G. Martinsson et. al, “Normalized power iterations for the computation of SVD,” In Proc. of NIPS Workshop on Low-Rank Methods for Large-Scale Machine Learning, 2011.

乱択アルゴリズムな特異値分解では、ランダム行列由来の直交行列との積で、元の行列のサイズを縮めるということをやっている。上記の論文では、（近似精度の高い?）直交行列を得るのに normalized power iterations というのを行う。が、これをしなくても、そんなに精度的には変わらなかった（けど実行速度は落ちる）ので、RandSVDでは、デフォルトでは切っている。 normalized power iterationsを試したければ、gesvd/gesddメソッドで、繰り返し回数を指定するとよい（デフォルトでは0回になっている）。

# 例えば5回おこなう場合:
t = 5
RandSVD.gesvd(mat, k, t)

おわりに

つまらないものですが、よろしくお願い致します。

洋食の日記

「だ・である」調ではなく「です・ます」調で書きはじめれば良かったなと後悔してる人のブログです