boost/range/adaptor 覚え書き

出力は以下の通り。

Convolutional Neural Network によるシーン認識

はじめに

　先のページで Convolutional Neural Network（CNN）を解説した。今回は、CNN をシーン認識に適用した結果を示す。

データセット

　このページに置いてある データ（LSP15）を使用した。このデータは15個のディレクトリから構成される。

1. bedroom
2. CALsuburb
3. industrial
4. kitchen
5. livingroom
6. MITcoast
7. MITforest
8. MIThighway
9. MITinsidecity
10. MITmountain
11. MITopencountry
12. MITstreet
13. MITtallbuilding
14. PARoffice
15. store

各ディレクトリ内にはその名前が表すシーンの画像が複数枚おさめられている。 今回のシーン識別では「MIT」の付くディレクトリだけを利用した。画像の大きさが統一されており扱い易かったからである。今回あつかう問題は8分類問題である。

データセットの増量

　各シーンの画像の枚数は以下の通りである。

シーン名	枚数
MITcoast	360(720)
MIThighway	260(520)
MITmountain	374(748)
MITstreet	292(584)
MITforest	328(656)
MITinsidecity	308(616)
MITopencountry	410(820)
MITtallbuilding	356(712)
総数	2688(5376)

画像の枚数を増やすため、各画像の左右反転画像を追加した。括弧内が追加後の枚数である。

データセットの分割

　データセットを train、validation、test の3つに分割する。枚数の比を train : validation : test = 5 : 1 : 1 とした。 以下に具体的な枚数を示す。
　

ラベル名	シーン名	train	validation	test
0	MITcoast	514	102	102
1	MIThighway	371	74	74
2	MITmountain	534	106	106
3	MITstreet	417	83	83
4	MITforest	468	93	93
5	MITinsidecity	440	88	88
6	MITopencountry	585	117	117
7	MITtallbuilding	508	101	101
	総数	3837	764	764

また、シーン認識時に使用するラベルを上図の第1列のように定義した。 左右反転画像は元の画像とは独立した画像として扱い、train / validation / test の3つに分ける際は乱数を用いた。

ネットワークの構造

　今回使用したネットワークの構造は以下の通りである。

convolutional 層と pooling 層の対を4つ、隠れ層を1つとした。 各種パラメータは以下の通りである。
　 　 　 　 　 　

層名	フィルタサイズ	フィルタの個数	pooling
0cp	5$\times$5	2	2$\times$2
1cp	6$\times$6	10	2$\times$2
2cp	6$\times$6	50	2$\times$2
3cp	5$\times$5	200	2$\times$2

フィルタサイズとは、convolution 時に使用する矩形画像のサイズのことである（フィルタの各要素の値は訓練時に決定される）。また、convolution時の stride は全て1とした。poolingはmax poolingを採用した。poolingのあと画像サイズは縦横ともに1/2になる。 隠れ層と出力層のユニット数は以下の通りである。
　 　 　 　 　

層名	ユニット数
隠れ層	500
出力層	8

各層における画素数の変遷は次のようになる。
　 　 　 　 　 　

層名	入力サイズ	conv後のサイズ	pooling後のサイズ
0cp	146$\times$146	142$\times$142	71$\times$71
1cp	71$\times$71	66$\times$66	33$\times$33
2cp	33$\times$33	28$\times$28	14$\times$14
3cp	14$\times$14	10$\times$10	5$\times$5

0cp層への入力画像に対しては、画素値を255で割ったあと、平均値が0、分散が1となるように正規化した。 3cp層から出力される画素の数は5000=(5$\times$5$\times$200)である。これが隠れ層への入力値となる。今考えているラベルの総数は8であるから出力層のユニット数は8となる。

計算方法

　Convolutional Neural Networks (LeNet) にあるサンプルソース（theanoによる実装）に手を加え計算を行なった。主な変更点は以下の通りである。

1. 上に示したネットワーク構造に変更した。
2. $L_1$正則化項と$L_2$正則化項を追加した。
3. 各epochごとに学習率が減衰していくようにした。
4. 隠れ層の活性化関数をsigmoidとした（tanhより精度が良くなった）。

正則化項についてはこちらで説明した。 パラメータは以下の通りである。
　 　

学習率	減衰率	$L_1$正則化項の重み	$L_2$正則化項の重み
0.01	0.995	0.0002	0.0002

