ロジスティック回帰

ディープラーニングをやってみたいので
http://deeplearning.net/tutorial/intro.html
を読んでます。

自分なりにまとめていきます。

今回は
http://deeplearning.net/tutorial/logreg.html#logreg

ロジスティック回帰ってのはクラス分類問題の一種。

$D_i$ 次元入力ベクトル $x$ が入力されたとき、それがどのクラスに所属するのかってのを識別する問題で、学習モデルに

$P(Y = i | x, W, b) = softmax_i (W x + b) = \frac{\exp\{W_i x + b_i\}}{\sum_j \exp \{W_j x + b_j \}}$

を使う奴。

これはパラメタが $W$ と $\theta$ で与えられたとき、入力 $x$ が $i$ 番目のクラスである確率。

クラスが $D_o$ 個なら、 $W$ は $D_o$ 行 $D_i$ 列の行列。

$x$ が入力されたとき、 $P(Y = i | x, W, b)$ を最大にする $i$ に $x$ を分類するとする。

回帰のよくあるパターンは学習データに対する尤度関数を最大化するとかそういうやつ。

今回のモデルの対数尤度関数は以下。

$\mathcal{L} (\theta=\{W,b\}, \mathcal{D}) = \sum_{i=0}^{|\mathcal{D}|} \log(P(Y=y^{(i)}|x^{(i)}, W,b))$

最小化問題にしたいから対数尤度関数に-1をかけてこの $\ell$ を最小化する事を考える。

$\ell (\theta=\{W,b\}, \mathcal{D}) = - \mathcal{L} (\theta=\{W,b\}, \mathcal{D})$

$\frac{\partial \ell}{\partial W} = 0$ と $\frac{\partial \ell}{\partial b} = 0$ が解ければそのときの $W$ と $b$ を使えばいいんだけど、今後のことを考えて反復してパラメタを更新しながら最適化する方法をとる。

最急降下法と確率的勾配降下法

$\theta$ でパラメタライズされる関数 $f(\theta)$ の値を最小化したいとき、以下の反復を適当な初期値から始めて繰り返すことで $f$ の最小値を与える $\theta$ を得ることが出来る。

$\theta^{i+1} = \theta^i - \eta \frac{\partial}{\partial \theta}f(\theta^i)$

この方法が最急降下法。

ところで、機械学習において、訓練データ $X = \{x_0, \dots , x_{N-1} \}$ に対して、
$F(\theta) = \sum_{n = 0}^{N-1}f_{\theta}(x_n)$ みたいな形の関数がよく出てきてこれを最小化するってパターンがよくある。

これを最急降下法でとくと
$\theta^{i+1} = \theta^i - \eta \sum_{n = 0}^{N-1} \frac{\partial}{\partial \theta}f_{\theta^i}(x_n)$
となって、一回のパラメタの更新につきN回の偏微分が必要になる。

Nが10とか100ならまだしも何十万、何千万ってなってくると厳しい。

ということで、確率的勾配降下法。

$\theta^{i+1} = \theta^i - \eta \frac{\partial}{\partial \theta}f_{\theta^i}(x_n)$

を、各 $i$ に対して $n$ を順々に変えながらパラメタを更新していく。

ロジスティック回帰で言えば
$\theta^{i+1} = \theta^i - \eta \frac{\partial}{\partial \theta}( - \log(P(Y=y^{(n)}|x^{(n)}, \theta^i)))$
となる。

Theanoを使ってロジスティック回帰を確率的勾配降下法で解く

データ

以下のデータを使用して手書き文字画像をクラス分類する。

http://deeplearning.net/data/mnist/mnist.pkl.gz

pythonでこのファイルを読むには以下を実行する

import cPickle
import gzip
f = gzip.open('mnist.pkl.gz')
train_set, valid_set, test_set = cPickle.load(f)
f.close()

手書き文字画像は28×28ピクセルのグレー画像で各ピクセルは0から1の値をとる(0は黒、1は白)。

