xiangze's sparse blog

機械学習、ベイズ統計、コンピュータビジョンと関連する数学について

Forecasting synchronizability of complex networks from dataを読んだ

Forecasting synchronizability of complex networks from data(PHYSICAL REVIEW E 85, 056220 (2012),PDF)という論文を読んだのでその概略をまとめます。

日本語キーワード

圧縮センシング 力学系 結合振動子 複雑ネットワーク

概略

本論文では同じような方程式に従う力学系(常微分方程式)が離散グラフ状に配置され、相互作用するような場合を想定しています。そして通常離散グラフをあらわす隣接行列は疎であることから圧縮センシングの手法でその結合の仕方(トポロジー)が推測できるということを主張しています。またそこから時系列を推測し、さらに推測できたネットワークの情報からその制御(同期化)が出来ることにも触れています。

この論文は以前本ブログで紹介した"Predicting Catastrophes in Nonlinear Dynamical Systems by Compressive Sensing"(arxiv)の著者の1人Wen-Xu Wangさんが名を連ねているものです。以前の論文では少数自由度の力学系の数式を圧縮センシングで推測することでその挙動を再現していましたが、本論文ではそのような力学系が結合したものの結合度合い(トポロジー)が対象となります。

モデルと例

結合力学系
f:id:xiangze:20140413122015p:plain
(自己結合を内部動作に取り込んだ形式ではf:id:xiangze:20140413122023p:plain)

を対象として、各素子はローレンツ方程式、Henon写像の場合を試行しています。
ローレンツ方程式のパラメータは
σ=10
ρ=28
β=2
のカオスの状態で固定としていますHenon写像も同様にカオス領域のものを使用しています。
結合の様式としてはErdos-Renyi型を想定しており、そのノードあたりのエッジの分布は指数αを用いて
 p(k) \sim k^{- ¥alpha}
と表されます。

観測によるノイズを考慮に入れた定式化として
(以前の論文と同様に)多変数多項式とその係数のベクトル
f:id:xiangze:20140413122041p:plain
(変数の数n=3,最大次数m=3の場合は12次元ベクトル)
各時点での係数と項のベクトルを
f:id:xiangze:20140413122212p:plain
f:id:xiangze:20140413122115p:plain
観測時系列のベクトルC_i(t)を
f:id:xiangze:20140413122230p:plain
とおくと問題は
f:id:xiangze:20140413122239p:plain
という形の方程式になります。モデルを差分方程式にしたものに対してB,Cが疎であるという仮定の下に
この方程式を圧縮センシングによって解くことでb_i,c_iを求め、ネットワークのトポロジーと各素子の時系列を推測することが出来るというのがこの論文の主張です。

さらにb_i,c_iが求まったあとにネットワークの同期状態を判定することが出来ます。
完全に同期した状態からの(微小)変分を
 ¥delta x_i(t)=x_i(t)-s(t)
として
 \frac{d\delta x_i}{dt}= DF(s)\delta x_i-\xi \sum_{j=0}^N G_{ij} DH(s) \delta x_j
を考え(DF,DHはそれぞれ差分方程式の内部の動き、相互作用のヤコビアン)、この線型方程式を対角化したもの
 \frac{d\delta y_i}{dt}=[DF(s)-K_i DH(s)\delta y_i]
*1固有値の分布が同期の状態を表していることになります。固有値がすべて負であれば同期状態は安定であり、そうでなければ不安定ということになります。

結果

推測されたネットワークの正確さを2つの指標を用いて計測しています。
ここでは結合がある素子間とない素子での誤差を区別しており、結合がある素子間に関しての誤差を対象としてます。
予測誤差E_nz=に対して
観測させた変数の比率

R_m=(観測の数)/(未知の係数の数)

のエラーE_nzとなるRの限界値

R_c=inf{R_m| E_nz(R_m)<=0.01}

を定義してそれをネットワークの密度

R_nz=(0でない係数の数)/(未知の係数の数)

とネットワーク指数αに対してプロットし、
R_nzが大、つまりネットワークが密になるにしたがってR_cが上がり誤差が大きくなる傾向にあること、αが大きくノード間の結合数の格差が大きい場合にはR_cがあがり、推定がうまくいくことが示されています(FIG 4.)

2つ目の指標として結合行列の固有値の分布を用いています。最大の固有値と最小の固有値の実際値、観測値の差の比率
f:id:xiangze:20140413122314p:plain
という少しわかりにくい量でそれを評価しています。
観測ノイズが小さくなるにつれ誤差も当然小さくなっていくように思われますが、FIG7では誤差が最小になるのはR_mの値が最小の場合ではないことがプロットで示されています。
計算時間についても触れられています。素子数N,観測データ数R_mの双方に対して線形に増加する傾向にあります。
同期の度合いはMSFの固有値の分布Kによって特徴付けられますが、FIG9では固有値の分布の推測値と実際の値もよく一定することが示されています。

またこのような推測をすることが出来るのは同期していない状態のみですが、ネットワークに擾乱を与えることで一時的に同期していない状態にしてやることで推測が可能になることも記述されています。

制御

ネットワークを同期した状態へと移行させるような制御についても触れられています(page 10)。

追試可能性

pythonpycompsenseなどの圧縮センシングのパッケージとnetworkxを用いれば実験することが可能と思われます。

感想

普通の状態では複雑ネットワークの結合は疎であることから圧縮センシングの手法が使えると考え、その構造を推定できることを立証したのが新規性だと思います。
確率的な見方は直接は現れていませんが含まれた観測値から圧縮センシングの過程で最尤な値を選んでいると考えられます。素子のばらつき(システムノイズ)に関しては考慮されていないように思えました。

ローレンツ方程式、Henon写像を使っているのはあくまで例示であり、カオスを示すような系では普遍的に適用できる(らしい)ということです*2
ただし応用分野は多くなさそうです。工学的には基本的には予測の難しいカオス的挙動は排するべきですし、自然現象(特に生物関連)においても素子の力学系が既知でトポロジーが未知という状況はあまり想像がつきません。神経細胞もモデルが最も現実的かもしれません。気象モデルなどの系ではその結合の仕方(トポロジー)はN次元格子などで既知の場合が多いと思います。
実際には多数、高次の変数からなる力学系である場合は圧縮センシングの対象となる次元が変数N、最大次数mに対してN*mで増えていってしまい、計算が難しくなっていくように思われますが、多項式の空間のなかでの力学系の式が疎でない場合はあまり想像がつかないので大丈夫ではないかと思っています。
実世界の応用では素子が従う方程式が正確にはわかっていない場合もありますが、そのような場合には情報量基準を用いた考慮が必要になるかもしれません。
また対象となる変数が多すぎてその性質がつかみづらい場合には 以前本ブログで紹介したIdentifying dynamical systems with bifurcations from noisy partial observation(arxiv)のような手法を用いて少数自由度力学系に落とし込むことが出来ると思います。

このような離散グラフと力学系を組み合わせた研究ではグラフの統計的な性質と挙動の関係を考察するものを多く見ていたのですが、本論文ではグラフの情報をほぼ推測してしまうという点で斬新でした。複雑ネットワークの制御(同期化)の文脈はほとんど追えていないのですが、そちらの方向からの意義が大きいのかもしれないです。
一部の変数が観測できない場合も線形の制御理論などでは普通ですが、本論文でも取り上げられていないようです。

なんか文章がすごい上から目線になってしまいましたが、大自由度の力学系に研究を進展させていてすごいと思いました。

*1:MSF(master stability function)という形式

*2:一方でこのような研究を実際の現象に即したモデルで行った例はあまり見たことがありません