《講義ノート》

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)

動的に変化するネットワークと結合力学系— 数理モデル構築からデータ解析まで —京都大学大学院・情報学研究科

青柳

富誌生

多数の要素が相互作用することで生じる非自明な現象に関して、従来の規則的な空間上の系に止まらず、最近では複雑なネットワーク上の相互作用を仮定した解析も行われています。そこで重要なことは、単に静的なネットワーク上の相互作用を考えるだけでは不十分であり、個々の要素の状態に応じてダイナミックにネットワーク構造が変化することを考慮すべき状況にしばしば直面することです（例えば、神経や感染症のネットワーク）。このような系を解析する場合、感染症など
個々の現象の特性に特化した数理モデルを用いるのが王道でしょう。一方、幅広い現象に対して限定的でも適用可能な汎用性の高い数理モデルを構築して解析する方法もあります。実際、イジング
スピン系などは磁性体に限らない脳や社会現象のモデルなど幅広い現象に拡張して応用されています。このアプローチでは、できるだけ普遍（あるいは不変）的な本質を捉えるため、理論的汎用性からリミットサイクル振動子や散逸と拡散が最初に考えるべきダイナミクスの候補となります。
たとえ素子がシンプルなダイナミクスであっても、相互作用のネットワークが動的に変化する場合
は、非自明な興味深い現象がしばしば見られます。この集中ゼミでは、動的ネットワーク上の力学系を解明するための第一歩として、以上の 2 つの解析例を紹介します。最後に、現実のデータから相互作用ネットワークを同定することも今後の重要な課題です。現代的なベイズ統計と力学系の理論を組み合わせれば、そのようなことも可能になりつつあり、特にリズム間の相互作用を推定する
いくつかの事例を概説します。

1

はじめに物理学を学び研究する皆さんは、一見複雑に見える物理現象が実はシンプルな幾つかの基本原理に従っており、たとえ複雑な振る舞いを示しても原理的には理解可能な事であると信じているでしょう。実際に物理学の理論は、日食などの天文現象の正確な予測から超電導や超流動などの一見不可解な現象の説明まで、幅広い自然現象の説明に成功しています。この輝かしい成功から、生命現象や社会現象までを同様に説明できるのではないかと期待が膨らむのは自然なことです。ただ、現状では生命現象や社会現象を物理現象のように説明できると断言できる段階には遥かにほど遠い状況であることは事実です1 。例えば、新型コロナの感染症の流行を予測する数理モデルの構築は非常に難しい問題であることは、皆さんも何となく感じていると思います。ただ、水晶玉による予言や星占いと本質的に異なるのは、基本的な感染症伝播のプロセスは数理モデルで記述できており、実験的にもそ1

物理現象に限ったとしても、非平衡・非線形の絡む現象の理論的解明は現在でも継続して研究されている
ホットなトピックです。

《講義ノート》

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)

れを示す事実がある点です。最近では、現代の人の社会ネットワーク構造を取り入れたス
ケールフリーネットワーク上の感染症流行の理論的研究もされており、通常の平均場的な場合（中世の村社会に対応する）に比較して、極めて低い感染率でも流行が起こることが知られています [1, 2]2 。
それでは、なぜ生命現象や社会現象の理論的な解析が難しいのでしょうか？例えば、上記感染症の例で考えてみましょう。感染症の広がりは人との接触であり、例えば友人関係などの社会的ネットワークを通じて広がると考えられます。そのネットワーク構造は、物理系でよく仮定される規則格子とは異なり、ランダム性もある複雑で非一様な構造です。
これに関係して、いわゆる複雑ネットワーク (complex network) の研究は近年盛んに行われ、スピン系や振動子などのネットワーク上での振る舞いが理論的に研究されています。
単純なランダムネットワークとは異なる、そのような点が理論解析を困難にしている側面はありますが、一定の理論的進展もあり克服が不可能とまでは言えないと思います [3]。それ以外に予測を困難にしている、今まであまり考慮されていない重要な要因があります。
それは、ネットワーク構造が動的に変化することです。社会や生命・工学のシステムは、しばしばネットワーク状の相互作用をしている結合力学系として見なすことができます。
ここで、力学系とはある時刻で系の状態がわかれば、その後の系の状態が固有の時間発展のルールにより基本的に決定されるような系の事をいいます3 。ネットワークの用語で言えば、ノードは個々の力学系（感染症なら人、神経系ならニューロン）に、リンクはノード間の相互作用（感染する機会、シナプス結合）に対応します。例えば、感染が広がると、その情報に依存して人々の行動に変容が起こり、実質的なネットワーク構造が変化すると考えられます。どう変化するかは、人の心理や個性も絡み数理モデルに取り込むのは難し
いと考えられます。しかし、ここで強調したいのは、そもそもネットワーク構造（相互作用の構造という意味です）が変化するシステムの研究は、シンプルなものも含め少なく、しばしば予想外の結果が生じて解析が困難になる点です。そのような系を理解する上で困難となる要因は、一般的に次の３つに整理されます。

