nykergoto’s blog

機械学習とpythonをメインに、雑多な内容をとりとめなく扱うブログです。

ニューラルネットにおける変数の初期化について

最近 keras をつかって色々とやることにはまってます。 keras だとネットワーク記述から fit までが完結に記述できてとても気に入っています。 そのときにドキュメントや実装を読んだりもしますが、ネットで自分がやりたいタスクで検索をかけて似たようなことをやっている人のコードを参考にしています。 その時、keras レイヤーの重みの初期化の引数に he_normal というのを指定しているのを見かけました。 ドキュメントを参照すればこの関数の定義式自体はわかりましたが、この関数が何をやってるのか、なぜこの関数にするとうれしいのか、平均0分散1のガウス分布ではだめなのか、ということがわからず困ってしまいました。

ニューラルネットワークにおける変数の初期化の方法について全く知らなかったので、ちょうどいい機会とおもい各種論文を読んでみたので、その内容をまとめていきます。

黎明期

黎明期、という言い方はあれかもしれませんが、まだ relu も登場していない時代の話から始めます。 この頃は、ヒューリスティックに以下のような分布が一般に用いられていたようです。

$$ W_{j, k}^{i} \sim U \left[ - \frac{1}{\sqrt{n_i}}, \frac{1}{\sqrt{n_i}} \right] $$

ここで $W_{j,k}$ は 第 $i$ レイヤーの重みの $(j, k)$ 成分、 $U[-a, a]$ は下限 $-a$, 上限 $a$ の一様分布, $n_i$ は第 $i$ 番目の隠れ層の数を表しています。またバイアス項はすべて 0 とします。

この関数をどのぐらい積極的に使っていたのかは、その当時に僕が機械学習をやっていなかったのでわかりませんが、んとなくこういう分布が経験的によいと知られていたのでしょうか。

この経験則に対し、活性化関数が linear なときという仮定のもとで議論を行い、新しい重みの分布の提案を行った論文が Understanding the difficulty of training deep feedforward neural networks で、これは今では一般に Glorot の一様分布と呼ばれているものです。

Glorot の一様分布 (Glorot Uniform)

Glorot の一様分布と呼ばれている初期化の方法では、第 $i$ 番目のレイヤーに対して以下のような分布から重みを選びます。

$$ W \sim U \left[ -\frac{ \sqrt{6}}{\sqrt{n_i + n_{i+1}} }, \frac{ \sqrt{6}}{\sqrt{n_i + n_{i+1}} } \right] $$

ここで $n_i$ は第 $i$ 番目の隠れ層の数を表します。 これだけ言われても「その6どこからきたの?」とかいろいろ聞きたいことがありますよね。 以下では提案論文を元に、なぜこの分布が導出されたかを追っていきましょう。

導出

ニューラルネットワーク上の第 $i$ 番目のレイヤーについて考えます。 第 $i$ レイヤーの活性化関数を $f$, レイヤーへの入力ベクトルを $x^i \in {\mathbb R^{n_i}}$, 活性化関数を通す前の出力ベクトルを $y^i \in \mathbb R^{n_{i+1}}$, レイヤーの重みを $W^i \in \mathbb R^{n_{i+1}, n_{i}}$ とします。ここで $n_i \in \mathbb N$ は第 $i$ レイヤーの入力次元を表しています。 以上の定義を用いると、以下の関係が成り立ちます。

$$ y^{i} = W^{i} x^i + b^i \\ x^{i+1} = f(y^i) $$

つぎにネットワークの学習について考えたいので、backpropagation で必要とされるネットワークの誤差 $C$ に対する偏微分を考えていきます。すると以下を得ます。

$$ \frac{\partial C}{\partial y_j^i} = f'(y_k^i) W_j^{i+1} \frac {\partial C}{\partial y^{i+1}} \\ \frac{\partial C}{\partial w_{j, k}^i} = x_k^i \frac{\partial C}{\partial y_j^i} $$

ここで、初期化された段階において, 活性化関数は線形性を持っていて, 重みは独立に初期化されていて, 入力変数はすべて同じ分散 ($Var[x]$) を持つ, という仮定を加えます。

すなわち活性化関数に関しては $f(x) = x, f'(0) = 1$ が成り立っている, ということを意味します(これは yigmoid や tanh など原点対称な活性化関数で成り立ちます)。 また重みに関しては独立に初期化されているので, 第$i$ 層の重みに対して $Var[w^i_{j, k}]$ が一定です。これを以下では $Var[W^i]$ と表記します。

このとき第 $i$ 層の出力の分散 $Var[x^i]$ を考えてみましょう. $$ \begin{align} Var[x^{i+1}] &= Var[f(y^i)]\\ &= Var[y^i] \\ &= Var[W^i x^i + b^i] \\ &= Var[W^i x^i] + Var[b^i] \\ &= n_i Var[W^i] Var [x^i] \end{align} $$

これを繰り返し用いると

$$ Var[x^i] = Var[x] \prod_{l=0}^{i-1} n_{l} Var[W^l] $$

を得ます。 したがって全てで $d$ 層あるネットワークにおいては, 第 $i$ 層の出力に対する微分について以下が成立します。

$$ \begin{align} Var \left[\frac{\partial C} {\partial y^i} \right] &= Var\left[f'(y_k^i) W_j^{i+1} \frac {\partial C}{\partial y^{i+1}} \right] \\ &= n_{i+2} Var\left[W^{i+1} \right] Var \left[\frac {\partial C}{\partial y^{i+1}} \right] \\ &= \left( \prod_{l=i}^{d} n_{l+1} Var[W^l] \right) Var\left[ \frac{\partial C}{\partial y^{d}} \right] \end{align} $$

また重みに対する微分に対して、以下が成立します。

$$ \begin{align} Var \left[ \frac{\partial C}{\partial W^i} \right] &= Var\left[x_k^i \right] \times Var \left[ \frac{\partial C}{\partial y_j^i} \right] \\ &= Var[x] \prod_{l=0}^{i-1} n_{l} Var[W^l] \times \left( \prod_{l=i}^{d} n_{l+1} Var[W^l] \right) Var\left[ \frac{\partial C}{\partial y^{d}} \right] \end{align} $$

今度は、ネットワークの出力の伝搬について考えて行きます。 出力の伝搬のとき、各レイヤーを通してもその出力は無限に発散したり、0に減衰したりするとただしくロスを計算できないため困ります。その為各レイヤーの出力 $x^i$ は同じ分散を持つことが好ましいです。つまり

$$ \forall (i, i'), Var[x^i] = Var[x^{i'}] $$

が成り立っていてほしい、ということになります。 また逆伝搬についても同様に

$$ \forall (i, i'), Var\left[ \frac{\partial C}{\partial y^i}\right] = Var\left[ \frac{\partial C}{\partial y^{i'}}\right] $$

がなりたっている必要があります。 これらを先の方程式に代入すると以下の条件に変形できます。

$$ \forall i, n_i Var[W^i] = 1 \\ \forall i, n_{i+1} Var[W^i] = 1 $$

この2つの条件の中間をとろう、というのが Glorot の提案する分布です。すなわち

$$ \forall i,\ Var\left[W^i \right] = \frac{2}{n_i + n_{i+1}} $$

