2018-07-12

NICE: NON-LINEAR INDEPENDENT COMPONENTS ESTIMATIONを読んだのでメモ

はじめに

NICE: NON-LINEAR INDEPENDENT COMPONENTS ESTIMATIONを読んだのでメモ．

基本アイディア

Normalizing flowにおいて，変数変換における関数 $f$ をニューラルネットを使って定義したという話（厳密には，議論の展開的にnormalizing flowは関係なく単純な変数変換の話のみ．ただ，最終的には変数変換を複数回行うため意味合いとしては変わらない）．前に読んだVariational Inference with Normalizing Flowの先駆け的なもの．

まずデータ $x$ を二つのブロック $(x_1,x_2)$ に分け，以下のような変換を定義する．

$\displaystyle y_1=x_1,\:y_2=x_2+m(x_1)$

$m$ が変換の中心で，ReLUとMLPによる非線型変換を表し，データの一方を使ってもう一方を変換することを表現している．変換をこのように定義してあげると各 $m$ においてunit jacobian deteminantを持ち，逆変換も以下のように定義できるためNormalizing Flowとして良く機能するでしょうという話．

$\displaystyle x_1=y_1,\:x_2=y_2-m(y_1)$

以下，この式が導かれるまでの話．

non-linear independent components estimation

ここでの問題設定は最尤推定によりデータ分布を単純な分布へ移す連続かつ微分可能な変換を学習すること．基本的にはnormalizing flowをベースとしているので以下の変数変換の関係性を用いる．

$\displaystyle \log p_X(x)=\log p_H(f(x))+\log\left|\mathrm{det}\left(\frac{\partial f(x)}{\partial x}\right)\right|$

$p_H(h)$ は事前分布で，扱いやすい等方性のガウス分布などを考えている（ただし，固定でなく学習してもよいとのこと）．もしそのような分布で事前分布が定義されていたとしたら，以下のnon-linear independent components estimation (NICE) criterionが得られる．

$\displaystyle \log p_X(x)=\sum_{d=1}^D\log p_{H_d}(f_d(x))+\log\left|\mathrm{det}\left(\frac{\partial f(x)}{\partial x}\right)\right|$

単純に第１項目の $f(x)$ が $f_d(x)$ に変わっただけで $f_d(x)$ は $f(x)=(f_d(x))_{d\leq D}$ を満たす関数．

今までの議論からNICEは対数尤度を増加させる逆変換可能なデータセットの前処理を学習していることがわかる．さらに変数変換の関係性が正則化として働いており任意に尤度が増加することを防いでいる．

モデルの構成

モデルはforward(encoder $f$ )とbackward(decoder $f^{-1}$ )どちらにおいてもヤコビアンが扱いやすくstraightforwardに計算ができなければならない．仮にニューラルネットのように多層構造（ $f=f_L\circ\dots\circ f_1$ ）として構成した場合，ヤコビアンの行列式は各レイヤーのヤコビアンの行列式の積になり嬉しい．そこでまずは $f$ を初歩的な構成で定義することから考える．

まず一つの手段としてアフィン変換（ニューラルネット）を考える．ここで満たさなければいけないことは，ヤコビアンや逆変換が容易に計算可能なことである．そのような行列として三角行列があげられる．特に正則な行列はLU分解によって上三角行列と下三角行列の積に分解できる．この恩恵を受ける一つの構成方法はニューラルネットの重みを三角行列にし，活性化関数が全単射であるものを選ぶことであるが，これは制約が強すぎてモデルのデザインが限られるのであまり良くない．

そこで別な手としてヤコビアンが三角行列になる全単射の変換を考える．これならヤコビアンの行列式が容易に計算できるため問題がない．そこで次のような変換を考える．

$\displaystyle y_{I_1}=x_{I_1},\:y_{I_2}=g(x_{I_2};m(x_{I_i}))$

基本アイディアのところで出てきた式と同じで， $D$ 個のデータ $x$ を $I_1,I_2(d=|I_1|)$ に分けており， $g:\mathbb{R}^{D-d}\times m(\mathbb{R}^d)\rightarrow\mathbb{R}^{D-d}$ は逆変換可能かつ， $x_{I_2}$ を $x_{I_1}$ を使って変換するような関数を表しており，ここではカップリング則(coupling law)と呼ぶ．このように二つのデータに分けることで $y_{I_1}$ を $x_{I_1}$ で微分すれば単位行列になり， $x_{I_2}$ で微分すれば0行列になる．そのためヤコビアンは

$\displaystyle \frac{\partial y}{\partial x}= \begin{bmatrix} I_d\:\:\:\:\:0 \\ \frac{\partial y_{I_2}}{\partial x_{I_1}}\:\:\:\frac{\partial y_{I_2}}{\partial x_{I_2}} \end{bmatrix}$

とうまいことブロック三角行列になる．よってヤコビアンの行列式は $\mathrm{det}\frac{\partial y}{\partial x}=\mathrm{det}\frac{\partial y_{I_2}}{\partial x_{I_2}}$ と簡単になる．さらに逆変換も

$\displaystyle x_{I_1}=y_{I_1},\:x_{I_2}=g^{-1}(y_{I_2};m(y_{I_i}))$

とかける． $m$ をカップリング関数（coupling function），変換をカップリング層（coupling layer）と呼ぶ．

カップリング則として加法カップリング則（additive coupling law）を選べば， $g(a;b)=a+b$ となり，

$\displaystyle y_{I_2}=x_{I_2}+m(x_{I_i}) \\ \displaystyle x_{I_2}=y_{I_2}-m(y_{I_i})$

となる．すると逆変換を陽に計算する必要はなく，カップリング関数にどのようなものを使用しても問題がなくなる．例えば， $|I_1|=d,\:|I_2|=D-d$ であることに注意をすれば， $m$ を $d$ 次元入力， $D-d$ 次元出力のニューラルネットとして構成することができる．さらに嬉しいことに， $y_{I_2}$ の $x_{I_2}$ による微分が単位行列になることからヤコビアンの行列式は常に1となる．

そのほかのカップリング則として乗法カップリング則（multiplicative coupling law） $g(a;b)=a\odot b,\: b\neq 0$ やアフィンカップリング則（affine coupling law） $g(a;b)=a\odot b_1+b_2,\: b_1\neq 0$ を選ぶことができる．ただしこの論文では計算の安定性等から加法カップリング則を使っている．

加法カップリング則を使えばヤコビアンの行列式が常に1で計算が楽になって嬉しいということだったが，変数変換の関係からカップリング層を導入しても分布自体は変化がない．そこでスケーリング行列と呼ばれる対角行列を導入する．これを使うことでヤコビアンの対角成分が1ではなくスケーリング行列 $S$ の対角成分に等しくなり，ヤコビアンの行列式は $\prod_i S_{ii}$ となる．結果としての加法カップリング則にスケーリング行列を導入した際の目的関数は