(i) 個々の系のダイナミクスの複雑さ（人の行動やニューロンの挙動）(ii) 相互作用ネットワークの複雑さ（友人関係、シナプス結合の非一様性など）(iii) 相互作用ネットワークが、個々の系の状態に依存して動的に変化すること（動的ネッ
トワーク構造）
ここでは、最後の (iii) 動的ネットワーク構造に着目して話をしたいと思います。なぜなら、
2
現代になって、新たな病気が次々に出現するのは、どこかの機関が細菌兵器を開発している陰謀論を仮定しなくても、社会構造のネットワークの変化から基本的には説明が可能だということです。
3「力学系」という用語は英語の dynamical systems を翻訳したものであり、狭い意味での「力学」（mechanics）より広い概念です。ノイズ項などのランダムネスを導入する場合は、確率的なルールも含めて考えます。

《講義ノート》

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)

ࣀ࣮ࢻᅛ᭷ࡢࢲ࢖ࢼ࣑ࢡࢫࢆࡶࡘຊᏛ⣔Dynamics of Network

X

0

Dynamics on Networkࣜࣥࢡྥࡁࡶྵࡵࡓ┦஫స⏝

図 1: ネットワーク上の結合力学系の状態変化とネットワーク自体の構造が相互に影響しあって時間発展する系：相互作用ネットワークがノード上の力学系の次の状態に影響し（Dynamics on Network)、そのノード上の力学系の状態全体がネットワーク構造に影響（Dynamics of Network) を与え、結果として非自明な興味深い振る舞いが見られます。
既に (i) や (ii) は歴史的にも経験済みであり、理論的にも克服した事例が幾つかあるからです。
まず (i) の個別の系の数理モデル構築の困難さは、常に問題となる点です。これに関して肝要なことは、現実の全ての側面を忠実に再現することを目指すのではなく、本質的に重要なダイナミクスを再現し、かつ、なるべくシンプルな数理モデルを構築することを目指すことです4 。例えば、Ising spin 系の 2 状態モデルは磁性の本質的な部分を単純化し、多くの理論的知見をもたらしました。更に、その拡張したモデルを活用して、神経系の連想記憶モデルや社会における投票者モデル等で興味深い理論的知見も得られています。
ニューロンや人のダイナミクスとして何を本質と見るかは、説明したい現象に依存しセンスが問われる部分でもあり、何にでも通用する処方せんが存在するわけではない点が悩ましいところです。しかし、本稿で力学系という観点から一つの方針は示唆される点に触れたいと思います。一方、(ii) の相互作用の複雑さに起因する難しさですが、これに関してはランダム系の物理の手法がしばしば有用です。もちろん、単純なランダム性だけでは本質を見失う場合もあり、例えばスケールフリーネットワークなどの特徴も取り入れた理論的解析が必要です。これに関しては一定の理論的進展もありますので、必要な知見に関して簡単に触れる程度としたいと思います。
現実の興味深い系は、しばしば動的ネットワーク上で相互作用する結合力学系として見4

これは物理の人は馴染み深い考え方でしょうが、生命系や社会系ではその是非も含め慎重に考える必要もあると思います。

《講義ノート》

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)

なすことができ、そのネットワーク構造とノード上の状態が相互に影響しあう事で、全体として変化・発展していきます (図 1)。繰り返しになりますが、本稿ではこの (iii) に注目しますが、その準備としてダイナミクスを記述する数理モデルに関して少し議論したいと思います。我々が興味を抱くシステムは、多数の要素が相互作用している結合力学系と見なせるケースが多いです。脳神経系や感染症における人のネットワーク5 などの例をあげましたが、そこではまず個々の要素の状態とそのダイナミクスを記述する必要があります。
その際、「筋の良いモデル」を構築できるかは、しばしば個人のセンスや時の運に依存していることも多いでしょう。しかし、もう少し系統的に「筋の良いモデル」を見つける方策は無いのでしょうか？一つの方策として、力学系の「構造安定性」6 に着目することで、限定的な状況ですが、ある種の「筋の良いモデル」を構築できる場合があることをお話したいと思います。

1.1

分岐理論と構造安定性

個々の要素の状態が、複数の変数 X で記述され、その時間発展が次の微分方程式

dX= F (X, μ)dt(1)で記述される場合を考えてみましょう。ここで μ は外部から設定するパラメーターとします。具体例で説明すると、ニューロンの場合は、膜電位や細胞内カルシウムイオン濃度、
イオンチャネルのゲートの開閉度等が状態変数 X （複数の変数をまとめて書いています）に対応し、外部から電極により加える電流が μ に対応します。このケースでは、電流を増加させると、時間に依存しない定常な状態から、突然周期的なスパイク（活動電位と呼ばれる）を生成する状態に変化します。このように系の外部からコントロール可能なパラメーターを変化させたとき、定常状態から周期的なリズムが生じる現象は、広く様々な系に見られます。その事実から一般化すると、微分方程式 (1) は様々な形をしていても、定常状態から周期状態に移り変わる変化のダイナミクスは、何か普遍的な記述があるのではないかと期待できます。一般的に、μ を変化させたとき解の性質が質的に異なる現象を分岐と呼び、その普遍的なパターンを研究するのが力学系の分岐理論です [4]。
解の構造全体を考えた大局的な分岐を調べるのは非常に難しいですが、ある一つの定常解に着目した局所的な分岐現象に限っても興味深いことが分かります。例えば、本質的に 1次元の力学系で生じる分岐（この場合は解の個数や安定性が変化するパターン）は genericにはサドルノード分岐と呼ばれる以下の形しか存在しません。

