Cav1.2離子通道馬爾科夫隨機過程模型的建立與應(yīng)用

宋璟波,張鎮(zhèn)西

(西安交通大學生物醫(yī)學信息工程教育部重點實驗室,710049,西安)

Cav1.2作為一類L型離子蛋白通道調(diào)控Ca2+流入細胞的跨膜行為,在心肌細胞中集中分布于T管與肌質(zhì)網(wǎng)共同形成的納米尺度的小空間中[1],在許多重要的細胞生理過程中,例如興奮收縮偶聯(lián)、動作電位的產(chǎn)生和傳播等都起著關(guān)鍵作用。因此,Cav1.2通道的門控行為受到了廣泛的關(guān)注和研究,由此所產(chǎn)生的量化模型作為實驗的輔助方法大量應(yīng)用于細胞和心臟的數(shù)字化工程中[2]。傳統(tǒng)的Hodgkin -Huxley形式方程作為一種唯象方法采用了“獨立門控因子”來表征通道的電壓相關(guān)的時變電流變化[3],然而這一范式無法反映近年來越來越受到重視的基因、藥物和離子通道的相互作用機制。另一方面,馬爾科夫隨機模型中的狀態(tài)可以對應(yīng)離子通道的無記憶構(gòu)象,狀態(tài)之間的轉(zhuǎn)換可以對應(yīng)通道構(gòu)象的變化,并且在心肌細胞的單通道電生理記錄實驗中,能夠較好地反應(yīng)實驗數(shù)據(jù)的統(tǒng)計特性[4],因此隨機馬爾科夫模型開始廣泛應(yīng)用于離子通道的建模之中[5]。

本文基于Cav1.2通道的結(jié)構(gòu)信息,通過改進Mahajan[6]提出的隨機模型參數(shù),建立了一個優(yōu)化的7態(tài)馬爾科夫隨機模型,它的特點有:①狀態(tài)遷移參數(shù)具有有界性,即使在異常生理條件下(膜電位范圍為-150～200 mV)依然處于合理的數(shù)值范圍之內(nèi);②具有廣泛的適用性,在單通道特性、I-V關(guān)系和準穩(wěn)態(tài)失活研究中,可以較好地重現(xiàn)實驗結(jié)果;③結(jié)合Ca2+擴散偏微分方程,可以準確地表征在小空間中的Cav1.2通道瞬態(tài)隨機門控行為;④通過相關(guān)結(jié)構(gòu)信息分析,簡單的調(diào)整狀態(tài)遷移參數(shù)便可以較好地預(yù)測G406R/G432N突變造成的通道電生理行為變化。

1 方 法

1.1 Cav1.2隨機馬爾科夫模型

Cav1.2通道作為一類多聚物復(fù)合體蛋白,由同源四聚體成孔亞基α1c、α2-δ、β和γ組成,分布于細胞的T小管區(qū)域,與肌質(zhì)網(wǎng)構(gòu)成了約15 nm的小空間。其中α1c的4個同源體(I～IV)都由6個跨膜部分S1～S6構(gòu)成(見圖1)。門控電流實驗數(shù)據(jù)[7]表明,Cav1.2通道的活化過程實際上是在膜電位影響下S4上的帶電殘基移動和旋轉(zhuǎn)造成的一系列構(gòu)象變化過程(見圖2a),因此它可以簡化成2個閉狀態(tài)(C1,C2)和一個開狀態(tài)(O)[6,8]。考慮到單通道實驗中的平均開放時間不受膜電位的影響[9],本文假定C2和O之間的遷移速率與電壓無關(guān)。

圖1 Cav1.2的微觀環(huán)境結(jié)構(gòu)簡圖

(a)Cav1.2處于開、閉狀態(tài)

(b)Cav1.2處于VDI狀態(tài)

Cav1.2通道具有兩類如圖1所示的調(diào)控機制,即電壓相關(guān)失活(VDI)和Ca2+相關(guān)失活(CDI)。在離子通道的分子水平上,VDI過程一般被認為主要是由I-II鏈接基團或α1亞基的羧基(C)末端的鉸鏈蓋效應(yīng)形成,即這些鏈接基團可以對通道形成阻斷,其中I-II鏈接基團受到C端影響,形成相對緩慢的VDI過程(見圖2b)。

(c)Cav1.2處于CDI狀態(tài)圖2 Cav1.2離子通道門控狀態(tài)示意圖

