バンドパスフィルタを用いた振動子集団の同期の解析1(指導教員：工藤和恵)

バンドパスフィルタを用いた同期の検出

単一振動子, すなわち N = 1 の場合を扱う. 周期外力の振動数 Ω = 0.7, 結合強度 κ = 0.02 のとき, 振動子の固有振動数が Ω−κ ≤ ω1 ≤ Ω+κ すなわち 0.68 ≤ ω1 ≤0.72 を満たすとき振動子と周期外力は位相ロックする.
そこで周期外力と非同期の場合の例として振動子の固有振動数 ω1 = 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.8, 0.9, 1.0,
周期外力と同期の場合の例として振動子の固有振動数
ω1 = 0.68, 0.7, 0.72 で数値計算を行った. 実効的な振動数が周期外力の振動数と一致する場合, 同期を検出できたと見なす. 観測位相から算出する実効的な振動数は以下で求める. 今回は T1 = 9000, T2 = 10000 とした.
⟩ Φ (T ) − Φ (T )⟨ωf2
ωf1
Φ̇ωf =(3)T2 − T1(a)

1.2

モデル

0.80.60.4

ϕ̇k = ωk + κ sin(Ωt − ϕk )01 ∑f (ϕk )
N

0.8
0.6
0.4
0.2

00.2

0.4

0.6

frequency

0.8

0100.20.4

0.6

frequency

0.8⟨⟩
図 1: 横軸は ωf , 縦軸は Φ̇ωf
通過帯域幅 β によらず, 周期外力と同期する固有振動数のときは観測位相の実効的な振動数は周期外力の振動数と一致する. 一方で, 非同期の場合は通過帯域幅によって結果が異なる. β = 0.22 のときは実効的な振動数は固有振動数と一致する. 一方 β = 0.01 のとき,
ω1 , |Ω − 2ω1 |, Ω で plateau が見られる. これを摂動論を用いて定性的に説明する. 式 (1) の解が定数 Ck を用いて(1)次に観測量を定義する. 振動子集団の平均的な振る舞いを観測できる場合について考えるので, 一般に観測量 F (t) は出力関数 f (ϕk ) を用いて以下のように与えられる. ただし, N は振動子の総数である.
NF (t) =κ=0.02 β=0.22

10.2

周期外力を受ける集団位相振動子を考える. 振動子
k(1, . . . , N ) の位相を ϕk , 固有振動数を ωk , 周期外力との結合強度を κ, 周期外力の振動数を Ω として, 各振動子のダイナミクスが以下で与えられるとする [2].

(b)

1.2

κ=0.02 β=0.01

1frequency2単一振動子の数値計算・解析計算

3

異なる固有振動数を持つ振動子集団が作用すると, 同期と呼ばれる秩序化現象が起きる. ‘Mooney‘ face とよ
ばれる, 情報量の少ない, 白黒の 2 値画像を見せ, 被験者が顔認識した場合には脳内の長距離離れた 2 点間でγ 波に同期が見られることを示す実験がある [1]. γ 波は脳波のパターンの 1 つで, 知覚や意識に関連づけられており, 神経細胞集団が電気信号を 40Hz 程度で放出した際に発生する. そこで, 電極から得られたデータに対し, (被験者の γ 波の周波数)±3Hz のみを通す
バンドパスフィルタに通して位相を得ることで上記の結果を得た.
電極によって測定される脳波は, 様々な振動成分を持つ 1 次元の時系列と考えられる. そこで, 本研究では, 様々な振動成分を持つ 1 次元の時系列を簡単な数理モデルで作る. 今回は周期外力の加わる系を扱う.
振動子集団の平均的な振る舞いのみを観測できる場合に, バンドパスフィルタを用いることで, 周期外力との同期を検出できるのかを検証する.

frequency

稲垣志帆

(2)∞
∑

k=1

ϕk (t) = ωk t + Ck +今回は簡単のため, f (ϕk ) = cos ϕk とした.
フィルタは BiQuad フィルタ [3] を用い, 通過帯域はβフィルタの通過振動数を ωf としたときに ωf 2− 2 <βω < ωf 2 2 とした. β は通過帯域幅を表すパラメータで, 単位は octave である. 今回は β = 0.01 と β = 0.22とした.
観測位相は以下の手順で得られる.

κn fn (t)(4)n=1で与えられるとする. ここで fn (t) は κ について O(1)とする. 式 (4) を式 (1) に代入することによりf1 (t) =1cos {(Ω − ωk )t − Ck }ωk − Ω(5)を得る. これをさらに式 (2) に代入すると

1. 観測量の時系列 F (t) に通過振動数 ωf の BiQuad
フィルタを施し, 得られた時系列を F̃ (t) とする.

F (t) = a1 cos(ωk t + Ck ) + a2 sin {(Ω − 2ωk )t − 2Ck }− a3 sin Ωt + O(κ2 )2. F̃ (t) をヒルベルト変換して得られた時系列を H̃(t)とする [2].

(6)を得る. ここで

3. 複素数 z(t) = F̃ (t) + iH̃(t) の偏角を Φωf (t) とし,これを観測位相とする.

a1 = 1,

ここで, ヒルベルト変換とは各フーリエ成分の位相をπ2 遅らせる操作である.

a2 = a3 =κω2(ωk − Ω)(7)である. 図 2 は, 式 (6) で O(κ2 ) を無視したものをフィルタに通すと得られる F̃ (t) の時系列を式 (8) のように

11

隔, すなわち ωk+1 − ωk により結果が大きく異なることがわかる. これらの結果の解析的なアプローチは今後の課題である.

