ホリケン's diary

趣味はでぃーぷらーにんぐ

時系列解析をゆるっと解説(2)

knto-h.hatenablog.com

時系列解析をゆるっと解説(1)の続きである。

時系列解析


状態空間モデルによる時系列の解析

  • 状態空間モデル
    時系列解析で用いられる様々なモデルを、状態空間モデルによって統一的に取り扱うことができる。一方、時系列解析の多くの問題が状態空間モデルの状態推定の問題として定式化できる。かなり重要
     y_nl変量の時系列とする。この時、この時系列を表現する次のようなモデルを状態空間モデル(state-space model)と呼ぶ。 $$ x_n = F_n x_{n-1} + G_n v_n (システムモデル) $$ $$ y_n = H_n x_n + w_n (観測モデル) $$ ここで、x_nは直接には観測できないk次元のベクトルで、状態(state)と呼ばれる。[texv_n]はシステムノイズあるいは状態ノイズと呼ばれ、平均ベクトル0、分散共分散行列Q_nに従うm次元の正規白色雑音である。一方、w_nは観測ノイズと呼ばれ、平均ベクトル0、分散共分散行列R_nに従うl次元の正規白色雑音とする。F_n,G_n,H_nはそれぞれk \times k, k \times m, l \times kの行列である。
    状態空間モデルは次のように2通りの解釈ができる。
    (1)観測モデルを時系列y_nが観測される仕組みを表現する回帰モデルと考えると、状態x_nはその回帰係数となる。この時システムモデルはその回帰係数の時間的な変化の様子を表現するモデルとなる。
    (2)状態ベクトルx_nを推定すべき信号と考えると、システムモデルは信号の発生メカニズムを表すモデル、観測モデルはその信号を実際に観測するとき信号が変換されノイズが加わる様子を表している。

  • カルマンフィルタによる状態の推定
    観測値Y_j = [y_1,\cdots,y_j]に基づいて時刻nにおける状態x_nの推定において、j < nの場合は予測、j=nの場合はフィルタ、j > nの場合は平滑化と呼ばれる。
    これらの状態推定の問題に答えるためには、観測値Y_jが与えられたもとでの状態x_nの条件付き分布p(x_n|Y_j)を求めれば良い。この時、求める必要のあるものは条件付き分布を規定する平均ベクトルと分散共分散行列である。通常この計算には莫大な計算量を要するが、状態空間モデルにおいては逐次的な計算アルゴリズムによって効率的に計算することができる。これがカルマンフィルタである。

  • 平滑化のアルゴリズム
    平滑化に関しては固定区間平滑化というアルゴリズムがある。フィルタが時刻nまでの観測値だけを用いてx_nを推定しているのに対し、平滑化のアルゴリズムは得られている全ての観測値を用いて推定を行なっているため、平滑化を行えば一般にフィルタよりも精度の良い状態推定が行えることになる。

  • 欠損値の補間
    時系列の状態空間モデルを用いると時系列のモデルがすでに与えられて入れば欠損値を含むデータに対しても厳密な尤度を計算でき、パラメータの最尤推定値を求めることができる。またそのモデルを利用して(欠損値の部分はフィルタリングのステップを省略しながら)カルマンフィルタにより予測分布、フィルタ分布を求め。その後平滑化のアルゴリズムを適用することにより欠損値の推定が行える。


