汪芃 李倩昀 黃志精 唐國寧

(廣西師范大學物理科學與技術(shù)學院,桂林 541004)(2018年3月21日收到;2018年5月22日收到修改稿)

1 引 言

大腦神經(jīng)系統(tǒng)中包含大量神經(jīng)元,一個神經(jīng)元可以與許多神經(jīng)元同時存在耦合.因此在大腦處理電信號的過程中,神經(jīng)元集體電活動的時空斑圖是極其復(fù)雜的.研究發(fā)現(xiàn):在大腦視覺皮層內(nèi)可以自發(fā)出現(xiàn)平行波[1],在睡眠狀態(tài)下或是藥物引起的振蕩中,哺乳動物大腦皮層中也會自發(fā)出現(xiàn)螺旋波[2,3];此外在癲癇發(fā)作時,在大腦皮層中也觀察到螺旋波[4].由于螺旋波是一種繞波頭旋轉(zhuǎn)的波,因此它可以在介觀尺度上組織和調(diào)節(jié)皮層內(nèi)神經(jīng)元的群體活動,從而對大腦的功能產(chǎn)生影響.螺旋波等有序波除了可以出現(xiàn)在腦神經(jīng)系統(tǒng)外,也可以出現(xiàn)在心臟系統(tǒng)中[5].我們知道心臟是可激發(fā)介質(zhì),需要起搏源才能使心臟工作.由于螺旋波是不需要波源就能夠自維持的波,而心臟中螺旋波的頻率通常大于心臟竇房結(jié)(心臟的起搏源)產(chǎn)生的靶波頻率(即心臟的跳動頻率),所以當心臟出現(xiàn)螺旋波電信號時,將引起心動過速,當螺旋波破碎成時空混沌時,就會導致心顫,危及生命.因此,了解大腦神經(jīng)系統(tǒng)和心臟系統(tǒng)中螺旋波的產(chǎn)生機制及其動力學行為,對了解大腦功能和治療癲癇疾病、心臟疾病具有重要意義.

最近,人們從理論和實驗上對神經(jīng)元網(wǎng)絡(luò)中螺旋波的產(chǎn)生機制和動力學行為進行了大量研究[6?12].許多研究結(jié)果表明:在神經(jīng)元網(wǎng)絡(luò)中,噪聲可以誘發(fā)螺旋波[7,9,13?15].Qin等[11]在神經(jīng)元網(wǎng)絡(luò)中觀察到自突觸誘發(fā)螺旋波和靶波;Ma等[16]發(fā)現(xiàn),增加長程連邊概率將導致螺旋波破碎;Wang等[17]發(fā)現(xiàn)調(diào)節(jié)細胞膜上的離子通道電導率也可以導致神經(jīng)元網(wǎng)絡(luò)自發(fā)產(chǎn)生螺旋波;在近鄰耦合Hindmarsh-Rose(HR)神經(jīng)元陣列中,當神經(jīng)元處于一周期或二周期態(tài)時,如果系統(tǒng)從混沌初相位分布態(tài)開始演化,適當選擇耦合強度,可以在這樣的系統(tǒng)中觀察到自發(fā)出現(xiàn)螺旋波[18];Fohlmeister等[19]在考慮興奮耦合、抑制耦合、噪聲和軸突延時的情況時,在神經(jīng)元網(wǎng)絡(luò)中觀察到自發(fā)出現(xiàn)螺旋波;Xiao等[20]在由類型I和類型II激發(fā)性神經(jīng)元組成的網(wǎng)絡(luò)中觀察到螺旋波的自發(fā)產(chǎn)生;Tao等[21]在由錐體神經(jīng)元和抑制性中間神經(jīng)元組成的神經(jīng)元網(wǎng)絡(luò)中也觀察到螺旋波的自發(fā)產(chǎn)生.這些結(jié)果表明:噪聲作用和與神經(jīng)元狀態(tài)相匹配的耦合方式可以使神經(jīng)元網(wǎng)絡(luò)自發(fā)產(chǎn)生螺旋波,興奮性神經(jīng)元與抑制性神經(jīng)元組合也能導致神經(jīng)元網(wǎng)絡(luò)自發(fā)產(chǎn)生螺旋波.目前對噪聲誘發(fā)螺旋波已經(jīng)研究得比較多,但是對于耦合方式誘發(fā)螺旋波的研究相對比較少,需要進一步的研究.

