劉 巖,孫 毅,張桂敏,夏健明,張洪明*

(1.昆明理工大學 力學系，云南 昆明 650500；2.云南阜外心血管病醫(yī)院，云南 昆明 650032)

KCa2.3離子通道蛋白是分布在人體心肌細胞隸屬于小電導鈣離子激活的鉀離子通道(small-conductance Ca2+-activated K+channel，SK Channel)，對鉀離子有高度選擇性[1]。當細胞質(zhì)內(nèi)游離的鈣離子濃度上升后會被激活，從而開放鉀離子通道，產(chǎn)生鉀離子外流，鉀離子外流是通過各類鉀離子通道實現(xiàn)的，鉀離子通道是生物體內(nèi)分布最廣泛、最復雜的一類離子通道[2]，SK離子通道就是鉀離子通道中的一類。鉀離子外流產(chǎn)生的動作電位會引起超極化后電位(afterhyperpolarization potential， AHP)現(xiàn)象， 超極化后電位對動作電位反復激動現(xiàn)象有抑制作用，以保證細胞的各項生命活動[3]。但是KCa2.3離子通道開放程度過大會使心房肌的動作電位時長明顯縮短[4]，從而會引發(fā)一些心血管疾病。

近幾年的研究表明，通過減小鉀離子外流可以有效減少某些心血管疾病尤其是房顫(valvular atrial fibrillation，VAF)的病發(fā)率[5]。我國的孫曼青[6]、苑磊[7]、LI[8]等的研究里表述了SK離子通道在人的心肌細胞里廣泛分布，且當血壓發(fā)生變化時，SK離子通道內(nèi)的電流也有較明顯波動。XU等在2003年發(fā)現(xiàn)SK通道分布心房優(yōu)于心室[9]。TAKAI等在2013年發(fā)現(xiàn)KCa2.3離子通道受心臟內(nèi)湍流切應力調(diào)控[10]。因此，對KCa2.3離子通道蛋白表達的進一步深入研究有科學應用價值。

本文基于超級計算機和分子動力學模擬軟件Gromacs，模擬計算了KCa2.3離子通道蛋白在不同壓力下對鉀離子外流的調(diào)控程度，分析了不同血壓下鉀離子外流的速度。

1 計算模型及研究方法

1.1 模擬原理

分子在整個環(huán)境中所受力場分為動能和勢能，勢能包括非鍵結(jié)勢能(Unb)、鍵的伸縮項(Ub)、鍵角彎曲項(Uθ)、二面角扭曲項(Uφ)、離平面振動項(Uχ)、庫倫作用項(Uel)，總勢能公式如下：

U=Unb+Ub+Uθ+Uφ+Uχ+Uel

(1)

依照經(jīng)典力學，系統(tǒng)中任一原子在環(huán)境中所受的力為:

(2)

式中，r為體系中原子的位置，通過數(shù)值求解牛頓方程模擬原子的運動。

(3)

最后得出每個原子的速度、位置信息。

1.2 KCa2.3通道工作機理

KCa2.3離子通道是小電導鈣激活鉀離子通道的一種，由KCNN3編碼，α亞單位和β單位組成，如圖1所示：

注：螺旋形表示的是α亞單位，固形箭頭表示的是β單位。

KCa2.3離子通道蛋白有1個單通電導，并受細胞內(nèi)Ca2+調(diào)控。KCa2.3離子通道作用機理如下：

a)當細胞質(zhì)內(nèi)鈣離子處于高濃度狀態(tài)時，鈣離子及鈣調(diào)蛋白(CaM)會結(jié)合KCa2.3離子通道的α亞單位，見圖2(a)；

b)KCa2.3離子通道緩慢激活、開放，見圖2(b)；

c)細胞質(zhì)內(nèi)鉀離子外流，見圖2(c)；從而產(chǎn)生動作電位。

圖2 KCa2.3離子通道激活過程

1.3 軟件及控制變量的數(shù)值選取

