岩波データサイエンスシリーズ「時系列解析」読書メモ

機械学習技術 人工知能技術 デジタルトランスフォーメーション技術 センサーデータ&IOT ICT技術 ストリームデータの処理 ICTインフラ技術 確率的生成モデル 関係データ学習 サポートベクトルマシン スパースモデリング 異常検知・変化検知技術 時系列データ解析 本ブログのナビ

サマリー

時系列データとは、株価や気温、トラフィック量などの時間の経過に応じて値が変化するデータのことを呼ぶ。この時系列データに対して機械学習を適用することで、大量のデータを学習し、未知のデータに対する予測を行うことで、ビジネスの意思決定やリスク管理などに役立てることができる。

それらの中で代表的な手法としては、ARIMA(自己回帰和移動平均モデル)、LSTM(長短期記憶モデル)、Prophet(Facebookが開発した時系列データの予測に特化したライブラリ)、状態空間モデルなどが用いられる。これらの手法は機械学習に基づいた予測手法で、過去の時系列データを学習して未来の予測を行うものとなる。

ここでは岩波データサイエンスシリーズ:「時系列解析−状態空間モデル・因果解析・ビジネス応用」をベースに、この時系列データ解析に対する状態空間モデルの適用について述べる。

状態空間モデルは、時系列データの解析によく使われる統計モデルの一つであり、汎用的な枠組みとなる。状態空間モデルでは、観測される時系列データが、何らかの確率過程によって生成されると考え、その過程を記述する数式モデルを構築する。また、その確率過程を支配する「状態変数」と呼ばれる非観測変数を仮定し、その状態変数の時間的変化を表す遷移モデルと、状態変数から観測される値を生成する観測モデルを考える。

今回は読書メモについて述べる。

岩波データサイエンスシリーズ「時系列解析」読書メモ

読書メモをいかに記載する。

状態空間モデル

1. 状態空間モデル

概要
状態空間モデルは多くの時系列モデルを統一的に取り扱うためのモデル
時系列の予測、補間、成分分解やバラメータ推定等の状態推定の問題を取り扱う

(1) 状態とは 時系列の予測と状態
気温や株価のように時間とともに変動する現象の記録を時系列とよぶ
時刻nの観測値yn
物理的なシステムでは、過去の情報は現在を経由して将来に伝えられる
過去の情報のうち将来の動きに関するものは現時点に集約されている
集約された情報のことを状態と呼ぶ
状態さえあれば過去の情報は捨ててしまってみ未来が予測される