考慮到大腦皮層神經(jīng)元同時受到大量不相關(guān)的興奮和抑制突觸電流的影響[22],在皮層組織培養(yǎng)中,觀察到抑制性神經(jīng)元與總神經(jīng)元之比達到20%,而在海馬組織培養(yǎng)中,觀察到抑制性神經(jīng)元與總神經(jīng)元之比達到30%[23],因此本文采用HR神經(jīng)元模型[24]研究混沌神經(jīng)元網(wǎng)絡(luò)是否可以自發(fā)出現(xiàn)有序波.該網(wǎng)絡(luò)由近鄰興奮性耦合和長程抑制耦合層組成.我們發(fā)現(xiàn),適當選擇兩種耦合的耦合強度,系統(tǒng)可以自發(fā)產(chǎn)生迷宮斑圖、螺旋波、平面波和向內(nèi)方形波,這些結(jié)果對了解大腦中螺旋波和平面波的形成機制具有積極意義.

2 模 型

研究中采用HR神經(jīng)元模型[24],網(wǎng)絡(luò)通過在二維方形點陣上增加長程連邊形成,其中近鄰耦合為興奮性耦合,長程耦合為抑制性耦合.該網(wǎng)絡(luò)動力學方程如下:

式中x代表神經(jīng)元的細胞膜電位;y代表與內(nèi)電流相關(guān)的恢復(fù)變量;z代表與鈣離子激活的鉀離子電流相關(guān)的慢變調(diào)節(jié)電流;a,b,c,d,r,s和x0為HR神經(jīng)元模型參數(shù);Iext為外部直流激勵;D為擴散系數(shù)(以下稱為興奮性耦合強度);g為抑制性耦合強度;xth為抑制性耦合的閾值;xav為長程耦合神經(jīng)元x變量的平均值;m為長程耦合距離;角標i,j=1,2,···,N代表N×N個HR神經(jīng)元的位置坐標.考慮到實際大腦網(wǎng)絡(luò)中抑制性神經(jīng)元占比達到20%—30%,我們每間隔一個格點安排一個抑制性格點,其位置角標為i,j=1,3,5,···,N.在這些格點上,?i,j=1,在非抑制性格點上,?i,j=0,因此抑制性格點與總格點數(shù)之比約為1/4.從(1)—(4)式可以看出,所有格點上的神經(jīng)元參數(shù)相同,但是在抑制性格點上還要考慮長程抑制耦合.

本文固定取a=1.0,b=3.0,c=1.0,d=5.0,r=0.006,s=4.0,x0=?1.6,Iext=3.0,N=201,m=8,xth=?1.5.由于增大m會得到類似結(jié)果,因此本文只將D和g設(shè)為可調(diào)參數(shù),即讓這兩個參數(shù)在0.1≤D≤2.0和0.1≤g≤0.9范圍內(nèi)變化.根據(jù)文獻[25]的結(jié)果,在此參數(shù)下單個神經(jīng)元處于混沌振蕩狀態(tài).在數(shù)值模擬中,解方程采用四階龍格-庫塔法,時間步長取Δt=0.01,近鄰擴散耦合使用無流邊界條件,抑制性耦合采用周期邊界條件,即當抑制格點坐標i′=i?m＜1時,取i′=N+i+1?m,當抑制格點坐標i′=i+m＞N時,取i′=i+m?N?1,每次數(shù)值模擬時間長度為10000單位時間.

為了描述系統(tǒng)的振蕩狀態(tài),引入系統(tǒng)的膜電位方差σ:

當σ小幅度不規(guī)則振蕩時,系統(tǒng)處于有序振蕩狀態(tài);σ大幅度無規(guī)則振蕩時,系統(tǒng)處于混沌振蕩狀態(tài).

3 數(shù)值模擬結(jié)果