dx= μ ± x2dt56

この場合の人は、未感染状態、感染状態、免疫がある状態の 3 状態に理想化して考えます。
数学的には難しい問題ですが、ここでは概念としておおらかにこの言葉を用いています。

(2)

《講義ノート》(a)

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)(b)

サドルノード分岐
安定解

リミットサイクル解
標準型 (normal form)

不安定解

0

相互作用

位相で見た相互作用0

相互作用

図 2: 分岐理論およびリズム力学系の位相縮約 (a)1 次元で generic に見られるサドルノード
分岐。コントロールパラメーター μ の増加にともない下図の様にペアの安定・不安定固定点が生じます。上図は μ と定常解の関係を表示した分岐図。主として dx/dt(≡ ẋ) = μ − x2の形を説明しましたが、符号が変わっても右下の枠内のように本質的に同じです。(b) 位
相縮約の概念図。高次元空間のリミットサイクルをもつ力学系の相互作用は結合位相振動子の形に理論的に縮約できます。

generic という意味は、微分方程式 (1) の関数形 F (X, μ) をわずかに変えたとしても, パ
ラメーターや平衡解が少しずれる場合はあるが、本質的には式（2）でダイナミクスが記述できるということです。具体的には図 2(a) に示すように、なにも解がない状況から二つの安定・不安定平衡解がペアで出現する現象です。分岐理論から、高次元においても本
質的に 1 次元の分岐現象として generic に見られるものは、サドルノード分岐であり、適切な座標変化を行えば式（2）の形（標準形）に表現できるダイナミクスである、ということです。これは奇妙に思えるかもしれません。なぜなら、一般の方程式を考えている以
上、様々な解の個数や不安定性の変化のパターンがあっても良いはずで、無限とは言わな
いまでも、一つしかないとは思えないでしょう。ここで、generic という言葉に注意してください。このようなダイナミクスの振る舞いは構造安定とも言われます。パラメーターの特殊な値や関数の特別な形の場合にしか見られない現象ではない点が重要です。ただし、もし系に厳密に成立する付加的な拘束条件があれば、別のタイプの分岐現象も generic に生じます。例えば、常に原点 x = 0 に解があれば（2）に加えてトランスクリティカル分
岐、x → −x に対して式が不変（反転対称性）であればピッチフォーク分岐が generic に生じます。逆に言えば拘束条件が少しでも崩れれば、サドルノード分岐が代わりに出現す
るということです。以上の事実から、「筋の良いモデル」として、局所的にはサドルノー
ド分岐の標準形を使うことが考えられます。これにより、分岐点近傍という限定された条件ですが、特定の詳細モデルに依存しない、幅広い現象を普遍的に捕らえることが期待できます。

《講義ノート》

1.2

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)

リズムを生み出す力学系の普遍性

分岐理論以外にも何か力学系の理論で、筋の良い普遍性のある数理モデル構築のヒントになる理論はないのでしょうか？一般的な現象から普遍的に起こるダイナミクスの変化として、定常状態の不安定化が起きた後によく見られる挙動として周期運動があります7 。特に散逸力学系（エネルギーなどの保存則が無い系）で見られる周期解をリミットサイクルとよび、非線形力学の世界では非常に良く研究されています。生物で見られる多くのリズ
ムは安定リミットサイクルとして数理モデルで表現できます。実は、安定リミットサイク
ルをもつ似たような力学系8 が多数相互作用している場合、シンプルな位相振動子モデル
で表現できることが理論的に知られています [5]。その理論を少し概観してみましょう。
リズムを生み出す系も、（1）と同様にdX= F (X)dt(3)と常微分方程式の形で一般的に記述できる点は変わりませんが、安定周期解が存在することが重要です。その安定周期解（安定リミットサイクル）の周期を T とし、解を X T (t) としましょう。なお、簡便さのために自由度を M として、各変数 Xi (i = 1, · · · , M ) をまとめて X などと表記しています。例えば、ニューロンの代表的な数理モデルである Hodgkin-

Huxley モデルの場合は M = 4 です [6]。一方、脳波などを生成している神経ネットワー
クは大自由度の力学系であり、M は数万以上になると考えられます。このリズムを生成している系の状態を、M 次元空間上の一点 X で表現しましょう。この抽象的な空間は、１点を指定すると系の状態が完全に指定できるので、状態空間と呼ばれます9 。X T (t) が状
態空間の中をダイナミクスの時間発展と共に動いていく状況では、1 周期ごとに同じ点に戻ってきます。すなわち、周期解 X T (t) の軌道を描くと、図 2(b) の左のように時間の経過とともに系の状態は一つの閉曲線上を周期 T で周回します。ここで、適当な位置を起点として一定の時間間隔で軌道上に印を刻み、それを位相 φ と名付けましょう。位相 φは起点で 0 から始まり、丁度 2π の値で再び起点に戻ってくるように定義します。すなわ
ち、位相 φ は dφ/dt = 2π/T ≡ ω を満たし、この座標で見れば、周期解は φ(t) = ωt とシンプルに表現されます。
実際の生命システムや工学システムで見られるリズム現象では、しばしば複数のリズムが相互作用していると見なせる場合があります。そのようなリズム間相互作用がある一般

