慶應義塾大学・株式会社Nospareの菅澤です．
今回はガウス過程を用いた空間データの(階層)ベイズモデリングについて紹介します．
空間データの分析
緯度・経度などの位置情報が付随したデータは一般的に空間データと呼ばれます．例えば不動産価格のデータは価格や物件の特徴量に加えて住所の情報から緯度・経度の位置情報を得ることができます．
実は空間データの中にも様々な型のデータがあるのですが，今回は$n$個の地点$s_i \ (i=1,\ldots,n)$において被説明変数$y_i$と説明変数$x_i$が観測されている状況を考えます．例えば不動産の例ですと，$y_i$が不動産価格，$x_i$が不動産の特徴量，$s_i$が物件の緯度・経度に相当します．
基本的な目的は
- 説明変数$x_i$が被説明変数$y_i$に与える影響を調べる
- 新しい地点$s_0$において説明変数$x_0$を用いて未観測の$y_0$を予測する
です．
上記の目的を達成するべく，このような空間データを分析する上で重要な考え方は以下の2点です．
- 空間相関: 近い位置のデータは似たような値になる (場所ごとの「相場感」を反映することが可能)
- 空間異質性: 場所ごとに$x_i$が$y_i$に与える影響が異なる (場所によるモデルの構造の違いを捉えることが可能)
不動産価格の例ですと，「相場感」を考慮することは都心と郊外でベースの価格が異なることを考慮することに相当しますし，「場所によるモデルの構造の違い」を考慮することは都心と郊外で説明変数が与える影響が異なる(例えば都心と郊外で物件の敷地面積が1平米増えることによる価格の増加度合いが違う)点を考慮することに相当します．このように，空間データに対しては空間相関・空間異質性ともにモデリング上考慮するべき重要な概念であることがわかります．
空間データに対する回帰モデル
上記の点を考慮した回帰モデルを考えてみます．まずは「相場感」のみを考慮した以下のようなモデルを考えてみます．
y_i=x_i^\top \beta + \omega_i + \varepsilon_i, \ \ \ \ \varepsilon_i\sim N(0, \sigma^2), \ \ \ i=1,\ldots,n.
ここで，$\omega_i$が各地点固有の切片項になりますので，各地点における相場感を表現するパラメータと解釈することができます．また，$\varepsilon_i$は誤差項です．通常の線形回帰モデルは$\omega_i=0$としたケース (相場感を全く考慮しないモデル) に相当します．
上記のモデルをどのように推定すれば良いか考えてみます．地点ごとに異なるパラメータ$\omega_i$を推定する必要がありますが，基本的に各地点にデータが1つしかない (同じ場所に複数の物件は存在しない) ので，何の構造も入れずに$\omega_i$を推定することはできません．
そこで，ベイズ的な考え方を採用し，$(\omega_1,\ldots,\omega_n)$に対して事前分布という形で構造 (確率モデル) を導入することを考えてみます．
では，$(\omega_1,\ldots,\omega_n)$に対してどのような性質を満たす確率モデルを用いるのが望ましいでしょうか．$\omega_i$が果たすべき役割が「相場感」の表現だったことを思い出すと，「近い位置の$\omega_i$は確率的に似た値になる」ことが重要な性質だと考えられます．
これは，すなわち $||s_i-s_j||$が小さければ$\omega_i$と$\omega_j$は近い値になる ということになります．
このような構造を与える代表例がガウス過程です．
ガウス過程
ガウス過程は位置情報$s_i$が存在する空間上の確率過程です．今，$\omega(s)$をある地点$s$上の効果 (相場感) として$s$の関数とみなすと，ガウス過程とはこの関数$\omega(s)$に対する確率分布を与えるものです．
「関数に対する確率分布」と聞くと難しい感じがしますが，ガウス過程に従うことの定義は有限個の関数値に対する同時分布によって与えられます．すなわち，$s$の候補として任意に$M$個の点$s_1,\ldots,s_M$を取ると，$\omega_1\equiv \omega(s_1),\ldots,\omega_M\equiv \omega(s_M)$の同時分布が
(\omega_1,\ldots,\omega_M)\sim N(0, \tau^2 C(\phi))
で与えられることが (平均0の) ガウス過程の定義です．ここで，$C(\phi)$は$M\times M$の相関行列です．このように，ガウス過程に従うとは 「任意の有限個の地点を適当に取り出したとき，その地点上の相場感の同時分布が多変量正規分布になる」 と解釈できます．
多変量正規分布の相関構造は$C(\phi)$によって決まっています．代表的な$C(\phi)$の構造として，適当な相関関数$\rho(\cdot;\phi)$を用いて$C(\phi)$の$(i,j)$成分を$\rho(||s_i-s_j||;\phi)$のようにする方法があります．相関関数としていくつか有名なものがありますが，代表例としては指数型相関関数$\rho(x;\phi)=\exp(-x/\phi)$があります．
議論を簡単にするため，以降は指数型相関関数を用いた例を考えていきます．その場合，$\omega_i$と$\omega_j$の共分散は
{\rm Cov}(\omega_i, \omega_j)=\tau^2 \exp(-||s_i-s_j||/\phi)
となります．これは，$s_i$と$s_j$の距離が離れるほど相関が小さくなり，距離に応じた相関の減衰の度合いをパラメータ$\phi$がコントロールしていることを表しています．すなわち，$s_i$と$s_j$の距離が近いと相関が高くなり，$\omega_i$と$\omega_j$が似た値を取りやすくなります．
(階層)空間回帰モデル
相場感を表すパラメータ(地点ごとに異なる切片項)$\omega_1,\ldots,\omega_n$を推定するため，ガウス過程に基づく多変量正規分布を事前分布として用いた推定を考えます．その結果，以下のような階層モデルが得られます．
y_i=x_i^\top \beta + \omega_i + \varepsilon_i, \ \ \ \ \varepsilon_i\sim N(0, \sigma^2), \ \ \ i=1,\ldots,n,\\
(\omega_1,\ldots,\omega_n)\sim N(0, \tau^2 C(\phi)).
このモデルは(階層)空間回帰モデルと呼ばれています．
このモデルでは，$(\omega_1,\ldots,\omega_n)$は単なるパラメータではなく潜在変数のような扱いになっています．
上記のモデルのもとで，$(\omega_1,\ldots,\omega_n)$を周辺化した$(y_1,\ldots,y_n)$の同時分布は$N(X\beta, \tau^2 C(\phi)+\tau^2 I_n)$となります．ここで，$X=(x_1,\ldots,x_n)^\top $です．したがって，$(y_1,\ldots,y_n)$は互いに相関していることがわかります．特に，$C(\phi)$が地点間の距離による空間相関を表していたことを思い出すと，誤差項が空間相関している回帰モデルを考えていることと同値であることがわかります．
上記の階層モデルを推定する際は，パラメータ$\beta, \sigma^2, \tau^2, \phi$に対して事前分布を導入し，MCMCによって潜在変数$(\omega_1,\ldots,\omega_n)$とパラメータの同時事後分布からのサンプリングを行うことで事後分布の計算を行うと便利です．詳細は省略しますが，パラメータを所与としたときの潜在変数$(\omega_1,\ldots,\omega_n)$の条件付き事後分布は多変量正規分布になります．他のパラメータについてもほとんどの完全条件付き分布が知られた分布になるので，ギブスサンプラーによるMCMCを用いることができます．
空間変化係数モデル
前節では，回帰モデルの切片だけが空間変化することを考えました．冒頭で記載した「空間異質性」を考慮するために，以下のような回帰係数も空間変化することを考慮した拡張モデル(空間変化係数モデル)を考えることもできます．
y_i=\sum_{k=1}^p x_{ik}\beta_{ik} + \omega_i + \varepsilon_i, \ \ \ \ \varepsilon_i\sim N(0, \sigma^2), \ \ \ i=1,\ldots,n,\\
(\omega_1,\ldots,\omega_n)\sim N(0, \tau_0^2 C(\phi_0)), \ \ \
(\beta_{1k},\ldots,\beta_{nk})\sim N(0, \tau_k^2 C(\phi_k)), \ \ \ k=1,\ldots,p.
このモデルにおいては，切片$\omega_1,\ldots,\omega_n$だけでなく，各説明変数に対する回帰係数$(\beta_{1k},\ldots,\beta_{nk})$も空間変化しており，全てガウス過程によってモデル化されています．
（階層）空間回帰モデルでは不動産の特徴量が不動産価格に与える影響が全地点共通だったのに対し，空間変化係数モデルでは地点によって影響が異なることを考慮したモデル化ができています．その分空間回帰モデルよりも複雑な構造のモデルになっていますが，ギブスサンプラーによって事後分布からのサンプリングを行うことができます．
おわりに
今回はガウス過程を用いた空間データ分析について紹介しました．今回は空間回帰モデルに特化してガウス過程を紹介しましたが，「空間的に異なるパラメータはガウス過程によってモデル化できる」 点が本質なので，他の文脈やモデルでも広く使うことができます．
階層ベイズに基づく空間データのモデリングの方法についてはBanerjee et al. (2014)に詳細が記述されています．
株式会社Nospareでは統計科学の様々な分野を専門とする研究者が所属しております．統計科学の技術アドバイザリーやビジネスデータの分析につきましては株式会社Nospareまでお問い合わせください．