在給定耦合強度D和g下,系統(tǒng)是否出現(xiàn)螺旋波與系統(tǒng)初態(tài)有關(guān).為了比較不同初態(tài)的演化效果,我們規(guī)定系統(tǒng)初態(tài)中各神經(jīng)元膜電位x＞0的神經(jīng)元數(shù)量與神經(jīng)元總數(shù)之比為ρ.本文通過以下方式產(chǎn)生系統(tǒng)演化的初態(tài):給每一個神經(jīng)元賦隨機初值,然后讓神經(jīng)元在無耦合下演化2000個時間單位.我們發(fā)現(xiàn),繼續(xù)延長系統(tǒng)演化時間,ρ的值一般在5.5%左右變化.這時我們可以找ρ滿足給定值的系統(tǒng)狀態(tài)作為D和g不為0時系統(tǒng)演化的初態(tài),該初態(tài)具有混沌的初相位分布,且系統(tǒng)處于混沌態(tài).在下面的數(shù)值模擬中,我只選擇ρ=5%與ρ=6%兩個初態(tài)進行數(shù)值模擬研究.

圖1給出了不同ρ值下g-D參數(shù)平面上的相圖,這是系統(tǒng)從兩個固定初態(tài)演化通過改變耦合強度得到的結(jié)果.從圖1可以看出:不論ρ取什么值,適當選擇興奮性和抑制性耦合強度,系統(tǒng)都會出現(xiàn)單螺旋波、多螺旋波、迷宮斑圖、迷宮斑圖與螺旋波共存現(xiàn)象,個別參數(shù)出現(xiàn)平面波.在耦合強度相同的情況下,不同的系統(tǒng)初態(tài)可能得到不同類型的斑圖,但是有如下相同的規(guī)律:1)當g=0.1時,系統(tǒng)一般不會出現(xiàn)有序結(jié)構(gòu),但是在適當?shù)南到y(tǒng)初態(tài)(如ρ=6%的初態(tài))和興奮性耦合強度下,系統(tǒng)可以偶爾出現(xiàn)單個螺旋波;2)當g=0.9時,系統(tǒng)很容易出現(xiàn)迷宮斑圖、多螺旋波和迷宮斑圖與螺旋波共存現(xiàn)象,適當選擇興奮性耦合強度可以出現(xiàn)平面波;3)當g=0.5—0.7時,適當選擇興奮性耦合強度,系統(tǒng)可出現(xiàn)從邊界向中心傳播的方形波,簡稱向內(nèi)方形波;4)在100組參數(shù)中,圖1(a)和圖1(b)分別有21組和22組參數(shù)出現(xiàn)迷宮斑圖,迷宮斑圖出現(xiàn)概率分別為21%和22%,同理可以得到螺旋波(包括多螺旋波)出現(xiàn)概率分別達到29%和26%,只是占據(jù)整個空間的大螺旋波(即單螺旋波)出現(xiàn)概率遠小于小螺旋波出現(xiàn)概率.向內(nèi)方形波出現(xiàn)概率分別達到9%和11%,平面波出現(xiàn)概率都為2%.可見在均勻長程抑制耦合下,螺旋波出現(xiàn)概率最大,迷宮斑出現(xiàn)概率次之,平面波出現(xiàn)概率最小.將我們的結(jié)果與文獻[21]的結(jié)果對比可以看出,我們使用相同的興奮性神經(jīng)元,引入了抑制耦合,系統(tǒng)就可以自發(fā)出現(xiàn)螺旋波,其效果等效在網(wǎng)絡(luò)系統(tǒng)引入抑制性神經(jīng)元.這些結(jié)果表明,抑制耦合有利于時空斑圖的穩(wěn)定和有序斑圖的產(chǎn)生.