を満たすような確率分布から重みを初期化します。 ではこのような分散を持つような一様分布とはどんな分布でしょうか、というのに答えるのが、この章の先頭に示した分布です。再掲すると以下の分布でした。

$$ W \sim U \left[ -\frac{ \sqrt{6}}{\sqrt{n_i + n_{i+1}} }, \frac{ \sqrt{6}}{\sqrt{n_i + n_{i+1}} } \right] $$

$[-a, a]$ の範囲の一様分布の分散は $ a^2/3 $ ですから

$$ \frac{1}{3} \left(\frac{ \sqrt{6}}{\sqrt{n_i + n_{i+1}}} \right)^2 = \frac{2}{n_i + n_{i+1}} $$

となって、ピッタリさきの分散に一致するように選ばれています。 一見どこから来たかわからない 6 という数字にも意味があったのですね。

一様分布以外にも, この分散を取り平均が 0 であるガウス分布 (Glorotのガウス分布) も提案されています。

まとめ

Glorot の分布の議論の仮定、特徴をまとめると以下のようになります。

  • 活性化関数が原点対称かつ原点周りで線形であることを仮定して議論
  • 順伝搬, 逆伝搬で各層の値がおなじになる、ということを制約としている
  • $i$ 番目のレイヤーの初期重みを $2 / (n_i + n_{i+1})$ を分散に持つような分布から選ぶ

問題点

きれいに片付いたかに見えた初期化分布問題ですが、Glorot は活性化関数に大きな仮定をおいていました。それは

  • 原点対称である
  • 原点付近で線形である

ということでした。 またネットワークの構造についても全結合のネットワークが対象であり CNN のようなカーネルを用いたレイヤーは想定されていません。 というわけで relu とか CNN で解析してみたよ、というのが次に紹介する He の正規分布と呼ばれているものです。*1

He の正規分布 (He normal)

He の正規分布は第 $i$レイヤーの初期重みを 平均0, 分散 $2/n_i$ の正規分布から設定します。 すなわち

$$ Var[W^i] = \frac{2}{n_i} $$

を満たします。 この分散を持つべき理由を導出していきましょう。

導出

まずネットワークの順伝搬 (Forward Propagation) について考えていきます。 第 $i$ 層の CNNの出力は

$$ y^i = W^i x^i + b^i $$

と書くことができます。 ここで $x \in \mathbb R^{c k^2}$ は $c$ チャネルの $(k, k)$のピクセルを持つ入力です (以下では簡易化のため $n = c k^2$ と記述します)。 また $W \in \mathbb R^{d, n}$ は各行がフィルタの重みとなるような $d$ のフィルタを表します。 $b \in \mathbb R^n$ はバイアス項, $y\in \mathbb R^d$ は出力マップの値です。 $f$ を活性化関数とすると

$$ x^{i+1} = f(y^i) $$

が成り立ちます。 ここでつぎの3つを仮定として加えます。

  • $x^i$ の各要素が独立同分布である
  • $x^i$ と $W^i$ の各要素は独立である。
  • $W^i$ の分布の期待値は 0 である。

すると

$$ \begin{align} Var[y^i] &= Var[W^i x^i + b^i] \\ &= Var[W^i x^i] \\ &= n_i Var[w^i x^i] \\ &= n_i \mathbb E[ (w^i x^i)^2] - n_i \mathbb E [w^i x^i]^2 \\ &= n_i \mathbb E[ (w^i)^2] \mathbb E[ (x^i)^2] - 0 \\ &= n_i Var[w^i] \mathbb E[ (x^i)^2] \end{align} $$

と変形できます。 ここで $\mathbb E [x^i]^2$ は $x^i$ の期待値が 0 でない限り分散に一致しません。 すなわち

$$ \mathbb E[ (x^i)^2] \ne Var[x^i] $$

です。たとえば活性化関数が relu のとき $x^i = \max (0, y^{i-1})$ ですからこれは明らかに期待値が 0 にはなりません.

ここで $w^{i-1}$ が 0 のまわりで対称な分布であったとし、また $b^{i-1} =0$とします。すると $y^{i-1} = W^{i-1} x^{i-1}$ と $W, x$ の独立性の仮定より $y^{i-1}$の各成分は平均 0 の対称な分布となります。 したがって

$$ \mathbb E (x^i)^2 = \mathbb E [\max (0, {y^{i-1}}^2)] = 0 + \frac{1}{2} \mathbb E[(y^{i-1})^2] = \frac{1}{2} Var[y^{i-1}] $$

と変形することができます。 これより

$$ Var[y^i] = n_i Var[w^i] Var[y^{i-1}] $$

を得ます。これを繰り返し用いると第 $L$ 層の出力 $y^L$ の分散は

$$ Var[y^L] = Var[y^1] \left( \prod_{i=2}^L \frac{1}{2} n_i Var[w^i] \right) $$

となります。 この後半の積の部分が 1 でなければ、層が増えていくにしたがって、指数的に無限大に発散もしくは0に収束してしまうことになります。 よってこの部分は 1 にっていることが望ましいです。 したがって

$$ Var[w^i] = \frac{2}{n_i}\ \forall i $$

を満たすような分布を用いて、重みの初期化を行うことが望ましいといえます。 今回は順伝搬についての議論でしたが、逆伝搬に対してもほぼ同様の議論ができ

$$ Var[w^i] = \frac{2}{n^{i+1}}\ \forall i $$

が得られます。 Glorotのときにも同じような結果となったので、まあそうなるよね、といったところでしょうか。 このような分散を満たす分布を keras では he_normal と呼んでいます。 Glorot のときと同様, これの一様分布 varsion の he_uniform もあったりします。

ちなみに、議論中 relu を通した期待値を計算するところで prelu を用いると以下が得られます。

$$ Var[w^i] = \frac{1}{1 + a^2} \frac{2}{n^i}\ \forall i $$

ここで $a$ は prelu の係数です。 prelu は $a=1$ のとき linear な関数に一致しますが, その時 Glorot のときに出てきた式と一致します。 また $a=0$ のとき prelu は relu と一致し, 初期重みの式も一致しています。

比較

さて大きくわけて2つの初期化分布を紹介しました。 これの違いは単純に分散が二倍ちがうだけです。

$$ \begin{align} Glorot&: \frac{2}{n_i + n_{i+1}}\\ He&: \frac{2}{n_i} \end{align} $$

対して違わないようですがこの初期化分布のズレは、層を重ねるごとに掛け算されていきますから、指数的に影響を及ぼします。ですから深いネットワークになればなるほど、この違いが顕著に現れるのです。 というのを実験している結果が以下になります。(参考文献 2 Figure 3 より引用)

f:id:dette:20171015120652p:plain:w600
n_layer = 22 のとき

赤いラインが He の正規化を行ったネットワークでの学習を表しています。 青いラインは凡例では Xavier となっていますがこれは Glorot の初期化の別名です。 これを見ると 22 layers のときは He のほうがよさそうですがまあ甲乙つけがたい感じとなっています。 しかしこれが 30 layers になると差が明らかになります。

f:id:dette:20171015120816p:plain:w600
n_layer = 30 のとき

