Morris-Lecar 系統(tǒng)中通道噪聲誘導(dǎo)的自發(fā)性動作電位研究1)

李揚 劉先斌

(南京航空航天大學(xué)機械結(jié)構(gòu)力學(xué)及控制國家重點實驗室,南京 210016)

引言

在生物物理學(xué)中,有越來越多的現(xiàn)象源于分段確定性動力系統(tǒng)和連續(xù)時間馬氏過程的耦合作用,其數(shù)學(xué)模型為隨機混合系統(tǒng)(stochastic hybrid system).例如遺傳開關(guān)現(xiàn)象[1-2],電力驅(qū)動的細(xì)胞內(nèi)運輸[3-4],隨機神經(jīng)元網(wǎng)絡(luò)[5-9]和隨機離子通道等等[10-12].此時,由于連續(xù)狀態(tài)和離散狀態(tài)之間的耦合特性,使得混合系統(tǒng)比隨機微分方程更難處理.雖然可以通過采用擴散近似方法得到相應(yīng)混合系統(tǒng)的擴散過程近似,但是眾所周知,如果考慮亞穩(wěn)態(tài)動力學(xué)行為,除了在接近分岔極限[13]的特殊情況之外,該方法通常會產(chǎn)生顯著的誤差[1,14-15].此外,恰恰是由于開通道數(shù)的離散特性,混合系統(tǒng)比連續(xù)系統(tǒng)更符合實際物理模型.因此,為了揭示實際模型的內(nèi)在動力學(xué)行為,直接處理原始的混合系統(tǒng)是必要且本質(zhì)的.

很多混合系統(tǒng)均可表現(xiàn)出可激性(excitability)特征,此時即便是很小的隨機擾動也可能對其動力學(xué)行為產(chǎn)生深遠(yuǎn)的影響,如自發(fā)性動作電位(spontaneous action potential)現(xiàn)象.一般意義下,如果偏離靜息態(tài)的一個小擾動能夠?qū)е缕潆娢辉诜祷仂o息態(tài)之前出現(xiàn)一個大的偏差,即稱產(chǎn)生了一個動作電位[16].通常,自發(fā)性動作電位可以分為兩個階段:初始階段和興奮階段.初始階段是指系統(tǒng)在噪聲驅(qū)動下從靜息態(tài)到運動到某個閾值(threshold)的過程,因此可視為噪聲誘導(dǎo)的離出(exit)事件或逃逸事件.這個階段,由于噪聲強度較弱,將花費相當(dāng)長的時間.一旦達到閾值,系統(tǒng)將沿著確定性流演化,直至返回靜息態(tài),此為興奮階段.與前者相比,興奮階段是暫態(tài)的,其運動軌跡僅受確定性動力學(xué)運動規(guī)律控制.因此,只要能夠準(zhǔn)確刻畫初始階段的離出行為,便可清晰地描述整個動作電位.

事實上,Freidlin 和Wentzell[17]所建立的大偏差理論(large deviation theory)已被廣泛應(yīng)用于分析弱噪聲對具有多種極限集—-不動點[18-19]、極限環(huán)[20-21]和奇怪吸引子[22-23]等系統(tǒng)的長期動力學(xué)行為的影響.依據(jù)大偏差理論中定義的作用(量)泛函(action functional)的概念,通過它在任一絕對連續(xù)函數(shù)處的值能夠得到任一隨機軌線通過該函數(shù)附近時其難易程度的量度.根據(jù)自吸引子出發(fā)所有路徑集合上作用泛函的駐值,便可以確定最優(yōu)的離出路徑,即系統(tǒng)以最大的概率沿著該路徑移動.

