吳望生 (長江大學物理與光電工程學院,湖北 荊州 434023)
唐國寧 (廣西師范大學物理科學與技術學院,廣西 桂林 541004)
?
兩耦合HR混沌神經元同步研究
吳望生(長江大學物理與光電工程學院,湖北 荊州 434023)
唐國寧(廣西師范大學物理科學與技術學院,廣西 桂林 541004)
[摘要]基于Hindmarsh-Rose(HR)神經元動力學模型,通過數值計算耦合神經元系統的同步差和峰峰間期(ISI)隨耦合強度的變化研究了2個耦合混沌HR神經元在近鄰反饋耦合下的同步過程。結果表明,初態混沌的HR神經元在耦合中存在著復雜的同步過程;伴隨耦合強度的增大,神經元會交替出現混沌發放、周期發放等放電現象;在達到精確的完全同步之前,會交替出現簇同步、近似同步和完全同步等同步現象;在弱耦合條件下,神經元會出現反同步和去同步。
[關鍵詞]Hindmarsh-Rose神經元;同步;反同步
自17世紀Huygens發現同步現象[1]后的很長一段時間里,有關同步的研究始終是建立在周期運動的基礎之上。20世紀90年代,Pecora和Carro[2,3]首次在實驗室的電子線路上實現混沌同步,從此混沌同步研究在許多學科領域深入展開。特別是在神經科學領域,同步成為當前研究的熱點問題。生理實驗證實,大腦中存在著神經元系統的局部同步模式。1989年Gray等[4]發現,當給予同一個物體的輪廓作為共同刺激時,分布在相互分離模塊的視覺皮層細胞常常會產生同步化放電;1997年Riehle等[5]在猴的運動皮質區發現不同神經元之間出現了峰同步,證實了大腦局部同步區的存在。結合生物神經系統中發現的許多同步現象,Pouget等[6]提出神經系統的大量信息是由神經元系統共同處理。神經元系統是一個由快慢變量組成的多尺度系統,慢變量是神經元的簇行為產生的機制[7],導致神經元會呈現豐富的發放行為,形成各種周期和混沌模式,因而耦合神經元系統的同步過程也非常復雜。研究表明,神經元在達到穩定的完全同步之前,存在著不同類型的同步過程。胡崗等[8]研究了廣義同步和相位同步之間的同步過程,指出相位同步可能優先也可能落后于廣義同步;Dhamala等[9]發現耦合神經元達到完全同步前,會出現神經元簇的同步和峰的同步。目前來看,對不同類型的同步過程仍有待進一步研究,分析耦合神經元系統的復雜同步過程對理解現實的神經元同步行為具有重要的指導意義。為此,筆者選用近鄰反饋耦合,研究了2個初態混沌的Hindmarsh-Rose(HR)神經元的同步過程。
1動力學模型
2個初態混沌的耦合HR神經元的動力學方程[10]如下:
(1)
(2)
(3)
(4)
(5)
(6)
其中,角標1(或2)分別代表HR神經元1(或2),統一簡寫為HR1和HR2; u代表神經細胞的膜電位;v是神經元內部的快速恢復電流;w是慢變調節電流(u和v是快變量,w是慢變量,共同構成一個雙尺度系統);a、b、c、d、γ、s及χ都是系統參數;Iext表示外界直流激勵; k是耦合強度。各參數取值如下:a=1,b=3,c=1,d=5,γ=0.006,s=4,χ=-1.6和Iext=3.0,去暫態后單個神經元進入混沌初態。
在數值模擬中,采用四階龍格庫塔法求解耦合方程組,積分步長Δt=0.001。神經元實現精確的完全同步需要時間依賴耦合強度,為保證2個神經元基本達到同步的漸進態,總積分時間取為10000時間單位。同時為了保證數值計算的準確性,需通過消除一定時間單位的暫態過程而得到穩態結果。筆者通過對2個神經元賦不同的隨機初值(-0.5,0.5),去除500時間單位的暫態,讓各神經元單獨演化進入混沌初態。圖1展示了去暫態時間段內HR神經元的動力學行為,圖1(a)顯示暫態末期2個神經元先后進入混沌發放模式,圖1(b)可看出2個神經元在(ui,wi)相平面上的軌道先后進入混沌吸引子。