元のソースでは、early-stoppingと呼ばれる手法を使って訓練を行っている。「訓練データ（train）で学習したあと検証データ（validation）で誤り率を見る」を、検証データで改善が見られなくなるまで繰り返す手法である。この部分はそのまま踏襲した。

計算環境

　アマゾンのEC2インスタンス（g2.2xlarge）を利用した。これは1536CUDA コアと 4GB のビデオメモリを搭載しており、今回の計算ではtheanoを通してGPUを利用した。

結果

　0/1損失（0/1 loss）の結果を以下に示す。
　 　

best validation score	test performance
11.9736842105%	10.6578947368%

ここで、best validation score とは、early-stopping の反復計算の過程で得られた「検証データ」に対する最適値である。 「検証データ」に対する最適値を実現したときのモデルを使って「テストデータ」を評価した結果が test performance である。 ただし、反復時の epoch 数の上限値を500とした。計算時間は 182.905974067 分であった。

考察

1. 左右反転画像を追加する前の誤り率は20%程度であり、これ以上減らすことができなかった。
2. 左右反転画像を追加した後の誤り率は10%程度となり改善された。
3. 白色化や隠れ層の追加などを試みたがこれ以上は改善しなかった。
4. なんとかして誤り率を一桁にしたいのだが。。。

正則化項について

はじめに

　これまでいくつかのNeural Networkを考察してきた（ 1次元ロジスティック回帰、 多次元ロジスティック回帰、 多層パーセプトロン、 Convolutional Neural Network）。 どの場合もパラメータを求める際はNegative Log Likelihood（NLL）を最小化した。過学習を防ぐため、正則化項（regularization term）と呼ばれる項をNLLに追加することがある。ここでは、この正則化項の役割についてまとめる。

正則化項

　多次元ロジスティック回帰を例として考える。全観測値の実現確率は次式で与えられた。 \begin{equation} P(W,\vec{b})=\prod_{n=1}^{N}\;p(c^{n}\;|\;\vec{x}^{\;n}\;;\;W,\vec{b}) \end{equation} ここで、$\vec{x}^{\;n}$ は観測値、$c^{n}$ はラベル、$W,\vec{b}$ はモデルを調節するパラメータである。上式におけるパラメータは確率変数とは無関係な量として定義したが、これも確率変数とみなし以下のように事前確率 $p(W,\vec{b})$ を導入することができる。 \begin{equation} P(W,\vec{b})=\prod_{n=1}^{N}\;p(c^{n}\;|\;\vec{x}^{\;n}\;;\;W,\vec{b})\;p(W,\vec{b}) \end{equation} すなわち、パラメータが確率 $p(W,\vec{b})$ で与えられたときの観測値の実現確率を考えるのである。 両辺の対数をとりマイナスを付ける。 \begin{eqnarray} L(W,\vec{b})&=&-\ln{P(W,\vec{b})}\nonumber \\ &=&-\sum_{n=1}^{N}\;\ln{p(c^{n}\;|\;\vec{x}^{\;n}\;;\;W,\vec{b})}-\ln{p(W,\vec{b})} \label{nll} \end{eqnarray} これまで考えてきたNLLに $-\ln{p(W,\vec{b})}$ が付加されることが分る。これが正則化項である。 $p(W,\vec{b})$ として以下の形の正規分布を考えることが多い。 \begin{equation} p(W,\vec{b}) = \exp{\left(-\eta \sum_{i,j}\;|w_{i,j}|^2\right)} \end{equation} ここで、$w_{i,j}$ は $W$ の成分、$\eta$ は正の量である。このとき $L(W,\vec{b})$ は以下のようになる。 \begin{equation} L(W,\vec{b})=-\sum_{n=1}^{N}\;\ln{p(c^{n}\;|\;\vec{x}^{\;n}\;;\;W,\vec{b})}+\eta \sum_{i,j}\;|w_{i,j}|^2 \end{equation} 上式を最小化する過程で、右辺第2項は $w_{i,j}$ を小さくする力となる。すなわち、正則化項はNeural Networkのユニット間の結びつきを 弱めるものであり、これが過学習を抑制するのである。上記のように $|w_{i,j}|$ の2乗のものを $L_2$ 正則化項、1乗のものを $L_1$ 正則化項と呼ぶ。