7

Hopf 分岐として知られていますが、今から述べる話は分岐点近傍だけで成立するのではありません。
似てるとは少し異なるリズムの集団を考えており、そのような周期の異なるリミットサイクル集団で引き込み転移が起きることが理論的に解析できることが知られています。非平衡・非線形で理論的に解けることは珍しく、平衡系でのスピン系のような立ち位置に対応する非線形の数理モデルと言えます。
9
位相空間 (phase space) ともいいます。
8

《講義ノート》

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)

的な状況を想定して、次のような力学系モデルを考えましょう。

dX i= F (X i ) + δF i (X i ) +V ij (X i , X j )
dt

(4)

j

ここで、各リズムを生み出すシステム i のダイナミクスは基本的にはほとんど同様の形 F (X)をしており、わずかな違いを δF i (X i ) で表しました。相互作用関数 V ij (X i , X j ) は、システ
ム i が j から受ける影響を表しています。なお、共通部分のダイナミクス dX i /dt = F (X i )では安定周期解 X T (t) をもつとし、わずかな差異 δF i (X i ) はリミットサイクル解 X T (t)の形をわずかに変化させるだけであると仮定します。同時に、相互作用は比較的弱く、互いにリミットサイクル軌道の形が微少に変化するかも知れないが、基本的には各システムの状態は元のリミットサイクル解近傍に沿って運動している状況を考えます。一般には各システムの状態は M 次元状態空間の座標 X i で指定するべきですが、上記状況では位相φi という１自由度により近似的に記述できることが期待できます。言い換えると、リミッ
トサイクル解近傍では状態空間の各点に位相を一意に割り当てることが可能であり、位相場 φi (X i ) が定義できるという意味です。それに基づけば、多変数の微分方程式の結合力
学系の代わりに、それぞれのシステムの状態を位相 φi という一つの変数で表現し、位相変数 φi (i = 1, · · · , M ) のダイナミクスとして閉じた形で表せることが期待できます。そんな都合の良いことが可能だろうかと思うかも知れませんが、位相縮約という理論的枠組みを使うと可能な事が知られています。すなわち、複数のリズムが相互作用している力学系(4) を位相縮約した方程式で記述すると
dφ= ωi +
Γij (φi − φj )
dt

(5)

j

と表せ、定数 ωi および関数 Γij は元の力学系 (4) から具体的に計算可能です。定数 ωi の i
依存性は力学系の差異 δF i (X i ) 由来の項で、単なるスカラー量（位相に依存しない）になることに注意してください。位相で見た相互作用関数 Γij (φ) は V ij (X i , X j ) 由来であり、2π 周期関数であることが示されます。詳細は専門書を参考にして頂くとして [5]、ここでのポイントは、安定リミットサイクル解により表現されるリズム現象に限定すると、そのダイナミクスは各リズムの状態を記述する位相 φi を用いて本質的に (5) の力学系モデルで記述出来るという点です。

2

変化するネットワーク結合力学系
現実の興味深い系は、しばしば動的ネットワーク上で相互作用する結合力学系として見なすことができ、そのネットワーク構造とノード上の状態が相互に影響しあう事で、全体として変化・発展していきます。これが (iii) の理解の困難の要因になっていると言いまし《講義ノート》

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)

た。実際、極めてシンプルなダイナミクスを仮定した場合でも、動的なネットワーク構造を導入することで、想像以上に複雑な挙動が生じることを次に示したいと思います。

2.1

拡散散逸の動的ネットワーク結合力学系

一番シンプルなダイナミクスとは何でしょうか？まず１自由度の力学系 dx/dt = F (x) のタイプを考えましょう。F (x) が 0 や定数10 など自明な例を除くと、線形微分方程式のもっとも基本的な形として F (x) = −γx が考えられます。これは「資源」x の一定割合を常に消費して減衰する散逸過程であり、解は x(t) = exp(−γt) となり、時間とともに指数関数
的に x は 0 に減衰します。放射性物質の崩壊過程やバクテリアの初期の増殖過程（γ の符号は異なりますが）など、自然現象を説明する意味でも普遍的に存在する基本的ダイナミクスです。そこでネットワークの個々のノード上のダイナミクスはこれを仮定しましょう。
一方、相互作用としては拡散過程は基本的な相互作用のダイナミクスの一つです11 。両方に資源がある場合に、その間のリンクを通じて拡散する、リンクの太さ wij が拡散の大きさに対応するような以下のモデルを考えて見ましょう。

xi (t + 1) − xi (t) = −γxi + D × (diﬀusion process via weighted links)(6)

