病みつきエンジニアブログ

機械学習、Python、Scala、JavaScript、などなど

Theanoを使ってPythonで行列演算とロジスティック回帰

TheanoというPython用のライブラリがあります。

ちょっと勉強したので、チュートリアルを日本語に翻訳しつつ、使い方とかを紹介します。

Theanoとは

まずはじめにTheanoとは、について。

Theanoはおそらく「テアノ」と読むのが多分正しいです。ピタゴラス(Pythagoras)の妻です。PythonとPythagorasをかけてるっぽいです。

何ができるかというと、簡単に言えば、行列演算ができます。以下の特徴を持っています(公式サイトより抜粋)

  • 実行スピードの最適化: Theano は g++などを使って式をコンパイルし、CPUやGPU操作に変換します(つまり、pure Pythonなコードよりも速い)
  • symbolicな微分: Theano は勾配を計算するために自動的にsymbolic graphsを組み立てます(訳注:つまり微分に便利だということ)
  • 安定な数値計算のための最適化: Theanoは不安定な数値計算を認識し、安定なアルゴリズムに変換します(訳注:例えばlogsumexpなどを言っていると思います)

なぜこれができるかというと、変数をシンボルとして扱うからだと思います。例えば、--xxにしたり、x * y / xyに変換することなどが明記されています。*1

この、変数を変数のまま保持しておくような考え方は、Theanoの設計を理解する上でちょっと重要な気もします。

また、Deep Learningとかニューラルネットの実装のためにpylearn2というライブラリが使われることもあるようですが、pylearn2はTheanoに依存しています。Theanoは、行列演算と微分の形になるNeural Networkと相性が良さそうですしね。

チュートリアル

(ロジスティック)シグモイドを造ります。

 s(x) = \frac{1}{1+e^{-x}}

で表せますから、これをそのまま数式にしてあげます。

まずは、必要な物をimportします。

import theano.tensor as T
from theano import function

次に、関数の引数に必要な変数を定義します。Theanoは型をちゃんと指定する必要があります(なんとなく、こうすることによって内部的に計算が高速になるメリットがあるんじゃないでしょうか)。

x = T.dscalar() # double型の数値
# x = T.dmatrix() # double型の行列

次に、シグモイド関数の計算式というか、関数の形を定義してあげます。

s = 1 / (1 + T.exp(-x))

このままでは、シグモイド関数sは呼び出し可能ではありません(s(0.5)みたいには呼び出せない)。ということで、呼び出し可能な関数を作ってあげます。

function 関数を使いますが、第1引数は、関数に必要な引数(ない場合は[])。第2引数は、「returnしてほしいもの」言い換えると関数の形とかです。

logistic = function([x], s)

ここで多分コンパイルが走ってる気がします。そしたらあとは呼び出すだけです。

logistic(0) # => array(0.5)

ここで、先ほどのxの型を「行列(dmatrix)」とかにしてあげると、呼び出し方も変わって、

logistic([[0]]) # => array([[0.5]])

となります。また、引数は行列なので、もちろん[[0, 1], [1.5, -1]]みたいな感じで行列を突っ込んで上げても良いです。

また、出力ももう少し拡張性が高いです。例えば、関数の形を配列に突っ込んであげたり。

logistic = function([x], [s, s])
logistic(0) # => [array(0.5), array(0.5)]

値の保持

値をどこかグロバール領域に保持しつつ計算することもできます(shared variables)。チュートリアルのサンプルでは、勾配法でロジスティック回帰を行っていますが、重みを保持し(つまりreturnしない)、ひたすら書き換えることによって実装しています。

このメリットは、グローバル変数と同じく共有しやすいことと(結果として表記上わかりやすくなる)、主に計算上のメリットだそうです。

shared variableを使うにはshared関数と、update引数を使います。

from theano import shared
state = shared(0) # 0は初期値
inc = T.iscalar('inc') # integer型のスカラー
accumulator = function([inc], state, updates=[(state, state+inc)])
# equals function([inc], state, updates={state: state+inc})

functionの引数のうち、2つめは必須ではありません。あくまで実行前の状態をreturnしてくれるというだけです(状態が確認しやすいだけ)。updateの引数は、各shared varibleに対する操作の詳細を示します。(書き込み先のshared variable, 関数の形)というタプルで指定します。まあ意味的には辞書みたいなものなので、dict型でもいいです。

乱数

乱数を使うこともできます。

from theano.tensor.shared_randomstreams import RandomStreams
srng = RandomStreams(seed=234)
rv_u = srng.uniform((2,2))
f = function([], rv_u)
f() # such as => array([[ 0.28179047,  0.23616647], [ 0.5958365 ,  0.1385743 ]])

例えば、box-muller変換を使った標準正規分布は、

rv_X = srng.uniform((1,2))
rv_Y = srng.uniform((1,2))
box_muller = T.sqrt(-2 * T.log(rv_X)) * T.cos(2 * math.pi * rv_Y)
normal = function([], box_muller)()
normal() # => array([[-0.53976774,  1.09561059]])

ただし、同じ乱数は式中で同じ値になるので注意が必要です。

function([], rv_X - rv_X)() #=> array([[ 0.,  0.]])

ロジスティック回帰の実装

説明に疲れてきたので詳細は割愛します、本家のチュートリアルを見てください。

こちらでは、L2正則化つきのロジスティック回帰が実装されています。具体的な微分式を与えず、

gw, gb = T.grad(cost, [w, b])

とすることで誤差の微分(勾配)を求めています。そして、

train = theano.function(
          inputs=[x,y],
          outputs=[prediction, xent],
          updates=((w, w - 0.1 * gw), (b, b - 0.1 * gb)))

という風に、updatesに((w, w - 0.1 * gw), (b, b - 0.1 * gb))としています。これは更新率ηを0.1として、shared variableのwをひたすら書き換える、という感じです。

また、このoutputs引数に指定された[prediction, xent]は全く使われていませんので、[]にしても動くと思います

参考文献

*1:実際にどれくらい認識して数値計算を安定化、もしくは省略してくれるのかはわからないので、結局自力で一番良い計算方法をコードに落としこむことが多いと思いますが。。。