基于大偏差理論,人們對噪聲誘導(dǎo)的逃逸問題進行了大量的研究,但主要結(jié)果限于閾值具有明確定義的情況,一般為系統(tǒng)鞍點的穩(wěn)定流形,此時最優(yōu)路徑則為連接吸引子和鞍點的異宿軌道[24].然而遺憾的是,很多系統(tǒng)在相空間中并沒有鞍點.因此為了研究這些離出問題,必須定義合適的閾值.近期,Chen等[25]和Khovanov 等[26]通過選擇虛擬分界線作為FitzHugh-Nagumo 系統(tǒng)的閾值,描述了其離出方案和有限噪聲效應(yīng).Carson 等[27]研究了通道漲落引起的Hodgkin-Huxley 模型的放電事件,并通過一維郎之萬(Langevin)近似計算了其平均放電頻率.James 等[28]將該模型簡化為只保留簡單的鈉離子通道噪聲的模型,并計算出平均首次離出時間的精確表達式.基于Hodgkin-Huxley 模型同時保留鈉離子和鉀離子通道噪聲的另一種簡化系統(tǒng),此時其相應(yīng)的確定性模型只有一個穩(wěn)定的不動點,Brooks 等[12]發(fā)現(xiàn)了擬環(huán)(quasi-cycle)的出現(xiàn).Newby[11]利用WKB 近似方法對同一系統(tǒng)進行了研究,得到了離出路徑滿足的Hamilton 表述.本文以他們的工作為基礎(chǔ),重點研究了該模型在自定義閾值下的自發(fā)性動作電位,并討論了通道數(shù)對可激系統(tǒng)動力學(xué)行為的影響.

迄今為止,用于計算擬勢的主要方法是WKB 近似,它被廣泛應(yīng)用于研究擴散過程、生滅控制方程和隨機混合系統(tǒng)中的逃逸事件.此外,本文采用的數(shù)值方法包括幾何最小作用量法(geometric minimum action method,GMAM)[29]和有序逆風(fēng)法(ordered upwind method,OUM)[30].它們均基于幾何作用概念,即將作用泛函的時間參數(shù)化轉(zhuǎn)化為弧長,從而可將關(guān)于無窮時間的積分簡化為有限長度的積分.GMAM 可用于求解具有固定端點的最優(yōu)路徑和它的作用泛函,而OUM 可用于計算離散化的相空間節(jié)點上的擬勢.

本文具體組織如下:第2 節(jié)介紹了隨機混合Morris-Lecar 模型,第3 節(jié)中給出了極值路徑滿足的輔助Hamilton 系統(tǒng)和改進的混合系統(tǒng)的Monte Carlo模擬算法.通過定義虛擬分界線為閾值,第4 節(jié)研究了系統(tǒng)的自發(fā)性放電現(xiàn)象.在第5 節(jié)中,通過改變鈉離子和鉀離子通道數(shù),我們詳細(xì)討論了不同噪聲組合對系統(tǒng)擬勢和最優(yōu)路徑的影響.結(jié)論見第6 節(jié).

1 隨機模型

本文所研究的模型最初由Hodgkin 和Huxley[31]提出并用于描述魷魚巨軸突的膜電位,Newby[11]將其簡化為只考慮單通道鈉離子和鉀離子的情形,出于論文的完整性考慮,我們簡述其模型結(jié)構(gòu).首先,膜電位v的演化受控于離子電流和外加刺激電流Iapp,即

其中,n=0,1,2,···,N和m=0,1,2,···,M分別表示鈉離子和鉀離子的開通道數(shù),函數(shù)fi(v)=gi(vi-v),i=Na,K,leak 分別決定了鈉離子和鉀離子電流和漏泄電流.鈉離子和鉀離子開閉狀態(tài)的切換取決于

其中,aNa(v)=exp[4(γNav+κNa)],bNa=1,aK(v)=exp(γKv+κK)和bK(v)=exp(?γKv-κK).由此可見,電位v的演化受控于方程(1).在兩次相鄰跳躍之間,其值是確定性的,但會通過影響轉(zhuǎn)移率反過來修正隨機的跳躍時間,因此它們之間是耦合的.

根據(jù)耦合過程的馬氏性,其概率密度函數(shù)p(v,m,n,t)滿足下面 Chapman-Kolmogorov (CK)方程[9]

其中的算子定義如下

而轉(zhuǎn)移率分別為

假設(shè)βNa和M均很大且滿足此外,可將w=m/M視為第二個狀態(tài)變量,在弱噪聲極限ε→0 時,系統(tǒng)恢復(fù)為確定性的Morris-Lecar 方程[32]

