🧠

PRML 第11章解答例

2022/05/15に公開

はじめに

PRML解答例まとめを参照

演習 11.1

\widehat{f}=\frac{1}{L} \sum_{l=1}^{L} f\left(\mathbf{z}^{(l)}\right) \tag{11.2}

で定義される有限のサンプルによる推定量fE[f]に等しい平均と

\operatorname{var}[\widehat{f}]=\frac{1}{L} \mathbb{E}\left[(f-\mathbb{E}[f])^{2}\right] \tag{11.3}

で与えられる分散を持つことを示せ.


「有限和で近似した期待値」の期待値が元の期待値と一致することを示す.定義から

\begin{aligned} \mathbb{E}\left[\widehat{f}\right]&=\mathbb{E}\left[\frac{1}{L} \sum_{l=1}^{L} f\left(\mathbf{z}^{(l)}\right)\right]\\ &=\frac{1}{L}\sum_{l=1}^{L}\mathbb{E}\left[ f\left(\mathbf{z}^{(l)}\right)\right]\\ &=\frac{1}{L}L\mathbb{E}\left[ f\left(\mathbf{z}^{(l)}\right)\right]\\ &=\mathbb{E}[ f] \end{aligned}

確率変数\mathbf{z}^{(l)}は確率分布p(\mathbf{z})からサンプリングされるため\mathbb{E}[f]=\mathbb{E}[f(\mathbf{z}^{(l)})]となることを用いた.

次に分散(11.3)を示す.分散と期待値の二乗を紐付ける式から

\begin{aligned} \operatorname{var}[\widehat{f}] &=\mathbb{E}\left[(\widehat{f}-\mathbb{E}[\widehat{f}])^{2}\right]=\mathbb{E}\left[\hat{f}^{2}\right]-\mathbb{E}[\widehat{f}]^{2}=\mathbb{E}\left[\hat{f}^{2}\right]-\mathbb{E}[f]^{2} \\ &=\mathbb{E}\left[\left(\frac{1}{L} \sum_{l=1}^{L} f\left(\mathbf{z}^{(l)}\right)\right)^{2}\right]-\mathbb{E}[f]^{2} \\ &=\frac{1}{L^{2}} \mathbb{E}\left[\left(\sum_{l=1}^{L} f\left(\mathbf{z}^{(l)}\right)\right)^{2}\right]-\mathbb{E}[f]^{2} \\ &=\frac{1}{L^{2}} \mathbb{E}\left[\sum_{l=1}^{L} f^{2}\left(\mathbf{z}^{(l)}\right)+\sum_{i, j=1, i \neq j}^{L} f\left(\mathbf{z}^{(i)}\right) f\left(\mathbf{z}^{(j)}\right)\right]-\mathbb{E}[f]^{2} \\ &=\frac{1}{L^{2}} \mathbb{E}\left[\sum_{l=1}^{L} f^{2}\left(\mathbf{z}^{(l)}\right)\right]+\frac{L^{2}-L}{L^{2}} \mathbb{E}[f]^{2}-\mathbb{E}[f]^{2} \\ &=\frac{1}{L^{2}} \sum_{l=1}^{L} \mathbb{E}\left[f^{2}\left(\mathbf{z}^{(l)}\right)\right]-\frac{1}{L} \mathbb{E}[f]^{2} \\ &=\frac{1}{L^{2}} \cdot L \cdot \mathbb{E}\left[f^{2}\right]-\frac{1}{L} \mathbb{E}[f]^{2} \\ &=\frac{1}{L} \mathbb{E}\left[f^{2}\right]-\frac{1}{L} \mathbb{E}[f]^{2}=\frac{1}{L} \mathbb{E}\left[(f-\mathbb{E}[f])^{2}\right] \end{aligned}

以上により示された.

演習 11.2

zは区間(0,1)上の一様分布を持つ確率変数で, zy=h^{-1}(z)で変換することを考えることで,h(y)

z=h(y) \equiv \int_{-\infty}^{y} p(\widehat{y}) \mathrm{d} \widehat{y} \tag{11.6}

で与えられるyが分布p(y)を持つことを示せ.


yの分布をp^\star(y)とおくと(11.5)式のp(y)=p(z)\left|\frac{\mathrm{d} z}{\mathrm{~d} y}\right|と(11.6)式およびp(z)=1を用いて

p^{\star}(y)=p(z) \cdot\left|\frac{d z}{d y}\right|=1 \cdot h^{\prime}(y)=\frac{d}{d y} \int_{-\infty}^{y} p(\widehat{y}) d \widehat{y}=p(y)

となりyが(11.6)式で与えられる分布p(y)を持つことがわかる.

演習 11.3

区間(0,1)上で一様分布する確率変数zが与えられたとき,y

p(y)=\frac{1}{\pi} \frac{1}{1+y^{2}} \tag{11.8}

で与えられるコーシ一分布を持つようにする変換y=f(z)を求めよ.


(11.8)より

p(y) = \frac{1}{\pi}\frac{1}{1+y^2} (11.5)
\begin{aligned} z&=h(y)\\ &=\int_{-\infty}^yp(\widehat{y})d\widehat{y}\\ &=\int_{-\infty}^y\frac{1}{\pi}\frac{1}{1+\widehat{y}^2}d\widehat{y}\\ &=\left[\frac{1}{\pi}\tan^{-1}\widehat{y}\right]_{-\infty}^y\\ &=\frac{1}{\pi}\tan^{-1}y-\frac{1}{\pi}(-\frac{\pi}{2})\\ &=\frac{1}{\pi}\tan^{-1}y+\frac{1}{2} \end{aligned}