ARMAモデルの状態空間表現

  • モデル表現
    ARMAモデルに関してはpart1の方に書いてあります。(m,l)次のの定常なARMAモデルを状態空間モデルとして解く場合以下のようにモデリングします(k =max(m, l-1)。
    $$ F = \begin{bmatrix} a_1 & 1 & & \\ a_2 & & \ddots & \\ \vdots & & & 1\\ a_k & & &
    \end{bmatrix} , G = \begin{bmatrix} 1 \\ -b_1 \\ \vdots \\ -b_{k-1}
    \end{bmatrix} $$ $$ H = \begin{bmatrix} 1 & 0 & \cdots & 0 \end{bmatrix} $$ これにおいてi>mの時a_i=0, i>lの時b_i=0と定義しておけば良い。すなわち、ARMAモデルは、係数行列F,G,Hが時間によらず一定で、観測ノイズが0の状態空間モデルで表現できる。


トレンドの推定(時系列モデリング手法)

経済時系列などの長期間にわたっての変動傾向を含む時系列の解析法として多項式回帰モデル。

  • 多項式トレンドモデル
    トレンドを推定する最も簡単な方法。時系列y_nが多項式と残差の和によって表されるとしたもの。 $$ \begin{align} y_n = t_n + w_n \\ t_n = \sum_{ i= 0}^{m} a_i x_n^{i} \end{align} $$ これは最小二乗法のハウスホルダー変換によって最適化が行える。

  • トレンド成分モデル-構造の確率的変化のモデル-
    多項式トレンドモデルでは真のトレンドが実際に多項式に従いw_nが白色雑音とみなせる場合は有力だが、そうでない場合やデータに追従しすぎて偶然の変動を拾ってしまうことがある。それを改善する種王がトレンド成分モデルで、以下の制約式を加える。 $$ \Delta^{k} t_n = v_n $$ これは状態空間モデルを用いて表現することができる。 $$ F = \begin{bmatrix} c_1 & c_2 & \cdots & c_k \\ 1 & & & \\ & \ddots & & \\ & & 1 &
    \end{bmatrix} , G = \begin{bmatrix} 1 \\ 0 \\ \vdots \\ 0
    \end{bmatrix} $$ $$ H = \begin{bmatrix} 1 & 0 & \cdots & 0 \end{bmatrix} $$ $$ x_n = \begin{bmatrix} t_n & t_n-1 & \cdots & t_n-k+1 \end{bmatrix} $$

$$ x_n = F_n x_{n-1} + G_n v_n $$ $$ y_n = H_n x_n $$

  • トレンドモデル
    実際に観測される時系列にはトレンドにノイズが加わったものである。トレンド成分モデルに白色雑音w_nを加えたものがトレンドモデルである。 $$ x_n = F_n x_{n-1} + G_n v_n (システムモデル) $$ $$ y_n = H_n x_n + w_n (観測モデル) $$


季節調整モデル(時系列モデリング手法)

  • 季節成分モデル
    このモデリングの考えは、トレンドとは別に周期的な傾向のある時系列データを使用する際にトレンド成分t_n, 季節成分s_n, 白色雑音w_nの三つの成分に分解して y_n = t_n + s_n + w_nと表すことである。s_nはトレンド成分と同様に以下の状態空間モデルとして記述することができる。 $$ F = \begin{bmatrix} d_1 & d_2 & \cdots & d_{l(p-1)} \\ 1 & & & \\ & \ddots & & \\ & & 1 &
    \end{bmatrix} , G = \begin{bmatrix} 1 \\ 0 \\ \vdots \\ 0
    \end{bmatrix} $$ $$ H = \begin{bmatrix} 1 & 0 & \cdots & 0 \end{bmatrix} $$ $$ x_n = \begin{bmatrix} s_n & s_n-1 & \cdots & s_n-k+1 \end{bmatrix} $$

  • 標準的季節調整モデル
    上で定義した季節成分を取り入れてモデリングすると、以下のように記述できる。 $$ F = \begin{bmatrix} F_1 & 0 \\ 0 & F_2 \end{bmatrix} ,G = \begin{bmatrix} G_1 & 0 \\ 0 & G_2 \end{bmatrix} ,H = \begin{bmatrix} H_1 & H_2 \end{bmatrix} $$

  • 曜日効果成分を含むモデル
    トレンド成分と季節成分を組み合わせて状態空間モデルを作った要領で、曜日効果成分を加えるだけ。同様のやり方で別の成分を加えることができる。


非ガウス型モデル(時系列モデリング手法)

状態空間モデルのシステムノイズ、あるいは観測ノイズが非ガウス型の場合への拡張。パラメータの急激な変化や異常値に強い。

  • 非ガウス型状態空間モデルと状態推定
    通常の状態空間モデルでは、ノイズを正規分布に従うことを前提としていたが、非ガウス型では任意の密度関数の分布に従う白色雑音だとする。これを非ガウス型状態空間モデルと呼ぶ。この時、カルマンフィルタの更新式にはベイズ更新を用いて計算する。実際の計算の際には非ガウス型の密度関数を以下にあげる分布の近似を使って解く。
    1.線形ガウスモデル → カルマンフィルタで計算
    2.正規分布近似 → 拡張カルマンフィルタで計算
    3.数値的に近似 → 階段関数近似、折れ線近似、スプライン近似