其中,x∞(v)={1+tanh[2(γNav+κNa)]}/2,w∞(v)=[1+tanh(γKv+κK)]/2 和τw(v)=2/βKcosh(γKv+κK).βK和βNa決定了鉀離子和鈉離子通道開閉的快慢,通常有βK?βNa,意味著鈉離子通道的開閉遠(yuǎn)快于鉀離子.

在研究隨機模型之前,需要先考察由方程(6)定義的確定性系統(tǒng)的動力學(xué)行為,方程中的參數(shù)見附錄.如圖1(a)所示,棕色和粉紅色曲線分別代表v零傾線和w零傾線,它們相交于一個全局穩(wěn)定的平衡點SN,圖中的小箭頭指明了向量場的變化趨勢.此外依據(jù)四階龍格?庫塔法,圖中給出了5 條具有代表性的確定性路徑,即圖中的綠色曲線.通過觀察,可以發(fā)現(xiàn),雖然所有的路徑最終都收斂于SN,但從不同的初始點出發(fā),它們之間存在著本質(zhì)上的區(qū)別.其中一些迅速而直接地衰減到SN,而另外的則先向右移動,直至到達v零傾線的右支,然后向左轉(zhuǎn),這樣便產(chǎn)生了一個動作電位.圖中的紅色粗體曲線稱為虛擬分界線(ghost separatrix)[33],它是通過v零傾線的頂點P時間后向積分得到的.可以看出,初始點位于此分界線左側(cè)的路徑,運動的振幅較小,而右側(cè)對應(yīng)于大振幅運動.在圖1(a)中,右側(cè)三條路徑的初始位置之間的差異僅為0.01,很難區(qū)分,但在這些路徑的中間位置上則出現(xiàn)了宏觀上的有限距離.因此,這條曲線起到了擬閾值的作用[16].

圖1 向量場、確定性軌線和隨機軌線.棕色和粉紅色曲線分別表示v 零傾線和w 零傾線,SN 為不動點,綠色曲線是確定性路徑,紅色粗體線指代分界線,小箭頭代表向量場的方向.一條代表性隨機軌道表示為灰色鋸齒線,里面包含一次自發(fā)性動作電位Fig.1 Vector field,deterministic and stochastic trajectories.Brown and pink curves indicate the v-nullcline and w-nullcline respectively.SN denotes the fixed point,the green curves show the deterministic paths,the red bold line suggests the ghost separatrix,and the small arrows represent directions of vector field.A representative stochastic trajectory has been denoted by gray zigzag line,containing a spontaneous action potential

在圖1(b)中,灰色鋸齒形曲線是一條具有代表性的表示自發(fā)性動作電位的隨機路徑.由于SN 的吸引性和較弱的噪聲強度,系統(tǒng)將長期在SN 附近徘徊.隨后在某個偶然的時刻,它會躍過閾值,然后再次沿著確定的路徑返回到SN.這整個過程構(gòu)成一個自發(fā)性動作電位.

2 理論分析與數(shù)值模擬方法

2.1 WKB 近似

隨機混合Morris-Lecar 的WKB 近似的結(jié)果已由Newby 導(dǎo)出[11].出于論文完整性的考慮,我們簡述其主要結(jié)果以便于后面的理論分析及數(shù)值計算,并給出幾點重要的評述.

在弱噪聲條件下,假設(shè)平穩(wěn)概率密度具有下面近似的形式

這里C(v,w)是前因子項,本文不予討論.r(n|v,w)是已知v,w后n的條件概率分布.W(v,w)稱為擬勢,它度量了系統(tǒng)在隨機激勵下抵達相空間中某一位置的困難程度.將方程(7)代入平穩(wěn)CK 方程并合并所有ε 的最高階項系數(shù)得到Hamilton-Jacobi 方程

根據(jù)特征線法可知,這一方程的解由一個輔助Hamilton 系統(tǒng)生成

通過積分可以得到每條特征線上的擬勢,即

對于WKB 近似的形式及結(jié)果,我們認(rèn)為有必要給出以下三個重要的評述:

首先,相比于擴散過程,方程(7)的形式多了一項r(n|v,w),這是因為系統(tǒng)中存在離散變量.由此,概率密度的歸一化條件不能像擴散過程那樣進行直接積分且其積分值等于1,而必須先對離散變量求和然后再對連續(xù)變量進行積分.