これをyについて解くと

y = \tan\left\{\pi(z-\frac{1}{2})\right\}

が得られる

演習 11.4

図11.3に示すように,z_1z_2が単位円上で一様分布し,

y_{1}=z_{1}\left(\frac{-2 \ln r^{2}}{r^{2}}\right)^{1 / 2} \tag{11.10}

および

y_{2}=z_{2}\left(\frac{-2 \ln r^{2}}{r^{2}}\right)^{1 / 2} \tag{11.11}

で与えられる変数変換を行うとする.(y_1,y_2)

\begin{aligned} p\left(y_{1}, y_{2}\right) &=p\left(z_{1}, z_{2}\right)\left|\frac{\partial\left(z_{1}, z_{2}\right)}{\partial\left(y_{1}, y_{2}\right)}\right| \\ &=\left[\frac{1}{\sqrt{2 \pi}} \exp \left(-y_{1}^{2} / 2\right)\right]\left[\frac{1}{\sqrt{2 \pi}} \exp \left(-y_{2}^{2} / 2\right)\right] \end{aligned} \tag{11.12}

に従って分布することを示せ.


方針:計算を楽にするため極座標で表してヤコビアンの計算を以下のように行う.

\begin{aligned} \left|\frac{\partial\left(z_{1}, z_{2}\right)}{\partial\left(y_{1}, y_{2}\right)}\right| &=\left|\frac{\partial\left(z_{1}, z_{2}\right)}{\partial(r, \theta)} \cdot \frac{\partial(r, \theta)}{\partial\left(y_{1}, y_{2}\right)}\right| \\ \end{aligned}

r^2=z_1^2+z_2^2であるから,

z_1=r\cos\theta
z_2=r\sin\theta

とすると

\frac{\partial\left(z_{1}, z_{2}\right)}{\partial(r, \theta)}=\left[\begin{array}{ll} \partial z_{1} / \partial r & \partial z_{1} / \partial \theta \\ \partial z_{2} / \partial r & \partial z_{2} / \partial \theta \end{array}\right]=\left[\begin{array}{cc} \cos \theta & -r \sin \theta \\ \sin \theta & r \cos \theta \end{array}\right]

であるから,行列式の値は

\left|\frac{\partial\left(z_{1}, z_{2}\right)}{\partial(r, \theta)}\right|=r\left(\cos ^{2} \theta+\sin ^{2} \theta\right)=r

となる.つぎに(y_1,y_2)について考えると

y_{1}=r \cos \theta\left(\frac{-2 \ln r^{2}}{r^{2}}\right)^{1 / 2}=\cos \theta\left(-2 \ln r^{2}\right)^{1 / 2}

同様に

y_{2}=\sin \theta\left(-2 \ln r^{2}\right)^{1 / 2}

が求まるので,

\frac{\partial\left(y_{1}, y_{2}\right)}{\partial(r, \theta)}=\left[\begin{array}{ll} \partial y_{1} / \partial r & \partial y_{1} / \partial \theta \\ \partial y_{2} / \partial r & \partial y_{2} / \partial \theta \end{array}\right]=\left[\begin{array}{cc} -2 \cos \theta\left(-2 \ln r^{2}\right)^{-1 / 2} \cdot r^{-1} & -\sin \theta\left(-2 \ln r^{2}\right)^{1 / 2} \\ -2 \sin \theta\left(-2 \ln r^{2}\right)^{-1 / 2} \cdot r^{-1} & \cos \theta\left(-2 \ln r^{2}\right)^{1 / 2} \end{array}\right]

となる.したがって行列式の値は

\left|\frac{\partial\left(y_{1}, y_{2}\right)}{\partial(r, \theta)}\right|=\left(-2 r^{-1}\left(\cos ^{2} \theta+\sin ^{2} \theta\right)\right)=-2 r^{-1}

となる.以上により求めたかったヤコビアンを計算することができて,

\begin{aligned} \left|\frac{\partial\left(z_{1}, z_{2}\right)}{\partial\left(y_{1}, y_{2}\right)}\right| &=\left|\frac{\partial\left(z_{1}, z_{2}\right)}{\partial(r, \theta)} \cdot \frac{\partial(r, \theta)}{\partial\left(y_{1}, y_{2}\right)}\right| \\ &=\left|\frac{\partial\left(z_{1}, z_{2}\right)}{\partial(r, \theta)}\right| \cdot\left|\frac{\partial(r, \theta)}{\partial\left(y_{1}, y_{2}\right)}\right| \\ &=\left|\frac{\partial\left(z_{1}, z_{2}\right)}{\partial(r, \theta)}\right| \cdot\left|\frac{\partial\left(y_{1}, y_{2}\right)}{\partial(r, \theta)}\right|^{-1} \\ &=r \cdot\left(-2 r^{-1}\right)^{-1}=-\frac{r^{2}}{2} \end{aligned}

が得られる.次にr^2(y_1,y_2)で表すと(11.10),(11.11) 及び r^2=z_1^2+z_2^2 から

y_1^2+y_2^2=-2\ln r^2
r^2=\exp\left\{-\frac{y_1^2+y_2^2}{2}\right\}

いまp(z_1, z_2)=\frac{1}{\pi}であるので