另一方面,新的實驗[10-11]表明,CDI可以被視為一個快速的VDI過程。它簡化的構(gòu)象變化機制是:C末端的鈣調(diào)蛋白CaM在無Ca2+環(huán)境下形成的構(gòu)象使得它位于通道孔附近,造成對I-II鏈接基團的“剎車”效應(yīng),即減緩了鏈接基團對通道的阻斷;在Ca2+環(huán)境下,Ca-CaM復(fù)合體基團在C末端的結(jié)合可能導(dǎo)致一個位移效應(yīng),使得它遠離了原來的“剎車”位點,從而加速了I-II鏈接基團的鉸鏈蓋效應(yīng)。因此,CDI過程可被視為Ca2+介導(dǎo)的如圖2c所示的快速VDI過程。

從上述結(jié)構(gòu)信息中可以自然地假定,VDI和CDI可以簡化為對稱形式,區(qū)別在于VDI的狀態(tài)遷移參數(shù)只受到膜電位的調(diào)控,而CDI的狀態(tài)遷移參數(shù)受到膜電位和Ca2+濃度的雙重調(diào)控。從計算性能和實際復(fù)雜性考慮,本文采用Mahajan[6]提出的模型框架(見圖3),并對其狀態(tài)遷移參數(shù)進行改進。

C1、C2:閉狀態(tài);O:開狀態(tài);Iv1、Iv2:電壓依賴失活

參數(shù)估計的方法如下。

(1)構(gòu)造7個狀態(tài)的馬爾科夫轉(zhuǎn)移矩陣

Q=

矩陣中各項參數(shù)代表狀態(tài)之間的遷移速率。例如,vC1C2代表從狀態(tài)C1到狀態(tài)C2的遷移速率,它們受到膜電位V的影響,而遷移速率中帶有ICa符號的參數(shù)還受到Ca2+濃度的影響。Σ代表矩陣各行的非0項之和。

(2)連續(xù)時間馬爾科夫鏈的狀態(tài)方程為

πt=π0eQ t

(1)

式中:π0代表初始時刻馬爾科夫狀態(tài)的概率分布向量;πt代表t時刻狀態(tài)的概率分布向量;eQ t代表矩陣的指數(shù)運算,Q是馬爾科夫轉(zhuǎn)移矩陣。由于π中元素表示相應(yīng)狀態(tài)在特定時刻的概率,因此π必須滿足πi≥0且Σiπi=1。矩陣的指數(shù)運算通過Matlab中的expm函數(shù)實現(xiàn)。

(3)C、O和Iv狀態(tài)間的遷移速率v是膜電位V的函數(shù),定義為

(2)

vIVC=

(3)

式中:α1、β1、α2、β2、ζ和λ是函數(shù)的待定控制參數(shù)。上述遷移情況可以滿足遷移參數(shù)有界性的要求,同時V在數(shù)值較大時遷移速率進入飽和,與通道蛋白構(gòu)象變化的趨勢一致。ICa狀態(tài)之間以及ICa和C、O之間的遷移速率還受到Ca2+濃度c的影響,其函數(shù)形式定義為

(4)

(5)

式中:ζ1、η1、ζ2、η2、κ和δ是函數(shù)的待定控制參數(shù)。

可見:在膜電位不變且Ca2+濃度升高時,遷移速率增大;在Ca2+濃度升高到較大值時,遷移速率進入飽和狀態(tài),與前述的通道CDI時構(gòu)象變化趨勢一致。

(4)定義ρ1=[α1,β1,α2,β2,ξ,λ],τ1=[ζ1,η1,ζ2,η2,κ,δ],則模型待定參數(shù)向量θ=[ρ1,τ1,ρ2,τ2……]。在貝葉斯參數(shù)估計框架中,模型參數(shù)的后驗概率分布由下式確定