ここで D は、拡散の全体的な大きさをコントロールする正のパラメーターです。特に面白
くないダイナミクスと思うかも知れませんが、事実、拡散係数がリンクごとにばらばらな値だったとしても、実は簡単に解けて平衡状態に緩和することがすぐ分かります。
この系に、リンクの太さ wij が変化するダイナミクスを導入するとどうなるのでしょうか？リンクのダイナミクスとして基本的なものは何か？すぐには指針が浮かばないかも知れません。そのヒントを得るため、少し強引ですが、xi を都市 i の人口に対応付けて、2つの都市間の道路の太さを考えましょう。道路の太さはネットワークのリンクの太さ wijに対応します。そうすれば、インフラ整備の観点から、その道路の交通量に応じて道路の太さを拡充・縮小するのが自然でしょう。基本的にはランダムに人が往き来すると仮定す
ると、交通量はリンク両端の都市の人口の積に比例すると考えるのが妥当です（化学反応の質量作用の法則）。そのように考えると、リンクに付随する拡散係数 wij と両端の資源の量の関係は wij ∝ xi xj に近づくように変化させるのが妥当です。式で書くと次のようになります。

wij (t + 1) − wij (t) =  [xi (t)xj (t) − wij (t)](7)

ここで、 は十分小さく、xi よりもゆっくり変化すると仮定します (図 3 上参照）。
10

位相振動子系の単一のダイナミクスは定数が振動数に対応するものと見なせますが、次に取り扱います。
相互作用とよぶには、受動的なプロセスに思えるかも知れませんが、それでも多様な振る舞いが見られます。
11

《講義ノート》

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)

0.10

0.0

PDF

3

0.08

2
1

-0.5

0
0

0.4 0.8 1.2 1.6

0.06

2

i
0.04

1

10

100

-1.0

-1

10

-1.99

x

-2

10

0.02

-3

10

Resource on node

-4

10

-5

10

-6

10

0

-7

10

-1

10

0

10

1

10

2

10

Resource: x

3

10

-1.5
0

0.1

0.2

0.3

0.4

Exponent of the cumulative
distribution of resources

Decay const. : κ4

0.5

Diffusion const.: D

図 3: 資源の拡散と散逸の動的ネットワーク上の力学系。結合の太さ（拡散の強さを表す）が時間的に変化すると、リンクを通じた拡散が支配的な場合はスケールフリー構造が出現
し、ノード上の散逸が支配的だと通常のランダムネットワーク構造になります。右下の相図は拡散と散逸のパラメーターに依存してスケールフリーのべき分布が出現する領域に色をつけてあります。色はべきの指数を表します。
以上のような、資源が散逸・拡散する動的ネットワーク上のダイナミクスを調べてみると、予想もしない興味深い振る舞いを示します [7]。結論から言えば、図 3 左下にしめされたように拡散係数と散逸のバランスを変えると、大別して二つのタイプの挙動が見られます。拡散より散逸が相対的に大きい、すなわちその場で資源が消費されるプロセスでは、最終的にノードの資源分布はガウス分布的になり、ネットワークはランダムネットワークのようなものが形成されます。その逆の拡散すなわち資源の移動が支配的な状況では、ノードの資源分布はべき分布となり、ネットワークもスケールフリーネットワークの構造になります。大胆に言えば、地産地消は皆同じような資源の分布になり、自由貿易は格差の差が激しくなる、と言えるかも知れません。また、図には示していませんが、ある時刻に資源が最も集まったノードでも時間とともに衰退し、栄枯盛衰が繰り返され、定常状態になりません（詳しく見るとリヤプノフ指数が正のカオス状態であることが分かります）。このようにネットワークを動的に変化させただけで（変化規則はシンプルなものです）、このようなスケールフリー構造の獲得やカオスが出現するのは、すこし驚かないでしょうか？

2.2

リズムの動的ネットワーク結合力学系
抽象化したネットワークの力学系を考える上で、他に汎用性の高いダイナミクスの候補はないでしょうか？ その一つは、物理現象だけでなく生命・社会現象に幅広く顕れるリズム現象です。概日リズムや歩行運動などは、一定の振幅とリズムを刻み外乱に対しても

《講義ノート》

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)

リミットサイクル振動子ネットワーク

２クラスタ状態
Hebbian-like−π

位相振動子の発展方程式

自己組織化

結合構造が動的に変化

0

π

2 値情報
を保持

コヒーレント状態
STDP-like
−π

0

因果関係
を保持

π

結合重みの発展方程式
カオス状態
Anti Hebbian-like

結合と状態が共に変化する場合

−π

0

π

非定常に
状態変化

典型的に 3 つの状態が構造安定に出現！