\begin{aligned} p\left(y_{1}, y_{2}\right) &=p\left(z_{1}, z_{2}\right)\left|\frac{\partial\left(z_{1}, z_{2}\right)}{\partial\left(y_{1}, y_{2}\right)}\right|\\ &=\frac{1}{\pi}\mid-\frac{r^{2}}{2}\mid\\ &=\frac{1}{2\pi}\exp\left\{-\frac{y_1^2+y_2^2}{2}\right\}\\ &=\left[\frac{1}{\sqrt{2 \pi}} \exp \left(-y_{1}^{2} / 2\right)\right]\left[\frac{1}{\sqrt{2 \pi}} \exp \left(-y_{2}^{2} / 2\right)\right] \end{aligned}

演習 11.5

zを平均ゼロと単位行列の共分散行列を持つD次元のガウス分布に従う確率変数とし正定値対称行列\mathbf{\Sigma}がコレスキー分解\mathbf{\Sigma} = \mathbf{LL}^{\mathrm T}を持つとする。ここで,\mathbf{L}は下三角行列(すなわち,対角成分より上側がゼロになる行列)である。変数\mathbf{y} = \boldsymbol{\mu}+\mathbf{Lz}が平均\boldsymbol{\mu},共分散行列\mathbf{\Sigma}であるガウス分布に従うことを示せ.これは,平均0分散11変数ガウス分布からのサンプルを用いて, 一般の多変量ガウス分布からのサンプルを生成する技術を提供する.


ガウス分布に従う確率変数を線形変換して得られる確率変数はガウス分布に従うので
\mathbf{y}の平均と分散を計算すれば良い.\mathbf{z}\backsim\mathcal{N}(\mathbf{0},\mathbf{I}) であるから

\begin{aligned} \mathbb{E}\left[\mathbf{y}\right]&=\mathbb{E}\left[\boldsymbol{\mu}+\mathbf{Lz}\right]\\ &=\boldsymbol{\mu}+\mathbf{L}\mathbb{E}\left[\mathbf{z}\right]\\ &=\boldsymbol{\mu} \end{aligned}

次に分散共分散行列について\operatorname{cov}[\mathbf{z}]=\mathbb{E}\left[\mathbf{z z}^{\mathrm{T}}\right]-\mathbb{E}[\mathbf{z}] \mathbb{E}\left[\mathbf{z}^{\mathrm{T}}\right]=\mathbb{E}\left[\mathbf{z z}^{\mathrm{T}}\right]=\mathbf{I}, を用いて

\begin{aligned} \operatorname{cov}[\mathbf{y}] &=\mathbb{E}\left[\mathbf{y} \mathbf{y}^{\mathrm{T}}\right]-\mathbb{E}[\mathbf{y}] \mathbb{E}\left[\mathbf{y}^{\mathrm{T}}\right] \\ &=\mathbb{E}\left[(\boldsymbol{\mu}+\mathbf{L} \mathbf{z}) \cdot(\boldsymbol{\mu}+\mathbf{L z})^{\mathrm{T}}\right]-\boldsymbol{\mu} \boldsymbol{\mu}^{\mathrm{T}} \\ &=\mathbb{E}\left[\boldsymbol{\mu} \boldsymbol{\mu}^{\mathrm{T}}+2 \boldsymbol{\mu} \cdot(\mathbf{L z})^{\mathrm{T}}+(\mathbf{L z}) \cdot(\mathbf{L z})^{\mathrm{T}}\right]-\boldsymbol{\mu} \boldsymbol{\mu}^{\mathrm{T}} \\ &=2 \boldsymbol{\mu} \cdot \mathbb{E}\left[\mathbf{z}^{\mathrm{T}}\right] \cdot \mathbf{L}^{\mathrm{T}}+\mathbb{E}\left[\mathbf{L z z}^{\mathrm{T}} \mathbf{L}^{\mathrm{T}}\right] \\ &=\mathbf{L} \cdot \mathbb{E}\left[\mathbf{z z}^{\mathrm{T}}\right] \cdot \mathbf{L}^{\mathrm{T}}=\mathbf{L} \cdot \mathbf{I} \cdot \mathbf{L}^{\mathrm{T}} \\ &=\mathbf{\Sigma} \end{aligned}

以上から.\mathbf{y}\backsim\mathcal{N}(\boldsymbol{\mu},\mathbf{\Sigma})がわかる

演習 11.6

この練習問題では,棄却サンプリングが求めたい分布p(\mathbf{z})から実際にサンプルを抽出することをより注意深く示す。提案分布をq(\mathbf{z})とし,サンプル値\mathbf{z}が受理される確率が\widetilde{p}(\mathbf{z}) / k q(\mathbf{z})であることを示せ.ここで,\tilde{p}p(\mathbf{z})に比例する任意の正規化されていない分布であり,定数kk q(\mathbf{z}) \geqslant \widetilde{p}(\mathbf{z})をすべての\mathbf{z}の値に対して保証する最小値に設定される.値\mathbf{z}を抽出する確率はその値をq(\mathbf{z})から抽出する確率とその値が抽出されたときにそれが受理される確率の積であることに注意せよ.この事実と確率の和と積の規則を共に用いて,\mathbf{z}上の分布を正規化された形に書き下し,それがp(\mathbf{z})に等しいことを示せ.


まずP.242の流れをしっかりと確認する。

p(\mathbf{z})から直接サンプリングすることは困難であるが、任意の与えられた\mathbf{z}の値についてp(\mathbf{z})を求めることは、正規化定数を除いて容易だとする。すなわち、

p(\mathbf{z}) = \frac{1}{Z_{p}}\tilde{p}(\mathbf{z})

において\tilde{p}(\mathbf{z})はすぐに求まるが、Z_{p}はわからないとする。

