\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

興奮・抑制均衡入門

著者:新垣貴史
日付:2015 年 06 月 01 日
バージョン:0.1.1a0

注釈

未完ですが, 気合があれば読めるかもしれません.

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

まえがき

理工学系の学部三年生くらいから読める, 理論神経科学のトピックの紹介を目指して書いています. 生物・神経科学や大学院レベルの数学・物理の知識は前提としない予定です. 今の所, [vanVreeswijk1998] の解説+αです.

いくつかのバージョンで読めます:

HTML版(MathJax):

http://balance-tutorial-ja.readthedocs.org

常に最新版が読めます. ページの描画に少し時間がかかるかもしれません.

HTML版(pngmath):

http://balance-tutorial-ja.readthedocs.org/ja/pngmath/

HTML版(MathJax)より更新頻度は低いです. 数式は若干荒く見えますが, ページは素早く描写されます. Zipファイルでダウンロードしてオフラインで読むことも可能です.

PDF版:

http://tkf.github.io/balance-tutorial-ja/balance-tutorial-ja.pdf

HTML版より更新頻度は低いです.

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

はじめに

生体内 (in vivo) の大脳皮質のニューロンは典型的に, (1) 発火の時系列パターンが不規則であり (発火時間間隔 (interspike interval; ISI) の分布が指数分布に近い), (2) ニューロン同士はほぼ無相関に発火するという特徴を持つ. 一方, 生体外 (in vitro) の電気生理記録から, 大脳皮質のニューロンは決定論的にふるまうことが知られている. それでも, 大脳皮質のニューロンは相互に結合しているのだから, 入力として他のニューロンの発火の不規則なパターンを受けることにより不規則な発火パターンを出力するというメカニズムが考えられる. しかし, この素朴に考えれば「自然な」メカニズムは, 簡単な考察から破綻することが分かる.

N 個のニューロンが相互に結合した系を考える. それぞれのニューロンは平均して K 個のニューロンから入力を受けているとし, 行列 \bm J を用いて, i 番目のニューロンが j 番目のニューロンからの入力を受けているなら J_{ij} = J / K そうでなければ J_{ij} = 0 として結合を表現する. ここで結合力を J/K とスケールした理由は以下で明らかになる. 以下の議論は興奮性の集団 (J > 0) でも抑制性の集団 (J < 0) でも成り立つ. ニューロンの数やその結合の数は非常に大きいので, N, K \to \infty の極限で現れる普遍的な性質を調べる. ニューロン i への入力を x_i, その出力を y_i と書く. 入力 x_i がどのように出力 y_i に変換されるかは指定しないが, 何らかの決定論的な仕組みに依る, つまり出力の不規則性は入力の不規則性のみに依るという仮定をする. 簡単のために, すべてのニューロンの出力 y_1, \ldots, y_N は独立で, 同じ平均 m = \Avg{y_i} (「ネットワークの平均発火率」のようなもの) を持つとする. さらに, 入力 x_i が線形和

x_i = \sum_{j=1}^{N} J_{ij} y_j

で表されるとする. 右辺はたくさんの y_i の和なのだから, 大数の法則 (law of large numbers) により「その平均 (= m) × J_{ij} が非ゼロの時の値 (= J / K) × J_{ij} が非ゼロの個数 (\approx K)」で近似出来るはずである. つまり, x_i \approx J m となる. 結合力が J/K という選択は x_i が K \to \infty で発散しないために必要だったのである. しかし, この場合, 中心極限定理 (central limit theorem) により x_i の分散は \Var (x_i) = O(1/K) とスケールすることが期待される. つまり, x_i の不規則性は K \to \infty の極限で消えてしまう. よって, 出力 y_i の不規則性は入力 x_i の不規則性にのみに依るという仮定から, 出力 y_i の不規則性も消えてしまう. しかし, この結論は出力 y_i に不規則性があるという仮定のもとで導いたのだから, 矛盾している. つまり, このモデルは間違っている.

この現象を避けるには, 結合力を J/K から J / \sqrt K に変更すれば良い. この場合, 同じく 中心極限定理 (central limit theorem) により \Var (x_i) = O(1) となる. つまり, K \to \infty の極限でも不規則性が消えない. ところが, 平均は \Avg{x_i} = O(\sqrt K) となり, K \to \infty の極限で発散してしまう. 生物学的なモデルでは変数は常に一定の範囲内に収まっているべきであり, 平均が発散するモデルは使えない.

しかし, この困難は相互に結合した興奮性と抑制性の集団と, それらの集団への外部からの興奮性入力 [1] を考えると解決出来る. 興奮性の集団を上付添字 E, 抑制性の集団を上付添字 I, 外部の (興奮性) 集団を上付添字 0 で表し, 例えば i 番目の興奮性ニューロンへの入力は x^E_i で, 抑制性ニューロンから興奮性ニューロンへの結合の強さは J^{EI} / \sqrt K とする. 集団 E, 0 は興奮性なので J^{kE}, J^{k0} > 0 (k = E, I), 集団 I は抑制性なので J^{kI} < 0 (k = E, I) であることに注意. 集団 k (= E, I) への入力は

x^k_i = \sum_{j=1}^{N} J^{kE}_{ij} y^{E}_j + \sum_{j=1}^{N} J^{kI}_{ij} y^{I}_j + \sum_{j=1}^{N} J^{k0}_{ij} y^{0}_j \approx \sqrt{K} (J^{kE} m^E + J^{kI} m^I + J^{k0} m^0)

となる. ここで J^{kE} m^E + J^{k0} m^0 は正で J^{kI} m^I は負である. よって, これらの値がうまく釣り合えば, つまり, J^{kE} m^E - J^{kI} m^I + J^{k0} m^0 \approx 0 が成り立てば, 入力 x^k_i は発散しない. 本稿ではこの条件の成り立つネットワークを 興奮・抑制均衡ネットワーク と呼び, その特性を解析する. 素朴に考えるとそのためには各結合力 J^{kl} (k, l = E, I, 0) を注意深く選ぶ (fine tuning する) 必要があり, 恣意的なモデルに見える. よって, 興奮・抑制均衡ネットワークの重要なアイディアは, この条件が 自動的に [2] 達成されることにある. もちろんそのためには結合力 J^{kl} はある条件を満たす必要があるが, その条件は大変ゆるい.

[1]ここでの議論だけでは外部の集団を考える必然性は自明では無い. しかし, 興奮・抑制均衡ネットワークは外部の集団なしでは成り立たないことを後に述べる.
[2]力学系の言葉で言えば, この条件が達成された状態が漸近安定.

まず, 興奮・抑制均衡ネットワーク でネットワークの詳細には依らない特性を導く. そして, 二状態ニューロンのネットワーク で上記のアイディアが統計力学的な厳密な系で成立することを示す. 興奮・抑制均衡ネットワーク と 二状態ニューロンのネットワーク の大部分の議論はゆるやかにしか繋がっていないので, どちらから読み始めても良い. 二状態ニューロンのネットワーク で 興奮・抑制均衡ネットワーク の議論が必要になるのは, 安定性解析 の部分だけである.

課題

文献を紹介する

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

興奮・抑制均衡ネットワーク

興奮・抑制均衡ネットワークは興奮性と抑制性の入力が均衡した, 均衡固定点を持つ. はじめての興奮・抑制均衡 では単純化な興奮・抑制均衡ネットワークの場合についてこれを導き, 強フィードバック系 では, 興奮・抑制ネットワークよりさらに一般的な力学系についてこれを導く. 均衡固定点が存在するために必要な条件を 均衡条件不等式 で示し, 安定化可能性 でこの固定点は抑制性集団の時定数が十分小さければ安定であることを示す.

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

はじめての興奮・抑制均衡

はじめに の簡単な議論から, k = E, I について興奮性と抑制性の入力が均衡した状態

J_{kE} m_E - J_{kI} m_I + J_{k0} m_0 \approx 0

が成り立てば, ニューロンの活動の不規則性が説明出来ることが示唆された. ただし, 集団を表す添字を上付添字から下付添字に変更した. [1] [2]

[1]これ以降の議論では, 集団を表す添字は下付きで, ニューロンを表す添字を上付にする (例: J_{EI}^{ij}).
[2]記号 \approx の意味がきちんと定義されていないので条件の意味を把握するのは難しい. 実は, この節に現れる X \approx Y の形の式 (Y = 0 でも良い) を, すべて X = Y + O(1/\sqrt K) に書き換えれば厳密な議論となる. ただし, O(1/\sqrt K) は ランダウの記号 あるいは O-記法 (Big O notation) と呼ばれる漸近関係を表す記号のひとつである (漸近解析 (Asymptotic analysis) を参照).

外部入力の項を h_k = J_{k0} m_0 と書き,

\bm J = \begin{pmatrix} J_{EE} & J_{EI} \\ J_{IE} & J_{II} \\ \end{pmatrix}, \bm m = \begin{pmatrix} m_{E} \\ m_{I} \\ \end{pmatrix}, \bm h = \begin{pmatrix} h_{E} \\ h_{I} \\ \end{pmatrix}

なる行列・ベクトル表現を使えばこの式は

\bm J \bm m \approx - \bm h

と書ける. 行列 \bm J が可逆であれば,

(1)\bm m \approx \bm J^{-1} \bm h

である. つまり, ネットワークのパラメタ \bm J と外部からの入力 \bm h が決まれば, このネットワークの「グローバルな状態」 \bm m はただひとつに決まる. 興奮性と抑制性の入力が常に均衡していればこの状態は時間変化しないはずであるから, これは \bm m のダイナミクスの固定点である. よって式 (1) を 均衡固定点 (balanced fixed point) と呼ぶ. この式から \bm h = 0 なら \bm m \approx 0 なので, 外部入力がなければ興奮・抑制均衡ネットワークは活動出来ない (非ゼロの均衡固定点が存在しない) ことが分かる.

一般論へ

上記の議論を正確におこなうには, \bm m のダイナミクスが定義されていなくてはならない. しかも, それがひとつひとつは不規則なニューロンのダイナミクスから導出出来る必要がある. この導出は 二状態ニューロンのネットワーク の 平均場方程式の導出 でおこなう. しかし今は, \bm m のダイナミクスが

\bm \tau \frac{\D \bm m(t)}{\D t} = - \bm m(t) + \bm f(C \, \{\bm J \bm m(t) + \bm h\})

(ただし, C = \sqrt K, \bm m = (m_{E}, m_{I}), \bm \tau = (\tau_{E}, \tau_{I})^\intercal, \bm f(\bm m) = (f_{E}(m_{E}), f_{I}(m_{I}))^\intercal である) と書けることを仮定して議論を進めよう. 関数 f_E や f_I は集団レベルでの, 入力と出力の関係を定義している. 例えばロジスティクス関数 f_E(x) = f_I(x) = 1 / (1 + \exp(-x)) を使うことを想像しても良い. 実は関数 \bm f はあるゆるい条件を満たせば何でも良いことがわかる.

「フィードバック強度」 C は, はじめに の議論では J_{kl}^{ij} のスケーリングをどう定義するかに依っていたが, この C の起原は興奮・抑制均衡ネットワークの特性を議論するためにほとんど効いてこない. そこで, 強フィードバック系 ではこの均衡固定点を別の視点から「導出」する.

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

強フィードバック系

次の p 次元力学系を考えよう.

(2)\bm \tau \frac{\D \bm x(t)}{\D t} = - \bm x(t) + \bm f(C \, \{\bm J \bm x(t) + \bm h\}; \bm x(t))

ここで, \bm \tau = \diag(\tau_1, \ldots, \tau_p) \in \mathbb R^{p \times p} は時定数を決める行列, \bm J \in \mathbb R^{p \times p} は系の状態 \bm x(t) \in \mathbb R^p の各成分の相互作用を表す行列, ベクトル \bm h \in \mathbb R^p はこの系への外部入力を表す. 正数 C はこの系のフィードバックの強さを決めるパラメタである. 関数 \bm f の第一引数 \bm z = C \, \{\bm J \bm x(t) + \bm h\} を系への「全入力」と呼ぶことにする. 上述の記号, 特にパラメタ \bm \tau, \bm J, \bm h は C に依存しない. フィードバック強度 C への依存性は全入力 \bm z を経由する影響のみである. 本稿では, 極限 C \to \infty におけるこの系の振る舞いについて興味がある場合に, 式 (2) のような系を 強フィードバック系 (strong feedback system) と呼ぶ.

この節では, 強フィードバック系が均衡固定点と呼ばれる特別な固定点を持ち, さらにその性質が系の詳細 (\bm f の定義) に依存しないことを示す. 均衡固定点のアイディア自体は簡単である (しかし, 今ここでアイディアにピンと来なくても以降の議論は読めるはずである). ここでの「均衡」とは, 固定点 \bm x^0 において, フィードバック入力 \bm J \bm x^0 と 外部入力 \bm h が打ち消し合って (= 均衡して), 「ほぼゼロ」になる, という意味である. 正しくは, これらが O(1/C) のオーダーで打ち消しあう, つまり \bm J \bm x^0 + \bm h = O(1/C) となることを言う. これより \bm x^0 = - \bm J^{-1} \bm h + O(1/C) となるので, まず固定点 \bm x^0 が外部入力について \bm h 線型であるという性質が導かれる. さらに, 固定点 \bm x^0 からの小さい摂動 \bm x - \bm x^0 = O(1/C) が, 大きな全入力の変化 \bm z - \bm z^0 = C \bm J (\bm x - \bm x^0) = O(1) を生むので, 摂動への応答速度が速いだろうということも予想できる. 今, 均衡固定点が 存在すれば どのような性質を持つかを簡単に議論したが, この均衡固定点はいつも存在するのだろうか? もちろん力学系に何らかの条件が必要だろうが, それは何だろうか? また, 均衡固定点は他の種類の固定点に比べてどう特別なのだろうか?

表記を簡潔にするために, 状態 \bm x^0 からの 「摂動」 \bm y = C \, (\bm x - \bm x^0) に対して \bm f がどう変化するかを,

(3)\bm F(\bm y) = \bm F(\bm y; \bm x^0) := \bm f(\bm J \bm y + \bm z^0; \bm y / C + \bm x^0) - \bm f(\bm z^0; \bm x^0)

で定義する. 文脈から明らかな場合は, 第二引数は省略する.

正しく C \to \infty における系の振る舞いを表記するために, \sim, O, \Theta, \Omega などの漸近関係記号を用いる (漸近解析 (Asymptotic analysis) を参照). この節では, これらの記号はすべて C \to \infty についての漸近関係を表すとする (末尾の “as C \to \infty” を省略する).

関数 \bm f が満たすべき条件

\bm x^0 を力学系 (2) の任意の状態とし, \bm z^0 = C \, \{\bm J \bm x^0 + \bm h\} と定義する. (この節では \bm x^0 が固定点の場合しか考えないが, ここでの定義においては固定点でなくても良い.)

伝達性

\bm z^0 = O(1) ならば \bm F(\Omega(1); \bm x^0) = \Omega(1).

省略せずに書けば, すべての \bm z^0 = O(1) と \bm y = \Omega(1) について,

(4)\bm f(\bm J \bm y + \bm z^0; \bm y / C + \bm x^0) - \bm f(\bm z^0; \bm x^0) = \Omega(1).

飽和性

\bm z^0 = \omega(1) ならば \bm F(O(1); \bm x^0) = o(1).

省略せずに書けば, すべての \bm z^0 = \omega(1) と \bm y = O(1) について,

(5)\bm f(\bm J \bm y + \bm z^0; \bm y / C + \bm x^0) - \bm f(\bm z^0; \bm x^0) = o(1).

関数 \bm F は陰に C への依存性をもち, 状態 \bm x^0 も C に依存しても良いことに注意.

この関数 \bm f への条件は, 後述する均衡固定点の性質を導くための必要条件ではなく, さらに広いクラスでこれらの性質は成り立つ. しかし, この条件が応用上必要な力学系では成り立たないことは稀だろう. [1] 簡単のために 飽和性条件 の成り立たない力学系を扱うことも多いが, この場合, 非均衡固定点は発散してしまい極限 C \to \infty で存在しないのでその点まわりの振る舞いなどは議論する必要も無い. 非均衡固定点まわりでの振る舞い以外で 飽和性条件 は使わないので, この条件は本質的な条件ではない. 結局, 伝達率条件 さえ成り立てば, 強フィードバック系の重要な性質はすべて保証される.

[1]これが成り立たない場合の取り扱いについては, 発火率モデル (rate model) の区分的線型関数についてのコメント参照.

支配項均衡の方法による解析

力学系 (2) の固定点は

(6)\bm x^0 = \bm f(C \, \{\bm J \bm x^0 + \bm h\}; \bm x^0)

は満たす系の状態であり, それに対応する全入力を \bm z^0 = C \, \{\bm J \bm x^0 + \bm h\} と置く. 全入力の各項がキャンセルしあって C \to \infty で発散しない場合, つまり \bm z^0 = O(1), あるいは同値の条件

(7)\bm J \bm x^0 + \bm h = O(1/C)

を満たす固定点を 均衡固定点 (balanced fixed point) と呼ぶ. これが成り立たない場合, つまり固定点での全入力が発散して \bm z^0 = \Omega(C) となる, あるいは同値の条件

(8)\bm J \bm x^0 + \bm h = \Omega(1)

を満たす固定点を 非均衡固定点 (unbalanced fixed point) と呼ぶ. [2]

[2]ここでの均衡固定点や非均衡固定点に属さないクラスの固定点, つまり 1 \ll \bm z^0 \lesssim C (同値だが別の表記では, \bm z^0 = \omega(1) かつ \bm z^0 = o(C)) となる場合も考えられる. 例えば, \bm z^0 = \Theta(C^{1/2}) や \bm z^0 = \Theta(\log C) など. この場合, \bm J \bm x^0 + \bm h = o(1) なので, 「広義の均衡固定点」と呼んでも良いかもしれない. しかし, 以下で導く小さい摂動に対する応答速度のオーダーに関する性質は, (狭義の) 均衡固定点のみに成り立つ.

固定点まわりの座標系 \bm y = C \, (\bm x - \bm x^0) で微分方程式 (2) を書き直す. 座標変換 \bm x = \bm y / C + \bm x^0 を施すと,

