Processing math: 100%

2018/02/26

MCMCの2D histogram (corner plots):正しい 1-sigma contourの引き方 (how to estimate "correct" 1-sigma percentile areas of corner plots with MCMC samples)

参考(reference ...or answer?)
https://corner.readthedocs.io/en/latest/pages/sigmas.html

よくある間違いは一次元正規分布をそのまま適用してしまい、
規格化した事後分布の中央値から
[0.683, 0.955, 0.997] の範囲をとってしまうこと。

これだと実際の分布を大きく取りすぎてしまう。
実際は二次元正規分布を考える必要がある。

等高線をどの値に引けばいいのか(i.e. 中央値からの規格化された信頼区間)がわかればいいので、xyそれぞれ別の分散を持つような一般的な場合などは考える必要はなく、いずれも同じ分散を持つとして簡単にする。
このとき確率分布関数(PDF)は、
\mathrm{pdf(x, y)} = \frac{1}{2\pi \sigma^2} \exp \left( - \frac{\left( x^2 + y^2 \right)}{2\sigma^2}\right)

x = r\cos \theta, ~ y = r\sin \theta
として、

\mathrm{pdf(r, \theta)} = \frac{1}{2\pi \sigma^2} \exp \left( - \frac{r^2}{2\sigma^2}\right)

これを任意の r 座標 (ここでは z とする) まで面積分することでCDFを求める。
二次元極座標のヤコビアンが r であるから、

\mathrm{cdf(z)} = \frac{1}{2\pi \sigma^2} \int^{z}_{0} \int^{2\pi}_{0} \exp \left( - \frac{r^2}{2\sigma^2}\right) r \mathrm{d}r \mathrm{d}\theta

\theta についての積分はただちに 2\pi と求まり、
\mathrm{cdf(z)} = \frac{1}{\sigma^2} \int^{z}_{0} \exp \left( - \frac{r^2}{2\sigma^2}\right) r \mathrm{d}r

t = r^2 として、置換積分すると、 \frac{\mathrm{d}t}{\mathrm{d}r} = 2r なので、
\begin{eqnarray} \mathrm{cdf(z)} &=& \frac{1}{2\sigma^2} \int^{z^2}_{0} \exp \left( - \frac{t}{2\sigma^2}\right) \mathrm{d}t \nonumber \\ &=& -\left[\exp \left( -\frac{t}{2\sigma^2}\right) \right]^{z^2}_{0} \nonumber \\ &=& 1-\exp \left( -\frac{z^2}{2\sigma^2} \right) \end{eqnarray}

つまり(1)式によれば、正しい分布のpercentileの取り方は、
z = 1\sigma,~2\sigma,~3\sigma の場合となり、それぞれ、
\mathrm{cdf}(z) = 0.393,~0.865,~0.989
となるわけ。

事後分布を高さで規格化するなら、それぞれ1から引いた値をcontour levelにしてやればいい、気がする。

0 件のコメント: