sage/SageでTheanoのロジスティック回帰を試す
をテンプレートにして作成
[
トップ
] [
新規
|
一覧
|
単語検索
|
最終更新
|
ヘルプ
]
開始行:
[[FrontPage]]
#contents
2016/05/03からのアクセス回数 &counter;
ここで紹介したSageワークシートは、以下のURLからダウンロー...
http://www15191ue.sakura.ne.jp:8000/home/pub/71/
また、私の公開しているSageのサーバ（http://www15191ue.sak...
アップロードし、実行したり、変更していろいろ動きを試すこ...
** 参考サイト [#rd740127]
ここでは、
[[人工知能に関する断創録>http://aidiary.hatenablog.com/]]
のTheanoに関連する記事をSageのノートブックで実装し、Theno...
今回は、TheanoのTutorialからMNISTの手書き数字認識を以下の...
- [[Theanoによる2クラスロジスティック回帰の実装>http://ai...
** ２値分類 [#x27ba760]
２値分類とは、入力xの値によって2種類に区別する問題です。...
入力xからyの値を推定します。
入力xを指定したとき、y=1となる事後確率\(p(y=1|x)\)をニュ...
$$
h(x) = p(y=1|x) \approx g(x; w)
$$
$$
g(z) = \frac{1}{1 + e^{-z}}
$$
べき乗の性質を使ったトリックを使うと事後分布\(p(y|x; w)\)...
$$
p(y|x) = p(y=1|x)^y p(y=0|x)^{1-y}
$$
と表すことができます。ここで\(p(y=0|x) = 1 - h(x)\)である...
$$
L(w) = \prod_{n=1}^N p(y_n | x_n; w) = \prod_{n=1}^N \{ h...
$$
で与えられ、負の対数尤度は、以下のようになります。
$$
E(w) = - \sum_{n=1}^N [y_n log h(x_n) + (1-y_n) log\{1 - ...
$$
*** 必要なライブラリをインポート [#bc217829]
最初に使用するライブラリをインポートします。
sageへの入力:
#pre{{
import numpy as np
import matplotlib.pyplot as plt
import theano
import theano.tensor as T
}}
*** 訓練データのロード [#i9acc48c]
訓練データex2data1.txtをDATAディレクトリから読み込みます。
sageへの入力:
#pre{{
# 訓練データをロード
data = np.genfromtxt(DATA+"ex2data1.txt", delimiter=",")
data_x = data[:, (0, 1)]
data_y = data[:, 2]
# 訓練データ数
m = len(data_y)
}}
*** データの分布 [#e994bef7]
y=1を赤、y=0を青の●(point)でプロットし、データの分布をみ...
SageのGraphicsとpointを使うと簡単にデータの分布が可視化す...
sageへの入力:
#pre{{
# データのプロット
data_plt = Graphics()
for i in range(m):
if data_y[i] == 1:
data_plt += point(data_x[i], rgbcolor='red')
else:
data_plt += point(data_x[i], rgbcolor='blue')
show(data_plt, figsize=5)
}}
&ref(sage0.png);
*** コスト関数の自動微分 [#e7fe7e72]
バイアス項もwに追加するために、入力データxの1列目に1の列...
従って重みwは、\( w = [w_0, w_1, w_2]\)となります。
Theanoでコスト関数で定義します。
- X: 訓練データを表す共有変数
- y: 訓練データの正解ラベルを表す共有変数
- w: 重み（パラメータ）を表す共有変数
人工知能に関する断創録のコメントは鋭く、T.dot(X, w)のXとw...
これはXの行がサンプル、列が特徴をデータをまとめた行列なの...
g_wは、コスト関数をパラメータwで微分した数式を表している。
このようにとても簡単に微分を定義することができるのがThean...
sageへの入力:
#pre{{
# 訓練データの1列目に1を追加
data_x = np.hstack((np.ones((m, 1)), data_x))
# 訓練データを共有変数にする
X = theano.shared(np.asarray(data_x, dtype=theano.config....
y = theano.shared(np.asarray(data_y, dtype=theano.config....
}}
sageへの入力:
#pre{{
# パラメータを共有変数にし、0で初期化
# 訓練データに1を加えたのでバイアスもwに含めてしまう
w = theano.shared(np.zeros(3, dtype=theano.config.floatX)...
# コスト関数として負の対数尤度を定義
h = T.nnet.sigmoid(T.dot(X, w))
# 通常のコスト関数とは異なり1/mがファクターとして掛けてあ...
cost = - (1.0 / m) *T.sum(y * T.log(h) + (1 - y) * T.log(...
# コスト関数の微分
g_w = T.grad(cost=cost, wrt=w)
}}
** Theanoによる勾配降下法の実装 [#n35429a4]
学習率learning_rateと更新updatesにwの更新を定義し、theano...
これを繰り返し実行するだけで勾配降下法が実現できてしまう...
#pre{{
# 訓練用の関数を定義
train_model = theano.function(inputs=[], outputs=cos...
}}
functionで指定したupdatesは、関数の実行後に指定された更新...
updatesには（共有変数, 更新式のシンボル）のペアが並ぶリス...
train_modelの戻り値には、cost関数の戻り値が返されます、こ...
今回は収束判定はしないで、30万回パラメータの更新を繰り返...
sageへの入力:
#pre{{
# 勾配降下法
# パラメータ更新式
learning_rate = 0.001
updates = [(w, w - learning_rate * g_w)]
# 訓練用の関数を定義
train_model = theano.function(inputs=[], outputs=cost, up...
}}
sageへの入力:
#pre{{
# 高度な収束判定はせずにiterations回だけ繰り返す
iterations = 300000
for iter in range(iterations):
current_cost = train_model()
if iter%20000 == 0:
print iter, current_cost
}}
#pre{{
0 0.69314718056
20000 0.547295463694
40000 0.488295709188
60000 0.445118809496
80000 0.41259440171
100000 0.387389522711
120000 0.367353248365
140000 0.351070505109
160000 0.337586093857
180000 0.326237734121
200000 0.316554335209
220000 0.308192817261
240000 0.300897944917
260000 0.294476164256
280000 0.288778154733
}}
*** 結果をプロット [#i5307f4c]
求められた重みwからy=0とy=1の境界を求めます。
$$
w_0 + w_1 x_1 + w_2 x_2 = 0
$$
から
$$
x2 = -w_0/w_2 - w_1/w_2 x_1
$$
となります。これをSageのPlot関数を使って表示してみましょ...
sageへの入力:
#pre{{
# 更新されたパラメータを表示
t = w.get_value()
print "w:", t
}}
#pre{{
w: [-9.25573205 0.07960975 0.07329322]
}}
sageへの入力:
#pre{{
# 決定境界を描画
xmin, xmax = min(data_x[:,1]), max(data_x[:,1])
x1 = var('x1')
x2 = -t[0]/t[2] - t[1]/t[2]*x1
line_plt = plot(x2, [x1, xmin, xmax])
(data_plt + line_plt).show(figsize=5)
}}
&ref(sage1.png);
*** Theanoによる確率的勾配降下法の実装 [#sd266672]
上記の勾配降下法は、誤差関数に各サンプルの誤差\(E_n(w)\)...
$$
E(w) = \sum_{n=1}^N E_n(w)
$$
これに対してサンプルからランダムに1個を抽出し、パラメータ...
と呼び、SGDと略します。SGDでのwの更新は以下のように行われ...
$$
w^{(t+1)} = w^{(t)} - \epsilon \Delta E_n
$$
TheanoでSGDを実装するには、更新式を少し書き替えるだけで実...
#pre{{
# 訓練データのインデックスを表すシンボルを定義
index = T.lscalar()
# コスト関数の微分をindex番目のデータのみ使うように...
h = T.nnet.sigmoid(T.dot(w, X[index,:]))
cost = -y[index] * T.log(h) - (1 - y[index]) * T.log...
# 訓練用の関数の定義で、inputにindexを指定し、index...
train_model = theano.function(inputs=[index], output...
}}
人工知能に関する断創録では、5000で計算していましたが、解...
sageへの入力:
#pre{{
# 確率的勾配降下法
# 訓練データを再度設定
data_x = data[:, (0, 1)]
data_y = data[:, 2]
# 訓練データの1列目に1を追加
data_x = np.hstack((np.ones((m, 1)), data_x))
# データをシャッフル
p = np.random.permutation(m)
data_x = data_x[p, :]
data_y = data_y[p]
}}
sageへの入力:
#pre{{
# 訓練データを共有変数にする
X = theano.shared(np.asarray(data_x, dtype=theano.config....
y = theano.shared(np.asarray(data_y, dtype=theano.config....
# パラメータを共有変数にし、0で初期化
# 訓練データに1を加えたのでバイアスもwに含めてしまう
w = theano.shared(np.zeros(3, dtype=theano.config.floatX)...
# 訓練データのインデックスを表すシンボルを定義
index = T.lscalar()
}}
sageへの入力:
#pre{{
# コスト関数の微分を構築
# 確率的勾配降下法なので全データの和ではなく、index番目の...
h = T.nnet.sigmoid(T.dot(w, X[index,:]))
cost = -y[index] * T.log(h) - (1 - y[index]) * T.log(1 - h)
# コスト関数の微分
g_w = T.grad(cost=cost, wrt=w)
# 更新式
learning_rate = 0.0005
updates = [(w, w - learning_rate * g_w)]
# 訓練用の関数を定義
# index番目の訓練データを使ってパラメータ更新
train_model = theano.function(inputs=[index], outputs=cos...
}}
sageへの入力:
#pre{{
# 確率的勾配降下法
# 5000だと解がばらばらなので、10000に増やしてみた
max_epoch = 10000
for epoch in range(max_epoch):
for i in range(m):
current_cost = train_model(i)
if epoch%500 == 0:
print epoch, current_cost
# 更新されたパラメータを表示
t = w.get_value()
print "w:", t
}}
#pre{{
0 0.256765757093
500 0.199191597088
1000 0.159751041384
1500 0.136903825723
2000 0.121894838402
2500 0.110467015435
3000 0.100861036861
3500 0.0923233412348
4000 0.0845513113405
4500 0.0774410166041
5000 0.0709542020968
5500 0.0650615672185
6000 0.059727054437
6500 0.0549077671741
7000 0.0505575651269
7500 0.0466303292568
8000 0.0430820865081
8500 0.0398721475311
9000 0.036963576944
9500 0.034323261154
w: [-13.15095347 0.10151697 0.10169091]
}}
*** 決定境界の表示 [#n1cf4127]
解は勾配降下法の線の当たりをばらつきます。
sageへの入力:
#pre{{
# 決定境界を描画
xmin, xmax = min(data_x[:,1]), max(data_x[:,1])
x1 = var('x1')
x2 = -t[0]/t[2] - t[1]/t[2]*x1
line_plt = plot(x2, [x1, xmin, xmax])
(data_plt + line_plt).show(figsize=5)
}}
&ref(sage2.png);
** コメント [#d9b36c07]
#vote(おもしろかった,そうでもない,わかりずらい)
皆様のご意見、ご希望をお待ちしております。
#comment_kcaptcha
終了行:
[[FrontPage]]
#contents
2016/05/03からのアクセス回数 &counter;
ここで紹介したSageワークシートは、以下のURLからダウンロー...
http://www15191ue.sakura.ne.jp:8000/home/pub/71/
また、私の公開しているSageのサーバ（http://www15191ue.sak...
アップロードし、実行したり、変更していろいろ動きを試すこ...
** 参考サイト [#rd740127]
ここでは、
[[人工知能に関する断創録>http://aidiary.hatenablog.com/]]
のTheanoに関連する記事をSageのノートブックで実装し、Theno...
今回は、TheanoのTutorialからMNISTの手書き数字認識を以下の...
- [[Theanoによる2クラスロジスティック回帰の実装>http://ai...
** ２値分類 [#x27ba760]
２値分類とは、入力xの値によって2種類に区別する問題です。...
入力xからyの値を推定します。
入力xを指定したとき、y=1となる事後確率\(p(y=1|x)\)をニュ...
$$
h(x) = p(y=1|x) \approx g(x; w)
$$
$$
g(z) = \frac{1}{1 + e^{-z}}
$$
べき乗の性質を使ったトリックを使うと事後分布\(p(y|x; w)\)...
$$
p(y|x) = p(y=1|x)^y p(y=0|x)^{1-y}
$$
と表すことができます。ここで\(p(y=0|x) = 1 - h(x)\)である...
$$
L(w) = \prod_{n=1}^N p(y_n | x_n; w) = \prod_{n=1}^N \{ h...
$$
で与えられ、負の対数尤度は、以下のようになります。
$$
E(w) = - \sum_{n=1}^N [y_n log h(x_n) + (1-y_n) log\{1 - ...
$$
*** 必要なライブラリをインポート [#bc217829]
最初に使用するライブラリをインポートします。
sageへの入力:
#pre{{
import numpy as np
import matplotlib.pyplot as plt
import theano
import theano.tensor as T
}}
*** 訓練データのロード [#i9acc48c]
訓練データex2data1.txtをDATAディレクトリから読み込みます。
sageへの入力:
#pre{{
# 訓練データをロード
data = np.genfromtxt(DATA+"ex2data1.txt", delimiter=",")
data_x = data[:, (0, 1)]
data_y = data[:, 2]
# 訓練データ数
m = len(data_y)
}}
*** データの分布 [#e994bef7]
y=1を赤、y=0を青の●(point)でプロットし、データの分布をみ...
SageのGraphicsとpointを使うと簡単にデータの分布が可視化す...
sageへの入力:
#pre{{
# データのプロット
data_plt = Graphics()
for i in range(m):
if data_y[i] == 1:
data_plt += point(data_x[i], rgbcolor='red')
else:
data_plt += point(data_x[i], rgbcolor='blue')
show(data_plt, figsize=5)
}}
&ref(sage0.png);
*** コスト関数の自動微分 [#e7fe7e72]
バイアス項もwに追加するために、入力データxの1列目に1の列...
従って重みwは、\( w = [w_0, w_1, w_2]\)となります。
Theanoでコスト関数で定義します。
- X: 訓練データを表す共有変数
- y: 訓練データの正解ラベルを表す共有変数
- w: 重み（パラメータ）を表す共有変数
人工知能に関する断創録のコメントは鋭く、T.dot(X, w)のXとw...
これはXの行がサンプル、列が特徴をデータをまとめた行列なの...
g_wは、コスト関数をパラメータwで微分した数式を表している。
このようにとても簡単に微分を定義することができるのがThean...
sageへの入力:
#pre{{
# 訓練データの1列目に1を追加
data_x = np.hstack((np.ones((m, 1)), data_x))
# 訓練データを共有変数にする
X = theano.shared(np.asarray(data_x, dtype=theano.config....
y = theano.shared(np.asarray(data_y, dtype=theano.config....
}}
sageへの入力:
#pre{{
# パラメータを共有変数にし、0で初期化
# 訓練データに1を加えたのでバイアスもwに含めてしまう
w = theano.shared(np.zeros(3, dtype=theano.config.floatX)...
# コスト関数として負の対数尤度を定義
h = T.nnet.sigmoid(T.dot(X, w))
# 通常のコスト関数とは異なり1/mがファクターとして掛けてあ...
cost = - (1.0 / m) *T.sum(y * T.log(h) + (1 - y) * T.log(...
# コスト関数の微分
g_w = T.grad(cost=cost, wrt=w)
}}
** Theanoによる勾配降下法の実装 [#n35429a4]
学習率learning_rateと更新updatesにwの更新を定義し、theano...
これを繰り返し実行するだけで勾配降下法が実現できてしまう...
#pre{{
# 訓練用の関数を定義
train_model = theano.function(inputs=[], outputs=cos...
}}
functionで指定したupdatesは、関数の実行後に指定された更新...
updatesには（共有変数, 更新式のシンボル）のペアが並ぶリス...
train_modelの戻り値には、cost関数の戻り値が返されます、こ...
今回は収束判定はしないで、30万回パラメータの更新を繰り返...
sageへの入力:
#pre{{
# 勾配降下法
# パラメータ更新式
learning_rate = 0.001
updates = [(w, w - learning_rate * g_w)]
# 訓練用の関数を定義
train_model = theano.function(inputs=[], outputs=cost, up...
}}
sageへの入力:
#pre{{
# 高度な収束判定はせずにiterations回だけ繰り返す
iterations = 300000
for iter in range(iterations):
current_cost = train_model()
if iter%20000 == 0:
print iter, current_cost
}}
#pre{{
0 0.69314718056
20000 0.547295463694
40000 0.488295709188
60000 0.445118809496
80000 0.41259440171
100000 0.387389522711
120000 0.367353248365
140000 0.351070505109
160000 0.337586093857
180000 0.326237734121
200000 0.316554335209
220000 0.308192817261
240000 0.300897944917
260000 0.294476164256
280000 0.288778154733
}}
*** 結果をプロット [#i5307f4c]
求められた重みwからy=0とy=1の境界を求めます。
$$
w_0 + w_1 x_1 + w_2 x_2 = 0
$$
から
$$
x2 = -w_0/w_2 - w_1/w_2 x_1
$$
となります。これをSageのPlot関数を使って表示してみましょ...
sageへの入力:
#pre{{
# 更新されたパラメータを表示
t = w.get_value()
print "w:", t
}}
#pre{{
w: [-9.25573205 0.07960975 0.07329322]
}}
sageへの入力:
#pre{{
# 決定境界を描画
xmin, xmax = min(data_x[:,1]), max(data_x[:,1])
x1 = var('x1')
x2 = -t[0]/t[2] - t[1]/t[2]*x1
line_plt = plot(x2, [x1, xmin, xmax])
(data_plt + line_plt).show(figsize=5)
}}
&ref(sage1.png);
*** Theanoによる確率的勾配降下法の実装 [#sd266672]
上記の勾配降下法は、誤差関数に各サンプルの誤差\(E_n(w)\)...
$$
E(w) = \sum_{n=1}^N E_n(w)
$$
これに対してサンプルからランダムに1個を抽出し、パラメータ...
と呼び、SGDと略します。SGDでのwの更新は以下のように行われ...
$$
w^{(t+1)} = w^{(t)} - \epsilon \Delta E_n
$$
TheanoでSGDを実装するには、更新式を少し書き替えるだけで実...
#pre{{
# 訓練データのインデックスを表すシンボルを定義
index = T.lscalar()
# コスト関数の微分をindex番目のデータのみ使うように...
h = T.nnet.sigmoid(T.dot(w, X[index,:]))
cost = -y[index] * T.log(h) - (1 - y[index]) * T.log...
# 訓練用の関数の定義で、inputにindexを指定し、index...
train_model = theano.function(inputs=[index], output...
}}
人工知能に関する断創録では、5000で計算していましたが、解...
sageへの入力:
#pre{{
# 確率的勾配降下法
# 訓練データを再度設定
data_x = data[:, (0, 1)]
data_y = data[:, 2]
# 訓練データの1列目に1を追加
data_x = np.hstack((np.ones((m, 1)), data_x))
# データをシャッフル
p = np.random.permutation(m)
data_x = data_x[p, :]
data_y = data_y[p]
}}
sageへの入力:
#pre{{
# 訓練データを共有変数にする
X = theano.shared(np.asarray(data_x, dtype=theano.config....
y = theano.shared(np.asarray(data_y, dtype=theano.config....
# パラメータを共有変数にし、0で初期化
# 訓練データに1を加えたのでバイアスもwに含めてしまう
w = theano.shared(np.zeros(3, dtype=theano.config.floatX)...
# 訓練データのインデックスを表すシンボルを定義
index = T.lscalar()
}}
sageへの入力:
#pre{{
# コスト関数の微分を構築
# 確率的勾配降下法なので全データの和ではなく、index番目の...
h = T.nnet.sigmoid(T.dot(w, X[index,:]))
cost = -y[index] * T.log(h) - (1 - y[index]) * T.log(1 - h)
# コスト関数の微分
g_w = T.grad(cost=cost, wrt=w)
# 更新式
learning_rate = 0.0005
updates = [(w, w - learning_rate * g_w)]
# 訓練用の関数を定義
# index番目の訓練データを使ってパラメータ更新
train_model = theano.function(inputs=[index], outputs=cos...
}}
sageへの入力:
#pre{{
# 確率的勾配降下法
# 5000だと解がばらばらなので、10000に増やしてみた
max_epoch = 10000
for epoch in range(max_epoch):
for i in range(m):
current_cost = train_model(i)
if epoch%500 == 0:
print epoch, current_cost
# 更新されたパラメータを表示
t = w.get_value()
print "w:", t
}}
#pre{{
0 0.256765757093
500 0.199191597088
1000 0.159751041384
1500 0.136903825723
2000 0.121894838402
2500 0.110467015435
3000 0.100861036861
3500 0.0923233412348
4000 0.0845513113405
4500 0.0774410166041
5000 0.0709542020968
5500 0.0650615672185
6000 0.059727054437
6500 0.0549077671741
7000 0.0505575651269
7500 0.0466303292568
8000 0.0430820865081
8500 0.0398721475311
9000 0.036963576944
9500 0.034323261154
w: [-13.15095347 0.10151697 0.10169091]
}}
*** 決定境界の表示 [#n1cf4127]
解は勾配降下法の線の当たりをばらつきます。
sageへの入力:
#pre{{
# 決定境界を描画
xmin, xmax = min(data_x[:,1]), max(data_x[:,1])
x1 = var('x1')
x2 = -t[0]/t[2] - t[1]/t[2]*x1
line_plt = plot(x2, [x1, xmin, xmax])
(data_plt + line_plt).show(figsize=5)
}}
&ref(sage2.png);
** コメント [#d9b36c07]
#vote(おもしろかった,そうでもない,わかりずらい)
皆様のご意見、ご希望をお待ちしております。
#comment_kcaptcha
ページ名:
SmartDoc
![[PukiWiki]](image/pukiwiki.png "[PukiWiki]")