Convolutional Neural Network

はじめに

　先のページでは多層パーセプトロンを解説した。今回は、これを拡張することによりConvolutional Neural Networkを導く。

多層パーセプトロン

　隠れ層が一層の多層パーセプトロンは以下のような構造であった。

これに対応する式は以下のように与えられる。 \begin{eqnarray} y_k &=& \pi\left(\;f_k\left( \vec{z}; W^{(2)},\vec{b}^{\;(2)}\right)\;\right) \label{eprob-2}\\ \pi(a_k) &=& \frac{\exp{(a_k)}} {\sum_{j=1}^{K}\;\exp{(a_j)}} \label{softmax-1}\\ \vec{f}\left( \vec{z}; W^{(2)},\vec{b}^{\;(2)}\right) &=& W^{(2)\;T}\;\vec{z} + \vec{b}^{\;(2)}\label{projection-2} \end{eqnarray} \begin{eqnarray} z_h &=& \tanh\left(\;g_h\left( \vec{x}; W^{(1)},\vec{b}^{\;(1)}\right)\;\right) \label{eprob-1}\\ \vec{g}\left( \vec{x}; W^{(1)},\vec{b}^{\;(1)}\right) &=& W^{(1)\;T}\;\vec{x} + \vec{b}^{\;(1)}\label{projection-1} \end{eqnarray} ここで、$W^{(1)}$ は $D\times H$ 行列、$\vec{b}^{\;(1)}$ は $H$ 次元ベクトル、$W^{(2)}$ は $H\times K$ 行列、$\vec{b}^{\;(2)}$ は $K$ 次元ベクトルである。ただし、隠れ層の活性化関数として $\tanh{}$ を採用している。 最下層は入力層であり $D$ 個のユニットから構成される。隠れ層のユニット数は $H$ 個、最上層は出力層であり $K$ 個のユニットから構成される。これを適用する具体例として、$D$ 次元のベクトルとして与えられる観測値を $K$ 個のクラスのいずれかに分類する多分類問題を挙げることができる。
　上図における層間の接続に注目したものが下図である。各ユニットは、下層の全ユニットと接続している。

Convolutional Neural Network(CNN)

　いま、入力層(input layer)と隠れ層(hidden layer)の間に複数の層を挿入する。

上図は、convolution layerと pooling layerの対を2つ挿入したものである。 一般に、この対を任意の数だけ挿入したものをconvolutional neural networkと呼ぶ。以下各層について説明する。

① convolution layer

　多次元ロジスティック回帰や多層パーセプトロンでは、1次元のベクトルを入力値としてきたが、ここでは入力値として2次元ベクトルを考える。具体的には画像が入力値となる場合を考察する。まず最初に、$m$ 番目の入力画像の $l$ 番目のチャンネルを考え、位置 $(i,j)$ にある画素 $x_{m,l,i,j}$ について以下の計算を行なう。 \begin{equation} g_{m,k,i,j} = \sum_{p,q,l} \; w_{k,l,p,q}\;x_{m,l,i+p,j+q} \label{eq1} \end{equation} これは、矩形のフィルタ $w_{k,l}$ を使って convolution を行い、$l$ について和をとることと等価である。カラー画像なら $l$ の総数は3、グレースケールなら1である。convolution を $\circledast$ で表すことにすれば、式(\ref{eq1})は \begin{equation} g_{m,k} = \sum_{l} \; w_{k,l} \circledast x_{m,l} \label{eq2} \end{equation} と書くことができる。

② max pooling

　次に行なう処理は pooling である。pooling の方法にはいくつかあるが、ここでは max pooling を取り上げる。max pooling とは、 画像を $a\times a$ の領域に細分し、各領域内の最大値でその領域を代表させる処理である。$a=2$ なら $w\times h$の画像は$w/2\times h/2$ のサイズにダウンサンプリングされる。いま形式的に max pooling の処理を $p_{\rm{pooling}}(\cdot)$ と書くことにすれば、この処理の出力値は \begin{equation} q_{m,k} = p_{\rm{pooling}}\left(\sum_{l} \; w_{k,l} \circledast x_{m,l}\right) \label{eq3} \end{equation} となる。

③ 活性化関数