本文采用Gromacs軟件進行模擬，力場選取全原子(OPLS-AA/L all-atom)力場[11]，水分子模型采用SPC/E模型。體系壓力通過Parrinello-Rahman[12]控制，其數(shù)值分別設(shè)置為60、70、80、90、100、110、120、130、140、150、160、170、180 mmHg；溫度耦合由Nose-Hoover[13]控制，數(shù)值設(shè)置為310 K；靜電相互作用采用Particle-Mesh-Ewald算法進行處理；范德華力計算采用截斷Cut-off方式，計算截斷值設(shè)為0.12 nm。

1.4 構(gòu)建分子模型

構(gòu)建分子模型步驟如下：

a)在蛋白質(zhì)數(shù)據(jù)庫下載并修改KCa2.3離子通道的分子模型，修改后的模型如圖3所示。左圖為蛋白質(zhì)的全原子模型，右圖為蛋白質(zhì)泛素的模型。

(a)全原子模型 (b)泛系模型

b)利用分子建模軟件構(gòu)建鈣離子和鉀離子模型；正常生理狀態(tài)下，細胞內(nèi)、外鈣離子和鉀離子的濃度如表1所示。

表1 鈣離子和鉀離子濃度表

結(jié)合該數(shù)據(jù)，將一定數(shù)量的鈣離子、鉀離子分別固定在KCa2.3離子通道蛋白的兩側(cè)，整個模擬體系在6 nm×6 nm×52 nm的盒子中進行，如圖4所示。

圖4 分子動力學模擬體系

圖中大的離子表示鈣離子，小的離子表示鉀離子，體系中間為KCa2.3離子通道蛋白。

c)加入正/負離子(Na+/Cl-)以保證體系的電荷守恒。表2是加入Cl-之后，分子動力學模擬體系內(nèi)各單位的數(shù)量。

表2 各單位原子數(shù)量

d)進行能量最小化模擬，優(yōu)化初始模型。初步進行能量最小化模擬，采用最速下降法，將最大能量設(shè)置在1 000 kJ/mol，防止分子模型在初始階段受分子間勢能影響進行其他運動。

1.5 參數(shù)設(shè)置

本次模擬僅研究不同壓力下KCa2.3離子通道對鉀離子的調(diào)控作用，設(shè)置溫度耦合為恒溫環(huán)境，溫度場為36.85 ℃(310 K)，步數(shù)為1×105，步長為1fs(1fs=1×10-15s)，耦合總時長為0.1 ns。

壓力耦合設(shè)置選用60～180 mmHg共 13組壓力值，每組數(shù)據(jù)步數(shù)均為1×105，步長為1fs，耦合總時長為0.1 ns。

進行最終分子動力學模擬時，步數(shù)選用1.5×106步，步長為1fs，總時間為1.5 ns，每50 ps保存一次軌跡文件。

2 計算結(jié)果

2.1 擴散系數(shù)

通過觀察鉀離子的擴散系數(shù)，我們可以直觀的看到鉀離子通過KCa2.3離子通道蛋白流出細胞外的程度，擴散系數(shù)越大，鉀離子外流速率越快。而Gromacs軟件中并不能直接給出鉀離子的擴散系數(shù)。因此，我們需要根據(jù)其均方位移(mean-square displacement，MSD)函數(shù)計算出鉀離子在不同壓力體系下的擴散系數(shù)。

均方位移表示的是某一原子相對起始位置的位移量，其公式為：

MSD=〈|r(t)-r(0)|2〉

(4)

再結(jié)合愛因斯坦的擴散定律得到，

limt→∞〈|r(t)-r(0)|2〉=6Dt

(5)

計算得出鉀離子的擴散系數(shù)為，

(6)

式中,D為粒子擴散系數(shù)(diffusion constant)且為常數(shù)，所以模擬時間t與鉀離子的均方位移為一次函數(shù)關(guān)系。圖5為同一體系在不同壓力下鉀離子的MSD函數(shù)，其斜率為鉀離子的擴散系數(shù)。

