« 読了：Finley, Banerjee, Gelfand (2015) RのspBayesパッケージ | メイン | 読了：Jing & De Oliveira (2015) RのgeoCountパッケージ »
2017年10月26日 (木)
空間データのモデリングではバリオグラムというのを描くけれど、たいてい従属変数は量で、誤差に正規性が仮定されている。もし従属変数が二値とかカウントとかだったら、モデルはGLMとして組めるにしても、バリオグラムはいったいどうなるのか。
これ、ここんところの疑問の一つだったのだが、このたびやっと説明を見つけたのでメモしておく。すいません、純粋に自分のための覚え書きです。
$S(x)$を平均0, 分散$\sigma^2$の定常ガウス過程とする。観察$Y_i$は互いに独立で平均$\mu_i=g(\alpha+S(x_i))$, 分散$v_i = v(\mu_i)$とする。$g(\cdot)$はリンク関数$h(\cdot)$の逆関数である。
$Y$の素直なバリオグラムについて考えよう。横軸は$u=||x_i - x_j||$、縦軸は
$\tau_Y(u)=E[(1/2)(Y_i-Y_j)^2]$
これを展開するとこうなる。
$\tau_Y(u)$
$= (1/2)E_S[{g(\alpha+S(x_i)) - g(\alpha+S(x_j))}]$
$+ E_s[v(g(\alpha+S(x_i)))]$
第2項に$i$だけ出てきているのは、$S(x_i)$の周辺分布が場所を問わず一定だから。よくよくみればこの項は定数なので、以後$\bar{\tau}^2$と書く。これ、$S(\cdot)$の条件付き分布の分散を平均したもので、ガウシアンモデルで言うところのナゲットである。
さて。第一項に出てくる$g(\alpha+S(x_i))$が、$g(\alpha) + S(x_i) g'(\alpha)$とテイラー近似できるのはわかるね。[←わかんないよ！こっちは文系なんだから！]
すると結局こうなる。
$\tau_Y(u) \approx g'(\alpha)^2 \gamma_S(u) + \bar{\tau}^2$
つまり、バリオグラムの縦軸は、潜在ガウス過程$S(\cdot)$だけを取り出したバリオグラム$\gamma_S(u)$に比例し、そこにさきほどの平均ナゲット効果を乗せたものになる。
というわけで、上記のような通常のバリオグラムをナゲットのサイズ把握のために使うのはOK。しかし、GLMは$S(x)$と$Y$のあいだに非線形的関係を仮定しているわけだから、診断のために使うのはおかしい。[←「診断のため」ってどういうことだろう。バリオグラム・モデルの診断には使えないということだろうか]
バリオグラム・モデルをどうやって推定するか。
まずは一般的に、GLMの尤度関数について考えよう。ランダム効果$S$の下での$Y$の条件付き分布を$f_i(y_i; S, \theta)$とすると、パラメータ$\theta$の条件付き尤度は
$L(\theta | S) | \prod_i f_i(u_i|S, \theta)$
混合モデルの尤度は、$S$の同時分布を$g(S; \phi)$として
$L(\theta, \phi) = \int_S \prod_i f_i(u_i|S, \theta) g(s|\phi) ds$
さて、こんな風に多重積分が掛かっている尤度が計算できちゃうのは、ひとえに$S_i$が互いに独立だからである。しかし空間データでは$S_i = S(x_i)$は互いに独立でないので、こんなもん計算できっこない。実をいうとこの積分を近似する方法も提案されているんだけど、$S$の分散が大きいときには役に立たない。
というわけで... モンテカルロ最尤法、階層尤度、GEEなどの方法が提案されておる、という話になるのだが、力尽きたのでこの辺にしておこう。
以上, Diggle & Ribeiro (2010) "Model-based Geostatistics", Chapter 4-5より。勤務先の経費とはいえ2万円以上した本だったので、元を取らねば、と読んだ次第。
実をいうと、ポワソン対数線形モデルを使ったクリギングの際に、いい感じのバリオグラム・パラメータを視認で決める方法について知りたかったんだけど、それはよくわかんなかった。先の章に書いてあるのかも。
雑記：データ解析 - GLM空間モデルのバリオグラムとはなんぞや