18. モンテカルロ・MCMC・ブートストラップ

積分が難しい期待値、複雑な推定量の分布、解析的に解けない事後分布を数値で近似します。乱数を使った結果には計算上の誤差があり、データ由来の不確実性と区別します。乱数種、反復数、アルゴリズム、診断を残せば再現しやすくなります。

独立な $X_b$ を使うモンテカルロでは、$\hat I=B^{-1}\sum h(X_b)$ の標準誤差は $s_h/\sqrt B$。計算誤差を半分にするには反復数を4倍にします。

独立な $U,V\sim U(0,1)$ を発生し、$U^2+V^2\le1$ の割合を4倍します。真の確率は $\pi/4$ なので、$B=10000$ の理論的SEは

$$4\sqrt{\frac{(\pi/4)(1-\pi/4)}{10000}}\approx0.01642.$$

小数点以下の多数の桁を表示しても、その桁まで正確ではありません。逆変換法は $X=F^{-1}(U)$、棄却法は提案分布から受理・棄却して乱数を作る方法です。

元の $n$ 観測から、サイズ $n$ の標本を復元抽出し、毎回推定量を計算します。そのばらつきで標準誤差やバイアス、区間を近似します。パーセンタイル区間は再標本化推定値の分位点を使いますが、すべての推定問題で高精度とは限りません。標本最大値など非正則な問題では通常の方法が失敗する場合があります。

同一対象の前後はペアごと、集落標本は集落ごと、時系列はブロックごとに抽出します。観測構造を無視した再標本化は依存を壊します。ブートストラップは小さな標本から未知の尾や選択バイアスを自動的に作り出す方法ではありません。

目標分布を定常分布とするマルコフ連鎖を作ります。Metropolis–Hastingsの受理確率は

$$\alpha(x,y)=\min\left\{1,\frac{\pi(y)q(x\mid y)}{\pi(x)q(y\mid x)}\right\}.$$

正規化定数が比で消えることが利点です。Gibbsは各成分を完全条件付き分布から順次更新します。2変量標準正規で相関 $\rho$ なら、$X\mid Y\sim N(\rho Y,1-\rho^2)$ と $Y\mid X\sim N(\rho X,1-\rho^2)$ を交互に引きます。相関が強いと混合が遅くなります。

サンプルは独立ではありません。自己相関があると有効標本数は概ね $B/(1+2\sum_{k\ge1}\rho_k)$。複数連鎖のトレース、$\hat R$、有効標本数、MCSEを確認します。burn-inを捨てただけでは収束の証明になりません。

MCMCを10万回実行したら独立な10万標本と同じ精度か？

解答・解説

自己相関によって有効標本数が小さくなるため、同じとは限りません。連鎖の混合と、推定対象ごとのMCSEを調べます。