圖5 體系在不同壓力下的MSD函數(shù)

從表中我們可以看出在壓力在90 mmHg時，MSD函數(shù)斜率最大；壓力在140 mmHg時，MSD函數(shù)斜率最小。因此,當血壓為90 mmHg時，鉀離子外流速率最大，血壓為140 mmHg時，鉀離子外流速率最低。根據(jù)圖5的MSD函數(shù)曲線，可以得到不同壓力對應鉀離子擴散系數(shù)，如表3所示。

表3 不同壓力對應的鉀離子擴散系數(shù)

根據(jù)表3繪制壓力和鉀離子擴散系數(shù)的樣條曲線如圖6。從圖中可以看出壓力與鉀離子擴散系數(shù)的關(guān)系呈“U”型，當壓力在120～150 mmHg之間時，鉀離子擴散系數(shù)最低，而當壓力小于120 mmHg或者大于150 mmHg時，鉀離子擴散系數(shù)變大，鉀離子外流程度增大，較120～150 mmHg時增大了3.19%～6.12%。該數(shù)據(jù)與邢愛君統(tǒng)計研究[14]發(fā)現(xiàn)的當收縮壓處于120～140 mmHg時發(fā)生心血管疾病死亡風險最低的結(jié)果范圍大致一致。

圖6 壓力與鉀離子擴散系數(shù)的關(guān)系

2.2 體系的自由能分析

自由能是影響系統(tǒng)變化的重要因素，其中分子非鍵作用勢能占絕對主導地位，其主要包括分子間的范德華力和分子間的靜電力兩部分。

(7)

式中，εij是第i個原子和第j個原子之間相互作用勢能的最小值；σij是第i個原子和第j個原子之間達到平衡的距離；rij是第i個原子和第j個原子之間的絕對距離；qi，qj是該原子的電荷量；ε0是介電常數(shù)，ε0=1/4πk。

圖7是模擬體系在進行壓力耦合前做的能量最小化數(shù)據(jù)，在模擬達到6 325步后，系統(tǒng)趨于平衡，其能量最小值(Epot)為-7.215×106kJ/mol。此時，體系達到最穩(wěn)定的狀態(tài)，與模擬過程中體系內(nèi)相互作用能進行對比(圖8)，在模擬過程中能量最大為-7 140 348.5 kJ/mol，能量最小為-7 224 514.5 kJ/mol，波動范圍在84 116 kJ/mol之間。因此，隨時間變化體系內(nèi)的自由能穩(wěn)定性能未發(fā)生改變，說明在KCa2.3離子通道工作時未發(fā)生化學反應，且各原子動能未發(fā)生較大變化，各原子的運動速度是一個定值。

圖7 能量最小化過程中勢能的變化

圖8 鉀離子與KCa2.3離子通道蛋白之間的相互作用能

2.3 KCa2.3離子通道在不同壓力下的穩(wěn)定性

均方根偏差(root mean square deviation，RMSD)是衡量模擬過程中蛋白質(zhì)穩(wěn)定性的重要指標，當模擬后的蛋白質(zhì)位置與初始狀態(tài)的位置偏差較大時，RMSD值就會越大，其計算公式為：

(8)

式中,rij表示j時刻對應i原子的坐標位置;ri0表示初始時刻原子所在的位置。

圖9、圖10分別為體系在壓力為60 mmHg和180 mmHg時的RMSD變化趨勢。從圖中我們可以看出，隨著模擬時間的增加，蛋白質(zhì)與其初始位置的距離也變大，且均在模擬0.8 ns后，其結(jié)構(gòu)趨于穩(wěn)定，其位移大致為0.45 nm和0.49 nm。體系在其他壓力狀態(tài)下KCa2.3離子通道蛋白的位移均在0.45～0.49 nm范圍內(nèi)，對比心肌細胞直徑的15 μm可以忽略不計，說明壓力對KCa2.3離子通道蛋白的穩(wěn)定性沒有影響或是影響較小。