$\displaystyle \log p_X(x)=\sum_i[\log p_{H_i}(f_i(x))+\log |S_{ii}|]$

となる．

最後に実験的なところとして，この論文では事前分布として以下のような因子分解可能な分布を使った．

$\displaystyle p_H(h)=\prod_{d=1}^Dp_{H_d}(h_d)$

基本的には問題に合わせてガウス分布かロジスティック分布を使えばいいとのこと．

$\displaystyle \log p_{H_d}=-\frac{1}{2}(h_d^2+\log 2\pi) \\ \displaystyle \log p_{H_d}=-\log(1+\exp(h_d))-\log(1+\exp(-h_d))$

以上から，NICEの良いところは変数変換に用いることのできる関数の自由度，つまり，活性化関数にReLUのような全単射でないものを使えるし，batch normやresnetのような構造を入れることも可能．しかもコスト関数にL2のような再構成誤差に依存しない良さがある．ただし，スケーリング行列の点は個人的には少し違和感。。。

まとめ

Variational inference with normalizing flowを読んだ後なのでスケーリング行列あたりでこれでいいのか感がでたけど面白かった．Variational inference with normalizing flowで急に出てきた変換の関数の形の意味がよくわかった点は勉強になった．

2018-07-12

ボルツマンマシンについて勉強したメモ　まとめ

ボルツマンマシンについて勉強したメモのまとめ．

その1 ボルツマンマシンのモデルとその学習方程式について

その2 学習に用いるMCMCについて

その3 学習に用いる平均場近似について

その4 制約ボルツマンマシンのモデルとその学習方程式について

その5 学習に用いるコントラスティブダイバージェンスについて

その他，ボルツマンマシンを勉強するきっかけになった論文．

Discrete Variational Autoencoders

GumBolt: Extending Gumbel trick to Boltzmann priors

(上二つの論文もそのうちまとめる予定)

2018-07-12

Variational Inference with Normalizing Flowsを読んだのでメモ

生成モデル変分推論 normalizing flow

はじめに

Variational Inference with Normalizing Flowsを読んだのでメモ．OpenAIからでたGlowという生成モデルの論文が話題になっていて面白そうだったので関連知識を一から勉強してみる．まず手始めにnormalizing flowを使った変分推論について．

Variational Inference

この論文で問題としているのは生成モデルにおいて変分推定を適用するときに，近似事後分布の表現力が低いというもの．より良い近似を使えば強い生成モデルが作れるはずなのでその方法としてNormalizing Flowを使った変分推定を提案しますとのこと．

変分推定は以下のようにJensenの不等式を使って変分下限（または変分下界）を求め，それを最大化しようというもの．

$\displaystyle \log p_\theta \\ \displaystyle = \log\int p_\theta(\mathbf{x}|\mathbf{z})p(\mathbf{z})d\mathbf{z} \\ \displaystyle \geq -\mathbb{D}_{KL}[ q_\phi(\mathbf{z}|\mathbf{x})||p(\mathbf{z})]+\mathbb{E}_q[\log p_\theta(\mathbf{x}|\mathbf{z})]$

ここで $\mathbf{x}$ は観測値で $\mathbf{z}$ は潜在変数． $\theta$ はモデルのパラメータを表し， $q_\phi(\mathbf{z}|\mathbf{x})$ は近似事後分布を表す．近似事後分布は因子分解可能であることなどを仮定して解析的に扱いやすくしているため，元の分布を完全に表現することは不可能である．一応，近似事後分布にすることでの表現力の低下以外の問題として第2項目の勾配の計算コストが高いことも問題としてあげられる(VAEのreparameterizationはそこをうまく扱ったもの)．

Deep Latent Gaussian Models（DLGM）

Deep Neural Networksにおいて各層が潜在変数 $z_l$ を出力するとする．ここで $z_l$ は $l$ 番目の層に対応する潜在変数．このとき各々の潜在変数は自分の上位の層以外とは独立と仮定している．層の数を $L$ とするとDLGMにおける同時分布は各潜在変数の独立性より

$\displaystyle p(\mathbf{x},\mathbf{z}_1,...,\mathbf{z}_L)=p(\mathbf{x}|f_0(\mathbf{z}_1))\prod_lp(\mathbf{z}_l|f_l(\mathbf{z}_{l+1}))$

と表される．各々の潜在変数は事前分布として標準ガウス分布 $p(\mathbf{z}_l)=\mathcal{N}(\mathbf{0},\mathbf{I})$ が仮定されている．DLGMはDeepな分表現力は高いが，結局は分散が対角成分しかないような単純化されたガウス分布を用いるため真の分布を表現するには無理がある．

Normalizing flow

Normalizing flowは近似事後分布の作り方の一つで表現力の高い近似事後分布を与えることができる．イメージとしては確率分布の変数変換を繰り返し用いることで単純な分布を複雑な分布へと変換できるというもの．

まず，確率変数の変数変換として基本的な規則について．逆変換を持つ連続な写像 $f:\mathbb{R}^d\rightarrow\mathbb{R}^d,\:f^{-1}=g$ を考える．分布 $q(\mathbf{z})$ を持つ確率変数 $\mathbf{z}$ に対し変換を行った確率変数 $\mathbf{z}'=f(\mathbf{z})$ を考える．このとき， $\mathbf{z}'$ に対する確率密度を $q'(\mathbf{z}')$ とすると