C^{-1} \bm \tau \frac{\D \bm y(t)}{\D t} & = - (\bm y(t) / C + \bm x^0) + \bm f(C \, \{\bm J (\bm y(t) / C + \bm x^0) + \bm h\}; \bm x) \\ & = - \bm y(t) / C + \underbrace{ \bm f(\bm J \bm y(t) + \bm z^0; \bm x) - \bm f(\bm z^0; \bm x^0) }_{= \bm F(\bm y)}

なので, これの両辺に C をかけて

(9)\bm \tau \frac{\D \bm y(t)}{\D t} & = - \bm y(t) + C \bm F(\bm y)

を得る (関数 \bm F の定義については式 (3) を参照).

支配項均衡の方法 を使ってこの系の小さな摂動 \bm y = \Theta(1) (元の座標系では \bm x - \bm x^0 = \Theta(1/C)) への応答を調べよう. もし, C \bm F(\bm y) のオーダーが O(1) より大きい (つまり C \bm F(\bm y) = \omega(1)) ならば, C \bm F(\bm y) と \D \bm y/\D t のオーダーが等しい (つまり, \D \bm y/\D t \asymp C \bm F(\bm y) となる) ことが, 式 (9) の成り立つ必要条件である. これは, 元の座標系では \D \bm x/\D t \asymp \bm F(\bm y) を意味する.

仮定 \bm y = O(1) から, \bm F(\bm y) のオーダーは