(2) 状態空間モデル
状態を使って時系列がどのように表現されるのか?
状態xn-1を知っているときに、次の時刻の状態xnがどのように表されるか?
状態空間モデル、観測モデル、システムモデル
システムモデル
xn=Fn(xn-1)+Gn(vn)
Gn(vn):ノイズ
観測モデル
yn=Hn(xn)+wn
状態空間モデル
観測モデルとシステムモデルを連結したもの
一般形
システムモデル
Xn ~ qn(・|xn-1)
観測モデル
Yn ~ rn(・|xn)
例1:線形・ガウス型状態空間モデル
関数Fn,Gn,Hnが線形関数
Fn(x)=Fnx, Gn(v)=Gnv, Hn(x)=Hnx
システムモデル
xn=Fnxn-1 + Gnvn
観測モデル
yn=Hnxn + wn
例2:AR(自己回帰)モデル
2次のAR(自己回帰)モデル
Yn = a1yn-1 + a2yn-2 + vn vn ~ N(0, σ2)
線形ガウス空間モデル
時系列の変動特性を捉えるのに周波数成分に分けるのがポピュラー
フーリエ変換をしてペリオドグラムを求めると多数のピークが現れて意味のあるピークの判断が難しい
ペリオドグラムの平滑化には経験的な要素が強く入る
ARモデルの次数mが求められると周波数fのパワースペクトルが自動で求められる
スペクトルがk個のピークを持つためにはAR次数が2k以上必要
例3:成分分解モデル
時系列ynをJ個の成分xn(j), j=1,…,J と観測ノイズwnに分解する
Yn = xn(1) +…+xn(J) +wn
xn(j)は状態空間モデル表される
システムモデル
xn(j) = F(j)xn-1 +G(j)vn
観測モデル
Yn = H(j)xn +wn
状態空間モデル
Xn =Fxn-1 + Gvn
F,G,H
例4:時変係数モデル
自己回帰モデルにおいて、 自己回帰係数ajを時間とともに変化する係数an,jで置き換えたもの
an,j = an-1,j +vn,j vn,j=N(0,τ2)
自己回帰モデルと係数の変化を表す”ランダムウォークの概要とアルゴリズム及び実装例“でも述べているランダムウォークモデルは状態空間モデルを用いて
Xn = Fxn-1 + Gvn
vn~N(0,σ2)
Yn = Hnxn +wn
wn~N(0,τ2)
例5:時変分散モデル
時系列ynは正規分布に従うが、その分散は時事刻々と変化し、yn=εn,εn~N(0,σn2)となる場合
分散σn2の時間変化がlogσn2 = logσn-12 +vn
状態をxn=logσn2とすると状態空間モデルは
xn=xn-1 + vn
Yn = exp(xn/2)wn
Wnは標準正規分布
定常時系列モデル
時間とともに平均、分散、共分散関数などの確率的構造が変化しない
非定常時系列モデル
定常性が成り立たないもの

2. フィルタリング-状態の逐次推定

(1) 状態推定の問題

時系列の観測値Yj={y1,…,yn}と状態空間モデルが与えられた時
状態xnを推定する問題は「状態推定」と呼ばれる
観測値の最終時点nと推定する状態の時刻jの関係により 3つの場合に分けられる
n<jの場合 :予測
n=jの場合: フィルタ
n>jの場合: 平滑化
状態の予測値だけではなく 分散または標準偏差も分からなければ 合理的な意思決定には役立たない
明日の気温の予測が同じ10度でも 標準偏差が1とぬわの場合では判断は異なる
尤度に関しても予測値と分散が必要

(2) 逐次フィルタと平滑化

観測値Yjが与えられた時、 状態推定の問題は Yjの元で状態xnの条件付き分布p(xn|Yj)を求める問題になる
予測分布はp(xn|xn-1,Yn-1)=p(xn|xn-1)
以下の計算式を使って逐次的に計算可能
一期先予測
時刻n-1のフィルタ分布が与えられると、時刻nの予測分布が計算できる
フィルタ
予測分布と最新の観測分布から求めることができる
予測とフィルタの操作を繰り返すことで逐次フィルタができる
平滑化のアルゴリズム
固定区間平滑化
固定ラグ平滑化
固定点平滑化

(3) 線形・ガウス型状態空間モデルとカルマンフィルター

関数FnとHnが線形で、白色雑音vnおよびwnが正規分布に従う場合に「線形・ガウス状態空間モデル」と呼ばれる
カルマンフィルタアルゴリズム

(4) 時系列モデルの尤度計算とパラメータ検定

パラメータθを持つ時系列モデルの尤度は、長さがNり時系列Yn={y1,…,yn}が与えられるときに同時確率密度関数を用いてL(Θ)=fN(YN|Θ)

(5) 状態空間モデルによる成分分解 季節調整

線形・ガウス型の状態空間モデルの応用例
経済時系列の成分分解
月次経済時系列には、 上昇や下降などの中長期傾向と、 毎年同様なパターンを繰り返す年周期変動が見られる
観測した時系列をトレンドと呼ばれる傾向成分と 季節成分および偶然変動に成分分解する方法
季節調整法
トレンドの代表的なモデル
2次のトレンドモデル:tn=2tn-1-tn-2+vn
時系列の変動のうち毎年繰り返して現れる変動パターンが「季節成分」
さらに「曜日効果項」を加える
状態空間モデルでの曜日効果
欠損値の処理
逐次フィルタではある観測値が欠測となった場合は、対応するステップを省略するだけで良い
フィルタリング後に平滑化を行って欠損値を補間