p(θ|D∝p(D|θ)p(θ)

(6)

式中:D代表實驗數(shù)據(jù);p(θ)代表模型參數(shù)的先驗概率分布,假定θ服從(0,100)的均勻分布;p(D|θ)代表在參數(shù)確定條件下模型計算得到實驗數(shù)據(jù)的概率,假定實驗數(shù)據(jù)服從以模型計算結(jié)果為中心的正態(tài)分布,即Di服從參數(shù)為(Modeli,σi)的正態(tài)分布,進一步假定σi為相同的常數(shù)。利用JAGS軟件包,通過馬爾科夫蒙特卡洛(MCMC)方法,對后驗概率分布進行Gibbs采樣,尋求后驗概率分布最大區(qū)域(MAP)。

1.2 微觀結(jié)構(gòu)下Cav1.2隨機模型

在微觀環(huán)境下,Cav1.2分布于肌膜和內(nèi)質(zhì)網(wǎng)膜構(gòu)成的納米尺度的小空間中。本文構(gòu)建了一個簡單圓柱體來模擬Cav1.2通道分布的微觀環(huán)境,利用擴散反應(yīng)偏微分方程和隨機開閉的Cav1.2通道邊界,對鈣離子的瞬態(tài)行為進行了計算

(7)

(8)

D?zc=-J|z=0

(9)

表1 微觀結(jié)構(gòu)下隨機偏微分方程的參數(shù)

利用式(6)計算得到估計的馬爾科夫模型參數(shù)后,將通道的模型代入偏微分方程的邊界條件式(9)中,利用Gillespie隨機算法確定通道的隨機開閉事件,進而在計算過程中動態(tài)地改變邊界條件。利用開源有限元軟件FENICS對上述隨機偏微分方程進行計算,并以PARAVIEW對結(jié)果進行可視化。隨機偏微分方程求解步驟如圖4所示。

圖4 隨機偏微分方程求解步驟

1.3 Ba2+和電壓相關(guān)失活(VDI)

盡管近年來有文獻報道Ba2+自身同樣會造成類似于Ca2+相關(guān)失活的現(xiàn)象[12],但是在大量的Cav1.2失活研究中,Ba2+依然被視為理想的隔離CDI效應(yīng)的2價陽離子[13]。本文把Ba2+為載荷子的實驗數(shù)據(jù)視為只受到VDI的影響,即在模型中忽略CDI路徑。

2 計算結(jié)果

2.1 I-V關(guān)系

I-V關(guān)系指Cav1.2通道在不同的電壓脈沖下的離子電流峰值變化。圖5顯示了在Ba2+和Ca2+作為載荷離子條件下的實驗數(shù)據(jù)[14]與模型計算結(jié)果的比較,數(shù)據(jù)針對最大值進行了歸一化。從圖中可以看到,在復(fù)極化電位(≤-40 mV)下,Cav1.2通道幾乎沒有電流,這是因為在負的膜電位下,通道孔徑的開放概率較低;在高電位下,通道電流減少是因為依據(jù)恒定場公式,通道最大電流與電位大小成反比。通過模型計算的I-V曲線的峰值在～10 mV之間,并且在-40～60 mV的脈沖電壓范圍內(nèi)與實驗結(jié)果吻合較好。此外,通過模型計算得到的Ba2+和Ca2+電流最大峰值約為1.7 nA,這與文獻[6]報道的2 nA相近。實驗中,10 mmol/L乙二醇二乙醚二胺四乙酸(EGTA)作為Ca2+緩沖液抑制了肌質(zhì)網(wǎng)RyR釋放Ca2+的能力,從而排除鈣誘導(dǎo)鈣釋放(CICR)過程對Cav1.2離子通道門控行為的額外影響。膜片鉗的刺激協(xié)議在圖5中顯示。

(a)Ba2+的I-V關(guān)系

(b)Ca2+的I-V關(guān)系圖5 歸一化的Ba2+和Ca2+的I-V關(guān)系曲線

2.2 失活

在長時間的去極化電壓下或在細胞內(nèi)Ca2+濃度提高的條件下,Cav1.2通道都會產(chǎn)生離子電流減少的現(xiàn)象,被稱為通道失活。通道失活作為負反饋調(diào)節(jié)機制是細胞內(nèi)正常Ca2+循環(huán)的關(guān)鍵環(huán)節(jié)之一。在馬爾科夫模型中,失活代表著狀態(tài)概率大部分由開放態(tài)O和關(guān)閉態(tài)C轉(zhuǎn)移到失活態(tài)I。圖6a顯示了在單通道實驗中不同保持電位下Cav1.2的電壓失活。實驗數(shù)據(jù)由平均Ba2+電流除以通道單位電流得到[15],并進行了歸一化。實驗和模型計算結(jié)果表明,隨著靜息電位的升高,單通道電流在減少,即通道的平均開放概率在減少。圖6b顯示全細胞下的Ba2+和Ca2+電流在200 ms的去極化電壓持續(xù)期下的失活情況[16]。雙電位脈沖協(xié)議如圖6所示,可以看出:Cav1.2通道在Ba2+作為載荷子時,它的失活程度與去極化電壓的強度和持續(xù)時間成正比,在去極化電壓大于一定閾值時失活進入平穩(wěn)期,這一閾值大致位于0～20 mV之間;在Ca2+條件下,即CDI和VDI共同存在的情況下,通道的失活曲線一方面在負值電壓區(qū)間表現(xiàn)出了更大的失活程度,另一方面則顯示出了U型曲線,即在去極化電壓進入正值后,失活程度有所衰減;本文建立的模型能夠準確地反應(yīng)單通道和全細胞兩種環(huán)境下Cav1.2VDI和CDI的變化趨勢。實驗中,利阿諾定(Ryanodine)作為抑制劑抑制了肌質(zhì)網(wǎng)RyR釋放Ca2+的能力,從而排除CICR過程對Cav1.2離子通道門控行為的額外影響。

(a)單通道Ba2+失活曲線

(b)Ba2+和Ca2+的失活曲線圖6 單通道和全細胞條件下的通道失活情況

2.3 在小空間中的Ca2+瞬態(tài)

將小空間抽象為一個半徑為160 nm、高為15 nm的圓柱體(見圖7a),考察Cav1.2在毫秒級時間內(nèi)引發(fā)CICR的能力。圓柱體下表面代表T小管的細胞膜,5個離子通道軸對稱分布;上表面代表肌質(zhì)網(wǎng)膜,參與CICR的利阿諾定受體RyR就分布于其上(此研究中RyR不參與Ca2+瞬態(tài)行為)。Cav1.2遵循馬爾科夫鏈的規(guī)律隨機開閉,其狀態(tài)遷移參數(shù)受到鄰近區(qū)域的Ca2+濃度和膜電位的雙重影響。隨機偏微分方程的計算結(jié)果(見圖7b)表明,單個Cav1.2通道的開放可以導(dǎo)致通道上部的肌質(zhì)網(wǎng)相應(yīng)區(qū)域(以Cav1.2正上方15 nm為圓心,徑向半徑為100 nm,見圖7a)在小于100 μs時間內(nèi)Ca2+濃度超過1 μmol/L,見圖7b。這一結(jié)果符合之前的文獻[17]報道。圖8反映了5個Cav1.2通道分別在1、4.6和12.4 ms時的開閉狀態(tài),以及空間中Ca2+濃度的分布,圖形由PARAVIEW軟件生成。

(a)小空間簡化模型

(b)Ca2+濃度徑向分布圖7 小空間簡化結(jié)構(gòu)和計算結(jié)果

(a)1 ms

(b)4.6 ms

(c)12.4 ms圖8 不同時刻的Ca2+瞬態(tài)分布情況

2.4 G406R/G432N突變

圖9 野生型和G406R突變型通道Ba2+失活計算結(jié)果和實驗結(jié)果的比較

(a)活化特征時間之比

(b)去活化特征時間之比圖10 G432N與野生型活化和去活化特征時間之比

3 結(jié) 論

本文基于Cav1.2通道的結(jié)構(gòu)信息提出了改進的7態(tài)馬爾科夫隨機模型,并以此對一些通道的門控特性進行了模擬。計算結(jié)果表明,本文模型能夠在保持計算輕量化的同時較好地重現(xiàn)實驗結(jié)果,可以作為Cav1.2通道生理實驗的有效輔助方法。本文進一步以此模型作為計算框架,基于現(xiàn)有的通道結(jié)構(gòu)信息,通過合理簡單地修改模型參數(shù),即能準確地重現(xiàn)G406R突變所造成的失活性狀的改變,以及G432N對活化和去活化過程的延長影響,從而一方面表明該模型與通道結(jié)構(gòu)緊密相關(guān),另一方面顯示出該模型在心肌細胞宏觀電生理、微環(huán)境結(jié)構(gòu)和Cav1.2通道病理條件等研究方面的準確性和靈活性。后續(xù)的工作是把建立的Cav1.2馬爾科夫模型納入心肌細胞模型中去,進而考察在不同病理生理條件下細胞層級電生理活動的特征,并輔助促進相應(yīng)實驗的展開。

[1] SOELLER C, CANNELL M B. Examination of the transverse tubular system in living cardiac rat myocytes by 2-photon microscopy and digital image processing techniques [J]. Circulation Research, 1999, 84(3): 266-275.

[2] NOBLE D. Modeling the heart-from genes to cells to the whole organ [J]. Science, 2002, 295(5560): 1678-1682.

[3] LUO C H, RUDY Y. A dynamic model of the cardiac ventricular action potential: I Simulations of ionic currents and concentration changes [J]. Circulation Research, 1994, 74(6): 1071-1096.

[4] COLQUHOUN D, HAWKES A G. On the stochastic properties of bursts of single ion channel openings and of clusters of bursts [J]. Philosophical Transactions of the Royal Society: B Biological Sciences, 1982, 300(1098): 1-59.

[5] MOROTTI S, GRANDI E, SUMMA A, et al. Theoretical study of L-type Ca2+current inactivation kinetics during action potential repolarization and early afterdepolarizations [J]. J Physiol, 2012, 590(18): 4465-4481.

[6] MAHAJAN A, SHIFERAW Y, SATO D, et al. A rabbit ventricular action potential model replicating cardiac dynamics at rapid heart rates [J]. Biophys J, 2008, 94(2): 392-410.

[7] BEZANILLA F. The voltage sensor in voltage-dependent ion channels [J]. Physiol Rev, 2000, 80(2): 555-592.

[8] GAUTHIER L D, GREENSTEIN J L, WINSLOW R L. Toward an integrative computational model of the Guinea pig cardiac myocyte [J]. Front Physiol, 2012, 3: 244.

[9] CAVALIE A, PELZER D, TRAUTWEIN W. Fast and slow gating behaviour of single calcium channels in cardiac cells: relation to activation and inactivation of calcium-channel current [J]. Pflugers Arch, 1986, 406(3): 241-258.

[10]CENS T, ROUSSET M, LEYRIS J P, et al. Voltage-and calcium-dependent inactivation in high voltage-gated Ca(2+) channels [J]. Prog Biophys Mol Biol, 2006, 90(1/2/3): 104-117.

[11]BENITAH J P, ALVAREZ J L, GMEZ A M. L-type Ca2+current in ventricular cardiomyocytes [J]. Journal of Molecular and Cellular Cardiology, 2010, 48(1): 26-36.

[12]FERREIRA G, YI J X, RIOS E, et al. Ion-dependent inactivation of barium current through L-type calcium channels [J]. Journal of General Physiology, 1997, 109(4): 449-461.

[13]GRANDI E, MOROTTI S, GINSBURG K S, et al. Interplay of voltage and Ca-dependent inactivation of L-type Ca current [J]. Progress in Biophysics and Molecular Biology, 2010, 103(1): 44-50.

[14]ZHONG J M, HWANG T C, ADAMS H R, et al. Reduced L-type calcium current in ventricular myocytes from endotoxemic guinea pigs [J]. American Journal of Physiology-Heart and Circulatory Physiology, 1997, 273(5): H2312-H2324.

[15]REUTER H, STEVENS C F, TSIEN R W, et al. Properties of single calcium channels in cardiac cell culture [J]. Nature, 1982, 297(5866): 501-504.

[16]SUN L, FAN J S, CLARK J W, et al. A model of the L-type Ca2+channel in rat ventricular myocytes: ion selectivity and inactivation mechanisms [J]. The Journal of Physiology, 2000, 529(1): 139-158.

[17]LANGER G A, PESKOFF A. Calcium concentration and movement in the diadic cleft space of the cardiac ventricular cell [J]. Biophys J, 1996, 70(3): 1169-1182.

[18]SPLAWSKI I, TIMOTHY K W, SHARPE L M, et al. Ca(V)1.2 calcium channel dysfunction causes a multisystem disorder including arrhythmia and autism [J]. Cell, 2004, 119(1): 19-31.

[19]DEPIL K, BEYL S, STARY-WEINZINGER A, et al. Timothy mutation disrupts the link between activation and inactivation in Ca(V)1.2 protein [J]. J Biol Chem, 2011, 286(36): 31557-31564.

[20]SPLAWSKI I, TIMOTHY K W, DECHER N, et al. Severe arrhythmia disorder caused by cardiac L-type calcium channel mutations [J]. Proc Natl Acad Sci USA, 2005, 102(23): 8089-8096.