図 4: 動的ネットワーク上の振動子集団の振る舞い。結合の変化ルールを変えると 3 タイプの挙動が観測されます。
ロバストであることに特徴があり、力学系ではリミットサイクル振動解に対応します。既に述べたように、背後の現象が複雑な大自由度非線形力学系でも、このようなリミットサ
イクルが相互作用している系は、位相振動子モデルで理論的に記述可能である事が示されています 。逆に言えば、位相振動子のネットワークを考えることで、幅広いリズム現象をカバーしつつ本質的な挙動を捕らえることが可能だということです。例えば、局所的な神経回路で生成されたリズムが相互作用する EEG （脳波 Electroencephalography）などのダイナミクスを解明するヒントになる可能性があります。そこで位相振動子を動的ネット
ワーク上で考えてみることは非常に有用かも知れません。しかし、相互作用のダイナミク
スをどう考えるかは悩ましいところです。位相の並進対称性の条件から、相互作用の時間変化は、両端の振動子の位相差の関数として仮定しましょう。その中で、神経系の知見も活用しながら相互作用の形を幅広くパラメーターで変化させて調べてみることが考えられます。我々はその指針に基づいて、位相振動子のネットワーク結合力学系で、相互作用の強さに様々な可塑性のルールを仮定して、その力学的挙動を調べてみました [8, 9]。 結果を端的にまとめると、A. 同 (反) 位相だと結合が強く (弱く) なる Hebb 的学習ルールの場
合は、最終状態は初期条件に依存して 2 値情報を保持する最終状態へ収束 B. 相手より位相が進んで (遅れて) いれば結合を強め (弱め) るルールでは、 順序 (因果関係) を保持する
最終状態へ収束 C. 反 Hebb 的学習ルールの場合は最終的にカオス状態へ収束、の 3 タイプが現れることが示されました（図 4）。事実、系の詳細に依らず、神経ネットワークでは

Hebb 則により発火・非発火状態で符号化される連想記憶モデルが実現できることが知られており、一方、最近発見された B に相当するシナプス可塑性 Spike-Timing Dependent

Plasticity では、発火の順序 (因果) 関係が学習可能です。C のカオス状態は今のところ現実世界で対応するものは見つかっていません。ここでは全相互作用が同一の変化ルールに従うと仮定しましたが、A や B と C も組み合わせてハイブリッドにすることで、外部の《講義ノート》

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)

入力に応じて構造を柔軟に変化させる事が可能かもしれません。この研究では素子や結合のダイナミクスも同一と仮定していました。神経系などを考えると、神経細胞も多種類存在し、シナプス結合の学習則も均一ではありません。そのような非一様性を取り入れて、
例えば階層性の動的な出現などが見られるのかを調べることは今後の課題です。

3

時系列データからの力学系の結合関数推定
これまでの話で、散逸力学系のリズム（安定リミットサイクル解）が相互作用している

力学系に関しては, 対象の系の詳細に依らずに位相振動子ネットワーク (5) でダイナミク
スを本質的に記述出来ることを繰り返し述べました。その際には、個別の系の詳細を記述したダイナミクス (4) の特性は自然振動数 ωi や相互作用関数 Γij (φ) などに反映されます。
それでは現実の系からリズムの相互作用すなわち因果関係を調べるために具体的な相互作用関数などを同定するためにはどうすればよいのでしょうか？これまでの研究では、各種の実験データから、まず個別の系の詳細なモデル (4) を構築し、それが正しいという前提条件下で位相縮約を適用してダイナミクス (5) を導出する方法が標準的でした。しかし、これには問題があります。もし詳細モデル (4) の正確さが十分保証できず、何らかの要素を見逃していたら、位相縮約したダイナミクス自体も信頼出来ないものとなります。例えば、神経系のモデルでは、関連するイオンチャネルすべての影響を正確にモデルに反映す
るには、膨大な電気生理学的実験が必要です。また、脳波などの場合はもっと深刻で、背後にある神経回路網の正確なモデル化は非常に困難であることは理解できると思います。
ここで、最初の枠組みの位相振動子モデル (5) に立ち戻って考えて見ましょう。計測するデータ自体はリズムを内在する時系列データ X です。神経ネットワークの場合は、各ニューロンの膜電位や細胞内イオン濃度などに対応する時系列データになるでしょう。これを直接精緻にモデル化したダイナミクス (4) は大変複雑かも知れないが、説明すべきことをリズム現象に限れば、適切な状況下において位相振動子ネットワーク (5) で記述出来ることを我々は既に知っています。そこで、

dφi= ωi +
Γij (φi − φj ) + ξi (t)
dt

(8)

j

という形を最初から仮定して、観測した時系列データを説明できるように定数 ωi および関
数 Γij (φ) を詳細モデルを介さない形で決定する方法が考えられます (図 5)。ここで、実用上の観点から、様々な要因によりモデルでは説明出来ない部分 (観測誤差等）を表すノイズ項 ξi (t) を導入します。リズム現象に限れば、実験でリズムデータを計測し、縮約した力
学系 (8) をそこから構築すれば、詳細モデルでは取り込めなかった未知の効果も織り込み済みの数理モデルを構築できます。また、リズム現象に限定したことにより、推定するモ
デルに力学系の知見に基づく拘束条件が課され、逆に少ないデータでも正しく力学系を構《講義ノート》

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)

実データからリズムの相互作用ネットワーク抽出㸽᫬⣔ิࢆ෌⌧ࡍࡿ┦஫స⏝᥎ᐃ