其次,若在方程(9)中令p=0,則輔助Hamilton方程恢復(fù)為確定性的Morris-Lecar 系統(tǒng),這意味著動量p度量了離出的極值路徑與確定性流之間的距離.就是說,它刻畫了噪聲在驅(qū)動系統(tǒng)離出的過程中所扮演的角色.因此,與隨機微分系統(tǒng)的情形類似,可將動量值視為噪聲影響的一個指示器.

第三,作者已嚴(yán)格證明了混合系統(tǒng)一定不滿足細(xì)致平衡條件的事實[14],這一結(jié)論與擴散過程的相關(guān)結(jié)論有著本質(zhì)的區(qū)別.對于混合系統(tǒng),離出路徑不會出現(xiàn)可以簡單地視為確定性流的時間反向的情況.于是,尖點和焦散線等奇異性[34]的出現(xiàn)會更為普遍,這源于由極值路徑所追溯的流形出現(xiàn)折疊的特性,并因此導(dǎo)致擬勢的多值性[35].

2.2 Monte Carlo 模擬

基于Gillespie 算法[36]的Monte Carlo 模擬的基本原理是通過計算隨機跳躍時間和狀態(tài)轉(zhuǎn)移來演化樣本軌跡.Bressloff[8]和Newby[11]分別描述了他們算法的主要步驟,并應(yīng)用于神經(jīng)網(wǎng)絡(luò)和Morris-Lecar等模型.然而,他們的方法得以實施均有一個前提條件,即要求方程˙v=Iion(v,m,n)對固定的m,n有顯式解,因此很難推廣到更復(fù)雜的系統(tǒng).在此我們提出一種適用于幾乎所有隨機混合系統(tǒng)的新算法,并且此方法能夠計算出足夠精確的跳躍時間,從而可以近似到任何期望的精度.

為了計算下一次跳躍時間,我們需要考慮四種可能的轉(zhuǎn)移:n→n±1 和m→m±1.在開通道數(shù)的兩次相鄰跳躍之間,電位的演化受控于確定性動力學(xué)此處求解的主要的困難來源于轉(zhuǎn)移率和電位的耦合關(guān)系.為了解決這個問題,定義

從而可以得到具有初始條件v(0)=v0和G(0)=0 的微分方程組

只需要提供時間步長h和兩個均勻分布隨機變量Y和r,當(dāng)下列條件滿足時

跳躍時間和電位的估計就可由τ≈kh和v(τ)≈v(kh)給出.為了確定所發(fā)生的轉(zhuǎn)移是4 種可能中的哪一種,變量a定義為滿足下列不等式的最小整數(shù)

其中a=1,2,3,4 分別對應(yīng)于以下4 種轉(zhuǎn)移行為:n→n+1,n→n-1,m→m+1,m→m-1.至此可以看出,通過將函數(shù)F的積分定義為一個新的變量,不需要像前人一樣解析計算出v的表達式,而是將其轉(zhuǎn)化為一個微分方程組來同時進行數(shù)值積分,因此這一方法是通用的.

3 不同閾值下的自發(fā)性動作電位

如前所述,自發(fā)性動作電位可分為兩個階段:初始階段和興奮階段.由于噪聲強度較弱,系統(tǒng)將在吸引子附近游蕩很長時間,隨后,它將離開靜息態(tài),并在某個隨機時間穿過閾值.一旦到達閾值,它將沿著確定性軌跡返回到吸引子.此時便意味著一個動作電位出現(xiàn)了.事實上,初始階段的本質(zhì)是:在給定閾值的條件下,由噪聲誘導(dǎo)的從吸引域中的逃逸事件(離出問題).