　ここまでの出力値 $q_{m,k}$ にスカラー量であるバイアス項 $b_{k}$ を加え、$\tanh$ で活性化する。 \begin{equation} r_{m,k} = \tanh{(q_{m,k}+b_k)} \label{eq4} \end{equation} 出力値 $r_{m,k}$ は、feature mapと呼ばれる。これを新たな入力値 $x_{m,l}$ とみなして ①②③ の処理を繰り返す。その際、フィルタの数 $k$ を変えることができる。たとえば、最初の処理①②③では20種類のフィルタを使い、次の処理①②③では50種類のフィルタを使うなどである。

　convolution layerとpooling layerの対が2つ、最初の画像サイズが28$\times$28、フィルタサイズが5$\times$5の場合を考える。最初のconvolution layerにより画像サイズは 24$\times$ 24となり、それに引き続く pooling layer により画像サイズは12$\times$12となる。次のconvolution layerへの入力画像のサイズは12$\times$ 12となり、フィルタサイズが最初と同じなら、convolution layerにより画像サイズは8$\times$8、pooling layerにより画像サイズは4$\times$4となる。2つ目のcovolution に使うフィルタの数が50なら、4$\times$4の画像が50個出来上がる。すなわち、1枚の画像（サイズは28$\times$28）から50組の16（=$4\times 4$）画素が作成される。この画素を一列に並べた800（$=50\times 16$）次元のベクトルを多層パーセプトロンのhidden layerの入力とみなす。あとの処理は先のページと同じである。

計算を行なう際の留意点

• CNNの計算量はMLPと比べ著しく多くなる。 私のiMac（3.06GHz Intel Core 2 Duo, 4GB）で、参考文献にあるtheanoコードを実行したところ、MLPで8.32時間、CNNで20.99時間となった。0-1lossの値は、MLPで1.68%、CNNで0.91%であった。
• 各層の計算時間を同じにするため、(feature mapの数$\times$feature mapを構成する画素の数)を各層でほぼ同じにすることがよく行なわれる。
• 考察するデータセットごとに最適なフィルタサイズは変わる。入力画像が28$\times$28程度であれば5$\times$5のフィルタがよく使われる。一方、数百画素$\times$数百画素の画像に対しては12$\times$12や15$\times$15程度のフィルタを用いることが多い。
• max poolingを行なわない選択肢もある。行なうとしても、そのサイズは2$\times$2が限度である。これ以上のdownsamplingを行なうと、多くの情報が失われてしまう。
• whitening（白色化）を導入することにより識別率は高くなる。
• 確率的勾配降下法の各イテレーションごとに学習率を小さくすることにより収束速度を上げることができる。

まとめ

　多次元ロジスティック回帰は入力層と出力層の2層からなるNeural Network、多層パーセプトロンはこの2層の間に隠れ層を複数個挿入したNeural Network、そして、Convolutional Neural Networkは、多層パーセプトロンにconvolution層とpooling層の対を複数個挿入したNeural Networkである。多分類問題に適用した際の識別率は、一般的には、入力層と出力層の間の中間層が増える程高くなる。
　

参考文献

Convolutional Neural Networks (LeNet)
Theanoを使った実装が掲載されている。

多層パーセプトロン

はじめに

　先のページでは多次元ロジスティック回帰について考察した。ここでは、これを拡張することにより、多層パーセプトロンを導入する。

多次元ロジスティック回帰

　先のページでは、観測値として $D$ 次元のベクトル $\vec{x}$ を考え、これが $K$ 個のクラスのいずれかに属する多分類問題を考えた。式を書き下ろすと、 \begin{eqnarray} y_k &=& \pi\left(\;f_k\left( \vec{x}; W,\vec{b}\right)\;\right) \label{eprob}\\ \pi(a_k) &=& \frac{\exp{(a_k)}} {\sum_{j=1}^{K}\;\exp{(a_j)}} \label{softmax}\\ \vec{f}\left( \vec{x}; W,\vec{b}\right) &=& W^{T}\;\vec{x} + \vec{b}\label{projection} \end{eqnarray} となる。ただし、$(W,\vec{b})$ はモデルを記述するパラメータであり、前者は $D\times K$ 行列、後者は $K$ 次元ベクトルである。 これを図で示したものが下図である。

この図から、多次元ロジスティック回帰は入力層と出力層の2層からなるNeural Networkであることが分る。