\bm F(\bm y) & = \left\{ \begin{array}{lll} \Omega(1) & (\text{if } \bm z^0 = O(1)) \\ o(1) & (\text{if } \bm z^0 = \omega(1)) \end{array} \right.

と見積もることが出来る. ここで, \bm z^0 = O(1) の場合は 伝達率条件 を, \bm z^0 = \omega(C) の場合は 飽和性条件 を用いた. \bm z^0 が素性の良い振る舞いをするなら (正確には \lim_{C \to \infty} \bm z^0 が有限値または無限大に収束すれば), \bm z^0 \neq O(1) と \bm z^0 = \omega(1) は同値なので, \bm F(\bm y) のオーダーはこのふたつしかありえない.

ゆえに, 小さな摂動 \bm y = C \, (\bm x - \bm x^0) = \Theta(1) に関して, 均衡固定点 (\bm z^0 = O(1)) 周辺での応答速度のオーダーは \D \bm y/\D t = \Omega(C) (\D \bm x/\D t = \Omega(1)) となり, 非均衡固定点など, それを満たさない固定点 (\bm z^0 = \omega(1)) 周辺での応答速度のオーダーは \D \bm y/\D t = o(C) (\D \bm x/\D t = o(1)) となる. もし固定点 \bm x^0 が安定ならば, 均衡固定点 (\bm z^0 = O(1)) への収束は速く (\D \bm x / \D t = \Omega(1) なので, C \to \infty の極限で消えない), 非均衡固定点 (\bm z^0 = \Omega(C)) への収束は遅い (\D \bm x / \D t = o(1) なので, C \to \infty の極限で限りなく小さくなる) ことが読み取れる. 以上の議論は抽象的な定義にのみ基づいているから, これは幅広いクラスの力学系の均衡固定点について成り立つ非常に強力な性質であるといえよう.

以上の議論では, 伝達率条件 と 飽和性条件 のみを仮定し, 固定点における \bm z^0 のオーダーについて (ほとんど相互排他的な) ふたつの場合にわけ, それぞれが固定点まわりの応答速度のオーダーが強フィードバック極限 C \to \infty で消えない場合と限りなく小さくなる場合に対応することを示した. さらに, 応答速度が速い (\D \bm x / \D t = \Omega(1)) 場合が実は均衡固定点 (\bm z^0 = O(1)) に対応することを示した. つまり, 固定点まわりのダイナミクスについて漸近解析をすれば, 固定点は均衡固定点と均衡固定点でないものに自然と場合分けされるのである.

均衡固定点の線型性

均衡固定点 \bm x^0 は, 式 (7) より,

\bm x^0 = - \bm J^{-1} \bm h + O(1/C)

と書ける (ただし, \bm J が可逆だと仮定した). つまり, この固定点 \bm x^0 は外部入力 \bm h に対し, O(1/C) の誤差を除けば線型の関係を持つ. もしこの固定点が安定ならば, この系は外部入力に対し, 「線型な応答」をしていることになり, しかもそれが系の詳細を決める関数 \bm f に依らないという著しい性質を持つ. 線型性から自明だが, 均衡固定点は存在すれば (極限 C \to \infty で) 一意である.

これが良く引用される「均衡固定点は線型性をもつ」という性質であるが, この「線型性」は 厳密な意味での線型性 ではない. 固定点の満たすべき式 (6) を思い出せば, - \bm J^{-1} \bm h が \bm f の 像 の 内部 にあるという条件が必要である. 神経科学での応用上 x_i は常に発火率であるとみなされるから, f_i(\bm z) \ge 0 なので, - \bm J^{-1} \bm h の成分に負の成分があってはならない. また, f_i(\bm z) には上限がある (大きい入力については飽和する) ことが応用上自然だから, 結局, 均衡固定点の存在する \bm h の領域は有界である. つまり, 「均衡固定点は線型性をもつ」は間違いで, 「均衡固定点は閾値を超えてから飽和するまでは線型」が正しい.

均衡固定点の安定性

支配項均衡の方法による解析 から, 均衡固定点まわりのダイナミクスは,

\bm \tau \frac{\D \bm y(t)}{\D t} & = - \bm y(t) + C \bm F(\bm y) \sim C \bm F(\bm y)

つまり,

\frac{\D \bm y(t)}{\D t} \sim C \bm \tau^{-1} \bm F(\bm y)

で記述されると分かる. よって, この系の安定性の解析は, 行列

\bm \tau^{-1} \frac{\partial \bm F(\bm y)}{\partial \bm y} = \bm \tau^{-1} \frac{\partial}{\partial \bm y} \bm f(\bm J \bm y + \bm z^0; \bm y / C + \bm x^0) = \bm \tau^{-1} ((D_1 \bm f) \bm J + (D_2 \bm f) C^{-1})

の固有値で決まる (- \bm y(t) の寄与が消えていることに注意). ただし,

D_1 \bm f = \left. \frac{\partial \bm f(\bm z; \bm x^0)}{\partial \bm z} \right|_{\bm z = \bm z^0}, \qquad D_2 \bm f = \left. \frac{\partial \bm f(\bm z^0; \bm x)}{\partial \bm x} \right|_{\bm x = \bm x^0}

である. もし D_2 \bm f = O(1) なら第二引数からの寄与が消え, {\partial \bm F(\bm y)}/{\partial \bm y} \sim (D_1 \bm f) \bm J となる. これは, この方向への摂動が \bm y / C = O(1/C) と小さいためである.

強フィードバック系の例

発火率モデル (rate model)

p 個の集団からなる発火率モデルのネットワークで, すべての集団の入出力関係 (input-output relationship あるいは transfer function) がシグモイド関数 g (例えば, ロジスティック関数 g(t) = 1/(1+\exp(-t))) [3] で与えられているとすれば, 関数 \bm f の i 番目 (i = 1, \ldots, p) の成分は

f_i(\bm z; \bm x) = g(z_i)

と書ける.

[3]\arctan, \tanh, Q関数 (の x 軸を反転したもの) などでも構わない.

シグモイド関数 g は z^0 = O(1) からの非ゼロの変化に対して, 必ず非ゼロの変化をうむ, つまり

g(\Theta(1) + z^0) - g(z^0) = \Theta(1)

が成り立つ. これは 伝達率条件 の十分条件である. また, ロジスティック関数の場合は極限 z \to \infty で g(z) \to 1, 極限 z \to -\infty で g(z) \to 0 だから, 飽和性条件 も成り立つ. 一般に, シグモイド関数のように入出力関係が極限 z \to \pm \infty で有限の値に収束すれば 飽和性条件 は成り立つ. 有限の値に収束しなくても 飽和性条件 が成り立つ (かなり人工的な) 例として, g = \log がある.

もし, 関数 g が

g(z) = \begin{cases} 1 & (1 < z) \\ z & (0 < z \le 1) \\ 0 & (z \le 0) \end{cases}

のような区分的線型関数の場合は, 伝達率条件 が成り立たない. 例えば, \Delta z := 1, z^0 := 3 とおくと \Delta z = O(1), z^0 = O(1) だが,

g(\Delta z + z^0) - g(z^0) = 0 \neq \Omega(1)

である. この場合, \bm x^0 が均衡固定点であるためには, 各 z^0_i について, 0 < z^0_i < 1 なる条件が必要である.

二状態ニューロンから成るネットワーク

二状態ニューロンから成るネットワーク の平均場方程式も強フィードバック系である (平均場方程式の導出 を参照). この場合は, p = 2 個の集団からなる力学系で, 状態は集団平均発火率 \bm x = (m_1, m_2)^\intercal で, 外部入力は \bm h = (J_{10} m_0, J_{20} m_0)^\intercal である. 関数 \bm f は,

f_k(\bm z; \bm x) = H(-(z_k - \theta_k) / \sqrt{\alpha_k(\bm x)}) \alpha_k(\bm x) = (J_{k1})^2 x_1 + (J_{k2})^2 x_2

で定義される. また, フィードバックの強さは C = \sqrt K で決まる. 上記のシグモイド関数の場合と同様に, 伝達率条件 と 飽和性条件 が成り立つ.

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

均衡条件不等式

ここでは, はじめての興奮・抑制均衡 と 強フィードバック系 の解析を興奮・抑制ネットワークが外部からの興奮性入力を受ける系にあてはめ, このネットワークは「発火率がゼロの固定点は持たない」という条件から, 均衡条件不等式 (11) を導く. さらに, 「発火率が飽和するような固定点を持たない」という仮定を加えれば不等式 (13), (15) などの条件も導くことが出来るが, これらは単純な条件では無いことと, 外部入力が十分弱ければつねに成立することから, ニューラルネットワークの研究ではこれらの条件は通常は無視されている. また, 「発火率が発散しない」という条件から, 不等式 (14) も導くことが出来, 原著 [vanVreeswijk1998] ではこれも含めて均衡条件不等式と呼んでいる. 無活動・均衡-固定点の非存在条件 で 均衡条件不等式 (11) を導くまでが重要で, それ以降の議論は読まなくても他の章の内容は理解できる.

強フィードバック系 の解析と同様, 強フィードバック極限 C \to \infty を考え, すべての漸近関係はすべてこの極限についてのものとするので “(as C \to \infty)” は省略する.

興奮・抑制ネットワーク

強フィードバック系 で扱った力学系が二次元 p=2 で, 系の状態 \bm x = (x_1, x_2)^\intercal = (m_1, m_2)^\intercal の第1, 第2成分はそれぞ興奮性と抑制性の集団の発火率を表すとする. つまり, つねに m_1, m_2 \ge 0 である. また, ベクトルや行列の添字は 1, 2 の代わりに E, I を適宜使い, 例えば J_{12} と J_{EI} は同義である. 強フィードバック系 と同様, 系への「全入力」を

\bm z = \begin{pmatrix} z_{E} \\ z_{I} \end{pmatrix} = C \, \{\bm J \bm x + \bm h\} = C \, \begin{pmatrix} J_{EE} \, m_E + J_{EI} \, m_I + h_E \\ J_{IE} \, m_E + J_{II} \, m_I + h_I \end{pmatrix}

で定義する. 関数 \bm f は, いわゆる入出力関係 (input-output relationship; transfer function) を表すので, k = E, I についてどんな入力 \bm z についても f_k(\bm z) \ge 0 であり, さらに f_k は強い興奮性の入力 (z_E, z_I \to \infty) で飽和し, 強い抑制性の入力 (z_E, z_I \to -\infty) で出力が無くなる, という条件を追加する. 形式的に書くと

(10)z_k = |\Omega(C)| \Leftrightarrow f_k(\bm z) \sim m^{\max} z_k = -|\Omega(C)| \Leftrightarrow f_k(\bm z) \ll 1

となる. ここで, m^{\max} はフィードバック強度 C に依存しない正の定数で, 両集団の発火率の上限を与える. 定数 m^{\max} は集団ごとに違う値を考えても良いが, m_k と J_{kl} をスケールすれば同じ値にすることができるので, この仮定は一般性を損なわない. さらに, 発火率の上限が存在しない場合 (m^{\max} = \infty) でも, 以下の議論は成り立つ. [1]

[1]片方の集団の発火率の上限が存在し, もう片方の発火率の上限が存在しない場合は扱っていないが, そんな変なモデルを考える必要は無いだろう.

興奮性の集団と抑制性のふたつの集団が興奮性の外部入力を受ける系について解析するので, パラメタには

\bm J &= \begin{pmatrix} J_{EE} & J_{EI} \\ J_{IE} & J_{II} \end{pmatrix}, \qquad (J_{EE}, J_{IE} > 0, \enskip J_{EI}, J_{II} < 0) \\ \bm h &= \begin{pmatrix} h_{E} \\ h_{I} \end{pmatrix}, \qquad (h_{E}, h_{I} > 0)

なる制限がかかる.

興奮・抑制ネットワークの均衡固定点

均衡固定点の線型性 より, 均衡固定点は

\bm x^0 = - \bm J^{-1} \bm h + O(1/C)

を満たすことが分かった. これは単なる二次の線型方程式なので, 均衡固定点は形式的に

\begin{pmatrix} m_{E} \\ m_{I} \end{pmatrix} = - \frac{1}{\det \bm J} \begin{pmatrix} J_{II} & -J_{EI} \\ -J_{IE} & J_{EE} \end{pmatrix} \begin{pmatrix} h_{E} \\ h_{I} \end{pmatrix} = \frac{1}{\det \bm J} \begin{pmatrix} - J_{II} h_{E} + J_{EI} h_{I} \\ J_{IE} h_{E} - J_{EE} h_{I} \\ \end{pmatrix}

と書ける. これ以降の数式では固定点を表す上付き添字 0 は省略する.

この均衡固定点以外に, 非均衡固定点が存在する可能性がある. 定性的に場合分けをすると, 次の5つの場合が想定できる.

固定点 定義
無活動-固定点 両方の集団の発火率がゼロ
無活動・均衡-固定点 片方の集団の発火率がゼロで, もう片方の集団の発火率がゼロでもなく飽和もしていない
飽和-固定点 両方の集団の発火率が飽和している
飽和・均衡-固定点 片方の集団の発火率が飽和していて, もう片方の集団の発火率がゼロでもなく飽和もしていない
飽和・無活動-固定点 片方の集団の発火率が飽和していて, もう片方の集団の発火率がゼロ

ただし, 発火率の上限が無い (m^{\max} = \infty) 場合, ここでの「飽和」は発火率の発散を意味する.

無活動-固定点 は, \bm h の成分が正である限り存在できない. これは, z_k|_{\bm x = \bm 0} = C h_k \neq - |\Omega(C)| より, 自明である.

条件: 均衡固定点の発火率は非負

この均衡固定点が存在するためには \bm J の行列式

\det \bm J = J_{EE} J_{II} - J_{EI} J_{IE}

が非ゼロであることが必要条件である. さらに, 発火率が非負であるという条件を考慮すると, \bm J と \bm h は

\det \bm J > 0, \quad - J_{II} h_{E} + J_{EI} h_{I} > 0, \quad J_{IE} h_{E} - J_{EE} h_{I} > 0

または, 上記の不等号をすべてひっくり返した

\det \bm J < 0, \quad - J_{II} h_{E} + J_{EI} h_{I} < 0, \quad J_{IE} h_{E} - J_{EE} h_{I} < 0

なる不等式系を満たす必要があることが分かる. これを変形すると,

& \det \bm J = J_{EE} J_{II} - J_{EI} J_{IE} \lessgtr 0, \quad - J_{II} h_{E} + J_{EI} h_{I} \lessgtr 0, \quad J_{IE} h_{E} - J_{EE} h_{I} \lessgtr 0 \\ \Leftrightarrow & J_{EE} J_{II} \lessgtr J_{EI} J_{IE}, \quad J_{II} h_{E} \gtrless J_{EI} h_{I}, \quad J_{IE} h_{E} \lessgtr J_{EE} h_{I} \\ \Leftrightarrow & \frac{J_{EE}}{J_{IE}} \gtrless \frac{J_{EI}}{J_{II}}, \quad \frac{h_{E}}{h_{I}} \lessgtr \frac{J_{EI}}{J_{II}}, \quad \frac{h_{E}}{h_{I}} \lessgtr \frac{J_{EE}}{J_{IE}} \\ \Leftrightarrow & \frac{h_{E}}{h_{I}} \lessgtr \frac{J_{EI}}{J_{II}} \lessgtr \frac{J_{EE}}{J_{IE}}

を得る.

(11)\frac{h_{E}}{h_{I}} > \frac{J_{EI}}{J_{II}} > \frac{J_{EE}}{J_{IE}}

(12)\frac{h_{E}}{h_{I}} < \frac{J_{EI}}{J_{II}} < \frac{J_{EE}}{J_{IE}}

無活動・均衡-固定点の非存在条件

片方の集団の発火率がゼロで, もう片方の集団の発火率がゼロでもなく飽和もしていない, 無活動・均衡-固定点 存在しない(十分)条件を求める.

k = E, l = I または k = I, l = E とする. 集団 k の発火率がゼロ (m_k \ll 1) だとすると, この状態が固定点となる必要十分条件は, 式 (10) より, z_k = -|\Omega(C)| つまり,

J_{kl} m_l + h_k = - |\Omega(1/C)|

である. l = E ならば左辺の項はすべて正なので実現不可能である. ゆえに, k = E, l = I でなければならない. さらに, 集団 l = I の発火率が非ゼロで有限 (m_I = \Theta(1)) である必要十分条件は式 (10) より, z_I = \Theta(1), つまり,

J_{II} \, m_I + h_I = \Theta(1/C)

である. これを m_I について解けば

m_I = - h_I / J_{II} + \Theta(1/C)

となる. これを z_E = -|\Omega(C)| となる条件にあわせると,

& J_{EI} (- h_I / J_{II} + \Theta(1/C)) + h_E = - |\Omega(1/C)| \\ & \Leftrightarrow - h_I J_{EI} / J_{II} + h_E = - |\Omega(1/C)| \\ & \Leftrightarrow h_I \left( \frac{J_{EI}}{J_{II}} - \frac{h_E}{h_I} \right) = |\Omega(1/C)| \\ & \Leftrightarrow \frac{J_{EI}}{J_{II}} > \frac{h_E}{h_I}

最後の式は, 式 (12) の必要条件であり, また, 式 (11) の下では偽である. ゆえに, 式 (11) はこの固定点が存在しない十分条件である. 以下では, すべて式 (11) が成り立つと仮定して議論を行う.

飽和-固定点の非存在条件

両方の集団の発火率が飽和した固定点 (飽和-固定点) m_E = m_I = m^{\max} を考える. 条件 (10) (z_k = |\Omega(C)| \Leftrightarrow f_k(\bm z) \sim m^{\max}) より, k = E, I について,

J_{kE} \, m^{\max}_E + J_{kI} \, m^{\max}_I + h_k = |\Omega(1/C)|

が成り立つことが, 飽和-固定点が存在する必要十分条件である. m^{\max} = \Omega(1) より, 左辺のオーダーは無条件に \Omega(1) なので, これは条件

J_{kE} \, m^{\max} + J_{kI} \, m^{\max} + h_k > 0

と同値である. これをさらに同値変形することにより,

\begin{align*} & k = E, I: \quad (J_{kE} + J_{kI}) m^{\max} > - h_k \\ \Leftrightarrow \enskip & k = E, I: \quad J_{kE} + J_{kI} \ge 0 \enskip \text{or} \enskip m^{\max} < \frac{- h_k}{J_{kE} + J_{kI}} \\ \Leftrightarrow \enskip & \left( J_{EE} + J_{EI} \ge 0 \enskip \& \enskip J_{IE} + J_{II} \ge 0 \right) \tag{c1} \\ & \quad \text{or} \quad \left( J_{EE} + J_{EI} \ge 0 \enskip \& \enskip m^{\max} < \frac{- h_I}{J_{IE} + J_{II}} \right) \tag{c2} \\ & \quad \text{or} \quad \left( m^{\max} < \frac{- h_E}{J_{EE} + J_{EI}} \enskip \& \enskip J_{IE} + J_{II} \ge 0 \right) \tag{c3} \\ & \quad \text{or} \quad m^{\max} < \min_{k=E,I} \frac{- h_k}{J_{kE} + J_{kI}} \tag{c4} \end{align*}

となる. 前節までで導いた均衡条件不等式 (11) は

J_{EE} - \frac{J_{IE}}{J_{II}} J_{EI} < 0

を導く. よって, 式 (11) の条件下で 1 \le - {J_{IE}}/{J_{II}} (つまり J_{IE} + J_{II} \ge 0) ならば J_{EE} + J_{EI} \ge 0 は満たされないことが分かる. ゆえに, 式 (11) の条件下で J_{EE} + J_{EI} \ge 0 と J_{IE} + J_{II} \ge 0 は同時に成立しない. ゆえに 式 (11) の条件下で外部入力に依存しない飽和-固定点の存在する条件 (c1) は成立し得ない. 以上の議論をあわせると, 式 (11) の条件下で 飽和-固定点 が存在しない必要十分条件は 条件 (c2)–(c3) である. 条件 (c3) の否定をさらに強化した条件

(13)m^{\max} \ge \max_{k=E,I} \frac{- h_k}{J_{kE} + J_{kI}}

は, 条件 (c2)–(c3) の否定を導くので, 飽和-固定点 が存在しない, 十分条件である (単なる条件 (c3) の否定は, 上式において max ではなく min を使うことに注意). これから, h_E と h_I が十分小さければ, 飽和-固定点 は存在しないことが分かる.

発火率が発散しない条件

発火率の上限が無い (m^{\max} = \infty) 場合, 条件 (c2)-(c3) は成立しようが無いから, (c1) が存在しなければ 飽和-固定点 (この場合は発火率が発散するような解) は存在しない. さらに, 式 (11) の条件下で J_{IE} + J_{II} < 0 は

J_{EI} < J_{II} \frac{J_{EE}}{J_{IE}} < - J_{IE} \frac{J_{EE}}{J_{IE}} = - J_{EE}

つまり J_{EE} + J_{EI} < 0 を導く (J_{IE} + J_{II} < 0 のほうが J_{EE} + J_{EI} < 0 より厳しい条件). よって, 飽和-固定点 が存在しない, 必要十分条件は

(14)J_{EE} + J_{EI} < 0

である. [2]

[2][vanVreeswijk1998] の式 (4.10) の条件 J_E > 1 は J_{EE} = J_{IE} = 1, J_{kI} = - J_k なる正規化(変数変換)の元で J_{EE} + J_{EI} < 0 と同値である. 発火率の上限が無い場合はこの条件のみで 飽和-固定点 が存在しないことが保証されるが, 発火率の上限が有限の場合は, 条件 (c4) はまだ成立し得るので, 式 (4.10) だけでは飽和-固定点の存在を否定していない.

飽和・均衡-固定点の非存在条件

飽和・均衡-固定点 が存在しない十分条件を求める. k = E, l = I または k = I, l = E とする. m_k \sim m^{\max} かつ集団 l の発火率が非ゼロで飽和していない (m_l = \Theta(1)) という必要十分条件は, z_k = |\Theta(C)| かつ z_l = \Theta(1), つまり,

J_{kk} \, m^{\max}_k + J_{kl} \, m_l + h_k = |\Theta(1)| J_{lk} \, m^{\max}_k + J_{ll} \, m_l + h_l = \Theta(1/C)

と同値である. 前者の式に後者の式を m_l について解いた結果

m_l = \frac{\Theta(1/C) - J_{lk} \, m^{\max}_k - h_l}{J_{ll}}

を代入して同値変形すると,

& J_{kk} \, m^{\max}_k - J_{kl} \, \frac{J_{lk} \, m^{\max}_k + h_l}{J_{ll}} + h_k = |\Theta(1)| \\ & \Leftrightarrow \left( \frac{J_{kk}}{J_{kl}} - \frac{J_{lk}}{J_{ll}} \right) J_{kl} \, m^{\max}_k - \frac{J_{kl}}{J_{ll}} h_l + h_k = |\Theta(1)| \\ & \Leftrightarrow \left( \frac{J_{kk}}{J_{kl}} - \frac{J_{lk}}{J_{ll}} \right) J_{kl} \, m^{\max}_k - \frac{J_{kl}}{J_{ll}} h_l + h_k > 0

を得る. これを2つの場合に分けてさらに同値変形する.

Case 1: k = E, l = I

この 飽和・均衡-固定点 の存在する必要十分条件は,

& \left( \frac{J_{EE}}{J_{EI}} - \frac{J_{IE}}{J_{II}} \right) J_{EI} \, m^{\max}_E - \frac{J_{EI}}{J_{II}} h_I + h_E > 0 \\ & \Leftrightarrow \left( \frac{J_{EE}}{J_{IE}} - \frac{J_{EI}}{J_{II}} \right) J_{IE} \, m^{\max}_E + h_I \left( \frac{h_E}{h_I} - \frac{J_{EI}}{J_{II}} \right) > 0 \\ & \Leftrightarrow m^{\max}_E < \frac{h_I}{J_{IE}} \left( \frac{h_E}{h_I} - \frac{J_{EI}}{J_{II}} \right) / \left( \frac{J_{EI}}{J_{II}} - \frac{J_{EE}}{J_{IE}} \right)

Case 2: k = I, l = E

この 飽和・均衡-固定点 の存在する必要十分条件は,

& \left( \frac{J_{II}}{J_{IE}} - \frac{J_{EI}}{J_{EE}} \right) J_{IE} \, m^{\max}_I - \frac{J_{IE}}{J_{EE}} h_E + h_I > 0 \\ & \Leftrightarrow \left( \frac{J_{EE}}{J_{IE}} - \frac{J_{EI}}{J_{II}} \right) \frac{J_{IE} J_{II}}{J_{EE}} m^{\max}_I + h_I \frac{J_{IE}}{J_{EE}} \left( \frac{J_{EE}}{J_{IE}} - \frac{h_E}{h_I} \right) > 0 \\ & \Leftrightarrow \left( \frac{J_{EI}}{J_{II}} - \frac{J_{EE}}{J_{IE}} \right) J_{II} \, m^{\max}_I + h_I \left( \frac{h_E}{h_I} - \frac{J_{EE}}{J_{IE}} \right) > 0 \\ & \Leftrightarrow m^{\max}_I < - \frac{h_I}{J_{II}} \left( \frac{h_E}{h_I} - \frac{J_{EE}}{J_{IE}} \right) / \left( \frac{J_{EI}}{J_{II}} - \frac{J_{EE}}{J_{IE}} \right)

2つの場合を統合

以上の結果を合わせると, 飽和・均衡-固定点 の2つの場合のうち少なくともいずれか1つが存在する必要十分条件は,

m^{\max} < \frac{h_I}{\min \left\{J_{IE}, - J_{II} \right\}} \left( \frac{h_E}{h_I} - \frac{J_{EE}}{J_{IE}} \right) / \left( \frac{J_{EI}}{J_{II}} - \frac{J_{EE}}{J_{IE}} \right)

である. ここで, \max \left\{J_{IE}^{-1}, - J_{II}^{-1} \right\} = 1 / \min \left\{J_{IE}, - J_{II} \right\} なる関係を用いた. よって, いかなる 飽和・均衡-固定点 も存在しない, 十分条件は,

(15)m^{\max} \ge \frac{h_I}{\min \left\{J_{IE}, - J_{II} \right\}} \left( \frac{h_E}{h_I} - \frac{J_{EE}}{J_{IE}} \right) / \left( \frac{J_{EI}}{J_{II}} - \frac{J_{EE}}{J_{IE}} \right)

となる.

飽和・無活動-固定点の非存在条件

飽和・無活動-固定点 は, m_E = m^{\max}, m_I = 0 の場合は z_I = C (J_{IE} \, m^{\max} + h_I) \neq - |\Omega(C)| より存在できない. 一方, m_E = 0, m_I = m^{\max} は z_E = C (J_{EI} \, m^{\max} + h_E) = - |\Omega(C)| かつ z_I = C (J_{II} \, m^{\max} + h_I) = |\Omega(C)|, つまり J_{EI} \, m^{\max} + h_E < 0 かつ J_{II} \, m^{\max} + h_I > 0 ならば存在する. これは,

- \frac{h_E}{J_{EI}} < m^{\max} < - \frac{h_I}{J_{II}}

とも書ける. これが成立する必要条件は

\frac{h_E}{h_I} < \frac{J_{EI}}{J_{II}}

だが, これは式 (11) の下では成立しない. ゆえに, いかなる 飽和・無活動-固定点 も式 (11) の下では無条件で存在しない.

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

安定化可能性

均衡条件不等式 で求めた興奮・抑制ネットワークの均衡固定点の安定性を解析する. 均衡固定点の安定性 で示したように, 均衡固定点からのずれ \bm \xi = \bm m - \bm m^0 のダイナミクスは

\frac{\D \bm \xi}{\D t} \sim C \bm S \bm J \bm \xi

に従う. ただし,

\bm S = \diag (S_E, S_I) = \bm \tau^{-1} D_1 \bm f

である. 均衡固定点の安定性を調べるためには行列 \bm S \bm J の固有値の実部を調べれば良い. 行列 \bm S \bm J は 2 \times 2 行列なので, その固有値は

\lambda_{\pm} = \frac{(S_E J_{EE} + S_I J_{II}) \pm \sqrt{ (S_E J_{EE} + S_I J_{II})^2 - 4 (S_E J_{EE} S_I J_{II} - J_{EI} J_{IE}) }}{2}

となる. \max \Re \lambda_{\pm} の負の成分は S_I J_{II} であるので, S_I の大きい極限 S_E \ll S_I で \max \Re \lambda_{\pm} が負に出来なければ, この系が安定となる可能性は無い. これは, 例えば \tau_E \gg \tau_I なる場合を考えることと同様である. この極限をとるためには S_E = 0 とすれば良い.

\lambda_{\pm} = \frac{S_I J_{II} \pm \sqrt{ (S_I J_{II})^2 + 4 J_{EI} J_{IE} }}{2}

(S_I J_{II})^2 + 4 J_{EI} J_{IE} \le 0 なら, \max \Re \lambda_{\pm} = S_I J_{II} / 2 < 0 である. (S_I J_{II})^2 + 4 J_{EI} J_{IE} > 0 なら, 4 J_{EI} J_{IE} < 0 より, S_I J_{II} + \sqrt{(S_I J_{II})^2 + 4 J_{EI} J_{IE}} < 0 なので \max \Re \lambda_{\pm} = \lambda_+ < 0 である. ゆえに, 興奮・抑制ネットワークではどんな入出力関係 \bm f や結合パラメタ \bm J や外部入力 \bm h に対しても, 抑制性集団が興奮性集団に比べて十分速い時定数をもてば, 均衡固定点を漸近安定に出来る.

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

二状態ニューロンのネットワーク

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

モデルの定義

N_E 個の興奮性のニューロンと N_I 個の抑制性のニューロンが相互に結合したネットワークの特性を N = N_E + N_I \to \infty の極限で調べる. ここで, 興奮性ニューロンと抑制性ニューロンの非 N_E/N_I は N に依らない数に固定されているとする. これらのニューロンは平均 K 個のニューロンとランダムに結合され, その確率分布は

\Prob \{ J_{kl}^{ij} = J_{kl}/\sqrt K \} = K/N_l \Prob \{ J_{kl}^{ij} = 0 \} = 1 - K/N_l

で定義される. [1]

[1]原著 [vanVreeswijk1998] には確率は K/N_k と書いてあるがこれは間違い.

それぞれのニューロンの状態は 0 か 1 の二値変数で表現され, \sigma_k^i(t) = 1 (= 0) は活動 (無活動) 状態に対応する. ニューロンの状態は「バラバラ」なタイミング (後述) で更新され, そのニューロンへの入力 u_k^i(t) に依って, 状態は

\sigma_k^i(t) = \Theta(u_k^i(t))

に更新される. ただし, \Theta は \Theta(x) = 0 (x \le 0), \Theta(x) = 1 (x > 0) で定義される ヘヴィサイド関数 である. ニューロンへの入力は

u_k^i(t) = \sum_{l=E,I} \sum_{j=1}^{N_l} J_{kl}^{ij} \sigma_l^j(t) + u_k^0 - \theta_k

で定義される. u_k^0 は外部入力を表し, \theta_k はニューロンの閾値を表す. u_k^0 が外部集団からの入力を表すならば, その大きさは O(\sqrt K) のオーダーになるはずである (はじめに 参照). よって, u_k^0 は外部集団との結合強度 J_{k0} > 0 と外部集団の「活動率」 m_0 を使い,

u_k^0 = \sqrt K J_{k0} m_0 \qquad (k = E, I)

と置く.

集団 k のニューロン i はそれぞれ平均 \tau_k の時間間隔で「バラバラ」に (非同期的に) 更新される. ここでは, 各ニューロンの 更新の時間間隔は平均 \tau_k の独立な指数分布に従う (つまり, 各ニューロンの更新は独立なポアソン過程である) とする. [2]

[2]更新のタイミングがランダムでないバージョンのネットワークも [vanVreeswijk1998] では解析していて, その統計的な振る舞いは同一であることが示されている.
\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

平均場方程式の導出

ニューロン i の初期値 \sigma_k^i(0) とその更新のランダム性と初期値 \bm \sigma に 関する平均を \AvgDyn{\bullet} [1] と書き, ニューロン i の(局所)活動率を

m_k^i (t) = \AvgDyn{\sigma_k^i(t)}

と定義する.

[1]初期値 \bm \sigma に関する平均とは, 時刻 0 での集団活動率 m_k(0) が \sigma_k^i(0) が 1 の確率で, それぞれの i について \sigma_k^i(0) が独立, という確率分布に関する平均である.

課題

他の場所では, \AvgDyn{\bullet} は使われていない. 使うべき? 例えば, 他の場所では m_k = \PAvg{\sigma_k^i} だけど, ここでは m_k = \PAvg{\AvgDyn{\sigma_k^i}} である.

ニューロン i への入力が閾値を超える確率 [2] は, 更新時間のランダム性に関する平均 \AvgDyn{\Theta(u_k^i (t))} で表すことが出来る. [3] さらに, 更新のタイミングはポアソン過程で表され, 入力が閾値を超えていた場合に状態 1 へ遷移する単位時間当たりの条件付き確率は, 1 / \tau_k である. こられらを合わせると, 単位時間あたりにニューロン i が状態を 1 に遷移する確率は \AvgDyn{\Theta(u_k^i (t))} / \tau_k で与えられていることが分かる. 期待値の時間発展 の関係式を用いれば,

\tau_k \frac{\D}{\D t} m_k^i (t) = - m_k^i (t) + \AvgDyn{\Theta(u_k^i (t))}

と書くことができる.

[2]正確には, 系の状態 \bm \sigma が与えられた時の条件付き確率, である.
[3]確率変数 X について事象 X \in A が起こる確率は指示関数 (indicator function) 1_A(X) を用いて \Prob \{ X \in A \} = \Expect \{1_A(X)\} と書けることを思い出そう.

この式の集団平均をとる (つまり両辺に \PAvg{\bullet}_i を施す) と, 集団活動率 m_k のダイナミクスを表す式

\tau_k \frac{\D}{\D t} m_k (t) = - m_k (t) + \PAvg{\AvgDyn{\Theta(u_k^i (t))}}_i

を得る. この節では, 右辺第二項

F_k(m_E, m_I) := \PAvg{\AvgDyn{\Theta(u_k^i (t))}}_i

を計算する. 素朴に考えれば右辺は系の微視的な状態 \bm \sigma に依存しているはずだが, 左辺は巨視的な状態, つまり集団活動率 m_k のみに依存することを主張している. この微視的な状態への非依存性は F_k(m_E, m_I) の計算の過程で自動的に出てくる結果である.

確率 F_k(m_E, m_I) は以下の仮定 [4] のもとで計算することが出来る.

仮定

すべてのニューロンの活動が無相関である.

形式的に書けば, いかなるふたつのニューロン (i, k) と (j, l) (k, l \in \{E, I\}, i = 1, \ldots, N_k, j = 1, \ldots, N_l) についても, それぞれの活動 X = \AvgDyn{\Theta(u_k^i (t))}, Y = \AvgDyn{\Theta(u_l^j (t))} はすべての時間 t について無相関, つまり,

\lim_{N \to \infty} \AvgJ{ \left( X - \AvgJ{X} \right) \left( Y - \AvgJ{Y} \right) } = 0

が, 成り立つ.

[4]原著 [vanVreeswijk1998] での仮定は「すべてのニューロンについて, それに結合しているすべてのニューロンの活動が無相関である」であり, 本稿で使っている仮定より若干弱い. しかし, 無相関性の「証明」 より本稿で使っている仮定は [vanVreeswijk1998] の仮定と同じ条件 K \ll \log N で成り立つことが分かる. さらに, すべてのニューロンが無相関でなければ, 無相関変数に対する大数の法則 が使えない (自己平均性 (self-averaging property) を参照).

これは, K \ll \log N が成り立てば成り立つ. 詳しい議論については, 無相関性の「証明」 を参照.

自己平均性 (self-averaging property) を [\Theta(u_k^i (t))]_i の計算に適用すれば, [\bullet]_i と \AvgJ{\bullet} を交換することが出来て,

F_k(m_E, m_I) \approx \AvgJ{\AvgDyn{\Theta(u_k^i (t))}}

を計算すれば良いことが分かる.

ニューロン i が n_E(t) 個の興奮性ニューロンと n_I(t) 個の抑制性ニューロンから入力を受けているとすれば, その全入力は

u_k^i (t) = \sqrt K J_{k0} m_0 + \frac{J_{kE}}{\sqrt K} n_E (t) + \frac{J_{kI}}{\sqrt K} n_I (t) - \theta_k

となる. 確率 F_k(m_E, m_I) はこの入力が正である確率であり,

F_k(m_E, m_I) &\approx \sum_{n_1, n_2 = 0}^\infty p_1(n_1 | m_1) \, p_2(n_2 | m_2) \, \Theta \left( \sqrt K J_{k0} m_0 + \sum_{l=1,2} \frac{J_{kl}}{\sqrt K} n_l - \theta_k \right) \\

となる. ただし, p_l (n_l | m_l) は集団 l \in \{E, I\} の活動率が m_l の時にニューロン i が集団 l から n_l 個の入力を受ける確率であり,

p_l (n | m_l) & \approx \sum_{s=n}^\infty \underbrace{ \frac{K^s}{s!} \E^{-K} }_{\text{(P1)}} \underbrace{ \begin{pmatrix} s \\ n \end{pmatrix} (m_l)^n (1 - m_l)^{s-n} }_{\text{(P2)}} \\ & = \frac{(m_l K)^n}{n!} \E^{-m_l K}

となる. ここで, (P1) は集団 k のニューロン (どのニューロンでも成立する) が集団 l の s 個のニューロンからの結合を持つ確率であり, (P2) はその s 個のニューロンのうち n 個のニューロンが活動している (\sigma_l^j = 1 である) 確率である. 最後の等式は, \exp の定義に基づけば、以下の計算で確認できる.

& \sum_{s=n}^\infty \frac{K^s}{s!} \E^{-K} \begin{pmatrix} s \\ n \end{pmatrix} (m_l)^n (1 - m_l)^{s-n} \\ & = \sum_{s=n}^\infty \frac{K^s}{s!} \E^{-K} \frac{s!}{n! (s-n)!} (m_l)^n (1 - m_l)^{s-n} \\ & = \frac{\E^{-K} (K m_l)^n}{n!} \sum_{s=n}^\infty \frac{K^{s-n}}{(s-n)!} (1 - m_l)^{s-n} \\ & = \frac{\E^{-K} (K m_l)^n}{n!} \E^{K \, (1-m_l)} \\ & = \frac{(K m_l)^n}{n!} \E^{-K m_l}

この確率分布は平均と分散が m_l K の ポアソン分布 (Poisson distribution) なので, 極限 K \to \infty, つまりこの平均と分散が大きな極限では ガウス分布 (Gaussian distribution)

\frac{1}{\sqrt{2 \pi m_l K}} \exp \left( - \frac{(n_l - m_l K)^2}{2 m_l K} \right)

で近似できる. この極限 K \to \infty で,

F_k(m_E, m_I) & = \sum_{n_1, n_2 = 0}^\infty p_1(n_1 | m_1) \, p_2(n_2 | m_2) \, \Theta \left( \sqrt K J_{k0} m_0 + \sum_{l=1,2} \frac{J_{kl}}{\sqrt K} n_l - \theta_k \right) \\ & \overset{(1)} \approx \int Dx_1 \, Dx_2 \, \Theta \left( \sqrt K J_{k0} m_0 + \sum_{l=1,2} \frac{J_{kl}}{\sqrt K} (m_l K + x_l \sqrt{m_l K}) - \theta_k \right) \\ & \overset{(2)} = \int Dx_1 \, Dx_2 \, \Theta \left( u_k + \sum_{l=1,2} x_l J_{kl} \sqrt{m_l} \right) \\ & \overset{(3)} = \int Dx \, \Theta (u_k + \sqrt{\alpha_k} x) \\ & \overset{(4)} = H \left( \frac{- u_k}{\sqrt{\alpha_k}} \right)

と計算できる. ここで,

Dx & := \D x \exp(-x^2) / \sqrt{2 \pi} \\ u_k & := (J_{k0} m_0 + J_{kE} m_E + J_{kI} m_I) \sqrt K - \theta_k, \\ \alpha_k & := (J_{kE})^2 m_E + (J_{kI})^2 m_I \\ H(z) & := \int_z^\infty \frac{\D x}{\sqrt{2 \pi}} \exp(- x^2 / 2)

である. 上記の u_k と \alpha_k はただ変数に名前をつけただけだが, これらの物理的意味については 入力のゆらぎ を参照せよ. Dx は ガウス測度 (Gaussian measure) と呼ばれるただの省略記号である. 関数 H(z) は Q関数 と呼ばれる関数である. 上の計算では, (1) n_l \approx m_l K + x_l \sqrt{m_l K} なる近似と ガウス確率変数の変数変換, (2) u_k の定義, (3) ヘヴィサイド関数の多重ガウス積分とQ関数 の関係, (4) H(z) の定義をそれぞれ用いた.

無相関性の「証明」

以下の議論は [Derrida1987] に依る.

今, 初期状態から n 回の更新が起こったとする. いかなるニューロンも, n 回の更新の前まで遡れば最大でも K^n 個 [5] のニューロンの初期状態に依存している. 2つのニューロンから伸びる「木」はそれぞれ平均で K^n の「枝」をもつ. この中で最低でも1つの枝が同じニューロンに繋がっている確率は, (1) 2つの木からそれぞれの1つの枝を選ぶ方法の総数と, (2) 1つのニューロンの選び方の総数と, (3) ある1つのニューロンを2回選ぶ確率の積なので,

p = \underbrace{(K^n)^2 \vphantom{\frac 1 2}}_{(1)} \underbrace{N \vphantom{\frac 1 2}}_{(2)} \underbrace{\frac 1 {N^2}}_{(3)} = O(K^{2n} / N)

となる. これが 0 に漸近する, つまり p \ll 1 (as N \to \infty) という条件から, K^n \ll \sqrt N が導かれる. いかなる自然数 n でもこれが成り立つには K \ll \log N であれば十分である.

[5]ただし, 各ニューロンの結合の数が平均 K 個のまわりでゆらいでいる効果は無視している.
\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

入力のゆらぎ

ここでは, 入力のゆらぎ [(\Devi u_k^i (t))^2]_i と集団平均活動率 m_l(t) を結ぶ式

(16)[(\Devi u_k^i (t))^2]_i \xrightarrow{N \to \infty} \sum_{l = E, I} ( J_{kl} )^2 \, m_l(t) =: \alpha_k(t)

を導出する. ここで [ \bullet ]_i は集団 k \in \{ E, I \} 内のニューロンに関する 集団平均 (population average; 添字 i に沿った平均) を表し,

[ \bullet ]_i = \frac{1}{N_k} \sum_{i=1}^{N_k} \bullet

で定義される. 文脈からどの添字に関する平均かが明らかであれば, 添字 i は省略する. また, \Devi は集団平均からの偏差 \Devi X_i := X_i - [X_i]_i である.

入力の集団平均

u_k(t) & = [u_k^i (t)]_i \\ & = \left[ \sum_{l = E, I} \sum_{j=1}^{N_l} J_{kl}^{ij} \sigma_l^j(t) + u_k^0 - \theta_k \right]_i \\ & = \sum_{l = E, I} \underbrace{ \left[ \sum_{j=1}^{N_l} J_{kl}^{ij} \sigma_l^j(t) \right]_i }_{(*)} + u_k^0 - \theta_k

この (*) は以下のように近似的できる.

(*) & \overset{(1)} \approx \sum_{j=1}^{N_l} \left[ J_{kl}^{ij} \right]_i \sigma_l^j(t) \\ & \overset{(2)} \approx \sum_{j=1}^{N_l} J_{kl} \frac{\sqrt K}{N_l} \sigma_l^j(t) \\ & \overset{(3)} = J_{kl} \frac{\sqrt K}{N_l} \sum_{j=1}^{N_l} \sigma_l^j(t) \\ & \overset{(4)} = J_{kl} \sqrt K \left[ \sigma_l^j(t) \right]_j \\ & \overset{(5)} = J_{kl} \sqrt K m_l(t)

ここで, (1) は, [ \bullet ]_i の線型性 (よって和 \sum_{j=1}^{N_l} と集団平均 [ \bullet ]_i は演算順序を入れ替えて良い) と \sigma_l^j(t) が i に依らないこと (これの正しい解説は ニューロンの状態と結合係数の相関 を参照), (2) 大数の法則 (law of large numbers) より算術平均は期待値に収束する, つまり N \to \infty の極限で \left[ J_{kl}^{ij} \right]_i \to \AvgJ{J_{kl}^{ij}} が成り立つことと, \Prob \{ J_{kl}^{ij} = {J_{kl}}/{\sqrt K}\} = {K}/{N_l} と \Prob \{ J_{kl}^{ij} = 0\} = 1 - {K}/{N_l} から期待値は \AvgJ{J_{kl}^{ij}} = ({J_{kl}}/{\sqrt K}) ({K}/{N_l}) = J_{kl} {\sqrt K}/{N_l} となること, (3) J_{kl} \sqrt K / N_l が j に依らない定数であること, (4) 集団平均の定義, (5) m_l(t) の定義を用いた.

これらの計算を合わせ, u_k^0 = \sqrt K E_k m_0 を思い出せば, 入力の集団平均 u_k(t) は

u_k(t) & = ... \\ & \approx \sum_{l = E, I} J_{kl} \sqrt K m_l(t) + u_k^0 - \theta_k \\ & = \sqrt K \left( \sum_{l = E, I} J_{kl} m_l(t) + E_k m_0 \right) - \theta_k

となる.

ニューロンの状態と結合係数の相関

上記の式変換(1)で「\sigma_l^j(t) が i に依らない」ことを用いたが, これは正しくは, J_{kl}^{ij} と \sigma_l^j(t) の相関が無い [1] と仮定することで, \sigma_l^j(t) は i, つまり J_{kl}^{ij} に依らずに決まるから集団平均 [ \bullet ]_i の演算にとっては定数として扱えることから言える. この J_{kl}^{ij} と \sigma_l^j(t) が無相関であるという仮定は,

  1. 事象 J_{kl}^{ij} \neq 0 と J_{lk}^{ji} \neq 0 が独立 (J_{kl}^{ij} の定義より)
  2. \sigma_l^j(t) と \sigma_l^{j'}(t) (j' \neq j) が無相関という仮定

から正当化される. この無相関の仮定は有限の N, K では正しくないので, この式変形は完全な等号では結ばれず, \approx と書いている.

[1]

原著 [vanVreeswijk1998] では,

Note that on the right-hand side (r.h.s.) of equation 3.11 we have neglected the correlations between the random fluctuations in the activity of a cell and the particular realization of its output connectivity. This is justified since such correlations are weak in the large N limit.

—p.1329 (p.9), van Vreeswijk, Sompolinsky (1998)

と説明されている. この “equation 3.11” はここで扱っている入力の集団平均 u_k(t) のことである.

課題

式変形 [J_{kl}^{ij} \, \sigma_l^j(t)]_i = [J_{kl}^{ij}]_i \, \sigma_l^j(t) を正当化する議論をもっと形式化する. 説明に自然言語つかいすぎ！

先に確率平均に行く方法もあるかも?: J_{kl}^{ij} と \sigma_l^j(t) が独立だという近似のもと, [J_{kl}^{ij} \sigma_l^j(t)]_i \approx \AvgJ{J_{kl}^{ij} \sigma_l^j(t)} = \AvgJ{J_{kl}^{ij}} \AvgJ{\sigma_l^j(t)} であることを用いる. この \AvgJ{\sigma_l^j(t)} は式変形 (4) にあるように, さらに集団平均 [\bullet]_j がかかるから, \left[ \AvgJ{\sigma_l^j(t)} \right]_j = \AvgJ{[\sigma_l^j(t)]_j} = \AvgJ{m_l(t)} となる. この系は self-averaging なので (とどこかで説明する必要があるけど,) \AvgJ{m_l(t)} = m_l(t) となる.

入力のゆらぎ

& [(\Devi u_k^i (t))^2] \\ & \overset{(1)} = \left[ \left( \Devi \left\{ \sum_{l = E, I} \sum_{j=1}^{N_l} J_{kl}^{ij} \sigma_l^j(t)) \right\} \right)^2 \right]_i \\ & \overset{(2)} = \left[ \left( \sum_{l = E, I} \sum_{j=1}^{N_l} J_{kl}^{ij} \sigma_l^j(t)) \right)^2 \right]_i - \left[ \sum_{l = E, I} \sum_{j=1}^{N_l} J_{kl}^{ij} \sigma_l^j(t)) \right]_i^2 \\ & \overset{(3)} = \left[ \left( \sum_{l = E, I} \sum_{j=1}^{N_l} J_{kl}^{ij} \sigma_l^j(t)) \right)^2 \right]_i - K \left(\sum_{l = E, I} J_{kl} m_l(t) \right)^2

