初めに

これは?　技術論文「Empirical dynamic modeling for beginners」を読んだまとめ
内容は?　詳細な理論ではなく概要と適用にあたってのフロー

かなり雑記なので恐縮です

Empirical dynamic modeling (EDM) とは

非線形時系列解析手法
生態系における種間の影響の分析，降水量から河川水位を予測，と用途は広い

前提条件
- 対象の時系列データは，確率過程ではなく力学系に基づいている
  ※ 真の系は不明でも良い (式の前提が不要 = Equation-free)
  ※ 系を構成する全ての変数が観測できなくても近似できる
- 使用する時系列データを標準化 (平均0, 分散1) しておく
- 使用する時系列データから線形トレンドや周期性を取り除いておく (定常である)
- 各変数は固定された同一のサンプリング周期に則って観測されている

EDMの適用

力学系の複雑さ (次元) の判定
非線形力学系 or 線形確率系の区別と非線形性の定量評価
因果性を持つ変数の判定
予測

1. 力学系の複雑さ (次元) の判定

複雑さ = アトラクタの再構成に必要な独立変数の数，と定義される
※ Takensの定理によると単一の時系列のラグから最高性が可能とされている
{ $ x_{t}, x_{t-1\tau}, x_{t-2\tau}, ..., x_{t-(E-1)\tau} $ }
$E$ : 埋め込み次元数
$\tau$ : ラグ

埋め込み次元数 $E$ については，真の力学系の次元を $D$ とすると
$ E < 2D+1 $
という上限があるが，実際のタスクにおいては $E$ も $D$ も未知であることがほとんどであるため $E$ の推定が必要

$E$ の推定: Simplex Projection法
埋め込み次元数 $E$ を変えてSimplex Projection法を適用することで予測能力を比較する
予測能力の指標としては予測値と真値の「相関係数 $\rho$」や「平均絶対誤差 (MAE)」が挙げられる
($\rho$ とMAEのどちらを考慮するかは諸説ある)

手持ちのデータを半分に分けて評価する事が多い: $X, Y$
→ Leave-one-out Cross-validationすることもあるのでこの限りではない
$X$ : ライブラリ (アトラクタの再構成に使用される "学習データ")
$Y$ : Out-of-sampleな予測計算用データ

2. 非線形力学系 or 線形確率系の区別と非線形性の定量評価

EDMは非線形な力学系から観測されるデータを対象としていることから，手持ちのデータが力学系に基づいているのか確率過程なのかを区別することは非常に重要

力学系における非線形性 = 状態依存性

非線形性 = 状態依存性の定量化: S-map法
Simplex Projectionが予測対象近傍 (状態空間において) の点のみを考慮するのに対して，S-mapは全ての点を考慮する (その代わり部分的に重みをつける)

重み関数: $w({\bf d})=\exp(-\theta{\bf d}/d_m)$
${\bf d}$: 予測点?予測対象の時系列?からライブラリの各点までの距離 (predicteeが指すものがなにか自信ないですここ)
$d_m$: 全ペアの平均距離
$\theta$: 状態依存性 ($\theta=0$のときは線形自己回帰モデルと同義)

$\theta=0$ と $\theta>0$ のS-mapモデルで比較することで，非線形力学系か確率過程かが区別できる．
$\theta=0$ のときに相関係数 $\rho$ が最大であれば対象は確率過程である．

また，相関係数 $\rho$ について，$\Delta\rho=\max(\rho_\theta-\rho_{\theta=0})$ となる $\theta$ が状態依存性 (= 非線形性) の定量評価と言える．

3. 因果性を持つ変数の判定

EDMは変数間の因果性を明らかにするためにも使用される．

因果性検定: Cross-mapping algorithm
単変数から再構成された系のマニフォールドはもともとの系への写像を得る，という考えに基づく手法
変数 $M_1, M_2$ について，Cross-mappingの結果が "収束" すれば同じ系に属する (因果的に結びついている) と捉えられる (Cross-mapping Convergent: CCM) ．ここで言う "収束" とは，ライブラリ長の増加に伴ってcross-mapping能力 ($\rho$) が増加することを指す．

代表的な方法としては，下記2つの検定がともに有意だった場合，収束が有意であると判定される．
ここで，モデリングに使用するライブラリ長を $L_i$ としたとき，「$L_0=$ 埋め込み次元数 $E$」と「$L_{max}=$ 全ライブラリデータ」とする．

Kenddallの $\tau$ 検定:
ライブラリ長の変化によるcross-mapping能力 $\rho(L)$ に単調増加トレンドがあるか
Fisherの $\Delta\rho Z$ 検定
$\rho(L_{max})$ が $\rho(L_{0})$ に対して増加しているか

最後に，
「$M_2$ から $M_1$ に対して収束が有意である」
ことは
「$M_1$ が $M_2$ を引きおこす」
という意味であることに留意する．

4. 予測

前述したSimplex ProjectionやS-mapでも力学系についての予測を行うことができる．
(だからこそ $E$ や $L$ を評価できる)

EDMにおいては，状態空間にける近傍のベクトルは似通った時間的発展を遂げる，という特性に基づいて将来の値を予測するため，状態空間をどのように再構成するかが重要な点になる．

予測のための手法は代表的に3つ
Empirical dynamic modeling for beginnersでは「Resource-Consumer-Predatorモデル」を例に下記3手法を比較している

univariate embedding
名前の通り予測対象の変数1つを使用する手法
最適な埋め込み次元 $E=3$ を例として，{ $ C_1(t), C_1(t-1), C_1(t-2) $ } から状態空間が再構成される．
multivariate embedding
これも名前の通り．予測対象の変数の他に因果性を持つ変数を加える (例では $R, P_1$ ) ．
$C_1$ の最適な埋め込み次元が $E=3$ であることから，因果変数を加えた { $ R(t), P_1(t), C_1(t) $ } によって状態空間を再構成する．変数が3つで次元が4 ( $>3$ ) のときは $C_1(t-1)$ などラグを加えていく様子．ラグを含めると組み合わせが色々と考えられるので探索が必要なのかもしれないし，同時刻での再構成がmultivariate embeddingなのかもしれない (多分後者) ．
multiview embedding
multivariate embeddingの最後で触れた変数とラグの組み合わせについて，考えられる組み合わせを全てトライして予測能力の上位 $k$ アトラクタによって予測する手法．アンサンブルっぽい．
組み合わせ総数 $m$ とすると，元論文では $k=\sqrt m$ と設定している．

一般的に精度は univariate < multivariate < multiview になる．

続き

(5) Tracking strength and sign of interactions
(6) Scenario exploration of external perturbation
と続くが本記事では割愛

データにおける諸問題

あるタスクにおいてすべての変数を観測しきれるわけではない (そもそも既知であることが少ない)
プロセス誤差として現れる
EDM適用の必要データ数
理論的な最小データ長は無い
35-40くらいは最小で欲しい?
欠損値の取り扱い
rEDMでは自動的に除去する用になっているが，性能には間違いなく影響する (悪い意味で)
データの収集方法はしっかり練りましょう

前処理について

各変数を標準化しておく (平均0, 分散1)
各変数を定常化しておく (1階差分や回帰)
なるべく生値を使用する (力学系の特性を壊したり信号の線形性を強調する)
周期性や季節成分に注意 (解決手法などは研究段階)

最後に

夏休みを使って勉強する予定でしたが結局勤務時間中に読みました．
前述した手法を提案しているSugiharaらによってPythonやRで実装されているそうなので，次は手持ちのデータに対して前述したフローを適用していきます．
精度が良ければこちらも記事で書く予定です．