考慮到不同初態(tài),系統(tǒng)出現(xiàn)的現(xiàn)象基本相同,下面就ρ=6%情況詳細介紹觀察到的各種現(xiàn)象.圖2給出了不同耦合強度對應(yīng)的x變量斑圖,所有斑圖都是從相同初態(tài)演化10000個時間單位后的結(jié)果,斑圖亮度與x變量大小成正比.從圖2可以看出,系統(tǒng)可自發(fā)出現(xiàn)各種不同有序斑圖,有3種不同類型的迷宮斑圖,第一種迷宮斑圖如圖2(b)所示,D=0.1對應(yīng)的迷宮斑圖都屬于這類斑圖,這類斑圖形成過程如下:首先系統(tǒng)自發(fā)形成孤子波(圖中十字狀集團),這些孤子波一開始呈規(guī)則分布,而且它們可以自發(fā)左右或上下運動,最后排成行和列運動,從而形成迷宮斑圖.由于每一個孤子波是可以自由運動,當部分孤子波運動方向發(fā)生改變時,孤子波成行運動可以轉(zhuǎn)變?yōu)槌闪羞\動,如圖3所示.當孤子波成行或成列后都朝一個方向運動時,就會形成平面波,如圖2(a)所示.當系統(tǒng)同時存在孤子波成行和列運動,且運動方向不變,就會形成圖2(c)這種形狀的斑圖,我們稱這種穩(wěn)定迷宮斑圖為第二類迷宮斑圖.當孤子波成行和成列后停止運動,就會形成線振源,產(chǎn)生向外傳播的平面波,平面波相碰就會消失,這時就會得到不穩(wěn)定迷宮斑圖,如圖2(d)所示.如果小的孤子波組成大的振源,迷宮斑圖中還會出現(xiàn)靶波,如圖2(e)所示.

圖1 不同ρ值下g-D參數(shù)平面上的相圖 (a)ρ=5%;(b)ρ=6%Fig.1.Phase diagram in the g-D parameter plane for different values of ρ:(a) ρ =5%;(b) ρ =6%.

圖2(f)和圖2(g)中顯示螺旋波與迷宮斑圖共存現(xiàn)象,其產(chǎn)生過程如下:在不穩(wěn)定迷宮斑圖下,小的孤子波成行后以一定概率形成螺旋波,導致螺旋波與迷宮斑圖共存現(xiàn)象出現(xiàn).在本應(yīng)該出現(xiàn)螺旋

波與迷宮斑圖共存的情況下,如果螺旋波的波頭形成比較早,螺旋波可將平面波抑制掉,在一些情況下就會出現(xiàn)圖2(h)和圖2(j)所示的多螺旋波現(xiàn)象.從這些圖中可以看出,系統(tǒng)出現(xiàn)了雙臂螺旋波(參見圖2(h))和螺旋波對(參見圖2(j)),螺旋波對產(chǎn)生了靶波.在另一些情況下,系統(tǒng)自發(fā)出現(xiàn)單螺旋波,如圖2(k)所示.單螺旋波也在系統(tǒng)處于混沌態(tài)下自發(fā)產(chǎn)生,例如在g=0.1情況下自發(fā)出現(xiàn)的螺旋波就屬于這種情形.圖2(l)顯示的斑圖為向內(nèi)方形波斑圖,這種波不斷從邊界產(chǎn)生,然后向中心傳播,波形狀像口字,因此稱為向內(nèi)方形波.

圖2 不同耦合強度下的x變量斑圖 (a)g=0.9,D=0.2;(b)g=0.9,D=0.1;(c)g=0.9,D=0.3;(d)g=0.3,D=0.3;(e)g=0.9,D=1.2;(f)g=0.9,D=1.5;(g)g=0.5,D=0.5;(h)g=0.9,D=1.7;(i)g=0.9,D=1.8;(j)g=0.7,D=1.2;(k)g=0.5,D=1.7;(l)g=0.7,D=1.9Fig.2.Pattern of the variable x for different values of coupling strength:(a)g=0.9,D=0.2;(b)g=0.9,D=0.1;(c)g=0.9,D=0.3;(d)g=0.3,D=0.3;(e)g=0.9,D=1.2;(f)g=0.9,D=1.5;(g)g=0.5,D=0.5;(h)g=0.9,D=1.7;(i)g=0.9,D=1.8;(j)g=0.7,D=1.2;(k)g=0.5,D=1.7;(l)g=0.7,D=1.9.

