トップ 差分 一覧 ソース 検索 ヘルプ RSS ログイン

BugTrack-状態空間モデリング/8

R備忘録 /状態空間モデリング/donlp2/その他のメモ

状態空間モデリング - 記事一覧

欠測値の扱い

  • 投稿者: みゅ
  • カテゴリ: なし
  • 優先度: 普通
  • 状態: 完了
  • 日時: 2008年06月20日 11時20分51秒

内容

R言語の状態空間を扱える【sspir】パッケージで、欠測値を扱う

概要

v = y - a
  • の箇所で、vを0とおけばよいだけだから

sspir

  • sspirでは欠測値の取り扱いが現状では出来ない
  • が、元の関数をちょっといじってやるだけですぐにできる
yがNAだったら、参考書で言うvを0にする
  • これだけ
  • 実際のコードと、「co2」データにおける例は後日・・・

コード

  • sspirパッケージの関数定義で【kfilter.SS】という定義がある.この中で時系列方向に、フィルタリングの逐次計算を行っているのが【filterstep】関数
> filterstep
function (y, Fmat, Gmat, Vt, Wt, mx, Cx)
{
    a <- Gmat %*% t(mx)
    R <- Gmat %*% Cx %*% t(Gmat) + Wt
    f <- t(Fmat) %*% a
    Q <- t(Fmat) %*% R %*% Fmat + Vt
    e <- y - f
    A <- R %*% Fmat %*% mysolve(Q)
    m <- a + A %*% e
    C <- R - A %*% Q %*% t(A)
    if (length(y) > 1)
        loglikterm <- log(dmvnorm(as.numeric(y), as.numeric(f),
            Q))
    else loglikterm <- -0.5 * (log(2 * pi) + log(Q) + (y - f)^2/Q)
    list(m = m, C = C, loglikterm = loglikterm)
}
<environment: namespace:sspir>
  • この関数の【e】変数が参考書で言う「v」.なので該当箇所を
      if( any(is.na(y)) ) {
        e <- f
        e[] <- 0
      }else {
        e <- y - f
      }
        A <- R %*% Fmat %*% mysolve_c(Q)
      if( is.na(y) ) {
        A[] <- 0                       --- 
      }
  • のように書き換えてあげればまったく問題なく.
  • この場合、欠測値をNAでyに与えておく
    • いま気がついたのだけれども、なぜ,里茲Δ糞述をしたのか思い出せない・・・
    • 【e】が0なんだから【A】がどんな値をとろうが
    m <- a + A %*% e
    • には、影響がないように思うのだが・・・.それとも、次の行の
    C <- R - A %*% Q %*% t(A)
    • に関係あるんだったか・・・.条件のところが
      if( any(is.na(y)) ) {
    • でなく
      if( is.na(y) ) {
    • になってるから、また別のケースを考慮してるのだろう、きっと・・・

例(「co2」データを使った状態空間表現の例より)

  • フィルタリングまで

↑状態変数が1階のランダムウォークなので、欠測値があっても直近の値で補完されるだけ

  • 平滑化

季節調整込み

  • フィルタリング
    • もとの系列
    • トレンド + 季節成分

↑季節調整成分を入れると、このように季節成分で補完される

コメント