ここから教科書の流れと少し変わる。値\mathbf{z}が与えられている時に、それが受理される確率はp(\textrm{acceptance}\mid \mathbf{z})と書けて、かつ図11.4のように、区間[0,kq(\mathbf{z})]の一様分布から得られたサンプルu0\le u \le \tilde{p}(\mathbf{z})となったときに受理(acceptance)され、サンプル値\mathbf{z}として保持されるので

p(\textrm{acceptance}\mid \mathbf{z}) = \int_{0}^{\tilde{p}(\mathbf{z})}\frac{1}{kq(\mathbf{z})}du = \frac{\tilde{p}(\mathbf{z})}{kq(\mathbf{z})}

これがサンプル値\mathbf{z}の受理確率である。

次に受理されたサンプル値\mathbf{z}の確率密度分布は数式上でp(\mathbf{z} \mid \textrm{acceptance})とかけ、これがp(\mathbf{z})となることを示せば良い。

ベイズの定理から

p(\mathbf{z} \mid \textrm{acceptance}) = \frac{p(\mathbf{z},\textrm{acceptance})}{p(\textrm{acceptance})}

となり、分子p(\mathbf{z},\textrm{acceptance})は任意の\mathbf{z}が受理されるときの同時確率を表すので、これはq(\mathbf{z})p(\textrm{acceptance}\mid \mathbf{z})の積になり、

\begin{aligned} p(\mathbf{z},\textrm{acceptance}) &= q(\mathbf{z}) p(\textrm{acceptance} \mid \mathbf{z}) \\ &= q(\mathbf{z}) \frac{\widetilde{p}(\mathbf{z})}{k q(\mathbf{z})} \\ &=\frac{\tilde{p}(\mathbf{z})}{k} \\ &=\frac{Z_{p}}{k}p(\mathbf{z}) \quad (\because (11.13)) \end{aligned}

分母p(\textrm{acceptance})はP.243の(11.14)の導出と同様で

\begin{aligned} p(\textrm{acceptance}) &=\int q(\mathbf{z}) \frac{\widetilde{p}(\mathbf{z})}{k q(\mathbf{z})} \mathrm{d} \mathbf{z} \quad (\because (11.14))\\ &=\frac{1}{k} \int \widetilde{p}(\mathbf{z}) \mathrm{d} \mathbf{z}\\ &=\frac{Z_{p}}{k}\int p(\mathbf{z}) \mathrm{d} \mathbf{z}\quad (\because (11.13)) \\ &= \frac{Z_{p}}{k} \quad \left(\because \int p(\mathbf{z}) \mathrm{d} \mathbf{z} = 1 \right) \end{aligned}

以上から、

p(\mathbf{z} \mid \textrm{acceptance}) = \frac{p(\mathbf{z},\textrm{acceptance})}{p(\textrm{acceptance})} = p(\mathbf{z})

となり、題意が示された。

演習 11.7

yが区間[0,1]上の一様分布に従うとせよ.変数z = b\tan y + c

q(z)=\frac{k}{1+(z-c)^{2} / b^{2}} \tag{11.16}

で与えられるコーシ一分布に従うことを示せ.


(11.5)式の変換を用いる。教科書11.1.1節における説明に使われている変数の組y,zの役割がこの問題では逆になっていることに注意する。

y[0,1]の一様分布とし、z=f(y)=b\tan y +cの関数zの値が求めたい特定の分布q(z)に従うようになっていることを示す。ここではp(y)=1となっていることに注意する。

\displaystyle y = \tan^{-1}\left( \frac{z-c}{b} \right)と書けるので、

q(z) = p(y)\left| \frac{dy}{dz} \right| \tag{11.5}

に当てはめるためにdy/dzの値を求めると、\displaystyle \frac{d}{dx} \tan^{-1} x = \frac{1}{1+x^2}であることを用いて

\begin{aligned} \frac{d y}{d z} &=\frac{1}{1+\left(\frac{z-c}{b}\right)^{2}} \frac{d}{d z}\left(\frac{z-c}{b}\right) \\ &=\frac{1}{b} \frac{1}{1+\left(\frac{z-c}{b}\right)^{2}} \end{aligned}

(11.5)式に当てはめて

\begin{aligned} q(z) &=p(y)\left|\frac{d y}{d z}\right| \\ &=1 \times \frac{1}{b} \frac{1}{1+\left(\frac{z-c}{b}\right)^{2}} \end{aligned}

k = 1/bとみなせば、これは(11.16)式と同型になる。

※しかし、そもそも教科書の(11.16)式は、提案分布q(z)ではなく比較関数kq(z)と呼ぶのが適切なのではないかという意見。。

演習 11.8

連続性と正規化の条件を用いて適応的棄却サンプリングの包絡分布

q(z)=k_{i} \lambda_{i} \exp \left\{-\lambda_{i}\left(z-z_{i}\right)\right\} \quad \widehat{z}_{i-1, i}<z \leqslant \widehat{z}_{i, i+1} \tag{11.17}

の係数k_iの式を決定せよ.


f(z)=l_{i} \lambda_{i} \exp \left\{-\lambda_{i}\left(z-z_{i}\right)\right\}, \hat{z}_{i-1, i}<z \leqq \hat{z}_{i, i+1} とする。

対数を取ると\ln f(z)=\ln l_{i}+\ln \lambda_{i}-\lambda_{i}\left(z-z_{i}\right)

これが、z_{i}において、\ln p(z_{i})と接する。傾きが同じになるので

\lambda_{i}=(\ln p(z_{i}))^{\prime}=\frac{p\left(z_{i}\right)^{\prime}}{p\left(z_{i}\right)}