圖3 在g=0.9和D=0.1情況下不同時刻的x變量斑圖 (a)t=0;(b)t=955;(c)t=9875;(d)t=9905;(e)t=9945;(f)t=10000Fig.3.Pattern of the variable x at different time moments for g=0.9 and D=0.1:(a)t=0;(b)t=955;(c)t=9875;(d)t=9905;(e)t=9945;(f)t=10000.

為了對迷宮斑圖和向內(nèi)方形波演化有直觀印象,我們做j=100這一行和i=100這一列格點x變量的時空斑圖,圖4(a)與圖4(c)給出了在圖2(c)迷宮斑圖參數(shù)下對應(yīng)行與列x變量的時空斑圖,圖4(b)與圖4(d)則是圖2(l)向內(nèi)方形波參數(shù)下對應(yīng)行列x變量的時空斑圖.從圖2(c)和圖4(a)可以看出:圖4(a)左邊有4個波向右運動,右邊也有4個波向左運動,形成粗斜線,因為不同波列對應(yīng)不同斜線,粗斜線朝右傾斜,波列就朝右運動,如果粗斜線朝左傾斜,波列就朝左運動,如果粗線為水平線,則表示波列不作橫向運動.從圖2(c)和圖4(c)可以看出,圖4(c)左邊有6個波向上(j增大方向)運動,右邊也有4個波向下(j減少方向)運動,形成粗斜線,因為粗斜線朝右傾斜,波列就向上運動,如果粗斜線朝左傾斜,波列就向下運動,如果粗線為水平線,則表示波列不在縱向運動.可見圖2(c)迷宮斑圖是由朝不同方向運動的平面波形成的.類似由圖4(b)和圖4(d)可以看出:1)圖2(l)是向內(nèi)方形波斑圖,相當在邊界上自發(fā)形成波源,它產(chǎn)生的波由外至內(nèi)傳播;2)圖中反映短暫出現(xiàn)由內(nèi)向外傳播的波,說明中心自發(fā)形成的振源產(chǎn)生了向外傳播的波,這些波與向內(nèi)傳播的波相遇而消失,因此這種看起來向內(nèi)傳播的波不能稱為反靶波,因為傳播方向不是指向波源.

圖4 不同耦合強度下一行和一列格點的x變量時空斑圖 (a)g=0.9,D=0.3;(b)g=0.7,D=1.9;(c)g=0.9,D=0.3;(d)g=0.7,D=1.9Fig.4.Spatiotemporal pattern of the variable x of a row and a column of grid points for different values of coupling strength:(a)g=0.9,D=0.3;(b)g=0.7,D=1.9;(c)g=0.9,D=0.3;(d)g=0.7,D=1.9.

上述結(jié)果表明,適當選擇耦合強度系統(tǒng)會自發(fā)出現(xiàn)有序結(jié)構(gòu),為了了解系統(tǒng)有序程度,圖5給出了不同參數(shù)下系統(tǒng)方差隨時間的變化,所使用參數(shù)與圖2對應(yīng).從圖5可以看出,系統(tǒng)演化初期,系統(tǒng)處于混沌態(tài),系統(tǒng)方差隨時間大幅度無規(guī)則變化,在系統(tǒng)出現(xiàn)有序結(jié)構(gòu)后,系統(tǒng)方差在某個值附近無規(guī)則振蕩.這表明整個系統(tǒng)雖然處于有序態(tài),但是單個神經(jīng)元振蕩仍處于不規(guī)則振蕩狀態(tài),表現(xiàn)在單個神經(jīng)元振蕩時峰的數(shù)量和振幅在隨時間變化,這可以從圖2斑圖中波臂出現(xiàn)明暗無規(guī)則分布可以看出.從圖5還可以看出,穩(wěn)定迷宮斑圖和平面波的系統(tǒng)方差變化幅度最小,系統(tǒng)方差變化第二小是螺旋波、多螺旋波和部分不穩(wěn)定迷宮斑圖,系統(tǒng)方差變化第三小是向內(nèi)方形波和螺旋波與迷宮斑圖共存態(tài),混沌態(tài)的系統(tǒng)方差變化一般都比較大.對比圖2和圖5可以看出,系統(tǒng)越有序,其方差隨時間變化幅度越小.