迄今為止,我們應(yīng)該能夠意識到閾值的重要性！因為它決定了什么樣的動力學(xué)行為意味著一次放電.如果相空間中存在鞍點,那么其穩(wěn)定流形自然起到閾值的作用,最優(yōu)路徑就是連接吸引子和鞍點的路徑,即異宿軌道.穿過鞍點后,系統(tǒng)會沿著鞍點的不穩(wěn)定流形運動.然而遺憾的是,對于不具有鞍點的系統(tǒng)(如本文中的模型),很難定義這樣一個閾值.盡管如此,本文的系統(tǒng)與具有鞍點的系統(tǒng)有很多相同的動力學(xué)特征.為了了解這一事實,我們需要首先選擇一個具有規(guī)范定義的閾值來觀察初始逃逸事件.考慮到第2 節(jié)中關(guān)于模型的確定性動力學(xué)的分析,虛擬分界線在自發(fā)性動作電位中起著重要作用.如果系統(tǒng)受噪聲的驅(qū)動下穿過分界線,那么系統(tǒng)之后將幾乎沿著確定性軌跡繞一個大圈回到不動點,形成一個電位的尖峰.從這一點來看,把此分界線看成閾值是合理的.

在考察離出路徑之前,我們先來考察系統(tǒng)的擬勢.依據(jù)OUM 可計算出相關(guān)區(qū)域的擬勢,其結(jié)果由圖2 所示的等高線和三維圖給出.值得注意的是,該數(shù)值方法是通過從不動點開始逐次更新相鄰節(jié)點的值來計算規(guī)則網(wǎng)格上的擬勢,其結(jié)果必然是擬勢的全局最小值.對于此系統(tǒng),用800×800 的網(wǎng)格覆蓋區(qū)域[?0.5,1.8]×[0,0.6]足夠精確.圖2(b)所示,SN 附近存在一個很深的勢阱,在勢阱右上方有一個很大的平坦區(qū)域,擬勢幾乎保持不變,隨后在遠(yuǎn)處迅速增長.我們注意到,圖2(a)中在w=0 處有一個擬勢的截痕,這是因為w表示鉀離子的開通道數(shù)的比例,不可能是負(fù)的.

圖2 利用OUM 在800×800 的網(wǎng)格上計算的擬勢Fig.2 Quasi-potential calculated by OUM on a 800×800 grid

由輔助Hamilton 系統(tǒng)(9)的解所刻畫的軌跡,在t0→?∞時起始于SN 的條件下,在四維擴展相空間中張成二維不穩(wěn)定拉格朗日流形.這些路徑在坐標(biāo)空間上的投影決定了極值路徑的模式.在擴展的相空間中,由于增加了兩個不穩(wěn)定維度,原來的穩(wěn)定不動點轉(zhuǎn)化為鞍點.為了描述極限路徑的模式,我們使用作用量圖(action plot)[23]方法求解輔助Hamilton系統(tǒng)的初值問題.在二維拉格朗日流形上選取均勻分布于以SN 為中心半徑為0.001 的圓周上的點為初始點,并且每條極值路徑可由初始點的角位置參數(shù)化,與擬勢的特征線方程(10)一起積分,直到閾值后終止.其中幾條具有代表性的極值路徑如圖3(a)所示,它們起始時極為接近,在到達分界線前散開.當(dāng)極值路徑觸碰到分界線時,該點的擬勢值便被記錄,因此我們能夠得到分界線上所有點的擬勢,其與w的函數(shù)關(guān)系由圖3(b)給出.此外,利用插值原理,前面的OUM 方法的結(jié)果也可以用于計算閾值處的擬勢,其結(jié)果如圖3(b)中的黑線所示.可以看出,兩種方法的結(jié)果吻合得相當(dāng)好.

圖3 選取分界線為閾值Fig.3 The separatrix chosen as threshold

GMAM 是另一種計算擬勢的方法.其主要思想是在兩個端點給定的條件下計算連接此兩點的最優(yōu)路徑及其上的作用泛函.我們將一個端點選在SN,另一個位于分界線,通過不斷改變分界線上的點以計算不同位置的擬勢值.為了提高計算效率,這一方法的誤差會比前兩個結(jié)果稍大一些,但仍在可接受范圍內(nèi).從圖3(b)中可以看出,分界線上的擬勢在w≈0.023 (對應(yīng)于左圖中的S0點)附近存在明顯的最小值.對于雙穩(wěn)態(tài)系統(tǒng),鞍點一般是邊界上擬勢的最小值點,因此離出行為總發(fā)生在鞍點.對于本文的模型,盡管不存在鞍點,但是S0有著類似于鞍點的作用.