(6) 非線形・非ガウス型状態空間モデルのフィルタリング

概要
線形・ガウス型の状態空間モデルでは、 状態の条件付き分布p(xn|Yj)は正規分布となり、 カルマンフィルタにより効率的に計算できる

1. 解析的近似による方法

モデルが非線形あるいはノイズの分布が非ガウス型の場合は別のモデルが必要
ガウス近似
混合ガウス近似
M個の正規分布φiの和で近似
カルマンフィルタを適用可能
計算量が莫大になる

2. 数値的近似による方法

階段関数近似
分点数nの階段関数で近似
分点数を数百以上にすると、任意の分布をかなり高い精度で近似可能
4次元以上の状態空間モデルに対しては適用困難

3. 粒子近似による方法

粒子近似
概念図
独立に生成されたみなせる多数の粒子を用いて分布を表現
密度関数の高いところには多くの粒子が集中している

(7) 粒子フィルタ

概要
粒子フィルタ
予測、フィルタ、平滑化の分布をm個の粒子を用いて近似する
各粒子点で1/mずつステップする階段関数を作成し「経験分布関数」とする
粒子数mが増大すると真の分布関数に収束する
粒子数mはモデルの複雑さや必要な制度によって、1000個から10万個程度が用いられる
予測およびフィルタの分布を近似する粒子を 逐次的に生成するための粒子フィルタのアルゴリズム
予測分布に従う粒子{pn(1),…,pn(m)}は一期前のフィルタの粒子{fn-1(1),…,fn-1(m)}から生成できる
フィルタの粒子{fn(1),…,fn(m)}は予測分布の粒子{pn(1),…,pn(m)}のリサンプリングで生成可能

一般予測
予測のステップでは、一期前の状態xn-1の条件付き分布p(xn-1|Yn-1)に従う
粒子{fn-1(1),…,fn-1(m)}とシステムノイズvnの粒子{vn(1),…,vn(m)}がそれぞれm個与えられているものとする
新しい粒子pn(j)をpn(j)=Fn(fn-1(j), vn(j))によって定義

フィルタ
フィルタのステップ
リサンプリングはなぜ必要なのか
リサンプリングを行わずに予測とフィルタのステップだけを繰り返すと、多くの粒子の重みはゼロになって、分布の縮退が発生する
リサンプリングの目的は、 重み{αn(1),…,αn(m)}を持つ粒子{pn(1),…,pn(m)}で定義された分布関数を、 重みが等しい経験分布関数の形で表現し直すこと
厳格なランダムサンプリングを行うことは必ずしも必要でない
層化抽出を行うことで単純サンプリングよりも精度よく再近似を行うことができる

粒子フィルタのアルゴリズム
粒子フィルタアルゴリズム

(8) 平滑化

平滑化への粒子フィルタの適用

(9) 例 地震波の変動分の推定