ここで, (1) \Devi(x + \text{const.}) = \Devi x, (2) [(\Devi x)^2] = [x^2] - [x]^2, (3) 上記の u_k(t) の計算 (特に (*) の部分) を用いた.

& \left[ \left( \sum_{l = E, I} \sum_{j=1}^{N_l} J_{kl}^{ij} \sigma_l^j(t)) \right)^2 \right]_i \\ & = \left[ \sum_{l, l' = E, I} \sum_{j=1}^{N_l} \sum_{j'=1}^{N_{l'}} J_{kl}^{ij} J_{kl'}^{ij'} \sigma_l^j(t)) \sigma_{l'}^{j'}(t)) \right]_i \\ & = \sum_{l, l' = E, I} \sum_{j=1}^{N_l} \sum_{j'=1}^{N_{l'}} \underbrace{ \left[ J_{kl}^{ij} J_{kl'}^{ij'} \right]_i \sigma_l^j(t) \sigma_{l'}^{j'}(t) }_{(*)}

上式の (*) の和は, 恒等式 1 = \delta_{ll'} (\delta_{jj'} + (1 - \delta_{jj'})) + (1 - \delta_{ll'}) を用いて [2]

\sum_{l, l' = E, I} \sum_{j=1}^{N_l} \sum_{j'=1}^{N_{l'}} \bullet_{l,l',j,j'} = \sum_{l = E, I} \sum_{j=1}^{N_l} \bullet_{l,l,j,j} + \sum_{l = E, I} \sum_{\substack{j,j'=1 \\ j \neq j'}}^{N_l} \bullet_{l,l,j,j} + \sum_{\substack{l, l' = E, I \\ l \neq l'}} \sum_{j=1}^{N_l} \sum_{j'=1}^{N_{l'}} \bullet_{l,l',j,j'}

のように分解できる.

[2]\sum_{j,j'=1} (1 - \delta_{j,j'}) \bullet = \sum_{\substack{j,j'=1 \\ j \neq j'}} \bullet
第一項の計算 (l = l', j = j')

& \sum_{l = E, I} \sum_{j=1}^{N_l} \left[ (J_{kl}^{ij})^2 \right]_i (\sigma_l^j(t))^2 \\ & \overset{(1)} \approx \sum_{l = E, I} \sum_{j=1}^{N_l} \AvgJ{(J_{kl}^{ij})^2} \, \sigma_l^j(t) \\ & \overset{(2)} = \sum_{l = E, I} \sum_{j=1}^{N_l} \left( \frac{J_{kl}}{\sqrt K} \right)^2 \frac{K}{N_l} \, \sigma_l^j(t) \\ & = \sum_{l = E, I} ( J_{kl} )^2 \frac{1}{N_l} \sum_{j=1}^{N_l} \sigma_l^j(t) \\ & \overset{(3)} = \sum_{l = E, I} ( J_{kl} )^2 \, m_l(t)

ここで, (1) 大数の法則 (law of large numbers) と \sigma_l^j(t) の取りうる値は 0 か 1 なので \left( \sigma_l^j(t) \right)^2 = \sigma_l^j(t) [3], (2) \Prob \{ J_{kl}^{ij} = J_{kl}/\sqrt K \} = 1 - \Prob \{ J_{kl}^{ij} = 0 \} = K/N_l (結合確率の定義 を参照), (3) m_l(t) = [\sigma_l^j(t)]_j = \sum_{j=1}^{N_l} \sigma_l^j(t) / N_l, を用いた.

[3]二値変数のからむ計算ではよく使われるテクニック.
第二項の計算 (l = l', j \neq j')

& \sum_{l = E, I} \sum_{\substack{j,j'=1 \\ j \neq j'}}^{N_l} \left[ J_{kl}^{ij} J_{kl}^{ij'} \right]_i \, \sigma_l^j(t) \, \sigma_{l}^{j'}(t) \\ & \overset{(1)} \approx \sum_{l = E, I} \sum_{\substack{j,j'=1 \\ j \neq j'}}^{N_l} \AvgJ{J_{kl}^{ij} J_{kl}^{ij'}} \, \sigma_l^j(t) \, \sigma_{l}^{j'}(t) \\ & \overset{(2)} = \sum_{l = E, I} \sum_{\substack{j,j'=1 \\ j \neq j'}}^{N_l} \left( \frac{J_{kl}}{\sqrt K} \right)^2 \frac{K}{N_l} \frac{K}{N_l} \, \sigma_l^j(t) \, \sigma_{l}^{j'}(t) \\ & = K \sum_{l = E, I} (J_{kl})^2 \frac{1}{N_l} \sum_{j=1}^{N_l} \sigma_l^j(t) \left( \sum_{j'=1}^{N_l} \frac{1}{N_l} \sigma_{l}^{j'}(t) - \frac{1}{N_l} \sigma_{l}^{j}(t) \right) \\ & = K \sum_{l = E, I} (J_{kl})^2 \left( \left\{ \frac{1}{N_l} \sum_{j=1}^{N_l} \sigma_l^j(t) \right\}^2 - \frac{1}{{N_l}^2} \sum_{j=1}^{N_l} (\sigma_{l}^{j}(t))^2 \right) \\ & = K \sum_{l = E, I} (J_{kl})^2 \left( (m_l(t))^2 - \frac{1}{N_l} m_l(t) \right)

ここで, (1) 大数の法則 (law of large numbers), (2) j \neq j' なので J_{kl}^{ij} と J_{kl}^{ij'} が独立であることと, J_{kl}^{ij} の確率分布 \Prob \{ J_{kl}^{ij} = J_{kl}/\sqrt K \} = 1 - \Prob \{ J_{kl}^{ij} = 0 \} = K/N_l (結合確率の定義 を参照), を用いた. 残りは単純な式変形である.

第三項の計算 (l \neq l')

& \sum_{\substack{l, l' = E, I \\ l \neq l'}} \sum_{j=1}^{N_l} \sum_{j'=1}^{N_{l'}} \left[ J_{kl}^{ij} J_{kl'}^{ij'} \right]_i \, \sigma_l^j(t) \, \sigma_{l'}^{j'}(t) \\ & \overset{(1)} \approx \sum_{\substack{l, l' = E, I \\ l \neq l'}} \sum_{j=1}^{N_l} \sum_{j'=1}^{N_{l'}} \AvgJ{J_{kl}^{ij} J_{kl'}^{ij'}} \, \sigma_l^j(t) \, \sigma_{l'}^{j'}(t) \\ & \overset{(2)} = \sum_{\substack{l, l' = E, I \\ l \neq l'}} \sum_{j=1}^{N_l} \sum_{j'=1}^{N_{l'}} \frac{J_{kl}}{\sqrt K} \frac{J_{kl'}}{\sqrt K} \frac{K}{N_l} \frac{K}{N_{l'}} \, \sigma_l^j(t) \, \sigma_{l'}^{j'}(t) \\ & = K \sum_{\substack{l, l' = E, I \\ l \neq l'}} J_{kl} J_{kl'} \frac{1}{N_l} \sum_{j=1}^{N_l} \sigma_l^j(t) \frac{1}{N_{l'}} \sum_{j'=1}^{N_{l'}} \sigma_{l'}^{j'}(t) \\ & = K \sum_{\substack{l, l' = E, I \\ l \neq l'}} J_{kl} J_{kl'} \, m_l(t) \, m_{l'}(t)

ここで, (1) 大数の法則 (law of large numbers), (2) l \neq l' なので J_{kl}^{ij} と J_{kl'}^{ij'} が独立であることと, J_{kl}^{ij} の確率分布 \Prob \{ J_{kl}^{ij} = J_{kl}/\sqrt K \} = 1 - \Prob \{ J_{kl}^{ij} = 0 \} = K/N_l (結合確率の定義 を参照), を用いた. 残りは単純な式変形である.

合計

[(\Devi u_k^i (t))^2]_i & = \sum_{l = E, I} ( J_{kl} )^2 \, m_l(t) \\ & \qquad + K \sum_{l = E, I} (J_{kl})^2 \left( (m_l(t))^2 - \frac{1}{N_l} m_l(t) \right) \\ & \qquad + K \sum_{\substack{l, l' = E, I \\ l \neq l'}} J_{kl} J_{kl'} \, m_l(t) \, m_{l'}(t) \\ & \qquad - K \left(\sum_{l = E, I} J_{kl} m_l(t) \right)^2 \\ & = \sum_{l = E, I} ( J_{kl} )^2 \, m_l(t) + \frac{K}{N_l} \sum_{l = E, I} (J_{kl})^2 m_l(t) \\ & \qquad + K \underbrace{ \left( \sum_{l, l' = E, I} J_{kl} J_{kl'} \, m_l(t) \, m_{l'}(t) - \left(\sum_{l = E, I} J_{kl} m_l(t) \right)^2 \right) }_{= 0} \\ & = \sum_{l = E, I} ( J_{kl} )^2 \, m_l(t) + O(1/N) \\ & \xrightarrow{N \to \infty} \alpha_k(t)

これで, 入力のゆらぎと集団平均活動率を結ぶ式 (16) が示された.

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

空間ゆらぎと時間ゆらぎ

ここでは, 以下の2つの問への答えを「同時に」導く:

  1. 空間ゆらぎ (クエンチされたゆらぎ): 活動率はニューロンごとにどう違っているのか?
  2. 時間ゆらぎ: ニューロンの活動率は, 時間的に一定なのだろうか? それとも平均値のまわりを時間的にゆらいでいるのだろうか?

実は, 入力のゆらぎ はこのふたつのゆらぎに分解することが出来, さらにその計算のためにはすでに求めた平均活動率 m_k の他に,「秩序変数」 q_k さえ分かれば良い. そして, この計算の結果からこれらのゆらぎの大きさだけでなく, 時間平均活動率の分布 も計算することが出来る.

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

クエンチされたゆらぎの計算

ここでは, クエンチされたゆらぎ (quenched fluctuations) [1] が

(17)\left[ \left( \Devi \Avg{u_k^i(t)}_t \right)^2 \right] \xrightarrow{N \to \infty} \sum_{l=1,2} J_{kl}^2 q_l =: \beta_k

となることを示す. ただし, \Avg{\bullet}_t は長い時間にわたる平均である. ここで, q_k は オーダーパラメター (order parameter) と呼ばれ, ニューロン i の活動率の時間平均 m_k^i を用いて,

m_k^i &:= \Avg{\AvgDyn{\sigma_k^i(t)}}_t \\ q_k &:= \PAvg{(m_k^i)^2}_i

と定義される.

[1]無理やり日本語にすると「焼入れされたゆらぎ」と言うのだろうか.
偏差の分解

まず, \Avg{u_k^i(t)}_t の集団平均 [\Avg{u_k^i(t)}_t] からのズレ具合 (偏差) を次のように, 2つの成分に分解できることを示す.

\Devi \Avg{u_k^i(t)}_t = \underbrace{ \sum_l \sum_j \Devi J_{kl}^{ij} \, m_l }_{\text{(d1)}} + \underbrace{ \sum_l \sum_j J_{kl}^{ij} \, \Devi m_l^j }_{\text{(d2)}}

ここで (d1) は 「結合数のゆらぎ」, (d2) [2] は 「時間平均活動率のゆらぎ」である. 結合数は時間によらないので, そのゆらぎが「クエンチされている」のは当然であるが, 活動率の時間平均 m_l^j も (平均操作のおかげで) 時間によらないので, そのゆらぎもクエンチされたゆらぎに含める必要がある. つまり, クエンチされたゆらぎのうち, 直接の影響である (d1) 「結合数のゆらぎ」と, それが引き起こす間接的な影響である (d2) 「時間平均活動率のゆらぎ」の2つを勘定すれば良い, という主張である.

[2]

ここでの (d1) は, 原著 [vanVreeswijk1998] の式(5.5)

\delta_1 \langle u_k^i \rangle = \sum_{l=1}^2 \sum_{j=1}^{N_l} \delta J_{kl}^{ij} [m_l^j]

と同値であることは, [m_l^j] = m_l より分かる. しかし, この表記では, [m_l^j] が何を意味する不鮮明である. 親切に書くのならば, [m_l^{j'}]_{j'} として, 集団平均が添字 j' についてとられることと添字 j への依存性が無いことを示すべきであるが, それならそもそも m_l と書く方が良い.

これは, 地道に入力の時間平均 \Avg{u_k^i(t)}_t の偏差を計算することによって示せる:

\Devi \Avg{u_k^i(t)}_t & \overset{(1)} = \Devi \Avg{ \sum_{l = E, I} \sum_{j=1}^{N_l} J_{kl}^{ij} \sigma_l^j(t) }_t \\ & \overset{(2)} = \Devi \left( \sum_{l = E, I} \sum_{j=1}^{N_l} J_{kl}^{ij} m_l^j \right) \\ & \overset{(3)} = \sum_{l = E, I} \sum_{j=1}^{N_l} J_{kl}^{ij} m_l^j - \left[ \sum_{l = E, I} \sum_{j=1}^{N_l} J_{kl}^{i'j} m_l^j \right]_{i'} \\ & \overset{(4)} = \sum_{l = E, I} \sum_{j=1}^{N_l} J_{kl}^{ij} m_l^j - \sum_{l = E, I} [J_{kl}^{i'j'}]_{i'} \sum_{j=1}^{N_l} m_l^j \qquad (\forall j') \\ & \overset{(5)} \approx \sum_{l = E, I} \sum_{j=1}^{N_l} J_{kl}^{ij} m_l^j - \underbrace{ \sum_{l = E, I} [J_{kl}^{i'j'}]_{i'} \sum_{j=1}^{N_l} m_l }_{\text{nothing depends on } j} \\ & \overset{(6)} = \sum_{l = E, I} \sum_{j=1}^{N_l} \left\{ J_{kl}^{ij} (m_l^j - m_l) - (J_{kl}^{ij} - [J_{kl}^{i'j'}]_{i'}) m_l \right\} \\ & = \text{(d1)} + \text{(d2)}

ここで, (1) \Devi(x + \text{const.}) = \Devi x, (2) 定義 m_l^j = \Avg{\sigma_l^j(t)}_t, (3) 偏差 \Devi の定義, (4) ニューロンの状態と結合係数の相関 の議論, (5) \sum_{j=1}^{N_l} m_l^j = N_l N_l^{-1} \sum_{j=1}^{N_l} m_l^j = N_l [m_l^{j''}]_{j''} = N_l m_l であり, N_l = \sum_{j=1}^{N_l} 1 なので, 結局 ... = \sum_{j=1}^{N_l} m_l, (6) - J_{kl}^{ij} m_l + J_{kl}^{ij} m_l = 0, を用いた. 式変形 (4) の右辺とそれ以降の式中に現れる j' は, 1 から N_l のどの値をとっても良い. これは 大数の法則 (law of large numbers) より [J_{kl}^{i'j'}]_{i'} が同じ値に収束するからである.

ふたつの偏差の相関

上記の計算より導かれた2つの偏差の二乗平均をとって, ゆらぎを

\left[ \left( \Devi \Avg{u_k^i(t)}_t \right)^2 \right] = \left[ \text{(d1)}^2 \right] + \left[ \text{(d2)}^2 \right]

のように求めたいが, そのためにはそれらの偏差が無相関 \PAvg{\text{(d1)}\text{(d2)}} = 0 でなければならない. これは簡単に示せる:

& \left[ \text{(d1)} \text{(d2)} \right] \\ & \overset{(1)} = \left[ \sum_{ll'jj'} \Devi J_{kl}^{ij} \, m_l \, J_{kl'}^{ij'} \, \Devi m_{l'}^{j'} \right]_i \\ & \overset{(2)} = \sum_{ll'jj'} \left[ \Devi J_{kl}^{ij} \, J_{kl'}^{ij'} \right]_i m_l \, \Devi m_{l'}^{j'} \\ & \overset{(3)} = \sum_{lj} \left( \left[(J_{kl}^{i*})^2 \right]_i - \left[J_{kl}^{i*} \right]_i^2 \right) m_l \, \Devi m_{l'}^{j} \\ & = \sum_{l} \left( \left[(J_{kl}^{i*})^2 \right]_i - \left[J_{kl}^{i*} \right]_i^2 \right) m_l \, \underbrace{\sum_j \Devi m_{l'}^{j}}_{=0} \\ & = 0

式変形 (1) では ニューロンの状態と結合係数の相関 の議論を用いた. 式変形 (2) では, \left[\Devi J_{kl}^{ij} \, J_{kl'}^{ij'} \right]_i は (l, j) \neq (l', j') だと

\left[\Devi J_{kl}^{ij} \, J_{kl'}^{ij'} \right]_i = \left[ \Devi J_{kl}^{ij} \right]_i \left[ J_{kl'}^{ij'} \right]_i = 0

なので, 非ゼロになるのは (l, j) = (l', j') の場合のみであることを用いた. 式変形 (3) は, 偏差 \Devi の定義に沿って

\left[\Devi J_{kl}^{ij} \, J_{kl}^{ij} \right]_i = \left[ (J_{kl}^{ij})^2 \right]_i - \left[ J_{kl}^{ij} \right]_i^2

という計算をすれば良い. 式変形 (3) の右辺以降に現れる添字の * は, この部分の添字が何でも良いことを表す.

結合数のゆらぎ

[\text{(d1)}^2] & = \left[ \left( \sum_l \sum_j \Devi J_{kl}^{ij} \, m_l \right)^2 \right]_i \\ & \overset{(1)} = \left[ \sum_{ll'jj'} \Devi J_{kl}^{ij} \, \Devi J_{kl'}^{ij'} \, m_l \, m_{l'} \right]_i \\ & \overset{(2)} = \sum_{ll'jj'} \left[ \Devi J_{kl}^{ij} \, \Devi J_{kl'}^{ij'} \right]_i m_l \, m_{l'} \\ & \overset{(3)} = \sum_j J_{kl}^2 \left(1 - \frac K N_l \right) \left( m_l \right)^2

ここで, (1) 和の積の計算のための添字テクニック と (2) ニューロンの状態と結合係数の相関 の議論を用いた. 最後の式変形 (3) では, (l, j) \neq (l', j') だと

\left[ \Devi J_{kl}^{ij} \, \Devi J_{kl'}^{ij'} \right]_i = \left[ \Devi J_{kl}^{ij} \right]_i \left[ \Devi J_{kl'}^{ij'} \right]_i = 0

となり, (l, j) = (l', j') だと

\left[\left( \Devi J_{kl}^{ij} \right)^2 \right]_i & \overset{(1)} = \left[\left( J_{kl}^{ij} \right)^2 \right]_i - \left( \left[ J_{kl}^{ij} \right]_i \right)^2 \\ & \overset{(2)} \approx \left( \frac{J_{kl}}{\sqrt K} \right)^2 \frac{K}{N_l} - \left( \frac{J_{kl}}{\sqrt K} \frac{K}{N_l} \right)^2 \\ & = \frac{J_{kl}^2}{N_l} \left( 1 - \frac{K}{N_l} \right)

となることを用いた. この計算では, (1) 偏差 \Devi の定義を使い, (2) 大数の法則 (law of large numbers) と 結合確率の定義 による期待値の計算をした.

時間平均活動率のゆらぎ

[\text{(d2)}^2] & = \left[ \left( \sum_l \sum_j J_{kl}^{ij} \, \Devi m_l^j \right)^2 \right]_i \\ & \overset{(1)} = \left[ \sum_{ll'jj'} J_{kl}^{ij} \, J_{kl'}^{ij'} \Devi m_l^j \, \Devi m_{l'}^{j'} \right]_i \\ & \overset{(2)} \approx \sum_{ll'jj'} \left[ J_{kl}^{ij} \, J_{kl'}^{ij'} \right]_i \Devi m_l^j \, \Devi m_{l'}^{j'} \\ & = \sum_{\substack{ll' \\ l \neq l'}} \sum_{jj'} \bullet + \sum_l \sum_{\substack{jj' \\ j \neq j'}} \bullet + \sum_l \sum_j \bullet

ここで, (1) 和の積の計算のための添字テクニック と (2) ニューロンの状態と結合係数の相関 の議論を用いた. 上記の3つの項は以下のように計算できる.

\sum_{\substack{ll' \\ l \neq l'}} \sum_{jj'} \left[ J_{kl}^{ij} \, J_{kl'}^{ij'} \right]_i \Devi m_l^j \, \Devi m_{l'}^{j'} & = \sum_{\substack{ll' \\ l \neq l'}} \left[J_{kl}^{i*} \, J_{kl'}^{i*} \right]_i \sum_j \Devi m_l^j \sum_{j'} \Devi m_{l'}^{j'} = 0

\sum_l \sum_{\substack{jj' \\ j \neq j'}} \left[ J_{kl}^{ij} \, J_{kl'}^{ij'} \right]_i \Devi m_l^j \, \Devi m_{l'}^{j'} & = \sum_l \left[J_{kl}^{i*} \right]_i^2 \sum_j \Devi m_l^j \sum_{\substack{j' \\ j \neq j'}} \Devi m_l^{j'} \\ & \overset{(1)} = \sum_l \left( \frac{J_{kl}}{\sqrt K} \frac{K}{N_l} \right)^2 \sum_j \Devi m_l^j \left( \sum_{j'} \Devi m_l^{j'} - \Devi m_l^j \right) \\ & = K \sum_l (J_{kl})^2 \frac{1}{N_l} \sum_j \Devi m_l^j \left( \frac{1}{N_l} \sum_{j'} \Devi m_l^{j'} - \frac{1}{N_l} \Devi m_l^j \right) \\ & \overset{(2)} = K \sum_l (J_{kl})^2 \left( [\Devi m_l^j]_j [\Devi m_l^{j'}]_{j'} - \frac{1}{N_l} [(\Devi m_l^j)^2]_j \right) \\ & \overset{(3)} = O(K/N)

ここで, (1) 結合確率の定義 による期待値の計算, (2) 集団平均の定義 \PAvg{\bullet}_j = \sum_j \bullet / N_l, (3) [\Devi m_l^j]_j = 0 を用いた

\sum_l \sum_j \left[ (J_{kl}^{ij})^2 \right]_i (\Devi m_l^j)^2 & \overset{(1)} = \sum_l \left[ (J_{kl}^{i*})^2 \right]_i N_l \left[ (\Devi m_l^j)^2 \right]_j \\ & \overset{(2)} = \sum_l \left( \frac{J_{kl}}{\sqrt K} \right)^2 \frac{K}{N_l} N_l \left[ (\Devi m_l^j)^2 \right]_j \\ & = \sum_l J_{kl}^2 \left[ (\Devi m_l^j)^2 \right]_j \\ & \overset{(3)} = \sum_l J_{kl}^2 \left( [(m_l^j)^2] - [m_l^j]^2 \right) \\ & \overset{(4)} = \sum_l J_{kl}^2 \left( q_l - m_l^2 \right)

ここで, (1) \left[(J_{kl}^{ij})^2 \right]_i が j に依存しないこと, (2) 大数の法則 (law of large numbers) と 結合確率の定義 による期待値の計算, (3) [(\Devi x)^2] = [x^2] - [x]^2, (4) q_k と m_k の定義を用いた.

合計

\left[ \left( \Devi \Avg{u_k^i(t)}_t \right)^2 \right] & \approx \left[ \text{(d1)}^2 \right] + \left[ \text{(d2)}^2 \right] \\ & \approx \sum_j J_{kl}^2 \left(1 - \frac K N_l \right) \left( m_l \right)^2 + \sum_l J_{kl}^2 \left( q_l - m_l^2 \right) \\ & = \sum_l J_{kl}^2 \, q_l + O(N_l^{-1}) \\ & \xrightarrow{N \to \infty} \sum_l J_{kl}^2 \, q_l

これで, クエンチされたゆらぎが式 (17) で表されることが示された.

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

時間ゆらぎの計算

ここでは 時間ゆらぎ (temporal fluctuations) が

(18)\Avg{ \left[\left( u_k^i(t) - \Avg{u_k^i(t)}_t \right)^2 \right]_i }_{t} \xrightarrow{N \to \infty} \sum_l J_{kl}^2 (m_l - q_l) = \alpha_k - \beta_k

となることを示す. [1]

[1]外側の時間平均 \Avg{\bullet}_t は [vanVreeswijk1998] には無いが,これがあるほうが定義としては自然だろう. 外側の時間平均 \Avg{\bullet}_t が無くても計算結果は変わらないことは, 時間平均なしでも時間ゆらぎは計算できる で示す.

まずは \Avg{\bullet}_t の中身を計算する.

\left[\left( u_k^i(t) - \Avg{u_k^i(t)}_t \right)^2 \right]_i & \overset{(1)} = \left[ \sum_{ll'jj'} J_{kl}^{ij} J_{kl'}^{ij'} (\sigma_l^j(t) - m_l^j) (\sigma_{l'}^{j'}(t) - m_{l'}^{j'}) \right]_i \\ & \overset{(2)} = \sum_{ll'jj'} \left[ J_{kl}^{ij} J_{kl'}^{ij'} \right]_i (\sigma_l^j(t) - m_l^j) (\sigma_{l'}^{j'}(t) - m_{l'}^{j'}) \\ & \overset{(3)} \approx \sum_{lj} \left[ J_{kl}^{ij} J_{kl}^{ij} \right]_i (\sigma_l^j(t) - m_l^j)^2 \\ & \overset{(4)} = \sum_l J_{kl}^2 \left(1 - \frac K N_l \right) \underbrace{ \frac 1 N_l \sum_j (\sigma_l^j(t) - m_l^j)^2 }_{(*)}

ここで, (1) 和の積の計算のための添字テクニック, (2) ニューロンの状態と結合係数の相関 の議論, (3) \sigma_l^j(t) と \sigma_{l'}^{j'}(t) が無相関なので, 時間平均活動率のゆらぎ と同様に, (j, l) = (j', l') の項以外は O(1/N) の大きさしかもたないこと, (4) 大数の法則 (law of large numbers) と 結合確率の定義 による期待値の計算, を用いた.

この式を時間平均 \Avg{\bullet}_t したものを計算したい. 時間に依存するのは (*) の部分のみなので, これの時間平均をとる.

\Avg{(*)}_t & = \Avg{ \left[ (\sigma_l^j(t) - m_l^j)^2 \right]_j }_t \\ & \overset{(1)} = \Avg{ \left[ \left( \sigma_l^j(t) \right)^2 \right]_j }_t - \left[ (m_l^j)^2 \right]_j \\ & \overset{(2)} = \Avg{\PAvg{ \sigma_l^j(t) }}_t - \left[ (m_l^j)^2 \right]_j \\ & \overset{(3)} = m_l - q_l

ここで, (1) [(\Devi x)^2] = [x^2] - [x]^2 と同種の計算, (2) \sigma_l^j(t) の取りうる値は 0 か 1 なので \left( \sigma_l^j(t) \right)^2 = \sigma_l^j(t), (3) m_k と q_k の定義, を用いた.

以上の計算を統合すると,

\Avg{ \left[\left( u_k^i(t) - \Avg{u_k^i(t)}_t \right)^2 \right]_i }_t & \approx \sum_l J_{kl}^2 \left(1 - \frac K N_l \right) (m_l - q_l)

が導かれる. これの左辺の極限 N \to \infty をとると式 (18) が導かれる.

時間平均なしでも時間ゆらぎは計算できる

時間に依存する部分 (*) を, 時間平均なしで計算すると

(*) & = \PAvg{ (\sigma_l^j(t) - m_l^j)^2 }_j \\ & = \PAvg{ \left( \sigma_l^j(t) \right)^2 }_j - 2 \PAvg{ \sigma_l^j(t) \, m_l^j }_j + \PAvg{ (m_l^j)^2 }_j \\ & = \PAvg{ \left( \sigma_l^j(t) \right)^2 }_j - \PAvg{ (m_l^j)^2 }_j - 2 \PAvg{ (\sigma_l^j(t) - m_l^j) \, m_l^j }_j \\ & = m_k(t) - q_l + 2 \PAvg{ \zeta_l^j(t) \, m_l^j }_j

となる. ただし, \zeta_l^j(t) = \sigma_l^j(t) - m_l^j とおいた.

この (*) が定常状態では \Avg{(*)}_t と等しいことが示せる. まず, 第一項は m_k(t) = m_k である. よって, 第三項がゼロになることが示せれば良い. 無相関性の「証明」 の議論から, \sigma_l^j(t) つまり \zeta_l^j(t) = \sigma_l^j(t) - m_l^j が違うニューロン同士では無相関であることに注意すれば, 無相関変数に対する大数の法則 より,

& \PAvg{ \zeta_l^j(t) \, m_l^j }_j \approx \AvgJ{\AvgDyn{ \zeta_l^j(t) \, m_l^j }} = (\star)

となる. これは,

(\star) & = \Avg{ \AvgJ{m_l^j(t) \, m_l^j} - \AvgJ{m_l^j(s) \, m_l^j} }_s

と書ける. \AvgJ{m_l^j(t) \, m_l^j} と \AvgJ{m_l^j(s) \, m_l^j} が同じ (集団レベルでの, あるいはグローバルな) 量を違う時間で評価したものであることを考慮すれば, 定常状態ではこの量は時間に依らないはずであり, (\star) = 0 が導ける.

もっと形式的にこれを導出するためには, 自己相関関数

C(t, r) & := \AvgJ{(m_l^j(t) - m_l(t)) (m_l^j(r) - m_l(r))} \\ & = \AvgJ{m_l^j(t) \, m_l^j(r)} - m_l(t) \, m_l(r).

を使う. これを用いれば,

& (\star) = ... \\ & = \Avg{ \AvgJ{m_l^j(t) \, m_l^j} - \AvgJ{m_l^j(s) \, m_l^j} }_s \\ & = \Avg{ \AvgJ{m_l^j(t) \, m_l^j(r)} - \AvgJ{m_l^j(s) \, m_l^j(r)} }_{s,r} \\ & = \Avg{ C(t, r) + m_l(t) \, m_l(r) - C(s, r) - m_l(s) \, m_l(r) }_{s,r} \\ & \overset{(1)} = \Avg{ C(t - r) - C(s - r) }_{s,r} \\ & = \Avg{\Avg{C(t - r)}_r - \Avg{C(s - r)}_r}_s \\ & \overset{(2)} = \Avg{\Avg{C(r)}_r - \Avg{C(r)}_r}_s = 0

となることが導ける. ただし, (1) で系が定常状態にあると自己相関関数が時間差のみに依存すること C(t ,r) = C(t - r) を用い, (2) で時間平均が時間シフトの元で不変であることを用いた

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

秩序変数の計算

これまでの計算で, クエンチされたゆらぎ(空間ゆらぎ)と時間ゆらぎを集団平均活動率 m_k と秩序変数 q_k で表すことが出来ると分かったが, 秩序変数 q_k の計算方法がまだ分からないので, これでは答えを得たとは言えない. ここでは, 秩序変数 q_k が満たすべき関係式 (自己無頓着 (self-consistent) 方程式) を導く.

平均場方程式の導出 と同様に, ニューロンの状態 \sigma_k^i(t) = \Theta(u_k^i(t)) を標準ガウス確率変数で書きなおそう. 今度は, 入力 u_k^i(t) のゆらぎをクエンチされたゆらぎ(空間ゆらぎ)と時間ゆらぎをに分け, それぞれ独立な標準ガウス確率変数 x_i と y_i(t) で表す. [1] つまり, 確率変数 x_i と y_i(t) は \PAvg{x_i}_i = \Avg{y_i(t)}_t = 0 と \PAvg{(x_i)^2}_i = \Avg{(y_i(t))^2}_t = 1 を満たすとする. これらの確率変数を用いて入力を

u_k^i(t) = u_k + c x_i + d y_i(t)

と書いて実数 c, d を求めよう. クエンチされたゆらぎの計算 の結果と

\left[ \left( \Devi \Avg{u_k^i(t)}_t \right)^2 \right] = \left[ \left( \Devi (c x_i) \right)^2 \right] = \left[ \left( c x_i \right)^2 \right] = c^2

より, c = \sqrt{\beta_k}, 時間ゆらぎの計算 の結果と

\Avg{ \left[\left( u_k^i(t) - \Avg{u_k^i(t)}_t \right)^2 \right]_i }_{t} = \Avg{ \left[\left( d y_i(t) \right)^2 \right]_i }_{t} = d^2

より, d = \sqrt{\alpha_k - \beta_k} が言える.

[1]確率変数 x_i とある時間 t における時間ゆらぎ成分を表す確率変数 y_i(t) は独立だが, 違う時間 t' \neq t と比べて, y_i(t) と y_i(t') が独立, という意味では ない. 確率変数 y_i(t) と y_i(t') はもちろん相関を持ち, その相関構造は u_k^i(t) の自己相関関数を計算することで理解できる.

よって, ニューロンの状態 \sigma_k^i(t) = \Theta(u_k^i(t)) は

\sigma_k^i(t) = \Theta \left( u_k + \sqrt{\beta_k} \, x_i + \sqrt{\alpha_k - \beta_k} \, y_i(t) \right)

と書けることが分かった. この表式では, 時間平均 \Avg{\bullet}_t は確率変数 y_i(t) に関する平均と同値 (つまり, \Avg{f(y_i(t))}_t = \int \D y \, f(y)) なので, m_k^i = \Avg{\AvgDyn{\sigma_k^i(t)}}_t は

(19)m_k^i = m_k(x_i) = H \left( \frac{-u_k - \sqrt{\beta_k} x_i}{\sqrt{\alpha_k - \beta_k}} \right)

と書ける. この表式を用いて q_k = \PAvg{(m_k^i)^2} を計算すると,

(20)q_k = \int Dx \left( m_k(x) \right)^2 = \int Dx \left\{ H \left( \frac{-u_k - \sqrt{\beta_k} x_i}{\sqrt{\alpha_k - \beta_k}} \right) \right\}^2

となる. クエンチされたゆらぎ \beta_k が q_k に依存していることを思い出せば, この式は秩序変数 q_k が満たすべき関係式であり, q_k を陰に定義していることが分かる.

  • 時間ゆらぎが無い場合 \alpha_k - \beta_k = 0:

    \sigma_k^i(t) = \Theta \left(u_k + \sqrt{\beta_k} \, x_i \right) は y_i(t) に依存しないので, m_k^i = \Avg{\AvgDyn{\sigma_k^i(t)}}_t = \Theta \left(u_k + \sqrt{\beta_k} \, x_i \right) [2] となるから, (m_k^i)^2 = m_k^i より,

    q_k = \int Dx \, \left( m_k(x) \right)^2 = \int Dx \, m_k(x) = m_k

    である. 時間ゆらぎの計算 で導いた関係式より, \alpha_k - \beta_k = \sum_{l=E,I} J_{kl}^2 (m_l - q_l) = 0 となり, 時間ゆらぎが無いという仮定と整合性があるので, q_k = m_k は式 (20) の解のひとつである. この解を, 原著 [vanVreeswijk1998] にならい 凍結解 (frozen solution) と呼ぶ. 時間ゆらぎは二乗の平均なので正, つまり \sum_{l=E,I} J_{kl}^2 (m_l - q_l) は正である.

    また, 0 \le m_k(x_i) \le 1 より, 各点 x_i で m_k(x_i)^2 \le m_k(x_i) だから, \PAvg{m_k(x_i)^2} \le \PAvg{m_k(x_i)} = m_k, つまり凍結解 q_k = m_k は秩序変数 q_k の上限を与えることが分かる.

    [2]

    この場合の時間平均活動率 m_k^i = \Theta \left(u_k + \sqrt{\beta_k} \, x_i \right) は, \alpha_k - \beta_k > 0 の場合の Q関数 を用いた表式 (19) の極限 \alpha_k - \beta_k \to 0^+ でもある.

  • クエンチされたゆらぎがない場合 \beta_k = 0:

    式 (20) に \beta_k = 0 を代入すると, m_k(x_i) = H({-u_k}/{\sqrt{\alpha_k}}) は x_i に依存しなくなり,

    q_k = \int Dx \left\{ H \left( \frac{-u_k}{\sqrt{\alpha_k}} \right) \right\}^2 = \underbrace{ \left\{ H \left( \frac{-u_k}{\sqrt{\alpha_k}} \right) \right\}^2 }_{= (m_k)^2} \underbrace{ \int Dx }_{= 1} = (m_k)^2

    となる. しかし, クエンチされたゆらぎの計算 で求めた関係式 \beta_k = \sum_{l=E,I} J_{kl}^2 q_l にこれをあてはめると, m_k = 0 または J_{kl} = 0 というトリビアルな状況を除けば, \beta_k > 0 となり, 仮定 \beta_k = 0 とは整合性がとれない. よって, q_k = (m_k)^2 は解ではない. 一方で, イェンゼンの不等式 (Jensen’s inequality) [3] を用いれば, \PAvg{(m_k^i)^2} \ge (\PAvg{m_k^i})^2 = (m_k)^2, つまりクエンチされたゆらぎがない場合 \beta_k = 0 が秩序変数 q_k の下限を与えることが分かる.

    [3]

    (下に)凸関数 f(x) と x に関する平均 \Avg{\bullet} について, f(\Avg{x}) \le \Avg{f(x)} が成り立つ. これを, イェンゼンの不等式 (Jensen’s inequality) という. 参考: イェンゼンの不等式 - Wikipedia / Jensen’s inequality - Wikipedia
\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

時間平均活動率の分布

集団 k の時間平均活動率の分布は

\rho_k(m) := \lim_{N \to \infty} \frac{1}{N_k} \sum_{i=1}^{N_k} \delta(m - m_k^i) = \lim_{N \to \infty} \PAvg{\delta(m - m_k^i)}_i

で定義される. 秩序変数の計算 の表式を用いれば, \PAvg{f(x_i)}_i \approx \int \D x \, f(x) だから,

\rho_k(m) = \int Dx \, \delta(m - m_k(x))

で計算できる.

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

安定性解析

課題

書く: 安定性解析

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

記号表

課題

記号表を埋める

記号 定義等 節 原著 [vanVreeswijk1998] 注
\PAvg{\bullet}_i := N_l^{-1} \sum_i \bullet / 集団平均 / population average 入力のゆらぎ Sec. 3, p.1327 (p.7)  
\PAvg{\bullet} = \PAvg{\bullet}_i (添字が自明な場合)      
\Devi \Devi X_i := X_i - [X_i]_i / 偏差 / deviation 入力のゆらぎ Sec. 3, p.1329 (p.9) *2
\AvgDyn{\bullet} 更新時間のランダム性に関する平均 平均場方程式の導出 A.1, p.1365 (p.45) *2
m_k(t) := \PAvg{\AvgDyn{\sigma_k^i(t)}}_i time dependent local rate variable   Eq. 3.2, p.1327 (p.7) *1
m_k := \Avg{m_k(t)}_t   Seq. 3, p.1328 (p.8) (See also Eq. 3.7) *3
m_k^i := \Avg{\AvgDyn{\sigma_k^i(t)}}_t クエンチされたゆらぎの計算 Eq. 5.3, p.1334 (p.14) *1
m_k^i (t) := \AvgDyn{\sigma_k^i(t)} 平均場方程式の導出 Eq. A.1, p.1365 (p.45)  
q_k := \PAvg{(m_k^i)^2}_i / オーダーパラメター / order parameter クエンチされたゆらぎの計算 Eq. 5.3, p.1334 (p.14)  

  • *1 原著 [vanVreeswijk1998] では, \AvgDyn{\bullet} をとならないものとして定義されていた.

    課題

    \AvgDyn{\bullet} をつけるほうが正しいのか検証

  • *2 原著 [vanVreeswijk1998] とは若干違う記法を用いている.

  • *3 原著 [vanVreeswijk1998] では, m_k は m_k(t) のダイナミクスの固定点 (fixed point) として導入されているが, (もし固定点が存在すれば) 時間平均 \Avg{m_k(t)}_t はそれと同値である.

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

付録

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

漸近解析 (Asymptotic analysis)

記号   定義 意味
f = O(g) f \lesssim g \limsup_{x \to \infty} |f(x)/g(x)| < \infty f の増加は g より遅いか同じ
f = o(g) f \ll g \lim_{x \to \infty} f(x)/g(x) = 0 f の増加は g より遅い
— f \sim g \lim_{x \to \infty} f(x)/g(x) = 1 f が g に漸近する
f = \Theta(g) f \asymp g f \lesssim g \lesssim f f と g は同速度で増加
f = \Omega(g) f \gtrsim g g = O(f) f の増加は g より速いか同じ
f = \omega(g) f \gg g g = o(f) f の増加は g より速い
  • どの変数 (例: x) がどう変化 (例: x \to \infty) するかを書かなければ厳密には漸近関係とは呼ばない. 例えば 「f(x) = O(g(x))」 は漸近関係ではなく, 正しくは「f(x) = O(g(x)) (as x \to \infty) 」などと書く必要がある.
  • 表の定義ではすべて “as x \to \infty” での漸近関係だが, 変数の収束先は何でも良い.
  • 記号 O, o, \Theta, \Omega, \omega は ランダウの記号 (Landau symbol) や O-記法 (big O notation) と呼ばれる. 記号 \lesssim, \ll, \sim は Vinogradov notation と呼ばれる.
  • f \ll g を f = O(g) と定義するか, f = o(g) と定義するかは文献に依るようだ. ここでは, f \ll g は f = o(g) で定義し, f = O(g) と同値の二項関係は f \lesssim g と書く. しかし, この記法は慣例的なものではない.
  • f \sim g を f \lesssim g \lesssim f で定義する流儀もある [1] ようだが, ここでの定義とは違う.
  • 記号 \approx や \lessapprox は漸近関係では無く荒い議論で発見的な説明するときなどに使う.
[1]254A, Notes 0: A review of probability theory | What’s new

漸近関係

漸近関係の例

極限 x \to \infty での漸近関係をいくつか示す. 末尾の “(as x \to \infty)” は省略する.

x = O(5 x^2)

x = o(4 x^2)

x = \Theta(3 x)

x^3 + x \sim x^3 + \log x

ランダウの記号

ランダウの記号が関数 q の中に現れ

f(x) = q(O(g(x))) \quad \text{as } x \to \infty

と書かれた漸近関係は, 「ある関数 z(x) が 存在して f(x) \sim q(z(g(x))) (as x \to \infty)」と定義される. また, ランダウの記号が左辺の関数 p の中に現れ

p(O(f(x))) = g(x) \quad \text{as } x \to \infty

と書かれた漸近関係は, 「いかなる関数 y(x) についても, p(y(x)) \sim g(x) (as x \to \infty)」と定義される. 他のランダウの記号についても同様である. ランダウの記号が左辺にあらわれた場合と右辺にあらわれた場合の意味が異なることに注意せよ. 例えば, x = O(x) (as x \to \infty) は正しいが, O(x) = x (as x \to \infty) は正しくない.

上記のふたつの記法を組み合わせて使えば, ランダウの記号の「強さ」の関係を完結に表すことが出来る [2]:

o(f(x)) = O(f(x)), \enskip O(f(x)) = \Omega(f(x)), \enskip \Omega(f(x)) = \omega(f(x)), \Theta(f(x)) = O(f(x)), \Theta(f(x)) = \Omega(f(x)) .

[2]さらに記法を乱用すると o(f) = O(f) = \Omega(f) = \omega(f) とも書ける.

支配項均衡の方法

支配項均衡の方法 (method of dominant balance) とは, 3項以上がからんだ方程式のある少数の項 (例えば2項だけ) が他の項に比べて非常に大きい場合に, その少数の項だけの方程式とみなして (漸近的な) 解を求める方法である.

二次方程式の解の公式を知らないふりをして

x^2 + C x + 1 = 0

の C \to \infty についての漸近的な解を求めよう. 各項に a = x^2, b = C x, c = 1 と名前をつける. まず項 b = C x がどのオーダーを持つか仮定し, それに基づいて計算した他の項のオーダーを使って一番大きなオーダーが釣り合うか, つまり同じ大きさのオーダーがふたつあるか, を確かめる. 釣り合いが示せたならば, 釣り合った項を漸近関係 \sim で結んだ式を解く. 例えば, 項 b と c が釣り合い, 項 a がそれらに比べて小さい, つまり a \ll b \asymp c なら, b \sim - c という式を解けばいい.

  1. 仮定: C x \asymp 1

    仮定より x \asymp 1 / C だから, x^2 \asymp (1 / C)^2 となる. 項 b と c が釣り合い, 項 a はそれらに比べて小さいので無視できる. 方程式は C x \sim -1 と書くことが出来る. 解は x \sim - 1 / C である.

  2. 仮定: C x \asymp C

    仮定より x \asymp 1 だから, x^2 \asymp 1 となる. 項 a, c はともに 1 のオーダーであり, 項 b \asymp C とは釣り合わない.

  3. 仮定: C x \asymp C^2

    仮定より x \asymp C だから, x^2 \asymp C^2 となる. 項 a, b はともに C^2 のオーダーで, 項 c はそれらに比べて小さいので無視できる. 方程式は x^2 \sim - C x と書くことが出来る. 解は 0 か -C が考えられるが, x = 0 は x \asymp C を満たさないので, x \sim -C が解である.

以上の考察から, 解は x \sim -1/C, -C であることが分かった. 二次方程式の解の公式にあてはめると, 厳密解は

x_{\pm} = \frac{-C \pm \sqrt{C^2 - 4}}{2}

となる. この厳密解が支配項均衡の方法で求めた漸近解に漸近し x_- \sim -C, x_+ \sim -1/C となることは,

\lim_{C \to \infty} \frac{x_-}{-C} = \lim_{C \to \infty} \frac{-1 - \sqrt{1 - 4 / C^2}}{-2} = 1

\lim_{C \to \infty} \frac{x_+}{-1/C} = \lim_{C \to \infty} \frac{-1 + \sqrt{1 - 4 / C^2}}{-2 / C^2} \overset{(1)} = \lim_{\epsilon \to 0} \frac{1 - \sqrt{1 - 4 \epsilon}}{2 \epsilon} \overset{(2)} = \lim_{\epsilon \to 0} \frac{1}{\sqrt{1 - 4 \epsilon}} = 1

より確かめられる. ただし, (1) \epsilon = 1/C^2 なる変数変換と, (2) ロピタルの定理 を使った. 後者の計算から分かるように, x_+ のふるまいは, 厳密解がわかっていても自明では無い. しかし, 支配項均衡の方法を用いれば簡単に漸近解を導ける.

参考

Asymptotic analysis - Wikipedia
このページの Method of dominant balance という項目では, 常微分方程式を漸近解析を使って解く例が紹介されている.
Perturbation & Asymptotic Series - YouTube
摂動法と漸近解析の講義ビデオ. 一部の応用例の紹介で量子力学などの物理学の知識が必要な部分があるが, メインの流れは学部1年生のレベルの数学さえ理解していれば追いかけられるだろう. 理論科学一般に興味のある人には大変おすすめである.
\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

確率変数の収束

参考

確率変数の収束 - Wikipedia

Convergence of random variables - Wikipedia

大数の法則 (law of large numbers)

課題

LNN について書く

大数の法則 (law of large numbers)

\lim_{N \to \infty} \frac 1 N \sum_{i=1}^{N} X_i = \mu

参考

大数の法則 - Wikipedia

Law of large numbers - Wikipedia

無相関変数に対する大数の法則

大数の法則が成り立つためには確率変数が同一の分布の従う条件や, 独立であるという条件は実は必要ではない. 確率変数の分散に「有界性」があり, (独立ではなく) 無相関であれば十分である. 本稿では, 以下の定理 ([Cinlar2011] の Proposition III.6.13 (p.121) より) を繰り返し用いる.

無相関変数に対する大数の法則

確率変数 X_n (n = 1, 2, \ldots) が無相関で, ある非ゼロで正の増加列 b_1, b_2, \ldots を用いて \sum_{n=1}^\infty \Var (X_n / b_n) < \infty とできるならば, (\sum_{m=1}^n X_n - \Expect \sum_{m=1}^n X_n) / b_n \to 0 が L^2 収束と確率収束の意味で成り立つ.

本稿では確率変数の収束の細かい違いについて神経質に使い分けないが, 前者は \Expect |(\sum_{m=1}^n X_n - \Expect \sum_{m=1}^n X_n) / b_n|^2 \to 0, 後者は \Prob \{ (\sum_{m=1}^n X_n - \Expect \sum_{m=1}^n X_n) / b_n \to 0 \} = 1 を意味する.

中心極限定理 (central limit theorem)

課題

CLT について書く

参考

中心極限定理 - Wikipedia

Central limit theorem - Wikipedia

自己平均性 (self-averaging property)

確率変数 \bm X = (X_1, \ldots, X_N)^\intercal がすべての j = 1, \ldots, N について \lim_{N \to \infty} N^{-1} \sum_{i=1}^N X_i = \Avg{X_j} を満たす時, 確率変数 \bm X は自己平均的である (\bm X is self-averaging) と呼ばれる. ただし, \Avg{X_j} は 確率変数 \bm X についての平均で, “average over realisations” や “average over samples” と呼ばれる.

確率変数 \bm X は例えばネットワークの状態 \bm \sigma = (\sigma_1, \ldots, \sigma_N)^\intercal (モデルの定義 参照) である. この場合, 各 X_n は有界なので, b_n = n で 無相関変数に対する大数の法則 の「有界性」の条件 \sum_{n=1}^\infty \Var (X_n / b_n) < \infty を満たす. さらに, すべての確率変数の平均が同一 \Avg{X_j} = \mu (j = 1, 2, \ldots) ならば, この定理より, \lim_{N \to \infty} N^{-1} \sum_{i=1}^N X_i = \mu が言える.

参考

Self-averaging - Wikipedia

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

ガウス分布 (Gaussian distribution)

参考

正規分布 - Wikipedia

ガウス確率変数の変数変換

確率変数 X が平均 \mu, 分散 \sigma^2 のガウス分布 に従い, 確率変数 Z が平均 0, 分散 1 のガウス分布 に従うとする. これらの確率変数に関する期待値を \Avg{\bullet} と書く. 統計量 \Avg{f(X)} (f は適当な関数) は Z を用いて,

\Avg{f(X)} = \Avg{f(\mu + \sigma Z)}

と計算出来る.

形式的な変数変換で示そう. 平均 \mu, 分散 \sigma^2 のガウス確率密度関数を

g_{\mu,\sigma}(x) = \frac{1}{\sqrt{2 \pi \sigma^2}} \exp \left( - \frac{(x - \mu)^2}{2 \sigma^2} \right)

と書く. 変数変換 x = \mu + \sigma z をすれば,

\Avg{f(X)} & = \int_{-\infty}^\infty \D x \, g_{\mu,\sigma}(x) \, f(x) \\ & = \int_{-\infty}^\infty \D z \, \sigma \, g_{\mu,\sigma}(\mu + \sigma z) \, f(\mu + \sigma z) \\ & = \int_{-\infty}^\infty \D z \, \sigma \, \frac{1}{\sqrt{2 \pi \sigma^2}} \exp \left( - \frac{(\mu + z \sigma - \mu)^2}{2 \sigma^2} \right) f(\mu + \sigma z) \\ & = \int_{-\infty}^\infty \D z \, g_{0,1}(z) \, f(\mu + \sigma z) \\ & = \Avg{f(\mu + \sigma Z)}

ガウス測度 (Gaussian measure)

物理の計算では, いくつものガウス積分が出てくることが多々あり, その度に \D z \, g_{0,1}(z) などと書くのは煩雑なので, ガウス測度 (Gaussian measure) と呼ばれる以下の記法を導入する.

Dx := \D x \exp(-x^2) / \sqrt{2 \pi}

これを用いれば, ガウス確率変数の変数変換 の公式は

\Avg{f(X)} = \int_{-\infty}^\infty Dz \, f(\mu + \sigma z)

とかける.

注釈

ただの省略のための記法なので, 数学の測度論 (measure theory) とは深い関わりなど無いし, そもそも測度論の測度は集合をとる関数なので別物である (測度論で似た記法を 別の 場面で使う流儀はあるが...). 物理の分野では \exp(-x^2) / \sqrt{2 \pi} と繰り返し書くのが面倒なので使われているだけである.

Q関数

標準ガウス分布の裾 (tail) は Q関数 (Q-function) と呼ばれ, ガウス分布を取り扱う計算では頻出である.

Q(z) := \int_z^\infty \frac{\D x}{\sqrt{2 \pi}} \exp(- x^2 / 2)

興奮・抑制均衡のオリジナルの論文 [vanVreeswijk1998] ではこのQ関数は H(z) と表記され, 相補誤差関数 (complementary error function, erfc) と呼ばれているが, 広く使われている定義

\text{erfc}(x) := \frac{2}{\sqrt \pi} \int_x^\infty \D t \, \exp(-t^2)

からは, 少し違う. erfc と Q関数 の間には,

Q(z) = \frac{1}{2} \text{erfc}(z/\sqrt 2)

の関係がある.

参考

Q-function - Wikipedia

Error function - Wikipedia

誤差関数 - Wikipedia

ヘヴィサイド関数の多重ガウス積分とQ関数

ヘヴィサイド関数 \Theta(x) = 0 (x < 0), \Theta(x) = 1 (x > 0), の多重ガウス積分とQ関数には次の関係が成り立つ.

\int \left( \prod_i Dx_i \right) \, \Theta \left( \sum_i c_i x_i + u \right) = Q \left( - \frac{u}{\sqrt{\sum_i c_i^2}} \right)

ただし, 積と和は i=1 から i=M までについておこない, 積分区間は \mathbb R^M 全体とする.

まずは M = 2 の場合について示す. a = c_1, b = c_2, \alpha = a^2 + b^2 とする.

二次元ベクトル (a, b)^\intercal に沿った単位ベクトルを \bm r = (a, b)^\intercal / \sqrt \alpha と定義する. [1] これに直交した2つある単位ベクトル \bm s = (s_1, s_2)^\intercal のうち, 行列式が \det (\bm r, \bm s) = 1 となるものを選ぶ. [2]

[1]x_1 と x_2 に関するガウス積分を (a, b)^\intercal の軸に沿ったガウス積分とそれに直交する軸に沿ったガウス積分に分ける, というアイディアである.
[2]行列式が \det (\bm r, \bm s) = -1 となるものを選んでも同様に計算できる.

これらの仮定, つまり,

\bm r \cdot \bm s \propto a s_1 + b s_2 = 0 \det (\bm r, \bm s) = (a s_2 - b s_1) / \sqrt \alpha = 1

より, [3] \bm s = (-b, a)^\intercal / \sqrt \alpha

[3]第一式より a s_1 = - b s_2, 第二式の両辺に a をかけて (a^2 s_2 - b (a s_1)) / \sqrt \alpha = a なので, s_2 = a / \sqrt \alpha. 第一式より, s_1 = - b / \sqrt \alpha である.

変数変換

r = \bm r \cdot (x_1, x_2)^\intercal = (a x_1 + b x_2) / \sqrt \alpha s = \bm s \cdot (x_1, x_2)^\intercal = (-b x_1 + a x_2) / \sqrt \alpha

を施す. ヤコビアンは,

\det \frac{\partial (x, y)}{\partial (r, s)} = \left( \det \frac{\partial (r, s)}{\partial (x, y)} \right)^{-1} = \frac{1}{\det (\bm r, \bm s)} = 1

なので, \D x \D y = \D r \D s である. さらに, r^2 + s^2 = x_1^2 + x_2^2 なので, 結局, Dx_1 Dx_2 = Dr Ds である.

& \int Dx_1 Dx_2 \, \Theta (a x_1 + b y_2 + u) \\ & = \int Dr \underbrace{\int Ds}_{=1} \, \Theta (\sqrt \alpha r + u) \\ & = \int Dr \, \Theta (r + u / \sqrt \alpha) = \int_{-u / \sqrt \alpha}^\infty Dr \\ & = Q \left( - \frac{u}{\sqrt \alpha} \right)

一般の M に関して成り立つことは, 帰納法で確かめられる.

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

ポアソン分布 (Poisson distribution)

参考

ポアソン分布 - Wikipedia

Poisson distribution - Wikipedia

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

マスター方程式

N 個の二状態をとるニューロン (binary neuron) が相互に作用している確率的な系を考える. ニューロン i の状態は 0 または 1 をとる状態変数 \sigma_i で表され, 系全体の状態はベクトル \bm \sigma = (\sigma_1, \ldots \sigma_N)^\intercal で表される (つまり, \bm \sigma \in \{0, 1\}^N). ここではこれらの状態変数の連続時間 t \in \mathbb R に沿った変化を考えるが, 時間に関する依存性を表記上は省略して \bm \sigma(t) を \bm \sigma と書く.

表記の準備をする. i 番目のニューロンの状態を反転した状態ベクトルを \hat{\bm{\sigma}}^{i}, i 番目の成分を除いた状態ベクトルを \bm{\sigma}^{\setminus i} (つまり, \bm{\sigma}^{\setminus i} \in \{0, 1\}^{N-1}) と書くことにする. 形式的にかけば,

\hat{\bm{\sigma}}^{i} = (\sigma_1, \ldots, 1 - \sigma_i, \ldots, \sigma_N)^\intercal \bm{\sigma}^{\setminus i} = (\sigma_1, \ldots, \sigma_{i-1}, \sigma_{i+1}, \ldots, \sigma_N)^\intercal

である. ニューロン i 以外の状態 \bm{\sigma}^{\setminus i} が与えられた時に, ニューロン i が状態を \sigma_i' から \sigma_i に遷移させる単位時間あたりの確率を w(\sigma_i | \sigma_i', \bm{\sigma}^{\setminus i}) で表す. つまり, \Delta t 時間の間にこの遷移を起こす確率は,

P_{\Delta t}(\sigma_i | \sigma_i', \bm{\sigma}^{\setminus i}) = w(\sigma_i | \sigma_i', \bm{\sigma}^{\setminus i}) \, \Delta t

である. ここでは, 系が時間 t に 状態 \bm \sigma をとる確率 P_t(\bm \sigma) の時間発展が, マスター方程式 (master equation) とよばれる次の微分方程式 [1] で与えられることを示す.

マスター方程式

(21)\frac{\D P_t(\bm \sigma)}{\D t} = - \sum_{i=1}^N w(1 - \sigma_i | \sigma_i, \bm{\sigma}^{\setminus i}) P_t(\bm \sigma) + \sum_{i=1}^N w(\sigma_i | 1 - \sigma_i, \bm{\sigma}^{\setminus i}) P_t(\hat{\bm{\sigma}}^{i})

[1]状態 \bm \sigma が 2^N (= \# \{0, 1\}^N) 個あることを思い出せば, これは 2^N 次元常微分方程式と考えても良い.

マスター方程式の導出

時間微分 {\D P_t(\bm \sigma)}/{\D t} = \lim_{\Delta t \to 0} (P_{t+\Delta t}(\bm \sigma) - P_t(\bm \sigma))/\Delta t を, 時間 t と t+\Delta t の確率をつなぐ関係式

P_{t+\Delta t}(\bm \sigma) = \sum_{\bm \sigma'} P_{\Delta t}(\bm \sigma | \bm \sigma') P_t(\bm \sigma')

を用いて評価する. ここで, P_{\Delta t}(\bm \sigma | \bm \sigma') は時間 \Delta t の間に系の状態が \bm \sigma' から \bm \sigma へ変化する確率で推移確率 (transition probability)と呼ばれる. これは, 以下のように計算出来る.

& P_{\Delta t}(\bm \sigma | \bm \sigma') \\ & = \prod_{i=1}^N P_{\Delta t}(\sigma_i | \bm{\sigma'}) = \prod_{i=1}^N P_{\Delta t}(\sigma_i | \sigma_i', \bm{\sigma'}^{\setminus i}) \\ & \overset{(1)} = \prod_{i=1}^N \left( P_{\Delta t}(\sigma_i | \sigma_i, \bm{\sigma'}^{\setminus i}) \, \delta_{\sigma_i, \sigma_i'} + P_{\Delta t}(\sigma_i | 1 - \sigma_i, \bm{\sigma'}^{\setminus i}) \, \delta_{1 - \sigma_i, \sigma_i'} \right) \\ & \overset{(2)} = \prod_{i=1}^N \left( (1 - w(1 - \sigma_i | \sigma_i, \bm{\sigma}^{\setminus i}) \, \Delta t) \, \delta_{\sigma_i, \sigma_i'} + w(\sigma_i | 1 - \sigma_i, \bm{\sigma}^{\setminus i}) \, \Delta t \, \delta_{1 - \sigma_i, \sigma_i'} \right) \\ & = \prod_{i=1}^N \left( \delta_{\sigma_i, \sigma_i'} + \left\{ - w(1 - \sigma_i | \sigma_i, \bm{\sigma}^{\setminus i}) \, \delta_{\sigma_i, \sigma_i'} + w(\sigma_i | 1 - \sigma_i, \bm{\sigma}^{\setminus i}) \, \delta_{1 - \sigma_i, \sigma_i'} \right\} \, \Delta t \right) \\ & \overset{(3)} = \underbrace{ \prod_{i=1}^N \delta_{\sigma_i, \sigma_i'} }_{O(1)} \\ & \qquad + \underbrace{ \Delta t \sum_{i=1}^N \left\{ - w(1 - \sigma_i | \sigma_i, \bm{\sigma}^{\setminus i}) \, \delta_{\sigma_i, \sigma_i'} + w(\sigma_i | 1 - \sigma_i, \bm{\sigma}^{\setminus i}) \, \delta_{1 - \sigma_i, \sigma_i'} \right\} \prod_{\substack{j=1 \\ j \neq i}}^N \delta_{\sigma_j, \sigma_j'} }_{O(\Delta t)} \\ & \qquad + O(\Delta t^2)

この式を用いて, P_{t+\Delta t}(\bm \sigma) を評価する:

P_{t+\Delta t}(\bm \sigma) & = \sum_{\bm \sigma'} P_{\Delta t}(\bm \sigma | \bm \sigma') P_t(\bm \sigma') \\ & = \underbrace{ \sum_{\bm \sigma'} P_t(\bm \sigma') \prod_{i=1}^N \delta_{\sigma_i, \sigma_i'} }_{P_t(\bm \sigma)} \\ & \qquad - \Delta t \sum_{i=1}^N w(1 - \sigma_i | \sigma_i, \bm{\sigma}^{\setminus i}) \underbrace{ \sum_{\bm \sigma'} P_t(\bm \sigma') \, \delta_{\sigma_i, \sigma_i'} \prod_{\substack{j=1 \\ j \neq i}}^N \delta_{\sigma_j, \sigma_j'} }_{P_t(\bm \sigma)} \\ & \qquad + \Delta t \sum_{i=1}^N w(\sigma_i | 1 - \sigma_i, \bm{\sigma}^{\setminus i}) \underbrace{ \sum_{\bm \sigma'} P_t(\bm \sigma') \, \delta_{1 - \sigma_i, \sigma_i'} \prod_{\substack{j=1 \\ j \neq i}}^N \delta_{\sigma_j, \sigma_j'} }_{P_t(\hat{\bm{\sigma}}^{i})} \\ & \qquad + O(\Delta t^2)

これで {\D P_t(\bm \sigma)}/{\D t} = \lim_{\Delta t \to 0} (P_{t+\Delta t}(\bm \sigma) - P_t(\bm \sigma))/\Delta t を評価する準備が出来た.

& \frac{P_{t+\Delta t}(\bm \sigma) - P_t(\bm \sigma)}{\Delta t} \\ & \qquad = \underbrace{ - \sum_{i=1}^N w(1 - \sigma_i | \sigma_i, \bm{\sigma}^{\setminus i}) P_t(\bm \sigma) + \sum_{i=1}^N w(\sigma_i | 1 - \sigma_i, \bm{\sigma}^{\setminus i}) P_t(\hat{\bm{\sigma}}^{i}) }_{= \D P_t(\bm \sigma)/\D t} + O(\Delta t)

この式の極限 \Delta t \to 0 で消えない第1項と第2項は確かに式 (21) の右辺である.

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

期待値の時間発展

参考

[Ginzburg1994]
この章の議論は, Ginzburg & Sompolinsky (1994) の Appendix A に基づいている.

期待値の時間発展

マスター方程式 と同じ設定の下, ニューロン i の状態の期待値 \Avg{\sigma_i(t)} = \sum_{\bm\sigma} \sigma_i P_t(\bm\sigma) は次の方程式に従う.

\tau \frac{\D}{\D t} \Avg{\sigma_i(t)} = - \Avg{\sigma_i(t)} + \Avg{g_i(\bm\sigma(t))}

ただし, g_i は状態 \bm\sigma に基づき遷移確率を与える関数で,

g_i(\bm\sigma(t)) = w(1 | \sigma_i', \bm{\sigma}^{\setminus i}) \, \tau \qquad (\sigma_i' = 0, 1)

と定義される. つまり, w(1 | \sigma_i', \bm{\sigma}^{\setminus i}) が \sigma_i' に依らない場合にのみ上式は成り立つ.

一般の場合の期待値の時間発展

期待値 \Avg{\sigma_i(t)} = \sum_{\bm\sigma} \sigma_i P_t(\bm\sigma) の時間微分をとると,

\frac{\D}{\D t} \Avg{\sigma_i(t)} & = \sum_{\bm\sigma} \sigma_i \frac{\D}{\D t} P_t(\bm\sigma) \\ & = \sum_{\bm\sigma} \sigma_i \left\{ - \sum_{j=1}^N w(1 - \sigma_j | \sigma_j, \bm{\sigma}^{\setminus j}) P_t(\bm \sigma) + \sum_{j=1}^N w(\sigma_j | 1 - \sigma_j, \bm{\sigma}^{\setminus j}) P_t(\hat{\bm{\sigma}}^{j}) \right\} \\ & = - \underbrace{ \sum_{\bm\sigma} \sum_{j=1}^N \sigma_i w(1 - \sigma_j | \sigma_j, \bm{\sigma}^{\setminus j}) P_t(\bm \sigma) }_{\text{(A)}} + \underbrace{ \sum_{\bm\sigma} \sum_{j=1}^N \sigma_i w(\sigma_j | 1 - \sigma_j, \bm{\sigma}^{\setminus j}) P_t(\hat{\bm{\sigma}}^{j}) }_{\text{(B)}}

\text{(B)} & = \sum_{\bm\sigma} \sum_{j=1}^N \overbrace{ (\delta_{ij} + 1 - \delta_{ij}) }^{= 1} \sigma_i w(\sigma_j | 1 - \sigma_j, \bm{\sigma}^{\setminus j}) P_t(\hat{\bm{\sigma}}^{j}) \\ & = \underbrace{ \sum_{\bm\sigma} \sum_{j=1}^N \delta_{ij} \sigma_i w(\sigma_j | 1 - \sigma_j, \bm{\sigma}^{\setminus j}) P_t(\hat{\bm{\sigma}}^{j}) }_{\text{(B1)}} \\ & \qquad + \underbrace{ \sum_{\bm\sigma} \sum_{j=1}^N (1 - \delta_{ij}) \sigma_i w(\sigma_j | 1 - \sigma_j, \bm{\sigma}^{\setminus j}) P_t(\hat{\bm{\sigma}}^{j}) }_{\text{(B2)}}

\text{(B1)} & \overset{(1)} = \sum_{\bm\sigma} \sum_{j=1}^N \delta_{ij} \sigma_j w(\sigma_j | 1 - \sigma_j, \bm{\sigma}^{\setminus j}) P_t(\hat{\bm{\sigma}}^{j}) \\ & \overset{(2)} = \sum_{\bm\sigma} \sum_{j=1}^N \delta_{ij} (1 - \sigma_j) w(1 - \sigma_j | \sigma_j, \bm{\sigma}^{\setminus j}) P_t(\bm{\sigma}) \\ & \overset{(3)} = \sum_{\bm\sigma} (1 - \sigma_i) w(1 - \sigma_i | \sigma_i, \bm{\sigma}^{\setminus i}) P_t(\bm{\sigma})

\text{(B2)} & = \sum_{\bm\sigma} \sum_{j=1}^N (1 - \delta_{ij}) \sigma_i w(1 - \sigma_j | \sigma_j, \bm{\sigma}^{\setminus j}) P_t(\bm{\sigma})

\text{(A)} + \text{(B2)} & = \sum_{\bm\sigma} \sum_{j=1}^N \sigma_i w(1 - \sigma_j | \sigma_j, \bm{\sigma}^{\setminus j}) P_t(\bm \sigma) \\ & \qquad + \sum_{\bm\sigma} \sum_{j=1}^N (1 - \delta_{ij}) \sigma_i w(1 - \sigma_j | \sigma_j, \bm{\sigma}^{\setminus j}) P_t(\bm{\sigma}) \\ & = - \sum_{\bm\sigma} \sum_{j=1}^N \delta_{ij} \sigma_i w(1 - \sigma_j | \sigma_j, \bm{\sigma}^{\setminus j}) P_t(\bm{\sigma}) \\ & = - \sum_{\bm\sigma} \sigma_i w(1 - \sigma_i | \sigma_i, \bm{\sigma}^{\setminus i}) P_t(\bm{\sigma})

\frac{\D}{\D t} \Avg{\sigma_i(t)} & = \text{(A)} + \text{(B1)} + \text{(B2)} \\ & = \sum_{\bm\sigma} (1 - 2 \sigma_i) w(1 - \sigma_i | \sigma_i, \bm{\sigma}^{\setminus i}) P_t(\bm{\sigma}) \\ & = \Avg{ (1 - 2 \sigma_i) w(1 - \sigma_i | \sigma_i, \bm{\sigma}^{\setminus i}) }(t)

遷移確率 w(1 | \sigma_i', \bm{\sigma}^{\setminus i}) が自己の状態 \sigma_i' に依らない場合

関数 g_i の定義 w(1 | \sigma_i', \bm{\sigma}^{\setminus i}) = g(\bm{\sigma}) / \tau (\sigma_i' = 0, 1) は

w(1 - \sigma_i | \sigma_i, \bm{\sigma}^{\setminus i}) = \frac{1}{2 \tau} \left\{ 1 - (2 \sigma_i - 1) [2 g_i(\bm{\sigma}) - 1] \right\}

と書くことが出来る. これは上式に \sigma_i = 0 と \sigma_i = 1 を代入することで確かめられる.

この式を用いて, \frac{\D}{\D t} \Avg{\sigma_i(t)} を計算すると,

\frac{\D}{\D t} \Avg{\sigma_i(t)} & = \Avg{ (1 - 2 \sigma_i) \frac{1}{2 \tau} \left\{ 1 - (2 \sigma_i - 1) [2 g_i(\bm{\sigma}) - 1] \right\} }(t) \\ & \overset{(1)} = \Avg{ \frac{1}{2 \tau} \left\{ (1 - 2 \sigma_i) + (1 - 2 \sigma_i)^2 [2 g_i(\bm{\sigma}) - 1] \right\} }(t) \\ & \overset{(2)} = \Avg{ \frac{1}{2 \tau} \left\{ (1 - 2 \sigma_i) + 2 g_i(\bm{\sigma}) + 1 \right\} }(t) \\ & = \Avg{ \frac{1}{\tau} \left\{ - \sigma_i + g_i(\bm{\sigma}) \right\} }(t) \\ & = \frac{1}{\tau} \left\{ - \Avg{\sigma_i}(t) + \Avg{g_i(\bm{\sigma})}(t) \right\}

相関関数

課題

相関関数の従う方程式 (“two-time” second moment の時間発展) を導出する.

\Avg{f(\bm \sigma(t), \bm \sigma(t+s))} := \sum_{\bm \sigma} P_t(\bm \sigma) \sum_{\bm \varsigma} P_{t+s, t}(\bm \varsigma | \bm \sigma) f(\bm \sigma, \bm \varsigma)

\CAvg{f(\bm \sigma(t+s))}{\bm \sigma(t)} := \sum_{\bm \varsigma} P_{t+s, t}(\bm \varsigma | \bm \sigma(t)) f(\bm \varsigma)

& \Avg{\sigma_i(t) \, \sigma_j(t+s)} \\ & = \sum_{\bm \sigma} P_t(\bm \sigma) \sum_{\bm \varsigma} P_{t+s, t}(\bm \varsigma | \bm \sigma) \, \sigma_i \, \varsigma_j \\ & = \sum_{\bm \sigma} P_t(\bm \sigma) \, \sigma_i \underbrace{ \sum_{\bm \varsigma} P_{t+s, t}(\bm \varsigma | \bm \sigma) \, \varsigma_j }_{= \CAvg{\sigma_j(t+s)}{\bm \sigma(t)}} \\ & = \Avg{\sigma_i(t) \CAvg{\sigma_j(t+s)}{\bm \sigma(t)}}

\tau \frac{\D}{\D s} \CAvg{\sigma_j(t+s)}{\bm \sigma(t)} = - \CAvg{\sigma_j(t+s)}{\bm \sigma(t)} + \CAvg{g_j(\bm \sigma(t+s))}{\bm \sigma(t)}

\Avg{\sigma_i(t) \cdot \text{(l.h.s)}} & = \Avg{\sigma_i(t) \, \tau \frac{\D}{\D s} \CAvg{\sigma_j(t+s)}{\bm \sigma(t)}} \\ & = \tau \frac{\D}{\D s} \Avg{\sigma_i(t) \CAvg{\sigma_j(t+s)}{\bm \sigma(t)}} \\ & = \tau \frac{\D}{\D s} \Avg{\sigma_i(t) \, \sigma_j(t+s)}

\Avg{\sigma_i(t) \cdot \text{(r.h.s)}} & = \Avg{\sigma_i(t) \left\{ - \CAvg{\sigma_j(t+s)}{\bm \sigma(t)} + \CAvg{g_j(\bm \sigma(t+s))}{\bm \sigma(t)} \right\}} \\ & = - \Avg{\sigma_i(t) \, \sigma_j(t+s)} + \Avg{\sigma_i(t) \, g_j(\bm \sigma(t+s))}

\tau \frac{\D}{\D s} \Avg{\sigma_i(t) \, \sigma_j(t+s)} & = - \Avg{\sigma_i(t) \, \sigma_j(t+s)} + \Avg{\sigma_i(t) \, g_j(\bm \sigma(t+s))}

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

計算テクニック

和の積の計算のための添字テクニック

課題

説明を追加

\left( \sum_i a_i \right) \left( \sum_i b_i \right) = \left( \sum_i a_i \right) \left( \sum_j b_j \right) = \sum_i \sum_j a_i b_j

\left( \sum_i a_i \right)^2 = \sum_i \sum_j a_i a_j

\(\newcommand{\D}{\text{d}} \newcommand{\I}{\text{i}} \newcommand{\E}{\text{e}} \newcommand{\Prob}{\mathbb{P}} \newcommand{\Expect}{\mathbb{E}} \newcommand{\Var}{\text{Var}} \newcommand{\PAvg}[1]{\left[#1\right]} \newcommand{\Avg}[1]{\left\langle#1\right\rangle} \newcommand{\AvgJ}[1]{\Avg{#1}_{\bm J}} \newcommand{\AvgDyn}[1]{\Avg{#1}_{\text{dyn.}}} \newcommand{\CAvg}[2]{\Avg{#1}_{\left|#2\right.}} \newcommand{\Devi}{\mathfrak{d}}\)

参考文献

[vanVreeswijk1998]van Vreeswijk, C., & Sompolinsky, H. (1998). Chaotic balanced state in a model of cortical circuits. Neural Computation, 10(6), 1321–71. http://dx.doi.org/10.1162/089976698300017214
[Ginzburg1994]Ginzburg, I., & Sompolinsky, H. (1994). Theory of correlations in stochastic neural networks. Physical Review E, 50(4), 3171–3191. http://dx.doi.org/10.1103/PhysRevE.50.3171
[Derrida1987]Derrida, B., Gardner, E., & Zippelius, A. (1987). An Exactly Solvable Asymmetric Neural Network Model. Europhysics Letters (EPL), 4(2), 167–173. http://dx.doi.org/10.1209/0295-5075/4/2/007
[Cinlar2011]Çınlar, E., (2011). Probability and Stochastics. Springer New York (Vol. 230) http://dx.doi.org/10.1007/978-0-387-87859-1