データがどう格納されているのか確認する。
train_setは要素を二つもつタプル型で、１つ目は画像データの集合、２つ目はそれらの画像につけられたラベルの集合でこれらはどちらもnumpy.ndarray型

>> type(train_set)
tuple
>> train_images = train_set[0]
>> type(train_images)
numpy.ndarray
>> train_images.size / (28*28)
50000
>> train_labels = train_set[1]
>> type(train_labels)
numpy.ndarray
>> train_labels.size
50000
>> train_labels[0]
5 # <= 下の画像は数字の'5'を書いたもの
>> image_data = train_images[0]
>> image_data.shape
(784,)
>> image = image_data.reshape(28, 28)
>> import matplotlib.pyplot as plt
>> plt.imshow(image)
>> plt.show()

一つの画像は784(=28*28)要素のベクトルとして表される。これを入力 $x_n$ とし、それに付けられたラベルを $y_n$ として学習して $W$ と $b$ を求める。

実装

まずは

$P(Y = i | x, W, b) = softmax_i (W x + b) = \frac{\exp\{W_i x + b_i\}}{\sum_j \exp \{W_j x + b_j \}}$

から

import numpy
import theano
import theano.tensor as T

n_in = 28 * 28 # ベクトル化した入力画像の次元
n_out = 10     # 出力ベクトルの次元
input = T.matrix('x') # 

W = theano.shared(value=numpy.zeros((n_in, n_out), dtype=theano.config.floatX),
                  name='W', borrow=True)
b = theano.shared(value=numpy.zeros((n_out,), dtype=theano.config.floatX),
                  name='b', borrow=True)
p_y_given_x = T.nnet.softmax(T.dot(input, W) + b)

p_y_given_x_func = theano.function([input], p_y_given_x)

# Wの要素もbの要素もすべて0の状態で図示したデータが各クラスに分類される確率を表示してみる
p_y_given_x_func(train_images[0].reshape(1,784))
# => array([[ 0.1,  0.1,  0.1,  0.1,  0.1,  0.1,  0.1,  0.1,  0.1,  0.1]])

つぎはこれ
$\ell (\theta=\{W,b\}, \mathcal{D}) = - \mathcal{L} (\theta=\{W,b\}, \mathcal{D})$

negative_log_likelihood = -T.log(p_y_given_x)

negative_log_likelihood_func = theano.function(inputs=[input, y], outputs=negative_log_likelihood[0,y])

# Wの要素もbの要素もすべて0の状態で図示したデータが入力され、それを数字の'3'だと分類したときの尤度を表示してみる
negative_log_likelihood_func(train_images[0].reshape(1,784), [3])
# => array([[ 2.30258509]])

最後にこれ

$\theta^{i+1} = \theta^i - \eta \frac{\partial}{\partial \theta}( - \log(P(Y=y^{(n)}|x^{(n)}, \theta^i)))$

cost = negative_log_likelihood[0,y]
g_W = T.grad(cost=cost, wrt=W)
g_b = T.grad(cost=cost, wrt=b)
learning_rate = 0.1
updates =[(W, W - learning_rate * g_W), (b, b - learning_rate * g_b)]

ここまでくれば学習データをつかってパラメタを更新していく事が出来る。

train_model = theano.function(inputs=[input, y], outputs=cost, updates=updates)
i = 0
train_model(train_images[i].reshape(1, 784), [train_labels[i]])
i = i + 1

最後の二行を繰り替えしていくと新たな画像が入力されて、それをもとにパラメタが更新されていく。

3000枚ほどで学習してからふたたび図示したデータが各クラスに分類される確率を表示してみる。

 p_y_given_x_func(train_images[0].reshape(1,784))
# => array([[  3.79957900e-03,   1.98228237e-07,   1.53438291e-04,
#              1.91115501e-02,   1.97851087e-10,   9.45596439e-01,
#              3.04181820e-05,   2.23148970e-03,   1.02881405e-02,
#              1.87887465e-02]])

'5'に分類される確率が0.945で他より十分高くなっており、学習が出来ているっぽいことが確認できた。