地震波の観測記録の例
P波、S波の到着により、その分さんは大きく変化する
変化する分散を推定する方法(変数変換 zn=log((y2n-12+y2n)/2)
原系列の分散の大きさが、系列の高さに対応する
原データでは分散が正規分布でも、変換後は二重指数分布に従う
トレンド分析

粒子mの影響
粒子数を1/000にした場合
m=100でも中央値の推定はできるが、分布は判定できなくなっている

3. 非線形状態空間モデルでどんなことができるか

(1) 自己組織型の状態空間モデル

状態とパラメータの同時推定

(2) 非線形平滑化

(3) 非ガウス型分布モデル

欠陥値の処理
リサンプリングは必要なのか
リサンプリングにおいて忠実なランダムサンプリングは必要か
ラグLは大きい方が良いか

Rによる状態空間モデル

1. dimパッケージ

概要
Giovanni Petrisが開発したRのパッケージ
線形・正規分布の動的線形モデル (dynamic linear model)を扱う
モデルの記法
yt = Ftθt + vt vt ~ N(0,Vt)
Ytは時点tにおける観測値
Θtは時点tにおける状態
vtは観測ノイズ
Ftは係数行列
Vtは観測ノイズの分散共分散行列
Θt =Gtθt-1 + wt wt ~ N(0,Wt)
wtはシステムノイズ
Gtは係数行列
Wtはシステムノイズの分散共分散行列
Θ0 N(m0, C0)
m0は初期状態の期待値のベクトル
C0は初期状態の分散共分散行列

(1) dimによるローカルレベルモデル

サンプルデータ
水準のみがランダムウォークするモデル (一階差分モデル)
yt,θt,Ft,Gt,Vt,Wtは全てスカラー
Ft,Gt,Vt,Wtは全て時間tに関わらず一定 (F,G,V,W)
dlmでローカルモデルを扱うときはdlmModPoly関数を使う
多項式で表現される動的線形モデルを構築する
Order引数に1を与えると、ローカルレベルモデルを指定する
dlmではモデルをオブジェクトとして構築し、それに対して操作を行う
観測ノイズの分散(dV)とシステムノイズの分散(dw)を与える必要がある
ひとまず1とする

(2) dimによる季節調整モデル

2. KFASパッケージ

(1) KFASによる季節調整モデル

(2) ポアソン分布モデル

3. 終わりに

状態空間モデルのパラメータ推定

粒子フィルタを実装してみる

R言語のプログラム
プログラムの説明
動かしてみる
追加の説明
コーシー分布を使ってみる
「データ同化」という考え方
この先のこと

推定のための計算方法の比較

状態空間モデルのマーケティングへの応用

1. 市場反応モデル 基本構造

2. モデルの拡張

(1) 進化と深化
(2) 回帰係数を事変にする
(3) 潜在構造を取り込む

3. 解析事例

4. まとめ

VARモデルによる因果関係の推論

1. データの前処理

欠損値の補間
定常性の確保

2. VARモデルの推定とラグ選択

VARモデルの推定
VARモデルのラグ選択

3. 制約付きVARモデルの推定

4. 因果推定とインパルス応答関数

因果推定
インパルス応答関数
まとめに代えて

隠れマルコフモデルと状態空間モデル

時系列データから天体の形を再構成する 時間軸天文学への誘い
降着円盤とは?
降着円盤の形を反映している振動現象
降着円盤の新しいトモグラフィー手法
アンドロメダ座V455の降着円盤

シミュレーションとデータサイエンス

シミュレーション、データ同化、そしてエミュレーション

1. 内挿と外挿
2. データ同化(データあシミュレーション)
3. 状態空間モデルと逐次フィルタ
4. エミュレーション
5. DNNとの結合

エミュレータの活用

天気予報とデータサイエンス

静かなる革命

1. 天気予報の仕組み
2. データ同化
3. 台風予報と機動的観測
天気予報士の仕事

揺らぐタンパク質と老いる私

タンパク質は折りたたまれる
タンパク質は揺らぐ
タンパク質はミスフォールドする
タンパク質のシミュレーション
シミュレーションの多変量解析
時系列から動的な情報を取り出す
一般化固有値問題

分子シミュレーションの「逆問題」とは

リアルなSimCityの夢

1. 問題解決ツールとしてのSimCity
2. リアルなSimCityの実現に向けた課題
3. リアルなSimCityを実現するための技術
4. おわりに

夢と脳と機械学習

夢の生理学
夢の機能とフロイトの夢理論
夢と「生成」モデル
夢内容と脳
ブレイン・デコーディング
夢の脳計測実験
夢のデコーディング
おわりに
計算機で作る面白いナンプレ

コメント

  1. […] 岩波データサイエンスシリーズ:「時系列解析−状態空間モデル・因果解析・ビジネス応用」読書メモ […]

タイトルとURLをコピーしました