為了驗證我們之前的理論結(jié)果,現(xiàn)采用歷程概率分布方法模擬漲落路徑的分布.這一概念由Dykman[37]首先提出,歷程概率分布Ph(x,t|xf,tf)被定義為三時刻轉(zhuǎn)移概率Ph(xf,tf;x,t|xi,ti)與兩時刻轉(zhuǎn)移概率Ph(xf,tf|xi,ti)之比,初始時刻ti被設(shè)為?∞,從而ti和xi均可從表達式中去掉.這一分布包含了系統(tǒng)到達xf之前關(guān)于系統(tǒng)演化的所有信息.根據(jù)大偏差理論,系統(tǒng)會以壓倒性的概率通過最優(yōu)路徑的鄰域,因此歷程概率分布將在最優(yōu)路徑附近達到峰值.為了觀察系統(tǒng)的動力學(xué)行為,我們固定終點時間為tf=0 并向前追溯.基于第2 節(jié)中改進的Monte Carlo算法,可以得到眾多從SN 出發(fā)直至閾值的隨機路徑.根據(jù)核密度估計,最終得到歷程概率分布的結(jié)果如圖4(a)所示.自這些軌跡的樣本集合中,提取每一時刻歷程概率分布的最大值位置,如圖4(b)所示,其結(jié)果與GMAM 計算的最優(yōu)路徑的理論值吻合得相當(dāng)好.此外,這一結(jié)果證實了我們之前的觀點,即離出位置集中在S0上.

前面的研究表明,將自SN 出發(fā)到分界線的最優(yōu)路徑作為動作電位的初始階段是合理的.然而,為了說明以分界線作為閾值的合理性,還應(yīng)該觀察系統(tǒng)在越過分界線后是否存在大振幅運動,并且其結(jié)果是否與模擬結(jié)果一致.為了觀察系統(tǒng)的放電特性,將離出路徑進一步延伸,允許最優(yōu)路徑觸及v零傾線的右分支,看看本文的分析是否成立?

圖4 理論計算與數(shù)值結(jié)果的最優(yōu)路徑對比Fig.4 The comparison for the optimal paths between the theoretical and numerical results

首先,通過對隨機軌道樣本進行積分直到v零傾線右支,經(jīng)統(tǒng)計可得其歷程概率分布如圖5(a)所示.提取該分布的峰頂線,所得結(jié)果如圖5(b)所示,曲線在穿過分界線的時候非常接近S0.另一方面,由GMAM 給出的極值路徑的理論解自SN 開始,至S0結(jié)束,進一步可延伸至v零傾線的N0點.從圖中可以看出,理論解與數(shù)值模擬結(jié)果非常吻合,這說明將分界線作為閾值不僅對初始階段是合理的,對興奮階段同樣合理.圖5(b)中所示的小箭頭表示確定性向量場的方向,可以看出,在穿過分界線后,最優(yōu)路徑與確定性流的方向很接近.我們換個角度觀察動量隨時間變化情況,即圖5(c).可以看出最優(yōu)路徑穿越閾值后,動量值迅速衰減.如前所述,動量本質(zhì)上是對漲落路徑偏離確定性流的程度的度量,動量的衰減這一現(xiàn)象正說明離出路徑接近于確定性軌線這一事實.此外,擬勢的值源于動量的積分,而擬勢在分界線右邊變化緩慢,這也與OUM 的結(jié)果一致.因此,當(dāng)系統(tǒng)到達N0后,幾乎沿確定性軌跡返回到SN,發(fā)生一次完整的放電現(xiàn)象.

圖5 初始階段和興奮階段的最優(yōu)路徑及其動量、作用量Fig.5 The optimal path as well as the momentum,the action for initial and excitation phase

圖6 另一種閾值選擇Fig.6 Another choice for threshold

