理化学研究所 脳科学総合研究センター 脳数理研究チーム 岡田真人 スパイク解析理論入門 理化学研究所　脳科学総合研究センター 脳数理研究チーム 岡田真人

背景(1/2) 多数のニューロンの タスク スパイクのデータ 行動，特徴 相互情報量 スパイクの統計量 平均発火率 高次統計量 順問題 逆問題 順問題 拘束 ニューラルネットワークモデル 解剖学的知見 予測

背景(2/2) どのようにして統計的に取り扱うか スパイク列 unit 多数ニューロンスパイクの 同時記録 1次統計量以上の情報 time unit スパイク列 多数ニューロンスパイクの 同時記録 1次統計量以上の情報 (平均発火率) 統計量 2次や3次以上の統計量 近年の電気生理の発展により、多数ニューロンの同時記録が可能となってきた。 従来のニューロサイエンスは、スパイク列の１次統計量である平均発火率に着目していた。 今後は、２次や3次の統計量も原理的には利用できる。 そこで我々は、スパイク列から、各種統計量を推定する方法について議論する。 そのために、詳しくは後に述べるが、スパイクは確率分布から生成されるとする。 スパイクを生成する確率分布をこのように仮定する。 ここで、σはスパイク列を表している。 各要素は発火、非発火に対応して、＋１、－１の値をとる。 確率分布のは、分布の統計量を反映したパラメータθをもっているとする。 本研究の目的の一つは、スパイク列が与えられたときに、スパイク生成の確率分布のパラメータθを推定すること。 どのようにして統計的に取り扱うか

スパイク列をビット列に変換する 各区間にスパイクが高々1個となるように時間軸を離散化する． ビン 時刻 t の i 番目のニューロンの発火状態 発火： 1, 非発火： 0 と対応させる

スパイク生成の同時確率分布 同時確率分布 からスパイクが生成されたとする． 自由度 現実的には決められない!!

取り扱える問題にする(1/3) 2変数の場合を考える 周辺化によって自由度を落とす 一部のデータは見なかったことにする 独立だと思う 1 1 一部のデータは見なかったことにする 独立だと思う 独立同分布だと思う 時空間同時確率分布 時間と空間に分解

取り扱える問題にする(2/3) 時間に注目：1個のニューロンのダイナミクスに注目 S1 0 1 0 0 0 1 0 1 1 1 0 0 ニューロン間の相互作用も繰り込まれている 普通は T は大きい 自由度 これでも大変 定常性を仮定する と の関係と と の関係は同じ

自己相関関数 定常性を仮定する と の関係と と の関係は同じ 時定数

取り扱える問題にする(3/3) 空間に注目：同時刻のニューロン間の相関に注目 独立だと思う サンプル数が足りない 定常性を仮定し，独立同分布だと思う S1 0 1 0 0 0 1 0 1 1 1 0 0 ニューロンのダイナミクスも繰り込まれている S2 0 0 0 1 0 0 0 1 0 0 0 1

対数線形モデル(2ニューロンの場合) 1 自由度 3 1 自由度 3 対数線形モデル

対数線形モデル(Nニューロンの場合) 自由度： 1対1対応 自由度： 対数線形モデル

いろいろな座標系 ニューロン数N=2の場合はとてもよい性質がある．

情報幾何的視点 対数線形モデル 直交性 局所座標の性質 KL分解 カルバックダイバージェンスに関する ピタゴラスの定理 相関成分 発火率成分 Let me summarize two neuron case. Log-linear expansion, theta’s, eta is firing rate. We can have different coordinate systems. (maybe explain JPSTH too) Local property and KL dicomposition. カルバックダイバージェンスに関する ピタゴラスの定理 相関成分 発火率成分

検定 Χ2検定 相互情報量の分解 相関の成分 平均発火率の成分 与えられた帰無仮説 推定された帰無仮設   相互情報量の分解 推定された帰無仮設 It’s possible to have a statistical test using asymptotic theory. Distance is estimated by xi square distributions. MI decomposition is also possible. Can relate each correlational component to behavior or stimulus. 相関の成分 平均発火率の成分

相互相関関数との関係(同時刻) S1 0 1 0 0 0 1 0 1 1 1 0 0 S2 0 0 0 1 0 0 0 1 0 0 0 1 と　　は単調増加な関係

相互相関関数との関係(時間差t) S1 0 1 0 0 0 1 0 1 1 1 0 0 S2 0 0 0 1 0 0 0 1 0 0 0 1 に関して独立同分布の対数線形モデルを導入する と　　　　　は単調増加な関係

