« 読了：Commandeur, Koopman, & Ooms (2011) 状態空間モデルのためのソフトウェア | メイン | Du & Kamakura (2012) メモ »
2014年10月 2日 (木)
仕事の都合で状態空間モデルのソフトを使っているのだけれど、参考書によってもソフトによっても、記法がバラバラで大変に混乱する。メモをつくったので載せておく。純粋に自分のための覚え書きであります。
話を簡単にするために、単変量のガウシアン状態空間モデルについてのみ考える。状態変数の数をmとする。
Commandeur & Koopman(2007): とりあえず参考書の例として。Section 8.1に基づく。
- 測定方程式: $y_t = z'_t \alpha_t + \epsilon_t,\ \ \epsilon_t \sim NID(0, \sigma^2_\epsilon)$
- 状態方程式: $\alpha_{t+1} = T_t \alpha_t + R_t \eta_t,\ \ \eta_t \sim NID(0, Q_t)$
Koopman, Shepard & Doornik (1998): SsfPackの基になっている論文なのだが、記法がすごくややこしい...
- まず$y_t = \theta_t + G_t \epsilon_t,\ \ \epsilon_t \sim NID(0, I)$ と考え ($\epsilon_t$ の分散が1であることに注意！)
- 信号をさらに$\theta_t = c_t + Z_t \alpha_t$ と分解して、
- $\alpha_{t+1} = d_t + T_t \alpha_t + H_t \epsilon_t$ (ここでも$\epsilon_t$が出てくる！)
- 簡略化して: $(\alpha'_{t+1}, y'_t)' = \delta_t \Phi_t \alpha_t + u_t,\ \ u_t \sim NID(0, \Omega_t)$. ここで$\delta'_t = (d'_t, c'_t)'$, $\Phi_t = (T'_t, Z'_t)'$, $u_t = (H'_t, G_t)'$, $\Omega_t$はブロック対角でないややこしい行列。
- 初期条件: $\alpha_1 ~ \sim N(0, P)$
SsfPack: Commandeur&Koopman(2007)のSection 11.2 より。
- 測定方程式: $y_t = Z'_t \alpha_t + \epsilon_t,\ \ \epsilon_t \sim NID(0, H_t)$
- 状態方程式: $\alpha_{t+1} = T_t \alpha_t + R_t \eta_t,\ \ \eta_t \sim NID(0, Q_t)$
- 簡略化して: $(\alpha'_{t+1}, y'_t)' = \Phi_t \alpha_t + u_t,\ \ u_t \sim NID(0, \Omega_t)$. ここで$\Phi_t = (T'_t, Z'_t)'$, $u_t = (\eta'_t, \epsilon_t)'$, $\Omega_t=\left(\begin{array}{cc}R_t Q_t R't & 0\\ 0 & H_t\end{array}\right)$となる。
- 信号: $\theta_t = Z'_t \alpha_t$
- 初期条件: $\alpha_1 ~ \sim NID(0, P_1)$。さらに $\Sigma = (P'_1, \alpha_1)'$ とする。
- というわけで、結局システム行列は $\Phi_t$ (サイズm+1, m), $\Omega_t$(m+1, m+1), 初期状態$\Sigma$(m+1, m)の3つである。
Rのdlmパッケージ: vignetteによれば、
- 測定方程式:$y_t = F_t \theta_t + \nu_t,\ \ \nu_t \sim NID(0, V_t)$
- 状態方程式:$\theta_t = G_t \theta_{t-1} + w_t, \ \ w_t \sim NID(0, W_t)$
- 初期条件:$\theta_0 \sim N(m_0, C_0)$
- というわけで、システム行列は$F_t$(1, m), $V_t$(1, 1), $G_t$(m, m), $W$(m, m), 初期状態$m_0$(m,1), $C_0$(m,1)となる。
RのKFAS パッケージ: reference manual によれば、
- 測定方程式:$y_t = Z_t \alpha_t + \epsilon_t,\ \ \epsilon_t \sim NID(0, H_t)$
- 状態方程式:$\alpha_{t+1} = T_t \alpha_t + R_t \eta_t,\ \ \eta_t \sim NID (0, Q_t)$
- 初期条件:$\alpha_1 ~ \sim NID(a_1, P_1)$
- というわけで、システム行列は$Z_t$(1, m), $H_t$(1, 1), $T_t$(m, m), $R_t$(m, k), $Q_t$(k, k), 初期状態$\alpha_1$(m,1), $P_1$(m, m)となる (kは$\eta_t$の長さ)。
よく使うdlmとKFASを比較すると、
- 測定方程式の係数行列は、dlmでは$F$, KFASでは$Z$
- 測定方程式の誤差分散行列は、dlmでは$V$, KFASでは$H$
- 状態方程式の遷移係数行列は、dlmでは$G$, KFASでは$T$
- 状態方程式の潜在変数(というか誤差項というか)の負荷行列は、dlmでは単位行列、KFASでは$R$
- 状態方程式の潜在変数(というか誤差項というか)の分散行列は、dlmでは$W$, KFASでは$Q$
- 状態変数の初期状態の平均は、dlmでは$m_0$, KFASでは$a_1$
- 状態変数の初期状態の分散は、dlmでは$C_0$, KFASでは$P_1$
と呼ばれているわけである。やれやれ。
追記: はじめてMathJaxを使ってみた。うわあ、これは便利だ...
雑記：データ解析 - 統計学者のみなさん、記法は統一してもらえないかしら (状態空間モデルの巻)