根據(jù)圖2 中擬勢的形狀可知,若系統(tǒng)能夠躍出勢阱,則其后續(xù)運動在平坦區(qū)域的能量消耗很小.因此,Newby 選擇勢阱與平坦區(qū)域的接合處作為閾值,將從勢阱中的逃逸視為起始階段[11].其具體操作如下:首先,如圖6 所示,相空間中存在一個尖點(cusp point,CP),延伸出兩條焦散線(caustics)和一條切換線(switching line,SL).這一現(xiàn)象的出現(xiàn)是由于拉格朗日流形在擴展相空間中的折疊結(jié)構(gòu),使得其在坐標(biāo)空間上的投影具有多值性[21].而擬勢取全局最小值,由此會導(dǎo)致擬勢在切換線上的不可微性,出現(xiàn)一條折痕.因此,到達SL 兩側(cè)兩點的最優(yōu)路徑在拓?fù)渖鲜峭耆煌?于是,可以選擇SL 作為閾值的一部分,另一部分由尖點的擬勢等高線(contour line,CL)提供,分別用粗體的棕色和綠色曲線表示.可以看出最終進入平坦區(qū)域的極值路徑在穿過CL 時幾乎都位于某個S 點附近.

在結(jié)束這一節(jié)之前,我們比較一下這兩種閾值選擇的差別.首先,以分界線為閾值更多的是關(guān)注確定性系統(tǒng)的行為,而第二種主要考慮的是擬勢的結(jié)構(gòu).另外,前者只承認(rèn)大振幅運動為一次放電現(xiàn)象,這實際上更符合動力學(xué)的觀點.更重要的是,由于鞍點的擬勢在整個邊界上取最小值,使得在弱噪聲極限情況下,最優(yōu)路徑必然會通過鞍點,因此鞍點狀態(tài)在逃逸問題中起著至關(guān)重要的作用.但是,由于選擇擬勢等高線作為閾值,CL 上的每個點具有相同的概率,這削弱了S 的特殊意義.綜合來看,我們認(rèn)為將分界線作為閾值更為合理.最后,需要強調(diào)的是,這兩種選擇的平均首次離出時間幾乎是相同的,因為它們都在擬勢的平坦區(qū)域,具有幾乎相等的擬勢.

4 不同噪聲組合下的影響

迄今為止,我們研究了隨機混合Morris-Lecar 系統(tǒng)對于不同閾值在M=N=40 的情況下自發(fā)性放電現(xiàn)象.本節(jié),我們將分析鈉離子和鉀離子噪聲的各種不同組合,即不同的M和N的值將如何影響離出方案和擬勢形狀.

首先為探索M和N量級的影響,令M=N.可以發(fā)現(xiàn),當(dāng)具有不同M和N值和相同初始點時,極值路徑完全一致.然而,當(dāng)使用GMAM 計算分界線上SN 到S0的最優(yōu)路徑時,擬勢卻發(fā)生了明顯的變化,如圖7(a)所示,擬勢與M,N的量級呈嚴(yán)格的正比關(guān)系.換句話說,M,N的量級線性地決定了勢阱的深度.這一事實很好理解,因為如果M和N的值很大,那么由m和n的加減所帶來的變化對v的導(dǎo)數(shù)影響很小,于是便導(dǎo)致了更大的逃逸難度.

圖7 M 和N 不同組合下擬勢的對比Fig.7 Comparison for quasi-potential with various combinations of M and N

另一方面,通過研究M和N值的不同組合,以確定何種漲落對整個逃逸過程的貢獻更大.首先觀察分界線上的擬勢與M,N的不同組合之間的關(guān)系.將N分別取為1,10,40,400,擬勢與M=1,10,40,400之間的4 條折線關(guān)系如圖7(b)所示.可以看出,如果固定N,則擬勢將隨著M的增大而顯著增大.當(dāng)然,也存在一個N=1 時的特殊情況,即擬勢自M=40到M=400 時增大得不明顯.不過,這一現(xiàn)象也是可以理解的,若M與N的比值過大,鉀離子的漲落便會足夠小.

另一個角度,通過固定M我們可以分析鈉離子通道噪聲對擬勢的影響.若M=1,則當(dāng)N取不同值時,擬勢值幾乎無法區(qū)分.此外,當(dāng)M的值較大時,在N