2次の相互作用と相互相関関数の比較

中間のまとめ 多数のニューロンの タスク スパイクのデータ 行動，特徴 相互情報量 スパイクの統計量 平均発火率 高次統計量 順問題 逆問題 拘束 ニューラルネットワークモデル 解剖学的知見 予測

統計量には繰り込んだ変数の影響が入っている S1 0 1 0 0 0 1 0 1 1 1 0 0 ニューロンのダイナミクスも繰り込まれている S2 0 0 0 1 0 0 0 1 0 0 0 1 相関関数だけで構造を議論するのは注意が必要 逆問題

対数線形モデル(Nニューロンの場合) 自由度： パラメータ 個数： 1対1対応 対数線形モデル

ニューロン数が多い場合 相互作用の一様性を仮定 相互作用の次数に上限を仮定 この条件がどのような性質を持つか？ 自由度： ニューロン数が２の場合について、相互作用の強さを推定してみる。 系はこのように書ける。 N=2なので、系の取る状態は４つ。 スパイク列から相互作用の強さを推定するには、系がそれぞれの状態をとる確率をスパイク列から求め、それをもとに、２＾２－１次の連立方程式をとく。 また、ｐ－－などを求める際には、Nかける１００位のデータ数が必要。 Nが３なら２＾３－１。 Nが多いときに、このような解析法は無理。 そこで、仮定を置くことが必要となる。 まず、同じ次数の相互作用は一様であると仮定。 次に、相互作用には上限があると仮定。 この仮定は、結構一般的に使われている。 その結果、Hはこのような形になる。 この形は、統計力学的解析にぴったり。 後の解析のため、変数変換を行う。Jrが相互作用の強さを示すパラメータとなる。 この条件がどのような性質を持つか？

統計力学的解析 この条件がどのような性質を持つか？ N→∞の極限を考える ：r次の相互作用の強さ ニューロン数が２の場合について、相互作用の強さを推定してみる。 系はこのように書ける。 N=2なので、系の取る状態は４つ。 スパイク列から相互作用の強さを推定するには、系がそれぞれの状態をとる確率をスパイク列から求め、それをもとに、２＾２－１次の連立方程式をとく。 また、ｐ－－などを求める際には、Nかける１００位のデータ数が必要。 Nが３なら２＾３－１。 Nが多いときに、このような解析法は無理。 そこで、仮定を置くことが必要となる。 まず、同じ次数の相互作用は一様であると仮定。 次に、相互作用には上限があると仮定。 この仮定は、結構一般的に使われている。 その結果、Hはこのような形になる。 この形は、統計力学的解析にぴったり。 後の解析のため、変数変換を行う。Jrが相互作用の強さを示すパラメータとなる。

モデル time unit が与えられている場合 高次相互作用 系の性質を議論 スパイク列からの高次相互作用 の推定 まとめると、このようになる。 今後我々は、２通りの方法で、このモデルを評価する。 本当にやりたいことは、スパイク列が与えられたときの、相互作用の決定。 まず、高次相互作用が与えられているときの、このモデルの振る舞い。 その次に、スパイク列からの高次相互作用の推定。 系の性質を議論 スパイク列からの高次相互作用　　　の推定

統計力学的解析(1/2) 理論解析の結果 統計力学では 例えば 系の巨視的性質 mのみで記述可能 分配関数が重要 分配関数が重要 さまざまな統計量が分配関数からもとまる。 ニューロン発火の期待値も分配関数からかける。 ニューロン発火の期待値は、平均発火率と等価。 計算の結果、分配関数はこのようになる。 Jrは与えられているので、系はｍのみで記述できる。 系の巨視的性質 mのみで記述可能

統計力学的解析(2/2) mの物理的な意味 ＝ ニューロンの発火の期待値 系の対称性 系のｒ次の分散は0 平均発火率＝ 系の巨視的性質 　　　　　＝ ニューロンの発火の期待値 系の巨視的性質 ニューロンの発火の平均発火率のみで記述可能 同様の計算により 系の対称性 系の巨視的な性質を記述する、ｍの物理的な意味を求めてみる。 先ほどの計算の続きを行う。 ニューロン発火の期待値はｍとなる。 Ｍの物理的な意味はニューロン発火の期待値である。 ニューロン発火の期待値ｍは、線形な変換で平均発火率にすることができる。 系の巨視的な性質は、平均発火率のみで記述することができる。 系のｒ次の相関を計算すると、ｍｒになる。 このことと、系の対称性をこうりょすると、系のｒ次の分散は０になることがわかる。 系のｒ次の分散は0 平均発火率＝