多層パーセプトロン

　上記の構造に、中間層（隠れ層と呼ぶ）を挿入する。

これに対応する式は以下のように与えられる。 \begin{eqnarray} y_k &=& \pi\left(\;f_k\left( \vec{z}; W^{(2)},\vec{b}^{\;(2)}\right)\;\right) \label{eprob-2}\\ \pi(a_k) &=& \frac{\exp{(a_k)}} {\sum_{j=1}^{K}\;\exp{(a_j)}} \label{softmax-1}\\ \vec{f}\left( \vec{z}; W^{(2)},\vec{b}^{\;(2)}\right) &=& W^{(2)\;T}\;\vec{z} + \vec{b}^{\;(2)}\label{projection-2} \end{eqnarray} \begin{eqnarray} z_h &=& \tanh\left(\;g_h\left( \vec{x}; W^{(1)},\vec{b}^{\;(1)}\right)\;\right) \label{eprob-1}\\ \vec{g}\left( \vec{x}; W^{(1)},\vec{b}^{\;(1)}\right) &=& W^{(1)\;T}\;\vec{x} + \vec{b}^{\;(1)}\label{projection-1} \end{eqnarray} ここで、$W^{(1)}$ は $D\times H$ 行列、$\vec{b}^{\;(1)}$ は $H$ 次元ベクトル、$W^{(2)}$ は $H\times K$ 行列、$\vec{b}^{\;(2)}$ は $K$ 次元ベクトルである。隠れ層の活性化関数として、$\tanh{}$ を採用したが、他にsigmoid関数を考えることもできる。一般に、隠れ層の数を任意にした構造を多層パーセプトロンと呼ぶ。多次元ロジスティック回帰の場合と同様に、確率的勾配降下法を用いて解くことができ、多次元ロジスティック回帰よりも精度の良い識別率を実現できる。

計算を行なう際の留意点

• 勾配降下法では、大域的最小値を求めることはできない。
• 活性化関数としてよく用いられる関数は、$\tanh$とsigmoidの2つである。原点について対称な関数が好まれる理由は、次の層への入力値の平均値を0に近い値にできるからである。経験的にsigmoidよりも$\tanh$の方が収束性が良い。
• 行列$W$の初期値は、0近傍の十分小さな値とする。0近傍であれば、活性化関数の線形部分が適用されるので、勾配を大きくすることができる。
• 活性化関数を$\tanh$とするときは$[-\frac{\sqrt{6}}{\sqrt{N_{\rm in}+N_{\rm out}}}, \frac{\sqrt{6}}{\sqrt{N_{\rm in}+N_{\rm out}}}]$の範囲で、sigmoidを活性化関数とするときは$[-4\frac{\sqrt{6}}{\sqrt{N_{\rm in}+N_{\rm out}}}, 4\frac{\sqrt{6}}{\sqrt{N_{\rm in}+N_{\rm out}}}]$の範囲で一様乱数を発生させ、$W$を初期化すればよい。ここで、$N_{\rm in}$は入力層にあるユニットの数、$N_{\rm out}$は隠れ層にあるユニットの数である。

参考文献

Multilayer Perceptron
Theanoを使った実装が掲載されている。

多次元ロジスティック回帰（Softmax 回帰）

はじめに

　先のページでは1次元ロジスティック回帰について考察した。ここでは、多次元ロジスティック回帰を取り上げる。

多次元ロジスティック回帰

