洋食の日記

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

非線形最適化手法L-BFGS-BのGemを作成した

はじめに

非線形最適化手法の一つであるL-BFGS-Bは、scipy.optimize.minimizeのデフォルト手法として採用されているなど、最適化手法ではde facto standardの様な存在である。Rではoptim関数がサポートしていたり、Juliaではbindingライブラリがあったりするが、探した限りRubyにはなかったので作成した。

lbfgsb | RubyGems.org | your community gem host

使い方

インストールは、native extensionをビルドできれば、gemコマンドでインストールできる。別で外部ライブラリをインストールする必要はない。Lbfgsb.rbでは、ベクトルをNumo::NArrayで表現するので、一緒にインストールされる。

gem install lbfgsb

使い方は簡単というか、minimizeというmodule functionが一つあるだけである。

require 'lbfgsb'

result = Lbfgsb.minimize(
  fnc: 目的関数, jcb: 目的関数を微分したもの, 
  x_init: 初期ベクトル, args: 目的関数などへの引数
)

p ressult[:x] # 最適化したベクトル

例として、二値分類なロジスティック回帰を実装してみる。ちなみに、現在のScikit-learnのLogisticRegressionも、デフォルトではL-BFGS-Bで最適化を行う。

class LogisticRegression
  def fit(x, y)
    # 切片を得るために特徴ベクトルに1を継ぎ足す.
    xb = Numo::DFloat.hstack([
      x, Numo::DFloat.ones(x.shape[0]).expand_dims(1)
    ])
    # Lbfgsb.rbで最適化する.
    result = Lbfgsb.minimize(
      fnc: method(:obj_fnc), jcb: method(:d_obj_fnc), 
      x_init: xb.mean(0), args: [xb, y]
    )
    # 重みベクトルを得る.
    @weight_vector = result[:x]
  end

  def predict_proba(x)
    # どちらのクラスに属するかの確率を推定する.
    xb = Numo::DFloat.hstack([
      x, Numo::DFloat.ones(x.shape[0]).expand_dims(1)
    ])
    1.0 / (Numo::NMath.exp(-xb.dot(@weight_vector)) + 1.0)
  end

  def predict(x)
    # クラス確率をもとにラベルを付与する.
    probs = predict_proba(x)
    predicted = Numo::Int32.zeros(x.shape[0])
    predicted[probs >= 0.5] = 1
    predicted[probs <  0.5] =-1
    predicted
  end

  private

  # 目的関数
  def obj_fnc(w, x, y)
    Numo::NMath.log(1 + Numo::NMath.exp(-y * x.dot(w))).sum
  end

  # 目的関数を微分したもの. 勾配ベクトルを得る.
  def d_obj_fnc(w, x, y)
    (y / (1 + Numo::NMath.exp(-y * x.dot(w))) - y).dot(x)
  end
end

このロジスティック回帰に、適当なデータを与えて分類を行ってみる。今回は、LIBSVM Dataからsvmguide3というデータを取得した。

$ wget https://www.csie.ntu.edu.tw/~cjlin/libsvmtools/datasets/binary/svmguide3
$ wget https://www.csie.ntu.edu.tw/~cjlin/libsvmtools/datasets/binary/svmguide3.t

LIBSVM形式の読み込みと、正確度の計算にはRumaleを使う。

$ gem install rumale
require 'rumale'

# 訓練・テストデータを読み込む.
x, y = Rumale::Dataset.load_libsvm_file('svmguide3')
x_test, y_test = Rumale::Dataset.load_libsvm_file('svmguide3.t')

# ロジスティック回帰による分類器を学習する.
logit = LogisticRegression.new
logit.fit(x, y)

# ラベルを推定する.
predicted = logit.predict(x_test)

# 正確度を計算し, 結果を出力する.
evaluator = Rumale::EvaluationMeasure::Accuracy.new
puts 'ground truth:'
pp y_test
puts 'predicted:'
pp predicted
puts("\nAccuracy: %.1f%%" % (100.0 * evaluator.score(predicted, y_test)))

実行すると以下のようになる。 同様のことを、scikit-learnのLogisticRegressionでも行うと、正確度が61.0%となった。 L-BFGS-Bにより最適化を行うロジスティック回帰を実装できた。 ※scikit-learnの場合、L2正則化項がつくので、ハイパーパラメータを C=1000, max_iter=1000 と調整した。

$ ruby logit.rb
ground truth:
Numo::DFloat#shape=[41]
[1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, ...]
predicted:
Numo::Int32#shape=[41]
[1, -1, 1, 1, 1, -1, 1, 1, -1, 1, 1, 1, 1, 1, -1, 1, -1, 1, -1, 1, -1, -1, ...]

Accuracy: 61.0%

苦労話

L-BFGS-B自体は、FORTRAN 77で書かれている。Pythonには、F2PYという、PythonFortranとをつなぐinterfaceを生成するものがあり、わりと簡単にFortranで書かれた科学計算ライブラリを扱える。Rubyにはそういったものがないので、L-BFGS-BをC言語に翻訳して、native extensionsで叩くことにした。C言語への翻訳には、まず定番のF2Cを使用した。F2Cで生成したC言語のコードは、当然ながらF2Cに依存する箇所がある。例えば、文字列の比較に、strncmp関数を使わずに、s_cmpという独自の関数を使用する。これらF2C独自の関数を、標準的なC言語で書き換えて、F2C依存を消していった。機械的に生成されたコードなため可読性が低く、なかなか骨の折れる作業であった。

おわりに

最適化手法のライブラリがあると、機械学習手法を実装するのが簡単になることが多い。Lbfgs.rbも、Rumaleの開発であるといいな、と思ったので作成した。以前に作成した最適化手法gemのMopti内のアルゴリズムの一種として実装しようかと思ったが、binding libraryとして単体であったほうが、取り回しが良いかなと思って、このような形にした。L-BFGS-B自体がバージョンアップすることはなさそうなので、今後の開発は、bugfixやC言語なコードの修正が、主となる予定である。

github.com