計算機実験1-1 系の巨視的性質 を用いてスパイクを生成 time unit (a) (b) (c) mのみで記述可能 同じｍを与える　　　の組は複数存在 (a)相互作用なし のみ値をもつ (b)２体相互作用のみ 　　　　　のみ値をもつ (c)４体相互作用のみ 同一の m=-0.96 系の巨視的な性質は、ｍで記述できる。 Ｍを与える式はこのようなもの。 この式から、同じｍを与えるJrの組が複数存在することがわかる。 同じ平均発火率、異なる相互作用を持つような系をいくつかつくって、スパイクを発生させた。 一つは、２次以上の相互作用を持たない系。 もう一つは、２次の相互作用のみを持つ系。 もう一つは、４次の相互作用のみを持つ系。 全て、同一の平均発火率0.02を持つようにした。 結果、このようなスパイク列が得られた。

計算機実験1-2 (a)相互作用なし のみ値をもつ (b)２体相互作用のみ のみ値をもつ (c)４体相互作用のみ 同一の m=-0.96 　　　　　のみ値をもつ (c)４体相互作用のみ 同一の m=-0.96 平均発火率２％ 本当に、相互作用の異なる各領域でｍが同じ値を持つか、先のデータから計算してみた。 結果はこのように、設計した平均発火率が得られた。 平均発火率＝

高次相互作用の推定(1/4) time unit が与えられている場合 高次相互作用 スパイク列からの高次相互作用 の推定 スパイク列からの高次相互作用　　　の推定 推定方法の性能評価 先までは、相互作用が与えられている場合の、系の性質についての議論。 これから当初の興味の対象の、スパイク列から高次相互作用Jrを推定する話。 モデルは今までと同じ物を使う。 このようなスパイク列が与えられたときに、スパイク間の高次相互作用を推定する。

高次相互作用の推定(2/4) 高次相互作用 観測されたスパイク列 モデル 推定 の観測される確率 観測時のノイズ ：真のスパイク列 一般的に、このようなシチュエーションを想定する。 まず、分布関数ｐ（σ）から、スパイク列が生成されるとする。 スパイク列σはそのまま観測されるのではなく、観測時にノイズが載ってしまい、一般にσとは異なるデータ、τが観測されるものとする。 この、データτが観測される確率を、ｐ（τ）とおく。 このようにして観測されたデータτから、高次相互作用Jrを推定しようというわけです。 ：真のスパイク列 ：観測されるスパイク列

高次相互作用の推定(3/4) つまり ということは まとめると の推定値 ここで 対数周辺尤度 対数周辺尤度を最大化 スパイク列　　が観測される確率 に依存 つまり を変更 が変更 ということは スパイク列　　が得られたとき を変更 最も　　が起こりやすいような　　　　　が存在 まとめると の推定値 の推定値：　　　　　を最大化 ここで 観測されたスパイク列τからの、具体的な相互作用の推定方法について述べる。 スパイク列τが観測される確率ｐ（τ）は、ｐ（σ）にノイズが載っているとはいえ、真のスパイク列生成の相互作用の構造Jrに依存することは明らか。 Jrを変更すると、ｐ（τ）も変更される。 スパイク列τが観測されたときに、Jrをいろいろ動かしてみると、最もτが起こりやすいようなｐ（τ）を得ることができる。 そのようなパラメータを、Jrの推定値として採用することにする。 まとめると、Jrの推定値は、ｐ（τ）を最大化するようなパラメータである、ということがいえる。 しかし、ここで、ｐ（τ）の形は、明らかではありません。 しかし、全確率の公式を使って、ｐ（τ）を、このように、書き下すことができる。 これは、周辺尤度関数と呼ばれる。 つまり、Jrを推定するためには、周辺尤度を最大化したい。 ｐ（σ）などは指数関数なので、対数をとったほうが、後々の解析が容易。 そこで、対数をとる。これは、対数周辺尤度と呼ばれる。 Logは単調関数なので、ｐ（τ）を最大化するパラメータは、Kも最大化します。 つまり、今回我々がとる戦略は、対数周辺尤度を最大化することで、スパイクの高次相互作用を推定するというもの。 対数周辺尤度 対数周辺尤度を最大化 スパイクの高次相互作用を推定