　多分類問題を考える。$N$ 個の $D$ 次元ベクトル $\vec{x}^{\;n},n=1,\cdots,N$ が与えられ、それぞれが ${1,\cdots,K}$ のいずれかのラベルを持つとする。このとき、観測値が従う確率分布モデル $p(c=1,\cdots,K\;|\;\vec{x}\;;\;W,\vec{b})$ を求めたい。ここで、$W,\vec{b}$ はモデルを調節するパラメータである。$W$ は $D \times K$ 行列、$\vec{b}$ は $K$ 次元ベクトルである。各観測値が他の観測値の影響を受けないと仮定すれば、全ての観測値の確率は次式で与えられる。 \begin{equation} L(W,\vec{b})=\prod_{n=1}^{N}\;p(c^{n}\;|\;\vec{x}^{\;n}\;;\;W,\vec{b}) \label{likelihood} \end{equation} これが最大となるように $(W,\vec{b})$ を調節すればよい。計算を容易にするため次式の最小化を考える。 \begin{equation} L(W,\vec{b})=-\sum_{n=1}^{N}\;\ln{p(c^{n}\;|\;\vec{x}^{\;n}\;;\;W,\vec{b})} \label{nll} \end{equation} 上式を負の対数尤度（Negative Log Likelihood:NLL）関数と呼ぶ。確率 $p$ は次式で定義される。 \begin{equation} p(c\;|\;\vec{x}\;;\;W,\vec{b})=\prod_{k=1}^{K}\;(y_k)^{t_k} \label{prob} \end{equation} ここで、$y_k$ は $K$ 次元ベクトル $\vec{y}$ の $k$ 番目の成分、$t_k$ は $K$ 次元ベクトル $\vec{t}$ の $k$ 番目の成分である。$y_k$ は次式で定義される。 \begin{eqnarray} y_k &=& \pi\left(\;f_k\left( \vec{x}; W,\vec{b}\right)\;\right) \label{eprob}\\ \pi(a_k) &=& \frac{\exp{(a_k)}} {\sum_{j=1}^{K}\;\exp{(a_j)}} \label{softmax}\\ \vec{f}\left( \vec{x}; W,\vec{b}\right) &=& W^{T}\;\vec{x} + \vec{b}\label{projection} \end{eqnarray} 行列やベクトルの肩にある $T$ は転置を表す。$\vec{f}$ は $K$ 次元ベクトルであり、$f_k$ はその $k$ 番目の成分を表す。 $\vec{t}$ は、観測値 $\vec{x}$ が持つラベルに相当する成分を1、それ以外の成分を0とする $K$ 次元ベクトルである。 たとえば、ラベル $3$ を持つなら \begin{equation} \vec{t} = (0,0,1,0,\cdots,0)^T \end{equation} となる。式(\ref{softmax})をsoftmax関数と呼ぶ。 $a_k={0,\cdots,9}$としたときのsoftmax関数を下に示す。

次式が成り立つことに注意する。 \begin{equation} \sum_{k=1}^{K}\;\pi(a_k)=1 \end{equation} すなわち、 \begin{equation} \sum_{k=1}^{K}\;y_k=1 \end{equation} ということであり、観測値 $\vec{x}$ はいずれかのラベルに必ず属するということである。 式(\ref{nll})に式(\ref{prob})を代入すると次式を得る。 \begin{equation} L(W,\vec{b})=-\sum_{n=1}^{N}\;\sum_{k=1}^{K} t_k^n\;\ln{(y_k^n)} \end{equation} ここまでの議論をまとめる。

1. 式(\ref{projection})を使って、$D$次元ベクトル（観測値）を $K$ 次元ベクトルに射影する。
2. $K$ 次元ベクトルの各成分を使って確率(\ref{softmax})を導入する。
3. 確率的勾配降下法を用いて対数尤度を最小化し、最適なモデルを決定する。

計算を容易にするため、 $x^n_d,y^n_k,t^n_k$ を $x_{dn},y_{kn},t_{kn}$ に置き換え、これを行列の成分とみなすことにする。 \begin{equation} L(W,\vec{b})=-\sum_{n=1}^{N}\;\sum_{k=1}^{K} t_{kn}\;\ln{(y_{kn})} \label{nllspecial} \end{equation}

更新式