圖1 去暫態時間段內HR神經元的動力學行為
為了研究2個神經元的同步程度,將同步差ei(t)及平均同步差δi定義為:
eu(t)=|u1(t)-u2(t)|
(7)
ew(t)=|w1(t)-w2(t)|
(8)
(9)
其中,eu(t)描述快變量的空間同步程度;ew(t)描述慢變量的空間同步程度。即eu(t)描述了t時刻2個神經元膜電位u的同步差,ew(t)描述了t時刻2個神經元慢變調節電流w的同步差。同步差越小,表明同步程度越高。通常同步差是振蕩式下降,根據其大小,一般不容易判斷不同參數下神經元的同步程度,所以筆者用δi描述去掉同步暫態后的一段時間內神經元之間的平均同步情況。不失一般性,選擇時間范圍[9000,10000]。當δi=0時,2個神經元實現精確的完全同步,然而有限的時間內δi一般是逐步趨于零,并且快變量δu和慢變量δw的同步速率不同。筆者設定,當快變量和慢變量均小于等于0.001時,即δi≤10-3時,則2個神經元達到完全同步;當快變量和慢變量其中之一小于等于0.001時,即δu≤10-3或δw≤10-3時,2個神經元達到近似同步。
2數值模擬
采用近鄰反饋耦合,對耦合混沌的神經元的峰峰間期(ISI)和同步差進行數值計算。圖2給出了HR神經元1的峰峰間期(ISI)隨耦合強度k變化的分岔圖,圖3給出了耦合混沌神經元快變量的平均同步差δi隨耦合強度k的變化關系。由圖2和圖3可以看出,隨著耦合強度k的逐步增大,HR神經元呈現豐富的放電活動。隨著耦合強度的增加,圖2被交替分割成多個帶狀密集區和線型軌跡區,表明混沌耦合的神經元在混沌發放和周期發放2個狀態反復切換,直至進入混沌發放達到完全同步。

圖2 HR1的ISI隨k變化的分岔圖 圖3 平均同步差δi隨k的變化
圖2的第1段線型軌跡區對應耦合強度區間0.008≤k≤0.028,表明神經元HR1由開始的混沌發放模式變成了規則的周期發放,神經元HR2-ISI分岔圖在此區間也有同類現象,這表明耦合神經元發生了簇同步。不失一般性,取耦合強度k=0.024來分析神經元的簇同步狀態,圖4展示了2個耦合神經元膜電位的時間歷程和在(u1,w1)和(u1,u2)平面上所對應的相圖。結合圖4(a)和(c)可以發現,2個耦合神經元的周期簇同時發生,出現了簇同步;而圖4(b)顯示峰錯位和圖4(d)中相位定位在二四象限,都表明在神經元簇里的峰是不相關的,未發生峰同步。

圖4 k=0.024時,耦合神經元簇同步狀態的分析圖
從圖2中可以清楚的看到,當耦合較弱且參數范圍較小時,耦合強度可以抑制混沌而變為周期的發放行為。但結合圖3和圖4,雖然神經元由于耦合強度增大出現倍周期軌道,并出現簇同步,但2個神經元之間并未達到峰同步,直到k大于0.42后,2耦合神經元才出現峰同步現象。這說明混沌發放被抑制變為周期發放,出現簇同步,只是峰同步前的一個中間過程,還需要繼續增大耦合強度,神經元之間才有可能出現峰同步,并最終達成完全同步。
隨著耦合強度逐步增強,周期發放的神經元又進入混沌狀態,且不同步,但當0.151≤k≤0.180時,2個耦合神經元隨耦合強度的變化交替出現倍周期軌道、混沌軌道,不同步和反同步等復雜的過渡模式。具體來說,當0.151≤k≤0.157時,神經元會出現周期8的倍周期軌道,且出現反相的周期性同步狀態。下面,取k=0.151來分析倍周期和反同步現象(見圖5)。

圖5 k=0.151時,耦合神經元反同步狀態的分析圖
從圖5(a)中可以看到,耦合神經元HR1出現了周期8的倍周期軌道,比圖4(c)中周期數增加了,雖然從混沌軌道變成了周期軌道,但2個耦合神經元并未出現如圖4(a)和(b)所示的簇同步,而是出現了反同步狀態;從圖5(b~d)中可以清楚看到神經元3個快慢變量均處于反相的周期性同步狀態。經過分析,反同步所對應的平均同步差δi均出現較大的變化,如膜電位同步差δu>1就會出現反同步,這正好解釋了圖3在0.151≤k≤0.180區間內出現的突然增大的同步差峰包,其δu值均大于1,說明在峰包對應的耦合強度范圍內,神經元出現了反相的周期性同步。
當耦合強度離開這個反同步區,神經元的反同步被破壞,神經元重新由周期軌道進入混沌軌道,直到耦合強度k≥0.422時,2個混沌耦合神經元又出現近似同步。隨著耦合強度的的進一步增加,當耦合強度位于0.422≤k≤0.434時,耦合神經元在這個區間存在著復雜的過渡模式。