また値が同じになるので

\ln p\left(z_{i}\right)=\ln \ln l_{i}+\ln \lambda_{i}=\ln \ln l_{i}+\ln \frac{p^{\prime}\left(z_{i}\right)}{p\left(z_{i}\right)}

これより

p\left(z_{i}\right)=l_{i} \frac{p^{\prime}\left(z_{i}\right)}{p\left(z_{i}\right)}
\therefore l_{i}=\frac{p\left(z_{i}\right)^{2}}{p\left(z_{i}\right)^{\prime}}

よって

f(z)=l_{i} \lambda_{i} \exp \left\{-\lambda_{i}\left(z-z_{i}\right)\right\}, \lambda_{i}=\frac{p^{\prime}\left(z_{i}\right)}{p\left(z_{i}\right)}, l_{i}=\frac{p\left(z_{i}\right)^{2}}{p\left(z_{i}\right)^{\prime}}

f(z)の正規化係数は

Z_{q}=\int f(z) d z=\sum_{i=1}^{N} \int_{\hat{z}_{i-1, i}}^{\hat{z}_{i, i+1}} f(z) d z

Nはグリッドの個数、\hat{z}_{i-1, i}は接線の交点のz座標。\hat{z}_{0,1}=-\infty, \hat{z}_{N, N_{2} 1}=\infty

ここで

\begin{aligned} \int_{\hat{z}_{i-1}, i}^{\hat{z}, i, i+1} f(z) d z &=\int_{\hat{z}_{i-1, i}}^{\hat{z}_{i, i+1}} l_{i} \lambda_{i} \exp \left\{-\lambda_{i}\left(z-z_{i}\right)\right\} d z \\ &=l_{i} \lambda_{i} \int_{\hat{z}_{i-1, i}}^{\hat{z}_{i, i+1}} \exp \left\{-\lambda_{i}\left(z-z_{i}\right)\right\} d z \\ &=l_{i} \lambda_{i}\left[-\frac{1}{\lambda_{i}} \exp \left\{-\lambda_{i}\left(z-z_{i}\right)\right\}\right]_{\hat{z}_{i-1, i}}^{\hat{z}_{i, i+1}} \\ &=l_{i}\left[\exp \left\{-\lambda_{i}\left(\hat{z}_{i-1, i}-z_{i}\right)\right\}-\exp \left\{-\lambda_{i}\left(\hat{z}_{i, i+1}-z\right)\right\}\right] \end{aligned}

よって正規化係数は以下になる。