上述結(jié)果是選擇兩個初態(tài)得到的結(jié)果,選擇其他不同初態(tài),也得到類似結(jié)果,結(jié)果表明,出現(xiàn)大螺旋波的概率比出現(xiàn)多螺旋波(或小螺旋波)的概率小,這與文獻[15]給出星形膠質(zhì)細胞合胞體的實驗結(jié)果一致,實驗表明:大螺旋波出現(xiàn)概率小,小螺旋波出現(xiàn)概率大,這表明引入抑制耦合的做法是合理的.

圖5 不同耦合強度下系統(tǒng)方差隨時間變化 (a)g=0.9,D=0.2;(b)g=0.9,D=0.1;(c)g=0.9,D=0.3;(d)g=0.3,D=0.3;(e)g=0.9,D=1.2;(f)g=0.9,D=1.5;(g)g=0.5,D=0.5;(h)g=0.9,D=1.7;(i)g=0.9,D=1.8;(j)g=0.7,D=1.2;(k)g=0.5,D=1.7;(l)g=0.7,D=1.9Fig.5.Time evolution of the variance of the system for different values of coupling strength:(a)g=0.9,D=0.2;(b)g=0.9,D=0.1;(c)g=0.9,D=0.3;(d)g=0.3,D=0.3;(e)g=0.9,D=1.2;(f)g=0.9,D=1.5;(g)g=0.5,D=0.5;(h)g=0.9,D=1.7;(i)g=0.9,D=1.8;(j)g=0.7,D=1.2;(k)g=0.5,D=1.7;(l)g=0.7,D=1.9.

4 結(jié) 論

采用HR神經(jīng)元模型,研究了在興奮和抑制耦合共同作用下二維混沌神經(jīng)元網(wǎng)絡(luò)系統(tǒng)從隨機相位分布的初態(tài)演化是否能自發(fā)產(chǎn)生螺旋波等有序結(jié)構(gòu)的問題.研究發(fā)現(xiàn),當抑制性耦合強度比較小時,一般系統(tǒng)不會出現(xiàn)有序結(jié)構(gòu),增加抑制性耦合強度,混沌神經(jīng)元網(wǎng)絡(luò)系統(tǒng)很容易產(chǎn)生迷宮斑圖和螺旋波.有兩種穩(wěn)定迷宮斑圖和一種不穩(wěn)定迷宮斑圖,它們形成的機制略有不同;在給定抑制性耦合強度g＞0.3情況下,逐漸增加興奮性耦合強度,系統(tǒng)首先出現(xiàn)迷宮斑圖,接著系統(tǒng)一般出現(xiàn)螺旋波與迷宮斑圖共存態(tài)或多螺旋波態(tài),當興奮性耦合強度足夠大時,系統(tǒng)還可出現(xiàn)單螺旋波態(tài)和向內(nèi)方形波態(tài).總之,在適當選擇耦合強度下,系統(tǒng)可以自發(fā)出現(xiàn)迷宮斑圖、單螺旋波、雙臂螺旋波、螺旋波對、平面波、靶波和向內(nèi)方形波等有序結(jié)構(gòu),迷宮斑圖出現(xiàn)概率平均為21.5%,螺旋波出現(xiàn)概率平均為27.5%,向內(nèi)方形波出現(xiàn)概率平均為10%,平面波出現(xiàn)概率為2%,這說明抑制性耦合有促進混沌系統(tǒng)產(chǎn)生螺旋波和迷宮斑圖的作用.

本文研究了由抑制性耦合導致的規(guī)則波,采用的是電突觸耦合.在大腦神經(jīng)系統(tǒng)中,神經(jīng)元還可以通過化學突觸耦合和電磁場耦合,研究發(fā)現(xiàn)[12]:在電磁輻射作用下,神經(jīng)元之間能實現(xiàn)同步.因此開展多層網(wǎng)絡(luò)混合耦合模式放電態(tài)的自組織和電磁場輻射對自組織行為的影響研究,有助于人們認識大腦皮層中螺旋波和平面波產(chǎn)生的機制,這也是我們今后要研究的內(nèi)容,通過這些研究有助于治療與螺旋波等有序波有關(guān)的大腦疾病.