圖6 對數化后平均同步差δi隨k的變化
具體的表現行為通過圖6可以分析清楚。將同步差δi對數化后可清晰看出,隨著耦合強度的增加,同步差震蕩式下降直至為0。藍色虛線表示δi=10-3,為完全同步分界線,線下為完全同步區。在δi與藍色虛線交界處,即耦合強度位于0.422≤k≤0.434時,有時δu<10-3而δw>10-3或δu>10-3而δw<10-3,說明出現了近似同步;有時快慢變量的同步平均差均小于0.001,即δi<10-3,發生了完全同步。上述情況表明,耦合強度較大時,神經元在完全同步之前,即可能出現峰同步(即近似同步),也可能出現簇和峰均同步(即完全同步),而且上述同步現象隨耦合強度的變化而交替出現,表明耦合混沌神經元同步之前的過渡狀態中包含著復雜的放電信息。
為了進一步分析這些復雜的同步過程,考察完全同步和近似同步的不同之處,圖7在不同的耦合強度下研究神經元快慢變量的各種相圖,其中黑色軌跡對應k=0.423,紅色軌跡對應k=0.433。圖7(a)和(b)顯示,黑色軌跡的相位與一三象限角平分線完全重合,表明k=0.423時,神經元出現了完全同步。紅色軌跡的相位則定位在一三象限角平分線附近的小鄰域內,這表明k=0.433時,神經元產生了近似同步且處于混沌狀態。需要說明的是,在耦合強度0.422≤k≤0.434的這個區間內,沒有發現之前弱耦合狀態下,神經元的峰不同步而簇同步的周期放電現象,神經元不論同步與否,均處于混沌放電狀態。

圖7 神經元在不同耦合強度下的平面相圖
從圖3和圖6中2條曲線完全歸零的位置可查出,當耦合強度k≥0.474時,平均同步差δi=0,這時2個耦合神經元達到了完全的峰同步和簇同步。圖8顯示,當k=0.474時,2個耦合神經元的混沌放電時間歷程圖8(a)和混沌吸引子圖8(b)完全重合,表明神經元達成完全同步后仍然保持了穩定的混沌狀態,結合圖1說明初態隨機的2個混沌神經元通過近鄰反饋耦合達成了精確的完全同步。