確率的勾配降下法については先のページで説明したので、ここでは更新式だけを示す。
　式(\ref{nllspecial})を $W$ の成分 $w_{dk}$ で微分する。 \begin{eqnarray} \frac{\partial L}{\partial w_{dk}} &=& -\sum_{n=1}^{N}\;\sum_{s=1}^{K}\;t_{sn}\frac{1}{y_{sn}}\frac{\partial y_{sn}}{\partial a_{sn}}\frac{\partial a_{sn}}{\partial w_{dk}} \end{eqnarray} ここで、 \begin{eqnarray} \frac{\partial y_{sn}}{\partial a_{sn}}&=&y_{sn}(1-y_{sn}) \\ \frac{\partial a_{sn}}{\partial w_{dk}}&=&x_{dn}\;\delta_{ks} \end{eqnarray} が成り立つので、これらを代入すると、 \begin{eqnarray} \frac{\partial L}{\partial w_{dk}} &=& \sum_{n=1}^{N}\;t_{kn}(y_{kn}-1)\;x_{dn} \\ &=& \sum_{n=1}^{N}\;(y_{kn}-t_{kn})\;x_{dn} \end{eqnarray} を得る。ただし、$t_{kn}\;y_{kn}=y_{kn}$ を用いた。同様に $b_k$ で微分すると次式を得る。 \begin{equation} \frac{\partial L}{\partial b_{k}} = \sum_{n=1}^{N}\;(y_{kn}-t_{kn}) \end{equation} いま、$e_{kn}=y_{kn}-t_{kn}$ と置くと \begin{eqnarray} \frac{\partial L}{\partial w_{dk}} &=& \sum_{n=1}^{N}\; x_{dn}\;e_{kn}\; \\ &=& \left(X\;E^T\right)_{dk} \label{update1}\\ \frac{\partial L}{\partial b_{k}} &=& \sum_{n=1}^{N}\;e_{kn} \label{update2} \end{eqnarray} となる。ここで $E$ は $e_{kn}$ を成分とする行列、$X$ は $x_{dn}$ を成分とする行列である。これらより更新式は以下のようになる。 \begin{eqnarray} w_{dk} &\leftarrow& w_{dk} - \eta \; \left(X\;E^T\right)_{dk} \\ b_{k} &\leftarrow& b_{k} - \eta \; \sum_{n=1}^{N}\;e_{kn} \end{eqnarray} この更新式を使って、確率的勾配降下法を実行する。以下に各行列を書き下ろしておく。 \begin{eqnarray} Y &=& \left(\begin{array}{cccc} y_{11} & y_{12} & \cdots & y_{1N} \\ y_{21} & y_{22} & \cdots & y_{2N} \\ \vdots & \vdots & \ddots & \vdots \\ y_{K1} & y_{K2} & \cdots & y_{KN} \end{array} \right) = {\bf \pi}(W^T\;X+B)\label{pi}\\ \nonumber\\[3pt] W &=& \left(\begin{array}{cccc} w_{11} & w_{12} & \cdots & w_{1K} \\ w_{21} & w_{22} & \cdots & w_{2K} \\ \vdots & \vdots & \ddots & \vdots \\ w_{D1} & w_{D2} & \cdots & w_{DK} \end{array} \right)\\ \nonumber\\[3pt] X &=& \left(\begin{array}{cccc} x_{11} & x_{12} & \cdots & x_{1N} \\ x_{21} & x_{22} & \cdots & x_{2N} \\ \vdots & \vdots & \ddots & \vdots \\ x_{D1} & x_{D2} & \cdots & x_{DN} \end{array} \right)\\ \nonumber\\[3pt] B &=& \left(\begin{array}{cccc} b_{1} & b_{1} & \cdots & b_{1} \\ b_{2} & b_{2} & \cdots & b_{2} \\ \vdots & \vdots & \ddots & \vdots \\ b_{K} & b_{K} & \cdots & b_{K} \end{array} \right) \\ \nonumber\\[3pt] T &=& \left(\begin{array}{cccc} t_{11} & t_{12} & \cdots & t_{1N} \\ t_{21} & t_{22} & \cdots & t_{2N} \\ \vdots & \vdots & \ddots & \vdots \\ t_{K1} & t_{K2} & \cdots & t_{KN} \end{array} \right) \\[3pt] E &=& Y-T \end{eqnarray} $B$ は同じ列から構成される $K\times N$ 行列である。 式(\ref{pi})の関数 $\pi$ は先に定義したsoftmax関数であり、行列の各成分に作用させることとする。

実装

pythonを使って実装を行なった。 $K=3$、$D=2$、$N=1500$とした。 すなわち、2次元平面上の3分類問題である。1500個の点を3等分し、各集合に1,2,3のラベルを割り振った。

1. 8-9行目：式(\ref{softmax})である。
2. 12-13行目：式(\ref{eprob})である。
3. 16-21行目：式(\ref{update1}),(\ref{update2})である。
4. 24-28行目：誤差を描画する関数である。
5. 30-45行目：境界線を描画する関数である。
6. 47-57行目：観測値を作る関数である。
7. 59-64行目：ラベルを作る関数だる。
8. 67-104行目：確率的勾配降下法を実行するコードである。