高次相互作用の推定(4/4) 今回は、ノイズはのらないものとした 最急上昇法 高次相互作用の推定： を観測 を動かしてKを最大化 を解く time unit 計算の結果 この式をもとに 今回は、観測データにノイズは乗らないものとした。 また、Kの最大値を求めるにあたっては、この式で表される最急上昇法を用いた。 計算の結果、対数周辺尤度Kの微分はこのような関数になる。 この式をもとに、最急上昇法、すなわち、この微分方程式をとくということをおこなって、このようなスパイク列から、スパイク間の高次相互作用の推定を行った。 高次相互作用を推定

計算機実験2 ：2次までの相互作用を持ち、 相互作用は 一意には決まらない 推定結果 系の巨視的性質 mのみで記述可能 同じｍを与える 結果を示す。 スパイクを生成した真のｐ（σ）は、２次までの相互作用をもち、J1=0.5、J2=0.2。 推定結果はこのようになった。 各点はスパイク列τが観測されたときの、J1,J2の相互作用の推定値を表している。 点は、全部で２００個ある。 横軸が１次の相互作用の係数、縦軸が２次の相互作用の係数で、例えばこの辺りの点は、J1=0.1, J2=0.5である。 本当のJrはここ。 この図からわかるように、相互作用の推定値は一意には決まらずに、ある不定性をもっている。 これはなぜかというと、これまでの議論で、系の巨視的な性質はｍのみで記述可能だといいったが、このｍを与える式からわかるように、同じｍを与えるJの組が複数存在するから。 この黒い実線が、この不定性の条件を示す理論線。 推定結果は、理論線の周りでばらついていることがわかる。 同じｍを与える　　　 の組は複数存在

推定方法の一般的性質 Kを、真の に関して平均をとる 一般的な の推定の性能を評価 計算すると の条件で不定 を解くことで、高次相互 一般的な　　　の推定の性能を評価 計算すると の条件で不定 先ほどの実験は、与えられたスパイク列から相互作用を推定するというもであった。 でもそれだけでは、この推定法のもつ性能がよくわからない。 そこで次に我々は、一般に周辺尤度最大化の基準での高次相互作用の推定がどのくらいの性能を持つか、評価を行った。 対数周辺尤度Kを、真の分布ｐ（τ）に関して平均をとり、得られた[K]について、高次相互作用の推定値の収束の様子などをみてみた。 計算の結果、[K]はこのような形をしていることがわかった。 先ほどと同様に、この[K]について、最急上昇法を用いることで、Jrの推定値を求めてみた。 [K]の微分はこのような形をしている。 この関数の形からも、Jrの推定値が不定であることが見て取れる。 を解くことで、高次相互 作用推定法の一般的な 性質をみる

計算機実験3 [K]の形 実線は最大値を結んだもの 推定結果の不定性 推定値の収束の様子 得られた推定値 結果をしめす。 対数周辺尤度の平均値は、このような形をしていることがわかった。 実線は、[K]の最大値を結んだ線。 先ほどの理論から予測されたように、最大値を与えるパラメータはやはり不定となっている。 最急上昇法の結果得られた推定値が、この図の各点。 見方は先ほどの図と同じ。 この系はNが無限大の理想的な状況ですから、推定値は不定の理論線の上にきちんとのる。 理論線上でばらついているのは、Jrの初期値が異なるから。 最後に、推定値の収束の様子をみてみた。 横軸は、先の微分方程式をオイラー方で解いたときのタイムステップ。 指数関数的に素直に収束していることがわかる。 推定値の収束の様子 得られた推定値

均一な相関構造を仮定した場合のまとめ もちいた仮定 推定方法 結果 ニューロン数無限 相関はニューロンに依らず均一 相関の次数は有限で打ち切る． 推定方法 周辺尤度最大化 結果 相関は一意には決まらない 正解 推定結果

わかったこと 多数のニューロンの タスク スパイクのデータ 行動，特徴 相互情報量 スパイクの統計量 平均発火率 高次統計量 順問題 逆問題 拘束 ニューラルネットワークモデル 解剖学的知見 予測

高次統計量の理論(1/2) 2値ニューロン ・揺動散逸定数（線形応答理論，久保公式） ・マスター方程式の方法 平衡点まわりのゆらぎが 高次統計量 ・マスター方程式の方法 これらの知見が定性的にスパイキング ニューロンでもどの程度成り立つか?

高次統計量の理論(2/2) E-E間の詳細 r: 双方向結合 u: 一方向結合 av: 平均 d: 結合なし E-E間の平均 E-I間の平均 (Ginzburg & Sompolinsky, 1994) I-I間の平均

まとめ