読了: Rue, Martino, Chopin (2009) INLAを使って楽しい楽しい潜在ガウシアン・モデル (後編)

 前回に引き続き、以下の論文のメモ。

Rue, H., Martino, s., Chopin, N. (2009) Approximate Bayesian Inference for Latent Gaussian Models by Using Integrated Nested Laplace Approximations. Journal of Royal Statistical Society, B. 71(2), 319-392.

 どうやら私の能力を超えている… 辛い…

3. 積分段階的ラプラス近似
 本章では潜在ガウシアン場の事後周辺密度\(\pi(x_i|\mathbf{y})\)をINLAアプローチで近似する方法を示す。次の3つのステップからなる。(1)\(\theta\)の事後周辺密度をラプラス近似によって近似する。(2)選んだ\(\theta\)の値について\(\pi(x_i | \mathbf{y}, \theta)\)を近似する。(3)数値積分を使って上の2つを結合する。これは1.7節の$$ \hat{\pi}(x_i | \mathbf{y}) = \sum_k \tilde{\pi}(x_i | \theta_k, \mathbf{y}) \tilde{\pi}(\theta_k | \mathbf{y}) \Delta_k $$ を使えばよい。

3.1 \( \tilde{\pi}(\theta | \mathbf{y}) \) を調べる
 INLAアプローチの最初のステップは、\(\theta\)の事後周辺密度の近似\( \tilde{\pi}(\theta | \mathbf{y})\)を計算することである。その主な用途は、あとで\(x_i\)の事後周辺密度を近似するとき、\(\theta\)に関する不確実性を積分消去するためである。そのためには、\( \tilde{\pi}(\theta | \mathbf{y})\)をパラメトリックに表現する必要はなくて、数値積分のためのよい評価点を選ぶことができるくらいにちゃんと分かっているだけでよい。
 話を簡単にするため、\(\theta = (\theta_1, \ldots, \theta_m) \in \mathbb{R}^m \)とする。そうなるように再パラメータ化するのはいつでもできる。

  • ステップ1. \(\tilde{\pi}(\theta|\mathbf{y})\)のモードを決める。そのためには、\( \log\{\tilde{\pi} (\theta | \mathbf{y}) \} \)を\(\theta\)に関して最適化すればよい[最大化ってことね]。具体的には、なんらかのquasi-Newton法で、連続する勾配ベクトルの差を用いて\( \log\{\tilde{\pi} (\theta | \mathbf{y}) \} \)の二階導関数を近似すれば良い。勾配は有限差分を用いて近似する。得られたモードの配置を\(\theta^*\)とする。
  • ステップ2. 有限差分を用いて、ヘッセ行列に-1を掛けた奴\(\mathbf{H} \gt 0\)を求める。\(\Sigma = \mathbf{H}^{-1}\)を求める(密度がガウシアンだったときに\(\theta\)の共分散行列になる行列)。これを\(\Sigma = \mathbf{V} \Lambda \mathbf{V}^\top\)と固有値分解して、$$ \theta(\mathbf{z}) = \theta^* + \mathbf{V} \Lambda^{1/2} \mathbf{z} $$ と定義する。なにをやっているのかというと、数値積分を単純にするために、\(\theta\)ではなく、それを標準化した変数\(\mathbf{z}\)を使いたいのである。もし\(\tilde{\pi}(\theta | \mathbf{y})\)がガウシアンなら、\(\mathbf{z}\)は\(N(\mathbf{0}, \mathbf{I})\)となる。
  • ステップ3. \(\mathbf{z}\)パラメータ化をつかって\(\log \{ \tilde{\pi}(\theta | \mathbf{y}) \}\)を調べる。まず、モード\(\mathbf{z} = \mathbf{0}\)からスタートし、\(\log(\pi\{\theta(\mathbf{0})|y\})\)を求める。\(z_1\)軸に一定のステップ幅で点を置いておき、その点における\(\log(\pi\{\theta(\mathbf{z})|y\})\)が、原点におけるそれと比べてある閾値を超えて小さくなるまで続ける。逆向きにもそうする。\(z_1\)軸を一定のステップ幅で平行移動し、その線上でも同じことをする。なお、等間隔なので最後に数値積分するとき\(\Delta_k\)は定数で良い。[要するに、\(\theta\)の空間というか\(\mathbf{z}\)の空間上で、\log(\pi\{\theta(\mathbf{z})|y\})\)がある程度高いエリアに、グリッド上に点を置いておくわけである]
  • ステップ4. \(\pi(\theta_j | \mathbf{y})\) を近似する。やろうと思ったら\(\tilde{\pi}(\theta | \mathbf{y}) \)を数値積分して得ることができるけれど、計算が大変なので、Step 3.で決めた点をつかって補完してから、それらをつかって数値積分する。