Glorot の初期化では学習が全く進んでいない一方、 He の初期化では epoch 数はかかっているものの、先と同様に学習ができている事がわかります。 このことから層が深くなると初期化の影響が大きくなるという主張がある程度裏付けられた形となっています。 初期化って適当に選んでいましたが、大事なんですね。

まとめ

ざっくりというと以下の事がわかりました。

  • 初期化アルゴリズム Glorot と He は活性化関数に対する仮定が違う
    • Glorot: 原点周りにおいて線形な活性化関数
    • He: relu
  • ちゃんと初期化アルゴリズムを選ばないと、特に深いネットワークの学習を行う際に痛い目をみる

適当に使うとだめってことがわかったので、調べたかいがあったかな?

参考文献

  1. Understanding the difficulty of training deep feedforward neural networks, Xavier Glorot, Yoshua Bengio ; Proceedings of the Thirteenth International Conference on Artificial Intelligence and Statistics, PMLR 9:249-256, 2010.
  2. Delving Deep into Rectifiers: Surpassing Human-Level Performance on ImageNet Classification, Kaiming He, Xiangyu Zhang, Shaoqing Ren, Jian Sun

*1:https://arxiv.org/pdf/1502.01852.pdf この論文は活性化関数 prelu の提案論文でもあったりします。

本のクラスタリングをやってみよう - 吾輩は猫であるに近い本は何なのか

最近理論よりなことばかりやっていたので今回は実際のデータを使った解析をやってみます。

今回使うデータは、読書メーターからクロールさせてもらって作成した、ユーザーに紐付いた読書履歴のデータです。ユーザーごとに [だれの, どんな] 本を読んだかがわかるようなデータになっています。一例は以下のような感じです。

アガサ・クリスティー   おづ まりこ    トマス・H・クック ムア・ラファティ    川口俊和    ジョナサン・オージエ  村田 沙耶香    岡崎 琢磨   米澤 穂信   ピエール・ルメートル  金内 朋子

この人はミステリーが好きなのかもしれませんね。 先の例は作者でしたが、これと同じように本のタイトルも取得しています。 取得した本の数(累積)は 100万冊, ユーザー数は 2500 と気づいたら案外大きいデータセットになっていました。 このうち今回は本のタイトルのデータを使って、本のクラスタリングをやってみたいと思います。

github.com

方法

さてどうやってクラスタリングをやっていくか、ということですが、基本的なアイディアは「読書の履歴にたいする各本は、文章に対する単語の関係と似ているのではないか」という仮説です。

ユーザーは一定の嗜好パターンを持っていて、それに基づいて読む本を決めているはずです。 これは文章がある意味合い(例えば新聞記事とかラブレターとか)を持っていて、その意味合いによって出現する単語の分布が異なる、ということによく似ています。 このことから、文章と単語の関係から単語の意味表現を得たり、クラスタリングを行うような手法は、読書の履歴から本の特徴を得る場合にも用いれるはずです。はず。

以上の仮説に基づいて文章のクラスタ分類でおなじみの LDA (Latent Dirichlet Allocation) と Word2Vec を使います。

LDA

まずは LDA をやってみます。 LDAは文章中に出てくる単語の回数をカウントして文章の特徴とする BoW (Bag of Words) を用いるモデルです。

例えば、僕が文章を書くことを考えてみます。 僕には文章の嗜好や知識の偏りがあるので、野球や数学の話題は出て来るでしょうが、裁縫に関する話題はあまり出てこないように思われます。 このような話題のことをトピックと呼んでいます。 トピックモデルでは「文章には何らかの話題の偏りがある」ということを仮定します。

そして文章中に現れる単語は、文章に紐付けられたトピックのどれかから生成されるものと考えます。すなわち、すべての単語はそれに紐付いたトピックを持っている、と考えます。 またトピックは、どの単語がよく現れるか、また現れにくいかという分布を持っていて、この分布はトピックごとに異なっています。

例えば仮に野球を表すであれば バット, ヒット といった単語の確率が高くなっていて、反対に シュート という単語は低い、といった感じです。 このトピックごとの単語のかたよりも学習していきます。

トピックモデルはグラフィカルモデルの有用さを示すいい例なので数式をたどると面白いかもしれません。僕はベイズとグラフィカルモデル大好きなので、LDA大好きです。 こんな雑な説明でもこれから勉強したい!となった方は、ぜひさとういっせい先生の本をおすすめします。とても丁寧に式を追っていて、数式の気持ちというか「なぜこれをやってるんだろう」を丁寧に説明してくれていて、数式の途中でまいごになることが少なく、おすすめです。

www.amazon.co.jp

結局 LDA をして嬉しいのは何かというか「文章がどのような偏りを持ったトピックで生成されたのか」という情報を知ることが出来るということです。 この情報は先の 野球トピック のようなものに相当しています。 すなわちこちらから何も情報を与えなくても(文章だけ与えれば) バットとヒットは同じ文脈で出てきやすいんだな、というようなクラスの情報を得ることが出来るわけで、それはとてもハッピーですね。

というわけで早速やってみましょう。 feature.py の関数を用いてデータをロードした後に, 10個のトピックに分類してみます。

from features import load_feature, compile_corpus_and_dictionary
from lda import LDATrainer

data = load_feature(root_path="./data/raw/", target="title")
corpus, dictionary = compile_corpus_and_dictionary(data, no_below=5, save_to="checkpoints/dictionary.txt")

lda_trainer = LDATrainer(output_dir="./checkpoints/")
lda_trainer.run(corpus, dictionary, topics=10)  # 学習の実行

# 上位のトピックを保存
lda_trainer.show_topics(save_to_dir="./checkpoints/")

出てきたトピックの分布は以下のようになりました。 10個のトピックをすべて載せてもよいのですが、ちょっと見にくいのと、あとでも述べますがあまりトピックごとに差分が見えなかったので2つだけ載せています。

topic 1

[('阪急電車(幻冬舎文庫)', 0.0013286450514788421),
 ('和菓子のアン(光文社文庫)', 0.0010029719922129846),
 ('永遠の0(講談社文庫)', 0.00095220966510543164),
 ('氷菓(角川文庫)', 0.000931927887614856),
 ('火花', 0.00088232395541623144),
 ('レインツリーの国(新潮文庫)', 0.00087065130728654241),
 ('ビブリア古書堂の事件手帖3~栞子さんと消え…', 0.00083625459565052221),
 ('陽だまりの彼女(新潮文庫)', 0.00078834353368666273),
 ('夜は短し歩けよ乙女(角川文庫)', 0.00076325430885853943),
 ('ぼくは明日、昨日のきみとデートする(宝島社…', 0.00076090024369739853)]

topic2

[('永遠の0(講談社文庫)', 0.0012030737598141335),
 ('阪急電車(幻冬舎文庫)', 0.0011685757490659541),
 ('イニシエーション・ラブ(文春文庫)', 0.0011032171029623371),
 ('舟を編む', 0.00093492775265262469),
 ('西の魔女が死んだ(新潮文庫)', 0.00085344844568903389),
 ('夜のピクニック(新潮文庫)', 0.00084074342930131642),
 ('図書館戦争図書館戦争シリーズ(1)(角…', 0.00079957001036084336),
 ('カラフル(文春文庫)', 0.00078677125722370774),
 ('死神の精度(文春文庫)', 0.00077641049565474061),
 ('ぼくは明日、昨日のきみとデートする(宝島社…', 0.00076590742659798081)]