$\displaystyle q'(\mathbf{z}')=q(\mathbf{z})\left|\mathrm{det}\frac{\partial f^{-1}(\mathbf{z}')}{\partial\mathbf{z}'}\right|$

と与えられる．詳しいことは「確率密度　変数変換」などでググるとたくさん出てくるので割愛．一応こちらのブログはnormalizing flowの内容含めわかりやすい．

ここで関数 $f$ が $\mathbf{z}_K=f_K\circ\dots\circ f_2\circ f_1(\mathbf{z}_0)$ のように $K$ 個の関数の合成関数になっていると考える．すると上記の変数変換を繰り返し用いれば，

$\displaystyle q_K(\mathbf{z}_K)=q_{K-1}(\mathbf{z}_{K-1})\left|\mathrm{det}\frac{f^{-1}_K(\mathbf{z}_K)}{\mathbf{z}_K}\right|=\dots=q_0(\mathbf{z}_0)\prod_k\left|\mathrm{det}\frac{f^{-1}_{k}(\mathbf{z}_{k})}{\mathbf{z}_k}\right|=q_0(\mathbf{z}_0)\prod_k\left|\mathrm{det}\frac{f_{k}(\mathbf{z}_{k})}{\mathbf{z}_k}\right|^{-1}$

となる．最後の等式は逆関数定理から成り立つ．さらにこれの対数をとると

$\displaystyle \log q_K(\mathbf{Z}_K)=\log q_0(\mathbf{z}_0)-\sum_k\log\mathrm{det}\left|\frac{f_k}{\mathbf{z}_k}\right|$

となる．ただし関数 $f$ の引数を省略したことに注意．また論文では $\log$ ではなく $\ln$ となっていたが，本質的には変わらないはずなのでここでは $\log$ を使って表現する．ここで注目すべきは $q_K$ を計算するのに $q_0$ しか必要がないこと．これによって適当な関数 $h(\mathbf{z})$ の $q_K$ に関する期待値の計算も $\mathbb{E}_{q_K}[ h(\mathbf{z})]=\mathbb{E}_{q_0}[h(f_K\circ\dots\circ f_1(\mathbf{z}_0))]$ のように行うことが可能．

初期の分布 $q_0$ は従来通りindependent Gaussianなど使っても，DLGMと違い各潜在変数は何か分布を仮定しているわけではないので表現力が高く，かつ $q_0$ を使って計算できるため解析的に扱いやすいという利点がある．

上記では有限のフローを扱ったが，有限があれば無限も考えたくなるのが性というもの．ただし，以下の内容は前提知識が足りず理解しきれていないのでほぼ論文の和訳メモになっていることに注意．

まずLangevin Flow（ランジュバンフロー）について．Langevin Flowでは以下のランジュバン方程式（Langevin stochastic differential equation）によって確率変数 $\mathbf{z}$ が記述される．

$\displaystyle d\mathbf{z}(t)=\mathbf{F}(\mathbf{z}(t),t)dt+\mathbf{G}(\mathbf{z}(t),t)d\mathbf{\xi}(t)$

ここで $d\mathbf{\xi}(t)$ は $\mathbb{E}[\mathbf{\xi}_i(t)]=0$ ， $\mathbb{E}[\mathbf{\xi}_i(t)\mathbf{\xi}_j(t')]=\delta_{i,j}\delta(t-t')$ で与えられるWiener過程， $\mathbf{F}$ はドリフトベクトルで $\mathbf{D}=\mathbf{G}\mathbf{G}'$ は拡散行列を表す．ここで初期の分布が $q_0(\mathbf{z})$ である確率変数 $\mathbf{z}$ が上記の方程式により変化していく場合，密度分布の変化の規則はFokker-Planck方程式で与えられる．すなわち，

$\displaystyle \frac{\partial}{\partial t}q_t(\mathbf{z})=-\sum_i\frac{\partial}{\partial z_i}[F_i(\mathbf{z},t)q_t]+\frac{1}{2}\sum_{i,j}\frac{\partial^2}{\partial z_i\partial z_j}[D_{ij}(\mathbf{z},t)q_t]$

として変化を記述できる．機械学習の文脈では $F(\mathbf{z},t)=-\nabla_z\mathcal{L}(\mathbf{z})$ ， $G(\mathbf{z},t)=\sqrt{2}\delta_{ij}$ としたLangevin flowを用いる．ただし， $\mathcal{L}(\mathbf{z})$ は正規化されていないモデルの対数密度関数．重要なこととして，定常状態における解 $q_\infty(\mathbf{z})$ はボルツマン分布になることが知られている．つまり，適当な分布 $q_0(\mathbf{z})$ を持つ確率変数 $\mathbf{z}_0$ をLangevin flowで変換していくと最終的には $q_\infty(\mathbf{z})\propto\exp(-\mathcal{L}(\mathbf{z}))$ となる．ボルツマンマシンが任意の分布を近似可能だったことを思い出せばnormalizing flowも任意の分布を近似可能と言える（のかな？）．

Inference with Normalizing Flows

ここでは有限のnormalizing flowを用いた推定について考える．基本的には変数変換における関数 $f$ は逆変換可能なものかつヤコビアンンが効率的に計算できるものでなければならない．そこで以下のような変換を考える．

$\displaystyle f(\mathbf{z})=\mathbf{z}+\mathbf{u}h(\mathbf{w}'\mathbf{z}+b)$

ここで， $\lambda=\{\mathbf{w}\in\mathbb{R}^D,\mathbf{u}\in\mathbb{R}^D,b\in\mathbb{R}\}$ は学習パラメータで， $h$ はsmoothでelement-wiseな非線形関数．上記の変換はヤコビアンの項を $\mathcal{O}(D)$ 時間で計算ができる．具体的には非線形関数の導関数を $h'$ とすると，右辺第二項目の $\mathbf{z}$ に関する微分は

$\displaystyle \psi(\mathbf{z})=h'(\mathbf{w}'\mathbf{z}+b)\mathbf{w}$

とかけるので，行列式は

$\displaystyle \mathrm{det}\left|\frac{\partial f}{\partial\mathbf{z}}\right| =|\mathrm{det}(\mathbf{I}+\mathbf{u}\psi(\mathbf{z})'|=|1+\mathbf{u}'\psi(\mathbf{z})|$

と計算可能．よって計算コストのオーダーが $\mathcal{O}(D)$ になる．

上記で定義された変換 $f$ を，元の有限なnormalizing flowの式に当てはめれば，

$\displaystyle \log q_K(\mathbf{z}_K)=\log q_0(\mathbf{z}_0)-\sum_k\log|1+\mathbf{u}_k'\psi_k(\mathbf{z}_k)|$

が導かれる．このflowは $\mathbf{w}'\mathbf{z}+b=0$ で与えられる超平面に垂直な方向での収縮と拡大を与えており，それにより初期密度 $q_0$ を変化させると考えられる．そのためこれをplanar flowと定義する．

また別なアプローチとして，初期密度を確率変数 $\mathbf{z}_0$ 周りで変化させるような変換群を考える．そのような変換群は

$\displaystyle f(\mathbf{z})=\mathbf{z}+\beta h(\alpha,r)(\mathbf{z}-\mathbf{z}_0) \\ \displaystyle \mathrm{det}\left|\frac{\partial f}{\partial\mathbf{z}}\right|=[1+\beta h(\alpha,r)]^{(d-1)}[1+\beta h(\alpha,r)+h'(\alpha,r)r)]$

として考えられる．ただし， $r=|\mathbf{z}-\mathbf{z}_0|$ ， $h(\alpha,r)=1/(\alpha+r)$ とし，パラメータとして $\lambda=\{\mathbf{z}_0\in\mathbb{R}^D,\alpha\in\mathbb{R},\beta\in\mathbb{R}\}$ ．この群による変換もまた線形時間で行列式を計算可能である．このflowは放射状の(radial)収縮と拡大を与えるのでradial flowと定義する．planer flowとradial flowで球面ガウス分布を変換した例が論文のFigure 1．より複雑な分布になっていることがわかる．

ここで注意としてplaner flow，radial flowの全てが逆変換可能な関数になっているわけではなく，一定の条件下のみで逆変換可能となる．詳細なことはAppendixに．

ここで変分推定における目的関数にnormalizing flowを適用すると，

$\displaystyle -\log p_\theta(\mathbf{x}) \\ \displaystyle \geq \mathbb{D}_{KL}[ q_\phi(\mathbf{z}|\mathbf{x})||p(\mathbf{z})]-\mathbb{E}_q[\log p_\theta(\mathbf{x}|\mathbf{z})] \\ \displaystyle = \mathbb{E}_{q_\phi(z|x)}[\log q_\phi(\mathbf{z}|\mathbf{x})-\log p(\mathbf{z})-\log p(\mathbf{x}|\mathbf{z})] \\ \displaystyle =\mathbb{E}_{q_0(z_0)}[\log q_K(\mathbf{z}_K)-\log p(\mathbf{x},\mathbf{z}_K)]\\ \displaystyle =\mathbb{E}_{q_0(z_0)}[\log q_0(\mathbf{z}_0)] -\mathbb{E}_{q_0(z_0)}[\log p(\mathbf{x},\mathbf{z}_K)]-\mathbb{E}_{q_0(z_0)}\left[\sum_k\log |1+\mathbf{u}'_k\psi_k(\mathbf{z}_k)|\right]$

となる．4行目はnormalizing flowでは $q_K$ の期待値が $q_0$ の期待値でかけることを利用しており，最後の行はplaner flowを代入した．

この論文ではVAEと同様にニューラルネットを使ったモデルを使っている．encoderとdecoderの間で関数 $f$ を適用していて， $k$ の数を変えて実験を行なっている．NICEと比較した結果の可視化が論文のFigure 3(a)(b)(c)．提案手法の方が真の分布をうまく表現できていることが分かる．

ちなみに学習に関して最終的な変分下限の式には $p(\mathbf{x},\mathbf{z}_K)$ という式が現れているが，目的関数としては $p(\mathbf{x}|\mathbf{z}_K),p(\mathbf{z}_K)$ で分かれるはず． $p(\mathbf{z}_K)$ はtable 1の分布を使ったと書いてあって， $p(\mathbf{x}|\mathbf{z}_K),p(\mathbf{z}_K)$ どちらも勾配は元のVAEと同様モンテカルロ推定で計算する．あと目的関数の $p(\mathbf{x},\mathbf{z}_K)$ 項には温度パラメータ $\beta$ を係数としてかけて $\beta_t=\min(1,0.01+t/10000)$ としてアニーリングしていくらしい．

まとめ

Normalizing flow非常に賢いなという思い．とりあえずGLOWの元になっていてこの論文でも比較に使われていたNICEを次は読む．

2018-07-10

ボルツマンマシンについて勉強した　その５

機械学習ボルツマンマシン CD法

はじめに

RBMの最適化におけるサンプリング方法について勉強したのでメモ．教科書はいつも通り．

機械学習スタートアップシリーズこれならわかる深層学習入門 (KS情報科学専門書)

作者:瀧雅人
出版社/メーカー: 講談社
発売日: 2017/10/21
メディア: 単行本（ソフトカバー）

これでひとまずボルツマンマシンは終了

ブロック化ギブスサンプリング

RBMの条件付き独立性を有効に活用して可視変数と隠れ変数を交互にギブスサンプリングしようというもの．なので基本的にはギブスサンプリングを行うだけ．具体的には可視変数 $\mathbf{v}$ をランダムに初期化して $h_j$ を $P(h_j|\mathbf{v},\theta)$ からサンプリング．全ての隠れ変数がサンプリングできたらそれらを用いて $v_i$ を $P(v_i|\mathbf{h},\theta)$ サンプリング．以下，この手順を繰り返して連鎖を走らせた後にサンプルを取得する．各々の分布からサンプリングは，ギブスサンプリングで説明した一様乱数を使う方法で行う．

コントラスティブ・ダイバージェンス法（CD法）

RBMになったからといって連鎖を走らせる時間が短くなるわけではないので相変わらずギブスサンプリングの計算コストは高い．実は嬉しいことにRBMではギブスサンプリングの改良版であるCD法が適用できるとのこと．

ブロック化ギブスサンプリングとの違いの一つは， $\mathbf{v}$ の初期値をランダムではなくデータの経験分布 $q(\mathbf{v})$ からのサンプル(つまりは適当に選んだ訓練サンプル)にするというもの．つまり，

$\displaystyle \mathbf{v}(0)\sim q(\mathbf{v}),\mathbf{h}(0)\sim P(\mathbf{h}|\mathbf{v}^{(n)},\theta),...\\ \displaystyle \mathbf{h}(T-1)\sim P(\mathbf{h}|\mathbf{v}(T-1),\theta),\mathbf{v}(T)\sim P(\mathbf{v}|\mathbf{h}(T-1),\theta)$

と0から $T$ までサンプリングする．ただし，上記の式ではサンプリングされた回数を括弧で表していることと，各変数一つ一つについてサンプリングを書くと煩雑なのでベクトル表記でまとめてることに注意．で，どうやらこうすると連鎖を長く走らせなくても実用上は問題ない精度でサンプリングができるらしい．なんなら $T=1$ でもいいとのこと．

ギブスサンプリングとのもう一つの違いとして，RBMの学習方程式のモデルに関する期待値の計算方法がある．CD法ではモデルの期待値を

$\displaystyle \langle v_ih_h\rangle_{model}\approx \mathbb{E}_{P(\mathbf{h}|\mathbf{v}(T),\theta)}[ v_ih_j]$

と計算する．さらにすごいことに勾配を求めるのにサンプリング平均 $\mathbb{E}_{P(\mathbf{v}(T)|\mathbf{v}^{(n)})}$ を計算する必要がなく一つのサンプルから計算しても良いらしい．以上から，CD法における各学習係数の勾配をまとめると

$\displaystyle \delta b_i\propto v_i^{(n)}-v_i(T) \\ \displaystyle \delta c_j\propto P(h_j=1|\mathbf{v}^{(n)},\theta)-P(h_j=1|\mathbf{v}(T),\theta) \\ \displaystyle \delta w_{ij}\propto v_i^{(n)}P(h_j=1|\mathbf{v}^{(n)},\theta)-v_i(T)P(h_j=1|\mathbf{v}(T),\theta)$

となる．CD法ではモデルに対する期待値が， $P(\mathbf{h}|\mathbf{v}(T),\theta)$ による期待値に変わっていることから条件付き独立性が成り立ち，見通しが良くなっている．具体的には $\delta b_i$ では $\sum_\mathbf{h}P(\mathbf{h}|\mathbf{v}(T),\theta)=\prod_j\sum_{h_j}P(h_j|\mathbf{v}(T),\theta)=1$ の関係を用いており，その他二つにおいては前の記事で解説したRBMの学習方程式の導出過程で出てきた変換を使った．あとは勾配上昇法を使って学習を行えばいい．

教科書に載っていたこととして，文献によって重みに関する勾配の右辺第２項における $v_i(T)$ を $P(v_i=1|h(T-1),\theta)$ として $v_i(T)$ を直接使わない場合もあるらしい．また，ミニバッチでやる場合にはミニバッチ分だけ初期値と連鎖を用意してサンプリングを行うとのこと．

まとめ

長かったボルツマンマシンの勉強も終わり．ちなみに，ブログにまとめるのは疲れたので割愛しているが．教科書にはこんなに思い切った近似しているCD法がなぜうまくいくかについても解説している．

2018-07-09

ボルツマンマシンについて勉強した　その４

はじめに

今回は制限付きボルツマンマシン（Restricted Boltzmann Machine, RBM）について．教科書はいつもの．

機械学習スタートアップシリーズこれならわかる深層学習入門 (KS情報科学専門書)

作者: 瀧雅人
出版社/メーカー: 講談社
発売日: 2017/10/21
メディア: 単行本（ソフトカバー）
この商品を含むブログ (1件) を見る

制限付きボルツマンマシン（RBM）

ボルツマンマシンの問題点として期待値計算が難しいという点があげられた．隠れ変数を導入するとモデルに対する期待値だけでなくデータに関する期待値の項も計算量が爆発してしまう．そこで隠れ変数を導入しつつ効率化をはかるために制限を設けたものが制限付きボルツマンマシン．制限付きボルツマンマシン隠れ変数同士，可視変数同士に結合を持たないようにグラフ構造に制限を入れたもの．

隠れ変数同士，可視変数同士に結合を持たないので隠れ変数と可視変数を陽に書き表したエネルギー関数は以下のようになる．

$\displaystyle \Phi(v,h,\theta) = -\mathbf{b}'\mathbf{v}-\mathbf{c}'\mathbf{h}-\mathbf{v}'\mathbf{W}\mathbf{h}$

通常の隠れ変数を持つボルツマンマシンに比べて隠れ変数同士，可視変数同士のエッジに関する重みの項が減ってエネルギー関数が簡潔になった．ちなみにかなり強い制約を与えたが，隠れ変数の数次第で任意の分布を任意の精度で近似できるらしい．

ここでRBMの独立性をもう少し見るために，隠れ変数が与えられた時の可視変数の分布を考える．同時分布を周辺分布で割れば良いので，

$\displaystyle P(\mathbf{v}|\mathbf{h},\theta) \\ \displaystyle =\frac{P(\mathbf{v},\mathbf{h}|\theta)}{\sum_\mathbf{v}P(\mathbf{v},\mathbf{h}|\theta)} \\ \displaystyle =\frac{\exp(\mathbf{b}'\mathbf{v}+\mathbf{c}'\mathbf{h}+\mathbf{v}'W\mathbf{h})}{\sum_\mathbf{v}\exp(\mathbf{b}'\mathbf{v}+\mathbf{c}'\mathbf{h}+\mathbf{v}'W\mathbf{h})} \\ \displaystyle =\prod_i\frac{\exp(b_iv_i+\mathbf{c}'\mathbf{h}+\sum_jv_iW_{ij}h_j)}{\sum_{v_i=0,1}\exp(b_iv_i+\mathbf{c}'\mathbf{h}+\sum_jv_iW_{ij}h_j)} \\ \displaystyle =\prod_i\frac{\exp(\mathbf{c}'\mathbf{h})\exp(b_iv_i+\sum_jv_iW_{ij}h_j)}{\exp(\mathbf{c}'\mathbf{h})\sum_{v_i=0,1}\exp(b_iv_i+\sum_jv_iW_{ij}h_j)} \\ \displaystyle =\prod_i\frac{\exp(b_iv_i+\sum_jv_iW_{ij}h_j)}{1+\exp(b_i+\sum_jW_{ij}h_j)}$

と計算ができ，積の形で表されていることから隠れ変数が固定されれば可視変数は条件付き独立になっていることがわかる．隠れ変数についても同様の計算から条件付き独立性を見ることができる．独立になることで計算が飛躍的に簡単になり，そこには近似も入らないため，非常に嬉しいというのがRBM．おかげでギブスサンプリングの実装もシンプルになる．

学習方法

通常のボルツマンマシンと同様，対数尤度の各変数に関する微分を求める．基本的にはバイアスの項が一つ増えただけと考えられる．

可視変数にかかるバイアス $b_i$ に関する微分は隠れ変数を周辺化した対数尤度を $L(\theta)$ とし，途中までの式の導出は前の記事に譲ると

$\displaystyle \frac{1}{N}\frac{\partial L(\theta)}{\partial b_i} \\ \displaystyle =\frac{1}{N}\sum_nv_i^{(n)}\sum_\mathbf{h}P(\mathbf{h}|\mathbf{v},\theta)-\langle v_i\rangle_{model} \\ \displaystyle =\frac{1}{N}\sum_nv_i^{(n)}-\langle v_i\rangle_{model}$

と計算できる．ただし，最後の行は $\sum P(\mathbf{h}|\mathbf{v},\theta)=1$ であることと， $q(\mathbf{v})$ が経験分布であることを利用した．隠れ変数にかかるバイアス $c_j$ に関しても途中の導出まで前の記事に計算を譲って計算をすると，

$\displaystyle \frac{1}{N}\frac{\partial L(\theta)}{\partial c_j} \\ \displaystyle =\frac{1}{N}\sum_n\sum_\mathbf{h}h_jP(\mathbf{h}|\mathbf{v},\theta)-\langle v_i\rangle_{model} \\ \displaystyle =\frac{1}{N}\sum_n\sum_\mathbf{h_j}h_j\prod_jP(\mathbf{h_j}|\mathbf{v},\theta)-\langle v_i\rangle_{model} \\ \displaystyle =\frac{1}{N}\sum_n\sum_{h_j}h_jP(h_j|\mathbf{v}^{(n)},\theta)\prod_{j'\neq j}\sum_{h_{j'}}P(h_{j'}|\mathbf{v}^{(n)},\theta)-\langle v_i\rangle_{model} \\ \displaystyle =\frac{1}{N}\sum_nP(h_j=1|\mathbf{v}^{(n)},\theta)-\langle v_i\rangle_{model}$

となる．ただし，RBMの条件付き独立性と $\sum_{h_{j'}}P(h_{j'}|\mathbf{v}^{(n)},\theta)=1$ を使った．最後に重みに関する計算に関しては，バイアス $b_i$ に関する微分における各項の係数が変わるだけなので，

$\displaystyle \frac{1}{N}\frac{\partial L(\theta)}{\partial w_{ij}} \\ \displaystyle =\frac{1}{N}\sum_nv_i^{(n)}\sum_\mathbf{h}h_jP(\mathbf{h}|\mathbf{v},\theta)-\langle v_ih_j\rangle_{model} \\ \displaystyle =\frac{1}{N}\sum_nv_i^{(n)}P(h_j=1|\mathbf{v}^{(n)},\theta)-\langle v_ih_j\rangle_{model}$

となる．ただし，最後の行は先ほどと同様に条件付き独立の性質を使った．

以上からRBMの学習方程式は

$\displaystyle \frac{1}{N}\sum_nv_i^{(n)}=\langle v_i\rangle_{model} \\ \displaystyle \frac{1}{N}\sum_nP(h_j=1|\mathbf{v}^{(n)},\theta)=\langle h_j\rangle_{model} \\ \displaystyle \frac{1}{N}\sum_nv_i^{(n)}P(h_j=1|\mathbf{v}^{(n)},\theta)=\langle v_ih_j\rangle_{model}$

となる．隠れ変数を持つ通常のボルツマンマシンと比較するとRBMの学習方程式の左辺はサンプルに関する和を計算するだけになっておりとても簡略化されていることがわかる．

まとめ

RBMの学習方程式の導出まで．あとは学習方程式の右辺を効率よく計算する方法に関してのみ．それが終わったらもともと読みたかった論文を読む予定．まさかこんなかかるとは思わなかった．

2018-07-07

ボルツマンマシンについて勉強した　その３

ボルツマンマシン機械学習平均場近似

はじめに

ボルツマンマシンについて勉強したのでメモその３．今回は本題と少し外れて平均場近似について．教科書は例によって以下．

機械学習スタートアップシリーズこれならわかる深層学習入門 (KS情報科学専門書)

作者: 瀧雅人
出版社/メーカー: 講談社
発売日: 2017/10/21
メディア: 単行本（ソフトカバー）
この商品を含むブログ (1件) を見る

平均場近似

ボルツマンマシンは各ノード同士を掛け合わせる $\sum_{i,j}w_{ij}x_ix_j$ があるため期待値の計算が難しい．逆にこの項が存在しなければボルツマンマシンの分布は

$\displaystyle P(x|\theta)=\frac{\exp(\sum_ib_ix_i)}{Z(\theta)}$

となり，簡単化する．平均場近似では各変数が独立となっている分布（テスト分布）の集合を用意し，その中で最も元の分布に近いものをボルツマンマシンの近似として採用しようというもの（ $w_{ij}=0$ よりかは表現力のあるもの）．

まず，テスト分布を $Q(x)$ とするとテスト分布の各変数は独立であると仮定していることから，

$\displaystyle Q(x)=\prod_iQ_i(x_i)$

とかける．このようにかける分布の中から，ボルツマンマシンの分布とのKLダイバージェンスが最小となる分布を探すというのが目的．つまり，

$\displaystyle D_{KL}(Q||P)=\sum_xQ(x)\log\frac{Q(x)}{P(x|\theta)} \\ \displaystyle s.t.\:\sum_{x_i=0,1}Q_i(x_i)=1$

の最小化を行う．等式制約条件つきなのでラグランジュの未定乗数法を使えば

$\displaystyle L(Q_i;\lambda_i)=D_{KL}(Q||P)+\sum_i\lambda_i\left(\sum_{x_i=0,1}Q_i(x_i)-1\right)$

という目的関数が得られる．ただし， $\lambda_i$ はラグランジュ未定乗数．この目的関数を $Q(x_i)$ と $\lambda_i$ で変分することで得られる2つの方程式の解を求めれば分布が決まる． $\lambda_i$ に関する変分は元の拘束条件を表すため， $Q_i(x_i)$ についての変分を考える．まず，変分を考えやすくするため， $Q(x)$ の独立性を使って，目的関数のKLダイバージェンスの項を以下のように変形する．

$\displaystyle D_{KL}(Q||P) \\ \displaystyle =\sum_xQ(x)\log\frac{Q(x)}{P(x|\theta)} \\ \displaystyle =\sum_x\prod_jQ_j(x_j)(\sum_k\log Q_k(x_k)-\log P(x|\theta))$

上記を $Q_i(x_i)$ について変分すると

$\displaystyle \frac{\partial D_{KL}}{\partial Q_i(x_i)} \\ \displaystyle =\sum_{\mathbf{x}_{-I}}\prod_{j\neq i}Q_j(\sum_k\log Q_k(x_k)-\log P(x|\theta)) + \sum_{\mathbf{x}_{-i}}\prod_jQ_j(x_j)\frac{1}{Q_i(x_i)} \\ \displaystyle =\sum_{\mathbf{x}_{-i}}\prod_{j\neq i}Q_j(x_j)(\sum_k\log Q_k(x_k)-\log P(x|\theta)+1)$

となる．よって，目的関数を $Q_i(x_i)$ について変分したものは

$\displaystyle 0=\frac{\partial L}{\partial Q_i(x_i)}=\sum_{\mathbf{x}_{-i}}\prod_{j\neq i}Q_j(x_j)(\sum_k\log Q_k(x_k)-\log P(x|\theta)+1)+\lambda_i$

ともとまる．ここで右辺の第１項の $k=i$ に対して， $\log Q_i(x_i)$ が $\sum_{\mathbf{x}_{-i}}$ に関係ないことと制約条件を使って

$\displaystyle \sum_{\mathbf{x}_{-i}}\prod_{j\neq i}Q_j(x_j)\log Q_i(x_i)=\log Q_i(x_i)\prod_{j\neq i}\sum_{\mathbf{x}_{-i}}Q_j(x_j)=\log Q_i(x_i)$

と変形できる．また $k\neq I$ の部分に関しては $\sum_ {\mathbf{x}_{-i}}$ で期待値を取っているため確率変数によらない定数になる．さらに第二項目に関してもボルツマンマシンの具体形を代入すると

$\displaystyle \sum_{\mathbf{x}_{-i}}\prod_{j\neq i}Q_j(x_j)\log P(x|\theta) \\ \displaystyle =\sum_{\mathbf{x}_{-i}}\prod_{j\neq i}Q_j(x_j)(\sum_mb_mx_m+\sum_{m,n}w_{mn}x_mx_n -\log Z(\theta)) \\ \displaystyle =(b_i+\sum_{\mathbf{x}_{-i}}\prod_{k\neq i}Q_k(x_k)\sum_kw_{ki}x_k)x_i+\sum_{\mathbf{x}_{-i}}\prod_{j\neq i}Q_j(x_j)\left(\sum_{m\neq i}b_mx_m+\sum_{m,n\neq i}w_{mn}x_mx_n\right)-\log Z(\theta) \\ \displaystyle =(b_i+\sum_kw_{ki}\mu_k)x_i+\sum_{j\neq i}b_j\mu_j+\sum_{j,n\neq i}w_{jn}\mu_j\mu_n-\log Z(\theta)$

と変形できる．ただし，2行目から3行目において右辺第１項と同様に $x_i$ の項に関する関係性を使った．さらに3行目で導入された $k$ は，ノード $i$ に隣接するノードのインデックスを示す．また，4行目で導入された $\mu$ は期待値計算によって出てくる定数で，以下で定義される．

$\displaystyle \mu_j=\sum_xx_jQ(x)=\sum_{x_j=0,1}x_jQ_j(x_j)$

これがいわゆる平均場で，期待値計算になっていることがわかる．よって今までの変形を組み合わせて確率変数 $x_i$ に注目すると

$\displaystyle \log Q_i(x_i)-\left(b_i+\sum_kw_{ki}\mu_k\right)x_i+\mathrm{const.}+\lambda_i=0$

と簡単にかける．よって求める分布は

$\displaystyle Q_i(x_i)\varpropto \exp(-\lambda_i)\exp(b_ix_i+\sum_kw_{ki}\mu_kx_i)$

となる．ここで上記分布は拘束条件 $\sum_iQ_i(x_i)=1$ を満たすので， $x_i=0,1$ を代入して計算することでラグランジュ未定乗数は簡単にもとまって，求める分布は

$\displaystyle Q_i(x_i)=\frac{\exp(b_ix_i+\sum_kw_{ki}\mu_kx_i)}{1+\exp(b_i+\sum_kw_{ki}\mu_k)}$

となる．結局のところ，注目している変数以外は平均値で置き換えてしまおうというのがこの近似法．ここで $\mu_i$ を書き下して見ると， $x_i=0,1$ から $x_i=1$ の項のみが残り

$\displaystyle \mu_i=\frac{\exp(b_i+\sum_kw_{ki}\mu_k}{1+\exp(b_i+\sum_kw_{ki}\mu_k}$

が得らる．これを条件に平均場 $\mu$ の値を決定する必要がある．ちなみにこれは自己無撞着方程式や平均場方程式と呼ばれるらしい．しかし問題は，これは非線形な連立方程式になるため解析的に解くのは困難ということ．実践的には期待値をランダムに初期化して，得られた期待値を使って実際に計算することで期待値の値を更新するという手続きを繰り返して求めるらしい．

最終的にボルツマンマシンにおいては $\langle x_i\rangle_{model}\approx\mu_i$ ， $\langle x_ix_j\rangle_{model}\approx\mu_i\mu_j$ と近似するだけの話．

まとめ

平均場近似について勉強した．教科書によると平均場近似は簡単化しすぎているせいで変数間の相関を失って精度が落ちるらしい．独立の要請を緩めたベーテ近似やクラスター変分法という方法もあるらしいがとりあえず教科書にはないので機会があれば勉強してみようかくらい．

2018-07-03

ボルツマンマシンについて勉強した　その２

機械学習ボルツマンマシン MCMC

はじめに

ボルツマンマシンについて勉強したのでメモその２．主にボルツマンマシンの学習におけるサンプリングについて．教科書は以下．機械学習プロフェッショナルシリーズの深層学習に比べてサンプリングに関わる説明が手厚くサンプリングを初めて学ぶ身としてはありがたい．

機械学習スタートアップシリーズこれならわかる深層学習入門 (KS情報科学専門書)

作者: 瀧雅人
出版社/メーカー: 講談社
発売日: 2017/10/21
メディア: 単行本（ソフトカバー）
この商品を含むブログ (1件) を見る

サンプリングについて

基本的にはサンプリングを使う意味として，ボルツマンマシンの学習時に期待値計算がしんどいので，分布からサンプリングしてその平均を期待値とすれば計算量が減って嬉しいというもの．サンプリングで得られたサンプルの平均を期待値として使うには各サンプルがi.i.dに所望の分布からサンプルされてないと真の期待値から大きく外れて役に立たないのでそのようなサンプルを生成するうまいやり方がないかというのがサンプリング手法のポイントとなるところ．基本的にサンプルがi.i.dに得られていたら大数の法則からサンプル数が無限大の時に真の期待値に収束する．

マルコフ連鎖

ここではモンテカルロ法の一種のマルコフ連鎖モンテカルロ法（Markov Chain Monte Carlo method, MCMC）について考える．

まずMCMCで重要なマルコフ連鎖について．マルコフ連鎖はマルコフ性を満たす確率過程のことで，現在の状態が1時刻前の状態だけから決定されるというもの．つまり $P(x_t|x_{t-1},...,x_1)=P(x_t|x_{t-1})$ というもの．マルコフ過程においてある時刻 $t$ の確率分布 $P(x_t)$ は周辺化を用いて $P(x_t)=\sum_{x_{t-1}}P(x_t|x_{t-1})P(x_{t-1})$ と表現できる．さらに同時分布はマルコフ性から $P(x_t,...x_1)=P(x_t|x_{t-1})...P(x_2|x_1)P(x_1)$ と表される．

教科書ではマルコフ連鎖の応用としてページランクが紹介されていたが今回の勉強目的とは関係ないのでここでは割愛．

ここで，確率変数 $x$ が取りうる状態を $\alpha\in\{1,2,..\}$ とすると， $\mathbf{\pi}_t=\left(P(x_t=1),P(x_t=2),...\right)$ のように時刻 $t$ における各状態の実現確率 $P(x_t=1),P(x_t=2),...$ を並べたベクトルを考える．このベクトルを状態確率分布と呼ぶ．同様に状態 $a$ から状態 $b$ への遷移確率 $P(x_t=a|x_{t-1}=b)$ を $ab$ 成分とする行列を遷移確率行列 $\mathbf{T}$ と呼ぶ．すると時刻 $t$ の状態確率分布は $t-1$ の状態分布を用いて，

$\displaystyle \mathbf{\pi}_t=\mathbf{T}\mathbf{\pi}_{t-1}$

と表現できる．ここで均一なマルコフ連鎖（遷移確率行列が時刻 $t$ によらない）を考えれば，ある時刻 $t$ における状態確率分布は $\mathbf{\pi}_t=\mathbf{T}^t\mathbf{\pi}_0$ となる．ここで $t\rightarrow\infty$ の時，もはや遷移確率行列により状態確率分布が変化しなくなった場合，その分布を平衡分布と呼ぶ．マルコフ連鎖は必ず平衡分布に収束するわけではなく規約性と非周期性を満たす時に収束する（詳細は教科書参照）．収束の条件式として以下の詳細つり合い条件が与えられる．

$\displaystyle T_{ba}P(x=a)=T_{ab}P(x=b) \\ \displaystyle P(x=a)=\sum_bT_{ab}P(x=b)$

2行目は1行目の両辺を $b$ について和を取ることで $\sum_bT_{ba}=1$ の関係を用いて得られる．

マルコフ連鎖モンテカルロ法

ようやく本題のMCMC．MCMCはサンプリングしたい分布 $P(x)$ が複雑でサンプリングが難しい時に， $P(x)$ に収束するマルコフ連鎖からサンプリングしようというもの．そうすればマルコフ連鎖が収束した時には，欲しい分布からサンプリングしているのと同じとみなせる．問題はマルコフ連鎖をどう設計するかで，収束が早いほど嬉しいのは言うまでもない．ここでは例としてギブスサンプリングを考える．

ギブスサンプリングは簡素な方法でマルコフ連鎖を設計するため，アルゴリズムも単純で汎用性が非常に高い．いまボルツマンマシンの確率変数が $\mathbf{x}=(x_1,...,x_M)$ と定義されている時，ギブスサンプリングでは以下の条件つき分布からサンプリングを行う．

$\displaystyle P(x_i|x_1,...,x_{i-1},x_{i+1},...,x_M)$

要は，全ての確率変数を一気にサンプリングするのではなく，他の確率変数を条件として一つの確率変数のみをサンプリングしようというもの．全ての変数を一気にサンプリングするのは難しいが，一つずつなら楽にできるというのがギブスサンプリングのいいところ．この時サンプリングする順番は任意．

今，確率変数 $i$ 以外の変数の値を $\mathbf{x}_{-i}=(x_1,..,x_{i-1},x_{i+1},...,x_M)$ とすると，ギブスサンプリングで必要な条件つき分布は

$\displaystyle P(x_i|\mathbf{x}_{-i},\theta) \\ \displaystyle =\frac{P(\mathbf{x},\theta)}{\sum_{x_i=0,1}P(\mathbf{x},\theta)} \\ \displaystyle =\frac{\exp(-\Phi(\mathbf{x},\theta))}{\sum_{x_i=0,1}\exp(-\Phi(\mathbf{x},\theta))}$

となる．ただし分母は $x_i$ についての周辺化を表している．この式を陽に書くと

$\displaystyle \frac{\exp(-\Phi(\mathbf{x},\theta))}{\sum_{x_i=0,1}\exp(-\Phi(\mathbf{x},\theta))}\\ \displaystyle =\frac{\exp(\sum b_ix_i + \sum w_{ij}x_ix_j)}{\sum_{x_i=0,1}\exp(\sum b_ix_i + \sum w_{ij}x_ix_j)} \\ \displaystyle =\frac{\prod_i\exp(b_ix_i)\prod_{i.j}\exp(w_{ij}x_ix_j)}{\sum_{x_i=0,1}\prod_i\exp(b_ix_i)\prod_{i,j}\exp(w_{ij}x_ix_j)} \\ \displaystyle =\frac{\exp(b_i)\prod_j\exp(w_{ij}x_ix_j)}{\sum_{x_i=0,1}\exp(b_ix_i)\prod_j\exp(w_{ij}x_ix_j)} \\ \displaystyle =\frac{\exp(b_ix_i+\sum_jw_{ij}x_ix_j)}{1+\exp(b_i+\sum_jw_{ij}x_j)}$

と変形ができる．3行目から４行目は $x_i$ と関係ない項を約分しており，最後の行は $x_i=0,1$ を代入することで求まる．ここで，両辺の $\sum_j$ はノード $i$ に隣接しているノードのみでの足し合わせになっているため，ノード $i$ に隣接するノードの状態から $x_i$ が決まることがわかる．すなわちボルツマンマシンの局所マルコフ性のおかげで， $\mathbf{x}_{-i}$ を知らなくても周辺のノードのみで $x_i$ を求めることができ，それにより効率的に計算が可能である．また，ここで上記の式をよく見ると， $x_i=0,1$ であることからシグモイド関数の形になっていることがわかる．そのため，サンプリングも非常に単純に行うことができ，区間 $[0,1]$ から発生させた一様乱数が $P(x_i=1|\mathbf{x}_{-i},\theta)$ より大きければ， $x_i=1$ を，それ以外には $x_j=0$ をサンプリング値とすればいい．

ギブスサンプリングの問題点

サンプリングのための条件つき分布とそのサンプリング方法がわかったので，あとはギブスサンプリングの枠組みに当てはめて全ての確率変数をサンプリングすれば良い．また，簡単な計算で連鎖が定常状態になれば元の分布になることを示せる（教科書に丁寧な証明があるので割愛）．ただし，ギブスサンプリングには一つ問題がある．当然，一番最初の変数をサンプリングする際には他の値はわからないため適当な値で初期化を行う．そこから平衡状態になるまで（初期値の影響が完全に消え去るまで）サンプリングを繰り返す必要があるが（この繰り返しをburn-inと呼ぶらしい），どれくらいサンプリングを繰り返せば平衡状態になるかを見積もる明確な方法がない．もう一つ問題として，burn-inをしたのち1つサンプルを得た直後に別なサンプルを得るとそのサンプルは直前のサンプルと強い相関を持つため独立なサンプルとして使うことができない．そのため，ある程度時間をあけてサンプルを取得する必要があるため，全てのサンプルが得られるまでに時間が必要となる．

実践的なサンプリング方法

上記の最後に挙げたように，ギブスサンプリングにはサンプルを連続して取得すると相関を持ってしまう問題がある．そこで連鎖を複数作ってそれぞれの連鎖からサンプリングする方法を考える．これは，サンプリングしたい数だけ連鎖を用意すればサンプル間の独立性を保証できる．ただし，連鎖の数だけburn-inが必要なので計算コストはかかる（処理を並列化すれば解決は可能？）．そのため実践的には，連鎖の本数と計算時間のトレードオフを考えて，複数の連鎖を使って各連鎖から複数のサンプルを取得する方法が取られる．

まとめ

ギブスサンプリングにより分布の期待値を近似する方法を勉強した．これでボルツマンマシンの学習に必要な期待値を近似的に求めることができる．ただし，隠れ変数が入ってくるとモデルだけでなくデータに関する期待値の計算が困難なってくるためサンプリングが必要になり嬉しくない．そこで構造に制約を加えて計算をしやすくしようということで制約ボルツマンマシンができたようなので次回はそこら辺を勉強予定．

ちなみに教科書にはギブスサンプリングではなく平均場近似を使って近似する方法についての解説があったのでそっちも勉強してまとめたい．

機械学習とかコンピュータビジョンとか

CVやMLに関する勉強のメモ書き。

NICE: NON-LINEAR INDEPENDENT COMPONENTS ESTIMATIONを読んだのでメモ

はじめに

基本アイディア

non-linear independent components estimation

モデルの構成

まとめ

ボルツマンマシンについて勉強したメモ　まとめ

Variational Inference with Normalizing Flowsを読んだのでメモ

はじめに

Variational Inference

Deep Latent Gaussian Models（DLGM）

Normalizing flow

Inference with Normalizing Flows

まとめ

ボルツマンマシンについて勉強した　その５

はじめに

ブロック化ギブスサンプリング

コントラスティブ・ダイバージェンス法（CD法）

まとめ

ボルツマンマシンについて勉強した　その４

はじめに

制限付きボルツマンマシン（RBM）

学習方法

まとめ

ボルツマンマシンについて勉強した　その３

はじめに

平均場近似

まとめ

ボルツマンマシンについて勉強した　その２

はじめに

サンプリングについて

マルコフ連鎖

マルコフ連鎖モンテカルロ法

ギブスサンプリングの問題点

実践的なサンプリング方法

まとめ