Z_{q}=\sum_{i=1}^{N} l_{i}\left[\exp \left\{-\lambda_{i}\left(\hat{z}_{i-1, i}-z_{i}\right)\right\}-\exp \left\{-\lambda_{i}\left(\hat{z}_{i, i+1}-z_{i}\right)\right\}\right.

ここで(11.17)は以下でも表現できる。

q(z)=\frac{1}{z_{z}} f(z)
q(z)=\frac{1}{Z_{q}} l_{i} \lambda_{i} \exp \left\{-\lambda_{i}\left(z-z_{i}\right)\right\}

演習 11.9

11.1.1節で述べた,単一の指数分布からサンプリングする技術を用いて

q(z)=k_{i} \lambda_{i} \exp \left\{-\lambda_{i}\left(z-z_{i}\right)\right\} \quad \widehat{z}_{i-1, i} \lt z \leqslant \widehat{z}_{i, i+1} \tag{11.17}

で定義される区分的な指数分布からサンプリングするアルゴリズムを考案せよ.


区分的な指数分布からのサンプリング確率を知る必要があるが、(11.17)のq(z)は正規化されていない。従って、まず初めに正規化定数Z_{q}を求める。

\begin{aligned} Z_{q} &=\int_{\tilde{z}_{0,1}}^{\tilde{z}_{N, N+1}} q(z) d z=\sum_{i=1}^{N} \int_{\tilde{z}_{i-1, i}}^{\tilde{z}_{i, i+1}} q_{i}\left(z_{i}\right) d z_{i} \\ &=\sum_{i=1}^{N} \int_{\tilde{z}_{i-1, i}}^{\tilde{z}_{i, i+1}} k_{i} \lambda_{i} \exp \left\{-\lambda_{i}\left(z-z_{i}\right)\right\} d z_{i} \\ &=\sum_{i=1}^{N}-\left.k_{i} \exp \left\{-\lambda_{i}\left(z-z_{i}\right)\right\}\right|_{\tilde{z}_{i-1, i}} ^{\tilde{z}_{i, i+1}} \\ &=\sum_{i=1}^{N}-k_{i}\left[\exp \left\{-\lambda_{i}\left(\widetilde{z}_{i, i+1}-z_{i}\right)\right\}-\exp \left\{-\lambda_{i}\left(\widetilde{z}_{i-1, i}-z_{i}\right)\right\}\right]=\sum_{i=1}^{N} \widehat{k}_{i} \end{aligned}

\widehat{k}_{i}は以下のように定義する。

\widehat{k}_{i}=-k_{i}\left[\exp \left\{-\lambda_{i}\left(\widetilde{z}_{i, i+1}-z_{i}\right)\right\}-\exp \left\{-\lambda_{i}\left(\widetilde{z}_{i-1, i}-z_{i}\right)\right\}\right] \tag{A}

この導出から、i番目の区分からサンプリングする確率は、Z_{q}=\sum_{i=1}^{N} \widehat{k}_{i}とすると、 \widehat{k}_{i} / Z_{q}で与えられる。 そこで今度は区間[0,1]で一様な補助的な確率変数\etaを定義する。

i=j \quad \text { if } \eta \in\left[\frac{1}{Z_{q}} \sum_{m=0}^{j-1} \widehat{k}_{m}, \frac{1}{Z_{q}} \sum_{m=0}^{j} \widehat{k}_{m}\right], \quad j=1,2, \ldots, N \tag{B}

ここで便宜上、\widehat{k}_{0}=0と定義し、ここまでで選択されたi番目の区分を決定した。
次に、11.1.1節の手法を用いて、i番目の指数分布からサンプリングする。(11.6)より、次のように書くことができる。

\begin{aligned} h_{i}(z) &=\int_{\tilde{z}_{i-1, i}}^{z} \frac{q_{i}\left(z_{i}\right)}{\widehat{k}_{i}} d z_{i} \\ &=\frac{1}{\widehat{k}_{i}} \int_{\tilde{z}_{i-1, i}}^{z} k_{i} \lambda_{i} \exp \left\{-\lambda_{i}\left(z-z_{i}\right)\right\} d z_{i} \\ &=\left.\frac{-k_{i}}{\widehat{k}_{i}} \exp \left\{-\lambda_{i}\left(z-z_{i}\right)\right\}\right|_{\tilde{z}_{i-1, i}} ^{z} \\ &=\frac{-k_{i}}{\widehat{k}_{i}}\left[\exp \left\{-\lambda_{i}\left(z-z_{i}\right)\right\}-\exp \left\{-\lambda_{i}\left(\widetilde{z}_{i-1, i}-z_{i}\right)\right\}\right] \\ &=\frac{k_{i}}{\widehat{k}_{i}}\exp \left(\lambda_{i} z_{i}\right)\left[\exp \left\{-\lambda_{i} \widetilde{z}_{i-1, i}\right\}-\exp \left\{-\lambda_{i} z\right\}\right] \end{aligned}

なお、q_{i}(z)は正しく正規化されておらず、q_{i}(z)/\widehat{k}_{i}が正しい正規化された形になる。
いくつか並べ替えを行うと次のようになる。

\begin{aligned} h_{i}^{-1}(\xi) &=\frac{1}{-\lambda_{i}} \ln \left[\exp \left\{-\lambda_{i} \tilde{z}_{i-1, i}\right\}-\frac{\xi}{\frac{k_{i}}{\hat{k}_{i}} \exp \left(\lambda_{i} z_{i}\right)}\right] \\ &=\frac{1}{-\lambda_{i}} \frac{\ln \left[\exp \left\{-\lambda_{i} \tilde{z}_{i-1, i}\right]\right.}{\ln \frac{\widehat{\hat{k}_{i} \xi}}{k_{i} \cdot \exp \left(\lambda_{i} z_{i}\right)}} \\ &=\frac{\tilde{z}_{i-1, i}}{\ln \xi+\ln \frac{\widehat{k}_{i}}{k_{i}}-\lambda_{i} z_{i}} \end{aligned}

結論として、まず、区間[0,1]に一様なランダム変数\etaを生成し、(B)に従って値iを求め、次に、同じく区間[0,1]に一様なランダム変数x_iを生成し、z=h_{i}^{-1}(\xi)を用いてzに変換している。

ここで,z_{1}, z_{2}, \ldots, z_{N}の格子点が与えられれば,\lambda_{i}, \widetilde{z}_{i, i+1}k_{i}が得られることに注意。詳細は前問を参照。これらの変数が得られた後、(A)を用いて\widehat{k}_{i}も決定することができるため、h_{i}^{-1}(\xi)を決定することができる。

演習 11.10

p\left(z^{(\tau+1)}=z^{(\tau)}\right)=0.5 \tag{11.34}
p\left(z^{(\tau+1)}=z^{(\tau)}+1\right)=0.25 \tag{11.35}

および

p\left(z^{(\tau+1)}=z^{(\tau)}-1\right)=0.25 \tag{11.36}

で定義される整数上の単純なランダムウォークが\mathbb{E}\left[\left(z^{(\tau)}\right)^{2}\right]=\mathbb{E}\left[\left(z^{(\tau-1)}\right)^{2}\right]+1 / 2という性質を持ち,よって帰納法により\mathbb{E}[(z^{(r)})^2] = \tau / 2であることを示せ.


仮定より、

\begin{aligned} \mathbb{E}[(z^{(r)})^2] &= \sum (z^{(\tau)})^2 p(z^{(\tau)}) \\ &= 0.5 \cdot \mathbb{E}[(z^{(\tau-1)})^2]+ 0.25 \cdot \mathbb{E}[(z^{(\tau-1)}-1)^2]+ 0.25 \cdot \mathbb{E}[(z^{(\tau-1)}+1)^2] \\ &= \mathbb{E}\left[\left(z^{(\tau-1)}\right)^{2}\right]+1 / 2 \end{aligned}

である。今、\mathbb{E}[(z^{(0)})^2] = 0で、\mathbb{E}\left[\left(z^{(\tau)}\right)^{2}\right] - \mathbb{E}\left[\left(z^{(\tau-1)}\right)^{2}\right] = 1/2の等差数列と見做せるので、\mathbb{E}[(z^{(r)})^2] = \tau / 2である

演習 11.11

11.3節で述べたギブスサンプリングアルゴリズムは,

p^{\star}(\mathbf{z}) T\left(\mathbf{z}, \mathbf{z}^{\prime}\right)=p^{\star}\left(\mathbf{z}^{\prime}\right) T\left(\mathbf{z}^{\prime}, \mathbf{z}\right) \tag{11.40}

で定義される詳細釣り合い条件を満たすことを示せ.


ギブスサンプリングではある時刻\tauにおいて\mathbf{z}^{(\tau)}の全変数のうちz_kを除く変数すべてを固定してz_k^{(\tau + 1)}をサンプリングする。したがってこの一操作ではkを除く残りすべての\{z_{i}\}_{i\ne k}は不変である。つまり\{{z_{i}}^{\prime}\}_{i\ne k} = \{z_{i}\}_{i\ne k}である。

P.254の説明から遷移確率T\left(\mathbf{z}, \mathbf{z}^{\prime}\right) \equiv p\left(\mathbf{z}^{\prime} \mid \mathbf{z}\right)であるが、ここでは一操作あたりT\left(\mathbf{z}, \mathbf{z}^{\prime}\right) \equiv p\left({z_{k}}^{\prime} \mid \{z_{i}\}_{i\ne k}\right)と書ける(下の式ではp^{\star}としているが同じである)。

以上を用いると

\begin{aligned} p^{\star}(\mathbf{z}) T\left(\mathbf{z}, \mathbf{z}^{\prime}\right) &=p^{\star}\left(z_{k},\left\{z_{i}\right\}_{i \neq k}\right) p^{\star}\left(z_{k}^{\prime} \mid\left\{z_{i}\right\}_{i \neq k}\right) \\ &=p^{\star}\left(z_{k} \mid\left\{z_{i}\right\}_{i \neq k}\right) p^{\star}\left(\left\{z_{i}\right\}_{i \neq k}\right) p^{\star}\left(z_{k}^{\prime} \mid\left\{z_{i}\right\}_{i \neq k}\right) \\ &=p^{\star}\left(z_{k} \mid\left\{z_{i}^{\prime}\right\}_{i \neq k}\right) p^{\star}\left(\left\{z_{i}^{\prime}\right\}_{i \neq k}\right) p^{\star}\left(z_{k}^{\prime} \mid\left\{z_{i}^{\prime}\right\}_{i \neq k}\right)\quad \left(\because \{{z_{i}}^{\prime}\}_{i\ne k} = \{z_{i}\}_{i\ne k} \right)\\ &=p^{\star}\left(z_{k} \mid\left\{z_{i}^{\prime}\right\}_{i \neq k}\right) p^{\star}\left(z_{k}^{\prime},\left\{z_{i}^{\prime}\right\}_{i \neq k}\right) \\ &=T\left(\mathbf{z}^{\prime}, \mathbf{z}\right)p^{\star}\left(\mathbf{z}^{\prime}\right) \\ &=p^{\star}\left(\mathbf{z}^{\prime}\right)T\left(\mathbf{z}^{\prime}, \mathbf{z}\right) \end{aligned}

以上で(11.40)式が示された。

演習 11.12

図11.15に示す分布を考えよ.この分布に対する標準的なギブスサンプリングの手続きがエルゴード的かどうか,よって,この分布から正しくサンプリングするかどうか,について論ぜよ.


エルゴード的でない。なぜなら、z_1軸とz_2軸、どちらに射影してみても、影がつけられた二つの領域は重ならず、それゆえ片方の領域から最初の点としてサンプリングを行うと、ギブスサンプリングにより他方の領域に到達することはなく、初期サンプリングがどちらの領域に所属するかによって分布の収束先が変わってしまうからである。

演習 11.13

図11.16に示す単純な3ノードのグラフで, 観測ノードxが平均\mu,精度\tauのガウス分布\mathcal{N}(x\mid \mu,\tau^{-1})で与えられるものを考えよ.平均と精度の周辺分布が\mathcal{N}(\mu \mid \mu_{0},s_{0})および\operatorname{Gam}(\tau \mid a,b)で与えられるとせよ.ここで,\operatorname{Gam}(\cdot\mid\cdot,\cdot)はガンマ分布を表す事後分布p(\mu,\tau \mid x)にギブスサンプリングを適用するために必要となる条件付き分布p(\mu\mid x,\tau)p(\tau \mid x,\mu)の式を書き下せ.


演習 11.14

z_iが平均\mu_i,分散\sigma_{i}^{2}を持ち,\nuが平均0,分散1を持つとき,過剰緩和の更新式

z_{i}^{\prime}=\mu_{i}+\alpha\left(z_{i}-\mu_{i}\right)+\sigma_{i}\left(1-\alpha^{2}\right)^{1 / 2} \nu \tag{11.50}

が平均\mu_i,分散\sigma_i^{2}である値z_{i}^{\prime}を与えることを示せ.


※期待値と分散の定義にしたがって計算していけば良い。

z_{i}^{\prime}の期待値は

\begin{aligned} \mathbb{E}\left[z_{i}^{\prime}\right] &=\mathbb{E}\left[\mu_{i}+\alpha\left(z_{i}-\mu_{i}\right)+\sigma_{i}\left(1-\alpha^{2}\right)^{1 / 2}{ }_{\nu}\right] \\ &=\mu_{i}+\alpha\left(\mathbb{E}\left[z_{i}\right]-\mu_{i}\right)+\sigma_{i}\left(1-\alpha^{2}\right)^{1 / 2} \mathbb{E}[\nu] \\ &=\mu_{i}+0+0 \\ &=\mu_{i} \end{aligned}

z_{i}^{\prime}の分散は

\begin{aligned} \operatorname{var}\left[z_{i}^{\prime}\right] &=\mathbb{E}\left[\left(z_{i}^{\prime}\right)^{2}\right]-\mathbb{E}\left[z_{i}^{\prime}\right]^{2} \\ &=\mathbb{E}\left[\left(\mu_{i}+\alpha\left(z_{i}-\mu_{i}\right)+\sigma_{i}\left(1-\alpha^{2}\right)^{1 / 2} \nu\right)^{2}\right]-\mu_{i}^{2} \\ &=\alpha^{2} \mathbb{E}\left[\left(z_{i}-\mu_{i}\right)^{2}\right]+\sigma_{i}^{2}\left(1-\alpha^{2}\right) \mathbb{E}\left[\nu^{2}\right] \\ &=\alpha^{2} \sigma_{i}^{2}+\sigma_{i}^{2}\left(1-\alpha^{2}\right)\\ &=\sigma_{i}^{2} \end{aligned}

演習 11.15

K(\mathbf{r})=\frac{1}{2}\|\mathbf{r}\|^{2}=\frac{1}{2} \sum_{i} r_{i}^{2} \tag{11.56}
H(\mathbf{z}, \mathbf{r})=E(\mathbf{z})+K(\mathbf{r}) \tag{11.57}

を用いて,ハミルトン方程式

\frac{\mathrm{d} z_{i}}{\mathrm{~d} \tau}=\frac{\partial H}{\partial r_{i}} \tag{11.58}

r_{i}=\frac{\mathrm{d} z_{i}}{\mathrm{~d} \tau} \tag{11.53}

と等価であることを示せ.同様に,(11.57)を用いて,

\frac{\mathrm{d} r_{i}}{\mathrm{~d} \tau}=-\frac{\partial H}{\partial z_{i}} \tag{11.59}

\frac{\mathrm{d} r_{i}}{\mathrm{~d} \tau}=-\frac{\partial E(\mathbf{z})}{\partial z_{i}} \tag{11.55}

に等価であることを示せ.


(11.57)をr_{i}で偏微分すると、

\frac{\partial H}{\partial r_{i}}=\frac{\partial K}{\partial r_{i}}=r_{i} (\because 11.56) (11.53)

同様に(11.57)をz_{i}で偏微分すると、

\frac{\partial H}{\partial z_{i}}=\frac{\partial E}{\partial z_{i}} (\because 11.56)

を得る。

(11.55)と比較して(11.59)を得る。

演習 11.16

K(\mathbf{r})=\frac{1}{2}\|\mathbf{r}\|^{2}=\frac{1}{2} \sum_{i} r_{i}^{2} \tag{11.56}
H(\mathbf{z}, \mathbf{r})=E(\mathbf{z})+K(\mathbf{r}) \tag{11.57}

および

p(\mathbf{z}, \mathbf{r})=\frac{1}{Z_{H}} \exp (-H(\mathbf{z}, \mathbf{r})) \tag{11.63}

を用いて,条件付き分布p(\mathbf{r}\mid \mathbf{z})がガウス分布であることを示せ.


ベイズの定理と(11.54),(11.63)によれば、

\begin{aligned} p(\mathbf{r} \mid \mathbf{z})&=\frac{p(\mathbf{z}, \mathbf{r})}{p(\mathbf{z})}\\ &= \frac{1 / Z_{H} \cdot \exp (-H(\mathbf{z}, \mathbf{r}))}{1 / Z_{p} \cdot \exp (-E(\mathbf{z}))}\\ &= \frac{Z_{p}}{Z_{H}} \cdot \exp (-K(\mathbf{r}))(\because 11.57)\\ &= \frac{Z_{p}}{Z_{H}} \cdot \exp (-\frac{1}{2} \sum_{i} r_{i}^{2})(\because 11.56) \end{aligned}

従って、p(\mathbf{r} \mid \mathbf{z})はガウス分布に従う。

演習 11.17

2つの確率

\frac{1}{Z_{H}} \exp (-H(\mathcal{R})) \delta V \frac{1}{2} \min \left\{1, \exp \left(H(\mathcal{R})-H\left(\mathcal{R}^{\prime}\right)\right)\right\} \tag{11.68}
\frac{1}{Z_{H}} \exp \left(-H\left(\mathcal{R}^{\prime}\right)\right) \delta V \frac{1}{2} \min \left\{1, \exp \left(H\left(\mathcal{R}^{\prime}\right)-H(\mathcal{R})\right)\right\} \tag{11.69}

が等しくよってハイブリッドモンテカルロアルゴリズムで詳細釣り合い条件が満たされることを確認せよ.


\frac{1}{Z_{H}} \exp (-H(R)) \delta V \frac{1}{2} \min \left\{1, \exp \left(H(R)-H\left(R^{\prime}\right)\right)\right\} \tag{11.68}
\frac{1}{Z_{H}} \exp \left(-H\left(R^{\prime}\right)\right) \delta V \frac{1}{2} \min \left\{1, \exp \left(H\left(R^{\prime}\right)-H(R)\right)\right\} \tag{11.69}

H(R)=H\left(R^{\prime}\right)の時、両者は明らかに等しい。
H(R)>H\left(R^{\prime}\right)の時、(11.68)は

\frac{1}{Z_{H}} \exp (-H(R)) \delta V \frac{1}{2}

に減少する。
この時、(11.69)

\frac{1}{Z_{H}} \exp \left(-H\left(R^{\prime}\right)\right) \delta V \frac{1}{2} \exp \left(H\left(R^{\prime}\right)-H(R)\right) =\frac{1}{Z_{H}} \exp (-H(R)) \delta V \frac{1}{2}

ゆえに両者は同一である。

H(R)<H\left(R^{\prime}\right)の時も同様である。

Discussion