圖8 k=0.474時,耦合神經元完全同步狀態的分析圖
3同步過程分析
HR神經元是一個具有快慢變量的雙尺度系統,其中快變量產生了神經元的峰,慢變量導致了神經元的簇的形成,在不同的參數條件下,呈現出豐富的動力學行為,如峰峰間期分岔、混沌吸引子等。下面,筆者來比較不同耦合條件下神經元在(u,w)相平面上的吸引子變化情況。
圖1(b)顯示了未耦合狀態下,神經元的混沌吸引子;圖4(c)展示耦合強度k=0.024時,神經元HR1的周期5軌道;圖5展現了k=0.151時,神經元HR1的周期8軌道;圖8則表明神經元達到精確同步后仍然處于混沌吸引子中。隨著耦合強度的逐步增大,耦合神經元吸引子軌道由混沌先后過渡到周期5、周期8,然后又進入混沌狀態,直至以混沌狀態進入精確的完全同步。
對于這種吸引子變化的原因可以用文獻[11]中耦合的Lorenz系統的動力變化的解釋:在弱耦合的情況下,混沌的耦合破壞了原來Lorenz系統雙圈吸引子來回切換所遵循的一致性,結果造成原來吸引子的破壞而產生了新的更復雜的吸引子。上述解釋對于Lorenz系統是成立的,對同樣具有多圈吸引子的其他耦合系統也是成立的,因此在筆者所研究的HR神經元模型中,神經元在給定的初始參數下呈現了如圖1(b)所示的多圈混沌吸引子,隨著耦合強度的增加,混沌軌道被耦合破壞而進入倍周期軌道(見圖4(c)和圖5(a)),繼續增大耦合強度,神經元又脫離周期軌道而進入混沌,直至最終達成精確的完全同步仍然保持混沌狀態。
耦合的多尺度系統的混沌同步與快慢子系統的動力行為是緊密相關的,耦合神經元達到完全同步前,會出現神經元簇的同步和峰的同步[9]。筆者研究發現,耦合混沌的HR神經元隨著耦合強度的增加,會在弱耦合區間出現簇同步,即慢變量的同步(見圖4(a)和(b)),慢變量的同步隨耦合強度繼續增加而被破壞,直至耦合強度大到能使神經元出現峰的同步和簇的同步,即完全同步,意味著快慢變量都達到同步。因此,2個耦合混沌的HR神經元采用近鄰反饋耦合,隨著耦合強度的增加,會先后出現簇同步、峰同步直至完全同步。耦合強度的增加除了可以改變神經元的吸引子,還可以影響神經元的放電節律形成反同步。如圖5所示,當在某些耦合強度區間內出現周期性吸引子時,2個神經元的快慢變量如膜電位和慢電流會出現反相的周期性同步,事實上在這個區間附近,2個耦合神經元隨耦合強度的變化交替出現倍周期軌道、混沌軌道,不同步和反同步等復雜的過渡模式。
4結論
研究了采用近鄰反饋耦合2個混沌HR神經元的同步問題。在數值模擬中,采用四階龍格庫塔法求解耦合方程組,通過計算耦合神經元系統的同步差和峰峰間期隨耦合強度的變化,發現在耦合混沌的HR神經元模型中,存在著復雜的同步過程。
1)隨著耦合強度的增大,神經元的吸引子會被改變,會交替出現混沌發放、周期發放等放電現象,表明混沌軌道被耦合破壞而進入倍周期軌道,繼續增大耦合強度,神經元又脫離周期軌道而進入混沌,直至最終達成精確的完全同步仍然保持混沌狀態。
2)在達到精確的完全同步之前,隨著耦合強度的增大,神經元會交替出現混沌發放、周期發放等放電現象,神經元系統會依次出現神經元簇的同步、峰的同步,最終達成完全同步。
3)耦合強度的增加除了可以改變神經元的吸引子,影響神經元的放電節律,甚至在某些弱耦合區間,會導致神經元系統出現反相的周期性同步,抑制同步的發生。
[參考文獻]
[1]Huygens C.Horologium oscillatorium [M]. Parisiis, France, 1673.
[2] Pecora L M, Carroll T L. Synchronization in chaotic systems [J]. Physical Review Letters, 1990, 64(8): 821~824.
[3]Pecora L M, Carroll T L. Driving systems with chaotic signals[J]. Physical Review A, 1991, 44(4): 2374~2383.
[4]Gray C M, K?nig P, Engel A K, et al. Oscillatory responses in cat visual cortex exhibit inter-columnar synchronization which reflects global stimulus properties[J]. Nature, 1989, 338(6213): 334~337.
[5]Riehle A, Grün S, Diesmann M, et al. Spike synchronization and rate modulation differentially involved in motor cortical function[J]. Science (New York, N.Y.), 1997, 278(5345): 1950~1953.
[6]Pouget A, Dayan P, Zemel R. Information processing with population codes[J]. Nature Reviews Neuroscience, 2000, 1(2): 125~132.
[7]王青云, 石霞, 陸啟韶. 神經元耦合系統的同步動力學[M]. 北京: 科學出版社, 2008: 38~50.
[8] Zheng Z G, Hu G. Generalized synchronization versus phase synchronization[J]. Physical Review E, 2000, 62(6): 7882~7885.
[9] Dhamala M, Jirsa V K, Ding M Z. Transitions to synchrony in coupled bursting neurons[J]. Physical Review Letters, 2004, 92(2): 02810.
[10]Hindmarsh J L, Rose R M. A model of neuronal bursting using three coupled first order differential equations[J]. Proceedings of the Royal Society of London Series B, 1984, 221(1222): 87~102.
[11]Liu Z H, Lai Y C, Matias M A. Universal scaling of Lyapunov exponents in coupled chaotic oscillators. [J]. Physical Review E, 2003, 67: 045203(R).
[編輯]張濤
[文獻標志碼]A
[文章編號]1673-1409(2016)07-0022-06
[中圖分類號]O415.5
[作者簡介]吳望生(1979-),男,碩士,講師,現主要從事統計物理方面的教學與研究工作;E-mail:wshwu@yangtzeu.edu.cn。
[基金項目]國家自然科學基金項目(11165004)。
[收稿日期]2015-11-27
[引著格式]吳望生,唐國寧.兩耦合HR混沌神經元同步研究[J].長江大學學報(自科版),2016,13(7):22~27.