なんだかどちらのトピックも大差ない内容になってしまっていますね。 本のタイトルをクラスタに分けたいというモチベーションだったので、その意味でこれは失敗になってしまいました。 その後も色々と試行錯誤はしてみた*1のですが、あまり良い結果が得られなかったので LDA は宿題に回して Word2Vec の方に移ります。

Word2Vec

次は word2vec をやってみます。

word2vec は単語の分散表現 (word embedding) を獲得手法と表現されることが多いように、単語が何の意味を持っているのかをベクトルで表したい、ということがモチベーションの手法です。 *2

word2vec には Continuous Bag of Words (CBOW) と skip-gram の二種類がありますが、発想はほとんど同じです。 CBOW では周辺の単語から中央の言葉を予測していて、 skip-gram では中心の単語から周辺の単語の予測をします。そしてどちらの手法もその予測器は隠れ層1の auto-encoder です。 数式にすると中心の単語が入力になる skip-gram の目的関数は以下のようなります。

$$ {\min_{W_1, W_2}} {\rm Loss}(W_1, W_2) = \sum_{i \in {\rm train}} \sum_{j \in {\rm surround(i)}} {\rm CrossEntropy} ( W_2^{T} ( W_1 x_i ), y_j) $$

このとき $x, y \in \mathbb{R}^n$ は単語を表す one-hot-vector で $n$ は文章中の語彙数を表しています。また $W_1, W_2 \in {\mathbb R}^{m \times n}$ は隠れ層の次元数が $m$ の auto-encoder のパラメータです。

CBOW の場合には入力が周辺の単語になりますから、上の数式で $x_i$ と $y_j$ が入れ替わったような形になります。

こうして得られた auto-encoder の隠れ層の部分、すなわちある単語 $x$ にたいして $W_1 x$ で計算される $m$ 次元ベクトルが、単語の情報をより低次元で表現できるベクトル(これを単語の分散表現と呼びます)になっていると解釈できます。

現在はそれにいろいろな改良やアイディアが付け加えられています。 例えば、ロス関数の計算において発生する softmax の計算が語彙の数 $n$ のオーダーで必要なことが計算のボトルネックになっているとして、すべての単語に対する計算をさぼって、近似的に softmax を計算してしまおうというネガティブサンプリングという手法もそのひとつです。 これにより最適化における大幅な速度向上を達成しました。 また、word2vecの入力に文章のidを加えて、文章も分類してしまおうという方法 (doc2vec) も提案されています。

python で word2vec をやりたいときは chainer や keras などのディープなフレームワークを使ってもできますが内部が C で実装されていて速度の早い gensim を用いるのが良いと思います。 特にGPU環境でなければ。(一回 word2vec を実装した時には cpu しかなく全く使い物になりませんでした)

このライブラリ関数やクラスの定義の仕方に若干の癖があるという難点はありますが cpu 環境でも十分実用的な時間に計算できるという点はとても魅力です。 では実際に計算をやってみます。 今回は僕が作成した訓練のための V2WTrainer クラスを用いていますが, 実装見てもらうとわかるように内部で大したことはやっていないので直接 gensim を使うときもこんなのりで出来るんだーという感じで見てもらえれば。

from word_to_vec import W2VTrainer

w2v_trainer = W2VTrainer(output_dir="./checkpoints/")
w2v_trainer.run(data)

結果

Word2Vec では単語のベクトルが得られるので、内積を取ることでその近さを知ることが可能です。よってある単語に最も近い単語、最も遠い単語を計算できます。 今回は単語ではなく本の名前なので、ある本と似た本、全く似ていない本を計算出来ることになります。

ここでようやくタイトルの「吾輩は猫である」が出てきます。吾輩は猫であるも色々とバージョンがある(ハードカバーとか文庫版とか)ので、まずは最も読まれていた新潮文庫版で試します。 吾輩は猫であるとの類似度の計算結果が以下になります。

我輩は猫である(新潮文庫) に近い本・遠い本

近い本

[('三四郎(新潮文庫)', 0.8102374076843262),
 ('こころ(新潮文庫)', 0.7891632914543152),
 ('門(新潮文庫)', 0.7769464254379272),
 ('草枕(新潮文庫)', 0.7562223672866821),
 ('虞美人草(新潮文庫)', 0.7479068636894226),
 ('阿部一族・舞姫(新潮文庫)', 0.7449460625648499),
 ('坊っちゃん(新潮文庫)', 0.7428607940673828),
 ('それから(新潮文庫)', 0.7405833005905151),
 ('行人(新潮文庫)', 0.7404835224151611),
 ('道草(新潮文庫)', 0.7297923564910889)]

遠い本

[('死神姫の再婚-五つの絆の幕間劇-(ビーズ…', 0.13127756118774414),
 ('絶園のテンペスト2(ガンガンコミックス)', 0.1260046362876892),
 ('死神姫の再婚-始まりの乙女と終わりの教師-…', 0.1172419935464859),
 ('信長協奏曲1(ゲッサン少年サンデーコミッ…', 0.09469656646251678),
 ('レベルE(上)(集英社文庫―コミック版)', 0.09304021298885345),
 ('絶園のテンペスト(5)(ガンガンコミックス)', 0.09189711511135101),
 ('信長協奏曲7(ゲッサン少年サンデーコミッ…', 0.09187307208776474),
 ('信長協奏曲2(少年サンデーコミックス)', 0.08920496702194214),
 ('天使1/2方程式1(花とゆめCOMICS)', 0.08579257875680923),
 ('闇の皇太子偽悪の革命家(ビーズログ文庫)', 0.08444569259881973)]

予想以上に綺麗に分離してくれました。 同じ夏目漱石の文庫本が多く並んでいます。 また森鴎外が入っているのも納得できます。 反対の遠い本ではコミックス系が多くなっていて、これも納得できる結果になりました。 LDA よりはうまくいってそう…

いい感じだったので森見登美彦の「四畳半神話大系(角川文庫)」でもやってみました

四畳半神話大系(角川文庫)に近い本・遠い本

近い本

[('夜は短し歩けよ乙女(角川文庫)', 0.8690416216850281),
 ('太陽の塔(新潮文庫)', 0.8665028810501099),
 ('美女と竹林(光文社文庫)', 0.8509087562561035),
 ('有頂天家族(幻冬舎文庫)', 0.8440086841583252),
 ('新釈走れメロス他四篇(祥伝社文庫も…', 0.8257927894592285),
 ('ペンギン・ハイウェイ(角川文庫)', 0.8094362020492554),
 ('([も]3-1)恋文の技術(ポプラ文庫)', 0.8090072274208069),
 ('きつねのはなし(新潮文庫)', 0.7985487580299377),
 ('四畳半王国見聞録(新潮文庫)', 0.7920351028442383),
 ('宵山万華鏡(集英社文庫)', 0.7382686734199524)]

遠い本