3.2 \( \pi(x_i | \theta, \mathbf{y}) \) を近似する
 3つの方法がある。ラプラス近似がもっともよいが、単純化ラプラス近似はコストが小さい(精度はわずかに落ちる)。

  • ガウシアン近似。すなわち、平均\(\mu_i(\theta)\), 周辺分散\(\sigma^2(\theta)\)を用いてガウシアン近似\(\tilde{\pi}_G(x_i|\theta, \mathbf{y})\)を得る。周辺分散2.1節の\(\Sigma_{ij}\)の式を繰り返して得ればよい。すでの3.1節で\(\tilde{\pi}(\mathbf{x}|\theta, y\)\)を得ているから、周辺分散だけが必要である[わからん…どういうこと?]。このガウシアン近似はあまり良い方法でない。
  • ラプラス近似。下式から求める。$$ \tilde{\pi}_{LA}(x_i | \theta, \mathbf{y}) \propto \left. \frac{\pi(\mathbf{x}, \theta, \mathbf{y}}{\tilde{\pi}_{GG}(\mathbf{x}_{-i} | x_i, \theta, \mathbf{y})} \right|_{\mathbf{x}_{-i} = \mathbf{x}^*_{-i}(x_i, \theta)} $$ 分母はガウシアン近似。\( \mathbf{x}^*_{-i}(x_i, \theta) \) はモードの布置である。分母は\(\tilde{\pi}_G(\mathbf{x}|\theta, \mathbf{y})\)とは異なることに注意。残念ながら、分母は\(x_i\)と\(\theta\)のそれぞれの値について再計算しないといけない。そこで、2つの修正案を提案しよう。
    • \(\tilde{\pi}_G(\mathbf{x}|\theta, \mathbf{y})\)から\(E_{\tilde{\pi}_G} (\mathbf{x}_{-i} | x_i))\)を求め、それを\(\mathbf{x}^*_{-i}(x_i, \theta)\)の近似として使う。[…中略…]
    • 直観的にいうと、\(x_i\)に「近い」\(x_j\)だけが\(x_i\)の周辺分布に効果を持つはずである。さきほどの$$ \mathbf{x}^*_{-1}(x_i , \theta) \approx E_{\tilde{\pi}_G}(\mathbf{x}_{-1}|x_i)$$から下式が導ける。なんらかの\(a_{ij}(\theta)\)について、$$ \frac{E_{\tilde{\pi}_G}(x_j | x_i) -\mu_j(\theta)}{\sigma_j(\theta)} = a_{ij}(\theta) \frac{x_i – \mu_i(\theta)}{\sigma_i(\theta)} $$ \(i\)のまわりの「関心領域」\(R_i(\theta)\)を決める単純なルールとして以下が考えられる。$$ R_i(\theta) = \{j:|a_{ij}(\theta)| \gt 0.001\}$$ [残念ながら、おそらく続きを読んでも私には理解できないと思うので後略]
  • 単純化したラプラス近似。[以下2pにわたって説明があるが、これはたぶん、読んでもわかんないな… パス]

[正直、この節の内容はほぼ理解できなかったが、まあ\( \pi(x_i | \theta, \mathbf{y}) \) をなんらか近似する方法があるということだけはわかった]

4. 近似誤差: 漸近的性質と実務的な諸問題
[パス]

5. 例
[パス。構成のみメモする]

5.1 シミュレーション
5.2 縦断データの一般化線形混合モデル
5.3 確率的ボラタリティモデル
5.4 がん生起データの疾病マッピング
5.5 対数ガウシアンCox過程

6. 拡張
[パス。構成のみメモする]

6.1 下位集合\(\mathbf{x}_S\)の事後周辺密度の近似
6.2 周辺尤度の近似
6.3 予測指標
6.4 デビアンス情報基準
6.5 適切なハイパーパラメータ数

7. 考察
[パス。構成のみメモする]

[このあとに付録、そして44組の識者による膨大なコメントと回答がついている。怖い… もちろんパス]
—————
反省: INLAについて学びたいからといっていきなり原典にあたったのは無謀であった。MCMCについて学びたいからと言って大昔の物理学者の論文を読むようなものではないか。もうちょっとユーザ向けの教科書を探そう…