事實上,根據(jù)表1 中的數(shù)據(jù),我們從另一個方面也可以發(fā)現(xiàn)這些特性.對角線顯示了前文中的線性關(guān)系.此外,當(dāng)我們觀察另一條對角線就會發(fā)現(xiàn)另一個更有趣的性質(zhì)—M和N的乘積保持不變而擬勢存在明顯的最大值,即,存在一個能使?jié)q落強度最小的M與N的最優(yōu)比例.我們可以利用這一性質(zhì)控制噪聲強度的大小.另外可以發(fā)現(xiàn),在這個最優(yōu)值的條件下M是N的數(shù)倍,這也從另一方面證實了鉀離子通道噪聲的主導(dǎo)作用.

表1 M 和N 不同組合下分界線上S0的擬勢Table 1 Quasi-potential of S0on the separatrix versus various combinations of M and N

接下來,通過調(diào)整M和N的相對大小,我們觀察鈉離子和鉀離子通道噪聲如何分別影響離出路徑.對不同的M,將8 種情況分為N=1 和N=40 兩組,它們的最優(yōu)路徑分別如圖8 所示.可以看出,最優(yōu)路徑的動力學(xué)行為由M和N的比值決定,與量級無關(guān).如果M遠(yuǎn)大于N,則最優(yōu)路徑接近水平線.當(dāng)M大于或者幾乎等于N時,最優(yōu)路徑位于v零傾線以下.若初始時M=N,則增大M會顯著改變最優(yōu)路徑的形狀,而N的增大所引起的改變很小.類似于擬勢與通道噪聲強度比值的關(guān)系,最優(yōu)路徑也主要取決于鉀離子通道噪聲.最后我們指出,在具有快慢尺度的擴散過程里存在類似的現(xiàn)象,即調(diào)整快慢方向的噪聲強度比例可以控制最優(yōu)路徑水平或彎曲.實際上,快變方向的噪聲較強時,系統(tǒng)能夠直接克服水平勢壘從而沿水平方向離出,而快變方向噪聲較弱時無法克服,因此只能沿著快變方向零傾線離出,表現(xiàn)為彎曲的最優(yōu)路徑.在文獻[38]中對一相對簡單的快慢尺度擴散過程給出了理論解釋,本文研究的模型雖為混合系統(tǒng),但仍具有快慢特性,因此出現(xiàn)了類似的現(xiàn)象.

圖8 M 和N 不同組合下最優(yōu)路徑的對比Fig.8 Comparison for the optimal paths with various combinations of M and N

5 結(jié)論

本文研究了因為系統(tǒng)缺乏鞍點狀態(tài)而沒有明確閾值的隨機混合Morris-Lecar 系統(tǒng)的自發(fā)性動作電位.通過定義虛擬分界線為擬閾值,可以將一次放電行為分解為一個長時間的初始階段和一個暫態(tài)的興奮階段.前者可以看成是噪聲誘導(dǎo)的逃逸事件,而后者則幾乎會沿著確定性的路徑在大范圍運動后回到靜息態(tài).通過三種不同的數(shù)值方法,即GMAM、OUM和作用量圖方法,計算得出閾值上的擬勢,發(fā)現(xiàn)結(jié)果存在明顯的最小值.由此可以肯定離出位置幾乎集中于這個點附近,因此可以斷定該點發(fā)揮著類似鞍點的作用.通過改進的蒙特卡羅方法我們對歷程概率分布進行了模擬,驗證了我們前面使用理論方法得到的最優(yōu)路徑和離出位置的結(jié)果.另外,針對前人所給出的另一種閾值取法,即切換線和尖點擬勢等高線的組合,我們比較了兩種選擇的優(yōu)劣性,并認(rèn)為,從動力學(xué)的角度來講,前者更具有實際物理意義.

最后,研究表明,M,N的量級嚴(yán)重影響勢阱的深度,并且它們之間呈現(xiàn)嚴(yán)格的線性關(guān)系.此外,M和N的不同組合對最優(yōu)路徑的模式和擬勢有著顯著影響,鉀離子通道噪聲在自發(fā)性動作電位過程中起主導(dǎo)作用.研究還發(fā)現(xiàn),存在一個M和N的最優(yōu)比值,能使?jié)q落強度達到最小.

附錄

參數(shù)值

vNa=3.7,gNa=0.22,vK=?0.9,gK=0.4,vleak=?0.36,gleak=0.1,βK=0.04,Iapp=0.06,γNa=1.22,κNa=?1.188,γK=0.8,κK=?0.8,M=40,N=40.