[('血界戦線4―拳客のエデン―(ジャンプコ…', 0.07178327441215515),
 ('光とともに…13―自閉症児を抱えて', 0.06594447791576385),
 ('王家の紋章(17)(Princessc…', 0.0623263344168663),
 ('VIVO!1(マッグガーデンコミックス…', 0.05496706813573837),
 ('王家の紋章(15)(Princessc…', 0.05385550856590271),
 ('光とともに…14―自閉症児を抱えて', 0.052545785903930664),
 ('彼氏彼女の事情(15)(花とゆめCOMI…', 0.050349362194538116),
 ('王家の紋章(16)(Princessc…', 0.050084952265024185),
 ('Vassalord.6(マッグガーデンコ…', 0.04914539307355881),
 ('彼氏彼女の事情(16)(花とゆめCOMI…', 0.0489741712808609)]

こちらも上位がすべて森見登美彦氏で埋まっていますね。 同じ作家で固まることがよく起こるというのは、作家によらず、特定の作家の本を読むと同じ作家の本を再び読む、ということが言えるのかもしれません。まあそりゃ当たり前か。 一度気にいると同じ作家の本ばかり買ってしまう、という感覚が可視化された、ということでしょうか。

ここで「ハードカバーの方がコアな客が多いから、同じ作家で固まりやすいのでは(類似度がほぼ1になるとか」と思ったので検証してみます。 ハードカバーの四畳半神話大系の類似度は以下のようになりました。

[('きつねのはなし', 0.7952536344528198),
 ('リンダリンダラバーソール(新潮文庫)', 0.7498961091041565),
 ('新釈走れメロス他四篇', 0.7390938401222229),
 ('夜は短し歩けよ乙女', 0.7354095578193665),
 ('太陽の塔', 0.7339460849761963),
 ('有頂天家族', 0.7336863875389099),
 ('風に舞いあがるビニールシート', 0.7267143130302429),
 ('ぐるぐるまわるすべり台', 0.7245926260948181),
 ('天国はまだ遠く', 0.7150847911834717),
 ('二枚舌は極楽へ行く(FUTABA・NOVE…', 0.7144078016281128)]

「きつねのはなし」は森見登美彦氏ですが、類似度は先よりも落ちています。 二番目の本はしらなかったのでググったのですが大槻ケンヂ氏の本なんですね、知りませんでした。

www.amazon.co.jp

内容紹介
僕らのバンドが、メジャーデビューすることになった! その頃、日本はバンドブームに沸いていた。無名だった若者が、次々とスターになった。ライブ会場は熱狂に満ちた。でも、ブームはいつか終わるものだ。大人たちは、潮が引くように去ってゆく。誰もが時の流れと無縁ではいられないんだ。僕と愛すべきロック野郎たちの、熱くて馬鹿馬鹿しくて切なかった青春を、いま再生する。

普通に面白そう……。 あでもそう思うってことは、ちがう作家だけどその人が面白いと思いそうな本を拾えてるってことなのかもしれないですね。サンプル僕だけであれですが。

文庫版の森見登美彦氏の本も見当たらないのが不思議ですね。これはハードカバーを読む人と、文庫を読む人はそもそも文脈がだいぶ違っていて、あの本はハードカバーだけど今回は文庫までまとう、といった読み方はしないんだなあということが伺えます。

ちなみに吾輩は猫である岩波文庫版になるとだいぶ様相が違っていました

[('「超」文章法(中公新書)', 0.8472597002983093),
 ('無限論の教室(講談社現代新書)', 0.8468811511993408),
 ('天の瞳幼年編〈2〉(角川文庫)', 0.8414252400398254),
 ('日本語の教室(岩波新書)', 0.8320746421813965),
 ('不都合な真実', 0.8299296498298645),
 ('朝の少女(新潮文庫)', 0.8163926601409912),
 ('人生に意味はあるか(講談社現代新書)', 0.8158650398254395),
 ('けものたちは故郷をめざす(新潮文庫あ4…', 0.813538670539856),
 ('DeepLove―アユの物語完全版', 0.812150239944458),
 ('ブラック・ジョーク大全(講談社文庫)', 0.8112159967422485)]

夏目漱石が激減してしまったのが不思議です。岩波文庫でも夏目漱石の本はいくつか出ているのですが…

まとめと今後の課題

  • word2vec の応用の広さすごい
  • 著者のデータでも何かしたい
  • LDA かわいそうなのでなんとかチューニングして結果出してあげたい。ベイズだしね。

*1:トピック数, αの事前分布, 最適化方法などいろいろやったものの徒労に終わる悲しみ

*2:元論文: [1301.3781] Efficient Estimation of Word Representations in Vector Space

ニューラルネットへのベイズ推定 - Bayesian Neural Network

ニューラルネットワーク過学習防止としてもちいられる Dropout がネットワークの重みに対してベイズ推定を行っていることに相当する、ということを述べた論文について紹介します。

くわえてそれを実装して簡単な実験をやってみます。

Ref

ベイズ推定

はじめにベイズ推定について簡単に。
ベイズ推定で目指す目的とは何でしょうか。

仮定として, $N$ 個の入力データ $X=\{x_1, \cdots. x_N \}$ と、それに対応する出力 $Y=\{ y_1, \cdots, y_N \}$ をすでに観測しているとします。またデータは独立同分布から生成されている, とします。
この状態で、新しい入力 $x$ が入ってきたときに出力 $y$ を予測するという問題を考えます。
ベイズ推定で求めたいのは、データが与えられたときの、隠れ変数の事後分布 $p(w \mid X, Y)$ です。

重みの事後分布さえ入手できてしまえば、それを用いて新しいデータ $\hat{x}$ が与えられたとき、その予測値 $\hat{y}$ が従う分布は以下で計算することが可能です。

$$p(\hat{y} \mid \hat{x}, X, Y) = \int p(\hat{y} \mid w, \hat{x}) p(w \mid X, Y) dw$$

ですから、ベイズ推定の立場にたてば、事後分布 $p(w \mid X, Y)$ をいかにして求めていくか, ということが問題となります。

事後分布の計算

直接事後分布が計算できるような単純なモデルでは、先のアイディアに基づいて事後分布を計算すれば良いです。

しかし、複雑なモデルでは、陽に事後分布が計算できない場合が出てきます。 このような仮定のもとでは、何らかの関数で事後分布を近似する必要があります。
以下ではあるパラメータ $\theta$ を持つ関数 $p_\theta(w)$ を事後分布に近づいけていく, という方法を考えます。
このとき、2つの分布の近さを測る距離として, 2つの確率分布から実数空間へ射影する関数として KL-Divergence (以下では ${\rm KL}$ と表記) を用います。
KL距離は距離空間ではありませんが, 確率分布 $p,q$ に対して, ${\rm KL}(q, p) = 0$ のとき $q=p$ が成立します。よって KL 距離の意味で最小化を行えば, 分布 $q_\theta$ は $p(w\mid X, Y)$ に近づく事が期待されます.
最小化する関数を変形すると、以下のようになります。

$$\begin{align}{\rm KL}(q_{\theta}(w) | p(w | X, Y))&=\int q_{\theta}(w) \log \left( \frac{q_{\theta}(w)}{p(w | X,Y)} \right)dw \\&=\int q_{\theta}(w) \log \left( q_{\theta}(w) \frac{p(X,Y)}{p(X,Y|w) p(w)} \right) dw \\&=\int q_{\theta}(w) \left( \log \frac{q_{\theta}(w)}{ p(w)} + \log \frac{p(Y|X) p(X)}{p(Y|X, w) p(X|w)} \right) dw \\&=\int q_{\theta}(w) \left( \log \frac{q_{\theta}(w)}{ p(w)} + \log \frac{p(Y|X)}{p(Y|X, w)} \right) dw \\&= \int q_{\theta}(w) \log p(Y|X) dw - \int q_{\theta}(w) \log p(Y|X,w)dw + {\rm KL} (q_{\theta}(w) | p(w)) \\&= \log p(Y|X) - \int q_{\theta}(w) \log p(Y| X, w)dw + {\rm KL} (q_{\theta}(w) | p(w))\end{align}$$

ここで第4行目の変形で $p(X|w)=p(X)$ を用いました。今 Variance Inference を ${\rm VI}$ と表し

$${\rm VI}(\theta) = - \int q_{\theta}(w) \log p(Y|X,w)dw + {\rm KL} (q_{\theta}(w) | p(w))$$

と定義します。すると上記の式は以下のように変形することができます.

$${\rm KL}(q_{\theta}(w) | p(w | X, Y)) = {\rm VI}(\theta) + \log p(X,Y)$$

ここで右辺最終項 $\log p(X,Y)$ は定数ですから, 結局 KL 距離の最小化問題は以下のように書き換えられます。

$$\min_{\theta} {\rm KL}(q_{\theta}(w) | p(w | X, Y)) = \min_{\theta}{\rm VI}(\theta)$$

よって ${\rm VI}$ を $\theta$ に関して最小化することは, KL-Divergence の意味での最小化と一致する事がわかります。
また ${\rm VI}$ の第一項

$$-\int q_{\theta}(w) \log p(Y|X,w)dw$$

の部分は, $0 \le p(Y|X,w) \le 1$ であることを考慮すると必ず正の値です。また得られたデータ $Y$ が現れる確率 $p(Y|X,w)$ の値が高い $w$ に対して、分布 $q_\theta(w)$ も高い値を取るときに小さくなり、データへの当てはまりを表していることがわかります。

一方第二項の

$${\rm KL} (q_{\theta}(w) | p(w))$$

については, 事前分布 $p(w)$ と事後分布の推定密度関数 $q_\theta(w)$ との距離を表しています。これは事前分布からかけ離れた分布へのペナルティを与える項に相当しており、一般に正則化項と呼ばれるものです。(簡単な例で言えば互いにガウス分布であるときこのKLダイバージェンス部分は $L2$ ノルムとなり線形回帰であれば Ridge 回帰とよばれる枠組みになります)
このように自分が提案する分布 $q$ をもとめたい分布 $p$ に近づけていく方法を変分推論 (Variance Inference) とよびます.

ニューラルネットワークに対する既存の VI

つぎにニューラルネットワークに対して Variance Inference を行う既存の枠組みについて振り返ってみます。近年のニューラルネットワークへのベイズ推定のアプローチの一つである Hinston and Van Camp, 1993 では近似する分布 $q_\theta$ にたいして, すべての重みに対して独立したガウス分布を仮定しました。

すなわち $I$ 層の隠れ層を持ち, その各 $i$ 層で重み $w_{ijk}$ を持つネットワークに対して $m_{ijk}, \sigma_{ijk} \in \mathbb{R}$ とおいたとき

$$q_\theta (w) = \prod_{i, j ,k} N(w_{ijk}\mid m_{ikj}, \sigma_{ijk}^2)$$

と記述することになります。 この分布を仮定した場合の最適化はとてもむずかしく, 論文中では一つの隠れ層を持つニューラルネットワークに対しての実験にとどまっていました。 実際理論的には事後分布の推論にガウス分布を用いると一層の隠れ層をもつニューラルネットに対しては数式上綺麗に解析が行えるものの、実験上性能がよくありませんでした。理由としては、この方法が本来重要である重み同士の関係性を記述できないことなどが挙げられています。

これに対し Barver and Bishop, 1998 では重みに対して混合ガウスを仮定したモデルを提案しました。これによって同一層の重みに対して関連性を考慮することが可能となりました。しかし、その反面計算量が増大してしまうため、複雑なモデル等に対してうまく働きません。

Variance Inference の近似

そのままの形ではなかなか良い結果を得ることができていなかった変分推論ですが、今回は期待値の意味で一致する近似式を用いて最適化することを考えていきます。


まず $X,Y$ は独立であると仮定しているので $p(Y | w,X) = \prod_{i} p(y_i | w, x_i )$ が成り立ちます。これを用いて Variance Inference を変形すると以下のようになります。

$$\begin{align*}{\rm VI}(\theta) &= - \int q_{\theta}(w) \log p(X,Y|w)dw + {\rm KL} (q_{\theta}(w) | p(w)) \\&= -\sum_{i=1}^N \int q_{\theta}(w) \log p(y_i | x_i, w) dw + {\rm KL} (q_{\theta}(w) | p(w)) \\&= -\sum_{i=1}^N \int q_{\theta}(w) \log p(y_i | f^w(x_i)) dw + {\rm KL} (q_{\theta}(w) | p(w))\end{align*}$$

ここで $f^w(x_i)$ は重み $w$ を持つモデルに $x_i$ が入力されたときの出力を表します。

この形式を取り扱う上での難しさは, 主に以下の二点にあります

  1. $w$ に関する積分の部分が扱いやすい形式ではない
  2. データ $N$ に対する和を取らなくてはならないため, データの数が多くなると計算が困難になる

この内 2 に関してはすべてのデータのうちで \( M \) 個だけサンプルする (mini-batchを計算する) ことで回避することができます. この近似を行った $\hat{{\rm VI}}$ は以下のようになります。

$$\hat{{\rm VI}} = - \frac{N}{M} \sum_{i \in S} \int q_{\theta}(w) \log p(y_i | f^w(x_i)) dw + {\rm KL} (q_{\theta}(w) | p(w))$$

ここで, $S$ はサンプルされた添字の集合を表します。

しかし理由 1 によりこの計算は難しいままです。

$w$ の積分の近似

$w$ に関する積分を近似することを考えます。 積分の部分を見ると, $q_\theta$ と対数尤度との掛け算になっています。よって $q_\theta$ から $w$ をサンプリングすることができれば積分計算の近似を行うことができます。

しかし, $q_\theta$ の形式には仮定が置かれておらず任意の分布を取ることが可能です。したがってこの分布からサンプリングすることはできません。

そこで $q_\theta(w)$ がパラメータを持たない別の分布 $p(\epsilon)$ を用いて, $w = g(\theta, \epsilon)$ で表現できるという仮定を加えてみます。そうすると分布 $q$ に対する情報がわかっていなくても、分布 $p(\epsilon)$ からサンプルした値 $\hat{\epsilon}$ を用いて提案分布からのサンプル $\hat{w} = g(\theta, \hat{\epsilon})$ を生成することが可能となります。

平均 0 分散 1 のガウス分布 (正規分布) のサンプル値から 任意の分散と平均値をもつガウス分布のサンプルを作成することなどがこれに相当します。

これを用いて先の VI の式中の $q_\theta$ を $p(\epsilon)$ で表現することで, 分布の部分から $\theta$ を取り除きます。

$$ \begin{align} \hat{{\rm VI}} &= - \frac{N}{M} \sum_{i \in S} \int q_{\theta}(w) \log p(y_i | f^w(x_i)) dw + {\rm KL} (q_{\theta}(w) | p(w)) \\ &= - \frac{N}{M} \sum_{i \in S} \int p(\epsilon) \log p(y_i | f^{g(\theta, \epsilon)}(x_i)) d\epsilon + {\rm KL} (q_{\theta}(w) | p(w)) \end{align} $$

こうなると積分はパラメータを持たない分布 $p(\epsilon)$ からのサンプリングを行うことで効率的に近似をおこなうことが可能になります。 サンプリングをミニバッチのサンプルと同時に行うとし、サンプルされた実現値を $\epsilon_i \sim p(\epsilon) (i = 1, \cdots, M)$ とします。
するとこのモンテカルロ法によって積分が近似された $\hat{{\rm VI}}_{MC}$ は以下のようになります。

$$\hat{{\rm VI}}_{MC} = - \frac{N}{M} \sum_{i \in S} \log p(y_i | f^{g(\theta, \epsilon_i)}(x_i)) + {\rm KL} (q_{\theta}(w) | p(w))$$

このとき $\epsilon, S$ に対する期待値を取ると, $\mathbb{E}_{S, \epsilon} \left[\hat{{\rm VI}}_{MC} \right] = {\rm VI}$ が成立します。

この新しい $\hat{{\rm VI}}_{MC}$ を目的関数として例えば勾配法を用いて $\theta$ について最適化を行えば、期待値を取れば元の Variance Inference を最適化することと同値です。よって最適化の各 $t$ ステップにおいて

$$ \displaystyle \theta_{t} = \theta_{t-1} + \eta \frac{\partial}{\partial\theta} \hat{{\rm VI}}_{MC} $$

によりパラメータ $\theta$ を更新すれば良いことがわかります。

Dropout による学習

ここで一旦 Variance Inference のことをおいておいて, ニューラルネットワーク過学習を抑える手法の一つである Dorpout について考えます。
Dropout とは訓練データが与えられたとき, 各層においてすべての隠れノードを用いて出力を行わず, ランダムに選ばれたノードの値のみを用いて出力をし, backword においても出力に関わったノードの値のみを更新する, という方法です。

単純に一層の隠れ層を持つネットワークを考え, 入力から隠れ層, 隠れ層から出力層への重みをそれぞれ $M_1 \in \mathbb{R}^{n\times m}, M_2 \in \mathbb{R}^{m\times l} $とします.
また隠れ層の定数 $b \in \mathbb{R}^{m}$, 活性化関数 $\sigma$ とします。

今入力として $x \in \mathbb{R}^n$ のデータが与えられたとします。この時 $n$ 次元上の確率 $p_1 (0\le p_1\le 1)$ のベルヌーイ分布 $Q$ から 0-1 の $n$ 次元バイナリからなるベクトル $e_1 \in \mathbb{R}^n$ を生成します。この $e_1$ を用いて Dropoutを適用した隠れ層 $h$ は

$$ \begin{align*} h = \sigma ( (x \bullet e_1) M_1 + b) \end{align*} $$

となります。 ここで $\bullet$ は要素積 $x \bullet y = \sum x_i y_i$ を表します。

同様に隠れ層 $h$ に対しても確率 $p_2$ で要素を0にします。すなわち 先と同様に $m$ 次元のベルヌーイ分布から $e_2 \in \mathbb{R}^m$ をサンプリングして一定の隠れ層ノードの値を0にします。すなわち

$$\hat{h} = h \bullet e_2$$

をノードの値であるとします。出力はこの値を用いて

$$\hat{y} = \hat{h} M_2$$

となります。
この $\hat{y}$ は $\bullet e = {\rm diag}(e)$ と変形できることを用いれば以下のように変形できます。

$$ \begin{align}\hat{y} & = (h \bullet e_2) M_2 \\&= ( \sigma \left( (x \bullet e_1) M_1 + b \right) \bullet e_2) M_2 \\&= ( \sigma (x ({\rm diag}(e_1) M_1) + b)({\rm diag}(e_2) M_2 ) \\&= \sigma \left( x\hat{W}_1 + b \right) \hat{W}_2 \end{align} $$

ここで  {\rm diag}(e_1) M_1 = \hat{W}_1, {\rm diag}(e_2) M_2 = \hat{W}_2 と定義しました。

以上を用いるとニューラルネットワークの出力は確率変数 $\hat{\omega} = \{ \hat{W_1}, \hat{W_2}, b\}$ を用いて

$$\hat{y} = f^{\hat{W_1}, \hat{W_2}, b}(x)$$

と記述できます. 


ながながとゴニョゴニョしましたが, 以上からdropout による mask のかけられた出力も重みの確率変数をもつネットワークの出力として表現できる, ということが確認できました。

Dropout の目的関数

これらの記号を用いてニューラルネットワークの目的関数を記述していきましょう。入力値と正解ラベルから誤差を計算するロス関数を $E$ とおき、ニューラルネットワークが最小化する真の関数を記述すると以下の用になります。

$$L_{dropout} = \frac{1}{M} \sum_{n \in N} E \left( f(x_n), y_n \right) + \lambda_1 ||M_1|| + \lambda_2 || M_2 || + \lambda_3 || b ||$$

ここで $f(x_n)$ は $n$ 番目の入力データに対する出力値を表しています。また右辺第二項以降は重みに対する正則化 (weight decay)を表しています.

実際には $N$ 個すべてのデータを使うことは困難ですから、その中からある一定の大きさのサンプル $S$ を取得します。Dropout ではそれと同時にネットワークに対する mask をかけて出力を作ります。したがってネットワークの出力は確定値 $f(x)$ ではなく、確率変数 $\omega =\{ \hat{W_1}, \hat{W_2}, b\}$ によって決定する確率的な出力となります。よって Dropout を用いたニューラルネットワークのロス関数は

$$L_{dropout} = \frac{1}{M} \sum_{i \in S} E \left( f^{ \{ \hat{W_1^i}, \hat{W_2^i}, b \} }(x), y_i \right) + \lambda_1 ||M_1|| + \lambda_2 || M_2 || + \lambda_3 || b ||$$

となります. ここで $\hat{W_1^i}, \hat{W_2^i}$ はサンプルされた $i$ のデータに対する dropout のマスクがかけられた重みを表しています。

回帰問題においては, 関数 $E$ は定数部分を除いて負の対数尤度関数で書き換えることができます。すなわち

$$E (f(x),y) = \frac{1}{2} \| y - f(x) \|^2 = - \frac{1}{\tau} \log p(y | f(x)) + {\rm const}$$

と変形できます。 ここで尤度関数は $p(y | f(x)) = N(y;f(x), \tau^{-1}I)$ のガウス分布で, $\tau$ は精度パラメータです.
今回は回帰問題を考えましたが、分類問題においても同様に負の対数尤度を用いて定式化することが可能です. (その場合 $\tau = 1$ となります.)

これよりロス関数は

$$L_{dropout} =  - \frac{1}{M} \sum_{i \in S} \frac{1}{\tau} \log p( y_i | f^{ \{ \hat{W_1^i}, \hat{W_2^i}, b \} }(x)) + \lambda_1 ||M_1|| + \lambda_2 || M_2 || + \lambda_3 || b ||$$

つぎに確率変数 $\hat{\omega}$ について考えます。この集合はネットワークの重みという確定的な値とdropout による確率変数の部分をあわせたものでした。それを明示的に記述すると

$$\hat{\omega_i} = \{ \hat{W_1}, \hat{W_2}, b\} = \{{\rm diag}(e_1^i)M_1, {\rm diag}(e_2^i)M_2, b \} := g(\theta, \hat{e_i})$$

と書き換えることができます。ここで $\theta = \{M_1, M_2, b\}$ は確定的な値を要素に持つ集合と定義し、$\hat{e_i} = \{ e_1^i, e_2^i \}$ は $i$ 番目のミニバッチのサンプルによって作成される dropout のマスクを要素に持つ、確率変数の集合であると定義します。 また関数 $g$ は2つの集合 $\theta, \hat{e_i}$ から $\hat{\omega_i}$ をつくる射影であるとします。
また $1 \le i \le N$ にたいして $e_1^i \sim p(e_1)$, $e_2^i \sim p(e_2)$です. ここで, $p(e_j)\ (j=1,2)$ はそれぞれ確率 $p_j$ のベルヌーイ分布の積であるとします。
これらを用いると dropout のロス関数は確率変数 $\hat{e_i}$ を用いて

$$L_{dropout} = - \frac{1}{M\tau} \sum_{i \in S} \log p(y_i | f^{g(\theta, \hat{e_i})}(x_i)) + \lambda_1 ||M_1|| + \lambda_2 || M_2 || + \lambda_3 || b ||$$

となります。この目的関数の重み $\theta := \{M_1, M_2, b\}$ に関する勾配は

$$\frac{\partial}{\partial \theta} L_{dropout} = - \frac{1}{M\tau} \sum_{i \in S} \frac{\partial}{\partial \theta} \log p(y_i | f^{g(\theta, \hat{e_i})}(x_i)) + \frac{\partial}{\partial \theta} \left( \lambda_1 ||M_1|| + \lambda_2 || M_2 || + \lambda_3 || b || \right)$$

となり, この勾配を逆伝搬させてネットワークの重みを更新します.

ところでこれは先程の変分推論の式ととても良く似ています. 再掲すると

$$\hat{{\rm VI}}_{MC} = - \frac{N}{M} \sum_{i \in S} \log p(y_i | f^{g(\theta, \epsilon_i)}(x_i)) + {\rm KL} (q_{\theta}(w) | p(w))$$

であり, これは ${\rm KL} (q_{\theta}(w) | p(w)) = N\tau \left( \lambda_1 ||M_1|| + \lambda_2 || M_2 || + \lambda_3 || b || \right)$ とおけば

$$L_{dropout} = \frac{1}{M\tau} \hat{{\rm VI}}_{MC}$$

となりdropout と変分推論は厳密に一致します!!


厳密に一致するということは, dropout を用いて学習したネットワークは、変分推論による重みの事後分布となっているということです。

すなわち dropout を用いるのはベイズ学習をしていることにほかならない ということです。ですから学習によって得られている値は確定的な値ではなく、 dropout と組み合わせることで重みの事後分布を計算できます。

これによって。学習済みのネットワークにとあるdropoutをかけた出力は重みの事後分布から一つの重みをサンプルしているという風に解釈することが可能です。したがって、複数のdropoutの係数を用いて複数の出力を生成しそのばらつきを用いると, 事後分布の分散を推定することも可能です.

この枠組では, ニューラルネットワークの出力は確率過程からのサンプルにほかならないため、出力がどの程度信頼できるか(ばらつきをもっているか) を見積もることが可能になっています。すごい。

数値実験

簡単な人口データの回帰問題を解いてみます.

実験に用いたコードは以下においています。

github.com

条件

トレーニングデータは1次元の100個のデータで, ターゲット変数は $f(x) = x + \sin 5x$ に平均0, 分散 0.1 のガウスノイズを載せて作成します.

ネットワークは5層で各隠れ層が512次元としました。dropoutは確率 0.5 のベルヌーイ分布を用います。

relu

活性化関数を relu として 1000 epoch 計算したものの出力が下のグラフです。

f:id:dette:20170814172111p:plain

青い点がデータで, 橙色で事後分布からサンプルされたネットワークの出力を表しています。 これを見ると, データが存在する部分の分散は小さくなっていますがデータがない部分に行くに連れて分散が大きくなっていることがわかります. これは事後分布の分散が大きくなっているためです。このように従来のニューラルネットの予測値に加えて、その信頼度が出せていることが確認できます。

tanh

次に活性化関数を tanh に変えて実験したものが以下の図です。

f:id:dette:20170814174525p:plain

tanh に関してはreluのときのようにデータがない部分で極端に分散が大きくなりませんでした。またデータへの当てはまりも tanh の方は積極的に行っておらず、$x$ 軸性の領域ではほとんど 0 を予測しておりあたかも「データが足りないよ」と言っているようです。 

これは relu のほうが勾配消失に強く学習が素早く進むため、同じ epoch 数でもデータへ強くフィッティングしていくことが原因と考えられます。

また範囲外への予測も relu のようにそれまでの傾きを踏まえた予測ではなく、だんだんと 0 に向かうように予測しているように見受けられますが、これも活性化関数の形の影響 (reluは非ゼロのとき単調増加関数となる) を受けていると考えられます。

ちなみに tanh でも epoch さえ増やせばデータへしっかり当てはめてくれるようで、 epoch 数を更に増やして 4000 回としたものが以下になります。

f:id:dette:20170814175210p:plain

これを見るとデータがある部分に関しては平均値はほぼデータ通りで、分散もこちらが設定した値 0.1に近い値となっていることが伺えます。 1000 epoch では epoch が足りていなかったみたいですね。

感想

今まで確定値しか扱えなかったところに、Dropoutという既存の手法とベイズ推定の枠組みを結びつけて行くところは面白いなと思って読んでいました。CNNやRNNへ応用すれば「この画像はよくわかんないけど猫」とか「絶対犬です!」みたいな信頼度を同時にだすネットワークが自然に作れそうで今後発展すると面白そうだなと思います。

おまけ - 動画バージョン

tanh

 

f:id:dette:20170814180735g:plain

うにょうにょしながら学習していくの面白い。

機械学習のコード整理

タイトルのままですが、今まで適当に書いてローカルに放置してたりブログには書いたけどそのまま放置してたりしていたコードを一つのリポジトリにまとめる作業をしました。

まとめながら、昔書いたコードをいろいろと手直ししていましたが、ひどく密につながりすぎたひどいコードがたくさん出てきてよりきれいな洗練されたコードを書きたいなと思う次第でした。

github.com

個人的には エビデンス最大化 と

f:id:dette:20170712221500p:plain

RVM のグラフが好きです。

f:id:dette:20170712221506p:plain

matplotlibもこの二年でだいぶおしゃれになって、seaborn使わなくても最初からかっこよく描画してくれるようになりました。かっこいいグラフは出てくるだけでやる気が出るので良いです。