与えたときに, 各振動数成分の係数の絶対値をプロットしたグラフである. 図 2(a) では |ã1 | より |ã2 |, |ã3 |が大きくなる ωf が存在する. さらに |ã3 | が |ã1 | より大きくなる領域の方が, |ã2 | が |ã1 | より大きくなる領域よりも大きい. 一方で図 2(b) では |ã1 | より |ã2 |, |ã3 |が大きくなる ωf が存在しない. 以上から図 2(a) は図1(a) の, 図 2(b) は図 1(b) の定性的な説明と言える.
{}F̃ (t) = ã1 cos(ωk t + C̃k1 ) + ã2 sin (Ω − 2ωk )t + C̃k2

11(a)
κ=0.0

0.6

0.4

0.4

0.2

0.2

0+ ã3 sin(Ωt + C̃k3 )

0.4

0.6

0.8

frequency

κ=0.02

0.025

~ ~ ~
|a1|, |a2|, |a3|

β=0.01

0.003

0.002

β=0.22
0.02

0.01
0.005

000

0.2

0.4

ωf

0.6

0.8

1(c)
β=0.01

0.002
0.00150.001
0.0005

5
00.2

0.4

ωf0.003κ=0.020.00250.6

0.8ωf0.50.52(d)κ=0.020.0025

β=0.010.0020.0015
0.001

00.68

ωf0.70.72

図 2: 赤色は |ã1 |, 水色は |ã2 |, 紫色は |ã3 |, (c),(d) は(a) の拡大図

4

今後の課題：多振動子系多振動子の場合に同期が検出できるか数値計算を行った. まず 100 個の振動子の振動数分布は 0 以上 1 未満の一様分布とした. 通過帯域幅 β は β = 10−5 , 0.01, 0.22の場合を実行した. 周期外力と結合していないときは,通過帯域幅が広いと振る舞いが特殊になることがわかる (図 3(a)). 一方で, 周期外力と結合しているときは, 通過帯域幅が狭いと plateu が見られず, 同期を検出できないことがわかる (図 3(b)). 次に通過帯域幅を
β = 0.22 に固定し, 振動数分布を変える. 図 3(c), (d)で固有振動数は以下のように与えた.

red : ωk =k100

(0 ≤ k ≤ 99)

k(0 ≤ k ≤ 199)
100kgreen : ωk =(0 ≤ k ≤ 49)
25
2k(0 ≤ k ≤ 14)
light blue : ωk =15
purple : ωk =0.6

0.8(9)

1

1

0.5

00.5

11.5

2

0

0

0.5

1frequency⟨⟩
図 3: 横軸は ωf , 縦軸は Φ̇ωf )frequency

1.5

まとめと展望参考文献(10)[1] Rodriguez E, George N, Lachaux JP, MartinerieJ, Renault B, Varela FJ. Perceptions shadow: longdistance synchronization of human brain activity.
Nature Feb 4; 397(6718):430-3(1999).
[2] Arkady Pikovsky, Michael Rosenblum, JurgenKurths. 徳田功 (訳). 同期理論の基礎と応用: 数理科学、化学、生命科学から工学まで, 丸善 (2009).
[3] Robert Bristow-Johnson.
Cookbook for-(11)(12)

周期外力と結合していないとき, 振動数の間隔, すなわ
ち ωk+1 − ωk が大きいと特殊な振る舞いは見られなくなることがわかる (図 3(c)). また周期外力と結合しているとき, ωf ≈ Ω で plateau はいずれの場合でも見られるが, それ以外の中心振動数 ωf のときは振動数の間

mulae for audio EQ biquad filter coeﬃcient.
http://www.musicdsp.org/files/Audio-EQ-Cookbook.txt

21κ=0.1 β=0.22系を構成する振動子集団の平均的な振る舞いのみを観測できる場合に, バンドパスフィルタを用いることで, 周期外力との同期を検出できた. また, フィルタに通すことで, 同期・非同期だけでなく, 振動数解析が行えるが, フィルタの帯域幅により結果が異なる. これらの数値計算の結果を摂動論を用いて定性的に説明した.
また, 今回は周期外力の加わる系で解析を行ったが,振動子同士が相互作用する結合振動子系でも同様の解析を試みたい.
さらに展望としては, 実データに適用することを目指したい. 振動数分布によって, フィルタの適切な通過帯域幅は変わると考えられる. 今回のようにダイナミクスがわかる系では, 結果の照らし合わせによりフィルタの通過帯域幅が適切かの判断ができたが, ダイナミクスがわからない系では判定ができない. そこで, 実データが与えられた時に, 適切なフィルタの通過帯域を見積もる方法の提案が第一段階である. 次に振動数成分が幅広く分布する際や, どの振動数成分が存在するかわからない系では, 数多くのフィルタの振動数で数値計算を行うことは現実的ではない. そこで, 適切なフィルタの振動数を求める理論が第二段階である. これらの知見を生かして, ダイナミクスがわからない系での, 振動子の振動数, および, 振動子間の相互作用の強さや向きの推定を行いたい.

10.0005

00.48

0.4

frequency

0.015

0.0010.003~ ~ ~|a1|, |a2|, |a3|(b)κ=0.02~ ~ ~|a1|, |a2|, |a3|

~ ~ ~|a1|, |a2|, |a3|

0.004

(a)

0.2

1.5

0.5

0

0.03

0.005

0(d)

2

κ=0.0 β=0.22

1.5

C̃k1 , C̃k2 , C̃k3 はフィルタに通すことにより存在しうる位相遅れのため導入した.

01

frequency

frequency0.2(c)

2(8)

κ=0.1

frequency

0.6

0(b)

0.8

frequency

0.82