下図は、次式で定義される3つの量を、横軸に試行回数をとり、プロットしたものである。 \begin{eqnarray} \bar{e}_1&=&\frac{1}{\left|\Lambda\right|}\sum_{n\in\Lambda}\;e_{1n}\;\;(\mbox{赤線})\\ \bar{e}_2&=&\frac{1}{\left|\Lambda\right|}\sum_{n\in\Lambda}\;e_{2n}\;\;(\mbox{緑線})\\ \bar{e}_3&=&\frac{1}{\left|\Lambda\right|}\sum_{n\in\Lambda}\;e_{3n}\;\;(\mbox{青線})\\ \end{eqnarray} ここで、$\Lambda$ はミニバッチを表し、$|\Lambda|$ はそこに含まれる点の数を表す。

$K=3$の場合、$(W,\vec{b})$ は以下のようになる。 \begin{eqnarray} W &=& \left(\begin{array}{ccc} w_{11} & w_{12} & w_{13} \\ w_{21} & w_{22} & w_{23} \end{array} \right)\\ \nonumber\\[3pt] \vec{b} &=& \left(\begin{array}{c} b_{1} \\ b_{2} \\ b_{3} \end{array} \right) \end{eqnarray} また、計算される $\vec{f}$ は次式のようになる。 \begin{eqnarray} \vec{f} &=& \left(\begin{array}{c} f_{1} \\ f_{2} \\ f_{3} \end{array} \right) \\ f_1 &=& \vec{w}^{\;T}_1\cdot\vec{x}+\vec{b}_1 \\ f_2 &=& \vec{w}^{\;T}_2\cdot\vec{x}+\vec{b}_2 \\ f_3 &=& \vec{w}^{\;T}_3\cdot\vec{x}+\vec{b}_3 \\ \vec{w}_1 &=& \left(\begin{array}{c} w_{11} \\ w_{21} \end{array} \right) \\ \vec{w}_2 &=& \left(\begin{array}{c} w_{12} \\ w_{22} \end{array} \right) \\ \vec{w}_3 &=& \left(\begin{array}{c} w_{13} \\ w_{23} \end{array} \right) \end{eqnarray} いま、各点集合の境界線を考える。 たとえば、ラベル1を持つ点集合とラベル2を持つ点集合の境界線上では、ラベル1を持つ確率 $p_1$ とラベル2を持つ確率 $p_2$ が等しいはずである。 式(\ref{eprob})より \begin{eqnarray} p_1 &=& \frac{\exp{(f_1)}}{\exp{(f_1)}+\exp{(f_2)}+\exp{(f_3)}} \\ p_2 &=& \frac{\exp{(f_2)}}{\exp{(f_1)}+\exp{(f_2)}+\exp{(f_3)}} \\ \end{eqnarray} であるから、境界線上では、$f_1=f_2$ が成り立つことになる。 従って、この2つの点集合を分ける境界線は次式で記述される。 \begin{equation} (\vec{w}^{\;T}_1- \vec{w}^{\;T}_2)\cdot\vec{x}+(\vec{b}_1-\vec{b}_2)=0 \end{equation} 同様にラベル2と3を分ける境界線は \begin{equation} (\vec{w}^{\;T}_2- \vec{w}^{\;T}_3)\cdot\vec{x}+(\vec{b}_2-\vec{b}_3)=0 \end{equation} ラベル1と3を分ける境界線は \begin{equation} (\vec{w}^{\;T}_1- \vec{w}^{\;T}_3)\cdot\vec{x}+(\vec{b}_1-\vec{b}_3)=0 \end{equation} である。 これらを図示したものが下図である。ラベル1を持つ点を赤、ラベル2を持つ点を緑、ラベル3を持つ点を青で表している。 青線がラベル1と2の境界線、緑線がラベル2と3の境界線、赤線がラベル1と3の境界線である。

多次元ロジスティック回帰におけるNLLの最小化はconvex optimizationとなるので、初期値に関わらず大域的最小値に収束する。

pcl-1.7.2をmacでコンパイル

前回、pcl-1.7.1をコンパイルしたときは、 いろいろエラーが出ましたが、pcl-1.7.2は何の問題もなくコンパイルできました。

ここからソースをダウンロードして、解凍。 依存ライブラリは既にmacportsでインストールされている。