State Space Model

Summary 状態空間モデル(State Space Modeling)とは、モデル構造として状態過程と観察過程をもつ統計モデルの総称。JAGSとRスクリプトの備忘録。 State Space Modeling 生態学の時系列データ解析では、状態空間モデルがとても役立つ。なぜなら、「状態過程(State process)」と「観察過程(Observation process)」という二つプロセスを明示的に切り分けたモデル構造を持つからである。状態過程は、状態の時間変化をモデル化する(状態とは、集団サイズの増減のような実際には直接観察できないものを指す)。一方、観察過程では、状態に何らかの誤差が加わり、その結果としてデータが得られるという過程をモデル化する。例えば、状態過程では集団動態を表すモデルを当てはめ、そこに何らかの観察誤差が加わった結果として手元のデータが得られているという仮定を置く。これにより、集団増加率やその年変動など、将来予測に大きな影響を及ぼすパラメーターをより適切に推定できる(バイアスが少ない)。 JAGS script 次のような15年分の時系列データ(魚の個体数)があるとする。これに対し、状態空間モデルを当てはめてみる。 plot(dat$X, type = "o", ylab = "Fish abundance", xlab = "Year") State process:状態過程では、集団動態を指数増加関数としてモデル化することにする(対象に応じて変える)。各年の集団増加率\(r_t\)は、正規分布に従って年変動すると仮定し、その平均は\(\mu_r\)、分散は\(\sigma^2_r\)とする(対数スケールであることに注意)。 \[ ln~n_{t+1} = r_{t} + ln~n_{t}\\ r_t \sim Normal(\mu_{r},\sigma^2_r) \] Observation process:この集団動態モデルから予測される値に対し、観察誤差が加わった結果として観察値\(X_t\)が得られるとする。ここでは、観察誤差は過分散を考慮したポアソン分布(Overdispersed-Poisson distribution)に従うと仮定する。 \[ X_t \sim Poisson(n_{obs,t})\\ ln~n_{obs,t} \sim Normal(ln~n_t, \sigma^2_{obs})\\ \] JAGSのスクリプトで書き下す。 model{ # Prior mu.r ~ dnorm(0, … Continue reading State Space Model

Varying Intercept and Slope Model

Summary JAGSによるVarying Intercept and Slope Modelの書き方の備忘録。 Varying Intercept and Slope Model 切片と偏回帰係数がグループごとに変化するモデルを指す。書き下すと、次のような数式になる。 \[ y_i \sim Normal(\mu_i, \sigma^2)\\ \mu_i = \beta_{j(i),1} + \beta_{j(i),2}x_i\\ \beta_{k,j} \sim MVN(\mu_{\beta_k}, \sum) \] 切片および偏回帰係数である\(\beta_k\)がグループ単位\(j\)で変化すると仮定している。これにより、グループごとに切片および偏回帰係数の値が異なる(グループごとに応答変数\(y_i\)の平均水準と説明変数\(x_i\)の効果が異なる)ことを表現する。切片と偏回帰係数には多変量正規分布を事前分布とすることで、切片と偏回帰係数の相関も考慮する。ここでは次のテストデータに、このVarying Intercept and Slope Modelを当てはめてみる。 head(dat, 15) ## Y X groupID ## 1 1.6533548 3.715236 1 ## 2 1.9452796 1.334316 1 ## 3 1.0739305 7.079826 1 ## 4 … Continue reading Varying Intercept and Slope Model

Varying Intercept Model

Summary JAGSによるVarying Intercept Modelの書き方の備忘録。 Varying Intercept Model 切片がグループごとに変化するモデルを指す。書き下すと、次のような数式になる。 \[ y_i \sim Normal(\mu_i, \sigma^2)\\ \mu_i = \beta_1 + \beta_2x_i + R_{j(i)}\\ R_j \sim Normal(0, \sigma_R^2) \] 切片である\(\beta_1\)に対してグループ単位\(j\)で変化する変数\(R_j\)が加わった構造になっている。これにより、グループごとに切片の値が異なる(グループごとに応答変数\(y_i\)の平均水準が異なる)ことを表現する。意味としては同じだが、次のように書くこともできる。 \[ y_i \sim Normal(\mu_i, \sigma^2)\\ \mu_i = \beta_{1,j(i)} + \beta_2x_i\\ \beta_{1,j} \sim Normal(\mu_{\beta_1}, \sigma_{\beta_1}^2) \] ここで、切片のグループ間のばらつきを表す変数を、正規分布に従う確率変数とする。こうすることで、グループごとにそれぞれ切片を推定するのではなく、全体切片(\(\beta_1\)もしくは\(\mu_{\beta_1}\))と分散(\(\sigma_R^2\)もしくは\(\sigma_{\beta_1}^2\))という二つのパラメーターの推定で済ませることができる(グループ数が多くなるほどパラメーター数を節約できる)。ここでは次のテストデータに、このVarying Intercept Model(この場合GLMMと同じ)を当てはめてみる。 head(dat, 15) ## Y X groupID ## 1 0.7620270 5.1462547 1 … Continue reading Varying Intercept Model

Regression Example in JAGS

Summary 様々な確率分布を仮定した単回帰モデルのJAGSスクリプトの例。基本的には、リンク先のスクリプトのLiklihoodの箇所を変更し、各種パラメーターの事前分布を指定してやればよい。 Gaussian Yのとりうる値の範囲:実数 \[ y_i \sim Normal(\mu_i, \sigma^2)\\ \mu_i = \beta_1 + \beta_2x_i \] Y[i] ~ dnorm(mu[i], tau) # tau is precision mu[i] <- b[1] + b[2]*X[i] Poisson Yのとりうる値の範囲:0以上の整数 \[ y_i \sim Poisson(\lambda_i)\\ ln\lambda_i = \beta_1 + \beta_2x_i \] Y[i] ~ dpois(lambda[i]) log(lambda[i]) <- b[1] + b[2]*X[i] Binomial Yのとりうる値の範囲:0以上N以下の整数 \[ y_i \sim Binomial(N_i, p_i)\\ … Continue reading Regression Example in JAGS

How to use JAGS: Simple Linear Regression

Summary JAGSの基本的な使い方の備忘録。ここでは、もっとも単純な単回帰モデルをJAGSで実装する。実装の手順は大きく二つに分けられる。(1)JAGSモデルの記述:JAGS言語を使ってJAGSのモデルファイルを記述する。(2)Rによる実装:RからJAGSのモデルファイルを走らせる。 JAGSモデルの記述 サンプルデータ ここではもっとも単純な例として、(GLMでも分析できる)単回帰モデルの分析をJAGSでする。ここで使用するサンプルデータは次のようになっている。 head(dat1) ## Y X ## 1 11.859762 6.635249 ## 2 4.682198 2.290394 ## 3 14.809464 7.288998 ## 4 7.549894 3.507948 ## 5 8.841815 3.150085 ## 6 13.473219 8.793556 plot(Y ~ X, dat1) このデータに対し、応答変数\(y_i\)を説明変数\(x_i\)で説明する簡単な統計モデルを考える(添え字\(i\)はサンプルのID番号)。応答変数\(y_i\)の誤差は正規分布に従うと仮定し、以下の統計モデルを考える。 \[ y_i \sim Normal(\mu_i, \sigma^2)\\ \mu_i = \beta_1 + \beta_2 x_i \] ここで、\(\beta_1\)は切片(intercept)を、\(\beta_2\)は回帰係数(slope)を表しており、これらが推定対象となる主なパラメーターである。 JAGSの記法 上の統計モデルをJAGS言語で書き下すと以下のようになる。 # JAGS … Continue reading How to use JAGS: Simple Linear Regression