図 5: 力学系理論に立脚したデータ駆動型の位相振動子ネットワークの推定
築できることが期待出来るというわけです。やるべきことは、多変量時系列データ Xi (t)から力学系 (8) の ωi および関数 Γij 等を決定することです。
まず最初にやらなくてはいけないことは時系列データ X i から φi に変換することです。
実用上避けて通れないのが、一般に X は多自由度の変数ですが、現実にはその中の１変数しか計測できない場合が多い点です。例えば、神経細胞では膜電位しか計測していない
実験がほとんどです。その１変数を手掛かりに位相 φ を求めなくてはならず、一種の不良設定問題となります。これに対して、既に色々な手法が開発されており、例えばヒルベルト変換を用いて” 暫定位相” を求める方法や、もし各リズム力学系から２変数を計測可能
なら、２次元プロットから適当な原点を定め、その原点を中心とした極座標から 暫定位
相 を求める素朴な方法もあります。この” 暫定位相” は、とりあえず基本となる周期軌道
上に座標を刻んだもので、必ずしも時間に対して一定の刻み幅になるような座標（位相縮約の位相）になっていません。そこで、平均的な”暫定位相”の分布を用いて補正すること
で、相互作用や揺らぎによる位相の変化のみを取り出すことを行い、位相縮約の位相 φi の
時系列を統計的に得ることが可能になります。まとめると、各リズムを生成する系から少なくとも１つの変数を計測して、そこから既存の様々な手法の１つを用いて 暫定位相を求め、非結合時に一定の速度で増加する位相（位相縮約の位相）に変換する手順になります [10, 11, 12]。
さて、位相の時系列が得られた後、ωi および関数 Γij 等を決定するため、ベイズ統計の枠組みを用います。手法の骨子を説明するため、自然振動数 ωi を推定する状況を想定しましょう。ベイズ統計を用いる特徴として、推定すべきパラメーターを確率変数と見なすことがあげられます。すなわち、自然振動数 ωi はある確率分布 P (ωi ) にしたがう確率変数と仮定します。位相の時系列データ φi (tk )(k = 1, · · · , ND 、ND はデータ点の数) が得られたとし、この位相時系列データをまとめて Φ と表します。一方、ωi が与えられた条件下での位相時系列データ Φ は、ホワイトノイズを含むランジュバン方程式 (8) に従うの

《講義ノート》

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)

(a)

(b)

20

30

40

(c)

(d)

20 周期

2
0
-2

1.01

2
0
-2

0.99

50

t

160 周期1000 周期1000 周期

1

0.98

0

π
φ2-φ1

2π

0.99
0.98

1000 周期160 周期異なるタイプ

0.97の変数を計測

0.96

20 周期

0.95

0

真値

π
φ1-φ2

2π

図 6: ２つの van der Pol 振動子結合系への適用例 (a) 各 yi を計測した時系列データのグ
ラフ。(b) 推定した相互作用関数: 縦軸は ωi + j Γij (φi − φj ) であり、データ数と共に正しい関数に収束する様子がわかります。実線：位相縮約の理論、破線：推定結果、灰色：
関数の不確定性の範囲。(c) 計測データ数が増加するにつれて ω2 の確率分布関数が正しい値に収束する様子。縦の破線は正しい値を示します。(d) 異なる種類の観測量 (x1 と y2 ) を計測しても正しく推定可能です。
で、ある確率 P (Φ|ωi ) で実現すると考えられます12 。条件付き確率の関係式（ベイズの定理）から、次の式が導出できます。

P (ωi |Φ) =P (Φ|ωi )P (ωi )
P (Φ)

(9)

左辺の P (ωi |Φ) は、時系列データ Φ を知った条件下でのパラメーター ωi の分布を意味しており、データを知ったことにより ωi の分布 P (ωi ) が変化して P (ωi |Φ) となります。デー
タ取得後の P (ωi |Φ) は、正しく推定出来ていれば、ピークの位置が真の値に近づくでしょう。同時に、そのピークもより鋭くなり、推定の確信度が上がっていることも示されるはずです。実際は ωi 以外にも全てのパラメーターについて同様の計算を行います。その際、関数 Γij (φ) は原理上は無限の自由度をもちますが、シンプルには、周期関数である性質を利用し適当な次数までのフーリエ展開した有限個の係数を推定する問題に帰着できます。
どの次数まで採用するかに関しては、モデルの複雑さと実際のデータを説明する能力の高さの兼ね合いを評価して決定します（エビデンス近似と呼ばれるモデル選択の手法）。詳
しくは、ベイズ推定に関しては様々な良書があるのでそれを参照して下さい [13]。また、フーリエ展開した関数は暴れることがしばしばあるのですが、最近我々はガウス過程回帰（ランダムな周期関数のサンプルから推定する手法）を用いた方法を開発しました。
具体的な適用例として、図 6 のように非対称に相互作用している２つの van der Pol 方
12

実際には具体的な形を計算することが必要です。

《講義ノート》

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)