圖9 60 mmHg時分子模擬過程中的RMSD值

圖10 180 mmHg時分子模擬過程中的RMSD值

為了更直觀地觀察不同壓力下KCa2.3離子通道蛋白的變化情況，我們將模擬前后的KCa2.3離子通道蛋白構(gòu)型放到同一體系內(nèi)進行觀察對比，如圖11所示。

圖11 模擬前后KCa2.3離子通道蛋白構(gòu)型對比

圖11(a)、11(b)中較薄的構(gòu)型為體系模擬前的分子構(gòu)型，圖11(c)中較薄的構(gòu)型為體系在60 mmHg壓力下的分子構(gòu)型。通過對比發(fā)現(xiàn)KCa2.3離子通道蛋白在模擬前后，結(jié)構(gòu)上有細微差異；而60 mmHg壓力下與180 mmHg壓力下KCa2.3離子通道蛋白結(jié)構(gòu)沒有發(fā)生變化，僅在模擬后的位置上有細微差異，進一步驗證了壓力對KCa2.3離子通道蛋白的穩(wěn)定性沒有影響。

2.4 分子聚集特性

徑向分布函數(shù)(radial distribution function，RDF)其物理意義如下：

設(shè)中心距離為r到r+dr圓環(huán)內(nèi)的原子數(shù)目為dN，RDF函數(shù)g(r)為

ρg(r)4πr2=dN

(9)

解得,

(10)

當RDF值越接近1時，表示該區(qū)域范圍內(nèi)的密度越接近平均密度。我們可以根據(jù)徑向分布函數(shù)查看鉀離子之間最短距離，對比不同體系中鉀離子之間的徑向分布函數(shù)，發(fā)現(xiàn)其分布函數(shù)沒有較大差異。圖12是在100 mmHg壓力下體系以鉀離子為中心，距離中心r到dr的圓環(huán)內(nèi)鉀離子數(shù)目的徑向分布函數(shù)。

圖12 系統(tǒng)在分子模擬過程中的RDF函數(shù)

在0.338 nm處RDF函數(shù)值開始不為零，即在該系統(tǒng)下兩個鉀離子之間最短距離為0.338 nm，當距離為0.434 nm時，徑向分布函數(shù)值達到最大值為1.631。在其他壓力情況下，當中心距離越大時，RDF也都是趨近于1。說明鉀離子在外流的過程未出現(xiàn)聚集現(xiàn)象，鉀離子在整個模擬過程中都是均勻分布在細胞液內(nèi)。

3 結(jié)論

為了研究KCa2.3離子通道蛋白在不同壓力環(huán)境下各方面性能的差異，本文采用分子動力學模擬的方法，對同一體系在13種不同壓力下進行了分子動力學模擬。通過分析和對比KCa2.3離子通道蛋白在表達時不同壓力狀態(tài)下鉀離子的擴散系數(shù)、體系內(nèi)的勢能與KCa2.3離子通道蛋白的RMSD函數(shù)、體系能內(nèi)勢能的變化、鉀離子的徑向分布函數(shù)，得出以下結(jié)論：

1)KCa2.3離子通道蛋白受壓力影響，二者之間的關(guān)系呈“U”型，當壓力處于120～150 mmHg時，通過KCa2.3離子通道蛋白的鉀離子數(shù)量較其他壓力狀態(tài)下少3.19%～6.12%；

2)KCa2.3離子通道蛋白在表達時未發(fā)生化學反應，僅為各原子之間勢能的轉(zhuǎn)移；

3)壓力對KCa2.3離子通道的穩(wěn)定性沒有影響或是影響較小；

4)KCa2.3離子通道蛋白在表達時，鉀離子未出現(xiàn)聚集現(xiàn)象，其在整個模擬過程中都是均勻分布在細胞液內(nèi)。