程式ẋ1 = y1 + K(x2 − x1 ) + ξx,1 (t),ẏ1 = 1 (1 − x21 )y1 − x1 + Kx22 y2 + ξy,1 (t),ẋ2 = y2 − Kx21 y1 + ξx,2 (t),ẏ2 = 2 (1 − x22 )y2 − x2 + Kx1 y12 + ξy,2 (t),を考えましょう。まず、各振動子から (a) のように y1 ,y2 の時系列を計測している状況を考えます13 。(b) に本手法で推定したモデルの結果を示します。計測時間が増えるにつれ正しい関数形に漸近する様子が見て取れます。また、(c) ではデータ量に応じて P (ωi |Φ) の具体的形がどう変化するかを示しており、真値をピークとするシャープな分布に収束している様子がわかります。更に (d) では、x1 ,y2 の異なる変数の時系列を計測しても上手く推定できる事が示されています14 。ここで強調したい点は、各振動子の力学変数の１つが観測出来れば良く、違うタイプの変数でも縮約した位相振動子系を推定可能である点です。これは実用上応用範囲が広がる大きなメリットです。例えば、歩行運動や音声信号と脳活動、心臓や腸内細菌のリズム等の異なった種類の集団リズムの関係性も解析可能であり、新規性のある切り口で実験データの解釈に新たな視点を与えることが期待できます。特に、リズム間の相互作用が一方向的か双方向的か（リズム間の因果関係）も力学系として明確な意味で推定可能です。更に、長時間計測することで、マクロなレベルでネットワークの相互作用の変化を捕らえることで、普遍的な動的ネットワーク上の結合力学系構築のヒントも得られるかも知れません。

4

おわりに
動的ネットワーク上の結合力学系の研究は、物理現象以外の領域では特に重要なトピッ

クになると考えています。社会や生命のシステムを理解するためには、まず最初に「筋の良い」動的ネットワーク上の結合力学系の数理モデルを構築し、その知見を蓄えておくことが大事だと思います。一方で、社会や生命のシステムを精緻にモデル化するのは困難ですが、膨大なデータが取得可能な時代でもあります。その膨大な時系列データから、背後にある動的ネットワーク上の結合力学系の数理モデルを統計的に同定することも重要になると思います。皆さんが、将来そのような関係ある分野の研究者となったとき、少しでもヒントになれば幸いです。
13

リミットサイクル振動解が存在するパラメーター 1 = 0.3, 2 = 0.7, K = 0.01、σ = 0.03 を用いました。
この場合、互いの位相の不確定性は原理的に残ります。何故なら x1 と y2 の対応関係は計測データだけから知ることは出来ないからです。結果では、その不定性は補正して表示しています。なお、x1 ,x2 の場合も当然推定は成功します。
14

《講義ノート》

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)

謝辞
共同研究者の青木高明氏（香川大）および太田絵一郎氏（サイボウズ）に感謝します。

参考文献[1] R. Pastor-Satorras, C. Castellano, P. Van Mieghem, and A. Vespignani, “Epidemic
processes in complex networks,” Reviews of modern physics, vol. 87, no. 3, pp. 925–
979, 2015.
[2] R. Pastor-Satorras and A. Vespignani, “Epidemic Spreading in Scale-Free Networks,” Physical Review Letters, vol. 86, no. 14, pp. 3200–3203, 2001.
[3] A. Barrat, Dynamical Processes on Complex Networks. Cambridge: Cambridge
University Press, 2012.
[4] J. Guckenheimer and P. J. Holmes, Nonlinear Oscillations, Dynamical Systems, andBifurcations of Vector Fields. Applied Mathematical Sciences, New York: SpringerVerlag, 1983.
[5] Y. Kuramoto, Chemical Oscillations, Waves, and Turbulence. Springer Series in
Synergetics, Berlin Heidelberg: Springer-Verlag, 1984.
[6] A. L. Hodgkin and A. F. Huxley, “A quantitative description of membrane current and its application to conduction and excitation in nerve,” The Journal of
Physiology, vol. 117, pp. 500–544, 1952.
[7] T. Aoki and T. Aoyagi, “Scale-Free Structures Emerging from Co-evolution of a
Network and the Distribution of a Diﬀusive Resource on it,” Physical Review Letters, vol. 109, no. 20, p. 208702, 2012.
[8] T. Aoki and T. Aoyagi, “Co-evolution of phases and connection strengths in a
network of phase oscillators,” Physical Review Letters, vol. 102, no. 3, p. 034101,
2009.
[9] T. Aoki and T. Aoyagi, “Self-organized network of phase oscillators coupled by
activity-dependent interactions,” Physical Review E, vol. 84, no. 6, p. 10464, 2011.
[10] K. Ota and T. Aoyagi, “Direct extraction of phase dynamics from ﬂuctuating rhythmic data based on a Bayesian approach,” arXiv:1405.4126 [nlin, q-bio], 2014.

《講義ノート》

物性研究・電子版 Vol. 10, No. 1, 101213(2022年3月号)[11] K. Ota, I. Aihara, and T. Aoyagi, “Interaction mechanisms quantiﬁed from dynamical features of frog choruses,” Royal Society Open Science, vol. 7, p. 191693,
2020.
[12] T. Onojima, T. Goto, H. Mizuhara, and T. Aoyagi, “A dynamical systems approachfor estimating phase interactions between rhythms of diﬀerent frequencies fromexperimental data,” PLoS Computational Biology, vol. 14, no. 1, p. e1005928, 2018.
[13] C. M. Bishop, Pattern Recognition and Machine Learning. New York: Springer, 1st
ed. 2006. corr. 2nd printing 2011 版 ed., 2006.