專利名稱：應用角動量守恒原理模擬生物大分子運動的方法
技術領域：
本發(fā)明涉及一種基于角動量守恒的非平衡態(tài)分子動力學模擬方 法，適用于研究生物大分子的結構構象一功能的相互關系，以及區(qū)域 轉動對離子通道和結構構象的影響。
背景技術：
分子動力學模擬(Molecular Dynamics Simulation) (Rahman, A. and F. H. Stillinger (1971). "Molecular dynamics study of liquid water."』 Chem Phys 55(7): 3336-3359 )是在研究生物大分子結構和功能關系中應 用最為廣泛的一種計算模擬方法。分子動力學模擬是用計算機模擬由 原子核和電子所構成的多粒子體系中原子核的運動過程，其中每一個頓定律運動，得到體系隨時間變化的結構信息和計算體系的動力學或 熱力學性質。盡管計算機技術的發(fā)展使動力學模擬時間尺寸得到延伸， 最初的模擬只涉及真空里的單個分子，比如原始的牛胰島素抑制劑模 擬只包含大約500個原子(McCammon, J. A" B. R. Gelin， et al. (1977). "Dynamics of folded proteins." Nature 267(5612): 585-590 )，而今天模擬 包含有104—106個原子的體系是很平常的，并且，模擬的目標正從研究 分子和超分子體系向細胞水平發(fā)展。目前已經開始用分子動力學方法 模擬蛋白質聚合物等復雜生物體系，如對乙酰膽^5咸受體(Law，R.丄，R. H. Henchman, et al. (2005). "A gating mechanism proposed from a simulation of a human a7 nicotinic acetylcholine receptor." Proc Natl AcadSci U S A 102(19): 6813-6818 ) 、 ATP合酶(Dittrich， M.， S. Hayashi， et al. (2003). "On the mechanism of ATP hydrolysis in F廣ATPase." Biophvs 185(4): 2253-2266 )和核糖體等構象變化進行分子動力學模擬；另 一方 面，最初的牛胰島素模擬時間尺度也只有9.2皮秒(1 ps=l(T12 s),而 現代的模擬中僅例行的所謂平衡過程的預處理就遠遠超過這個時間， 一般的模擬時間尺度都是在10納秒(1 ns=l(T9 s )的數量級上(Xu， Y.， F. J. Barrantes, et al. (2005). "Conformational dynamics of the nicotinic acetylcholine receptor channel: a 35-ns molecular dynamics simulation study," J Am Chem Soc 127(4): 1291-1299 ),甚至在1998年4艮道了對蛋 白質折疊途徑的樣i秒級(1 jis-l(T6 s )的分子動力學模擬(Heymann， B. and H. Grubmiiller (2000). "Dynamic force spectroscopy of molecular adhesion bonds." Phvs Rev Lett 84(26): 6126-6129 )。然而目前常規(guī)分子 動力學模擬的主要挑戰(zhàn)仍然來自時間尺度的巨大差異，因為我們所感 興趣的生物反應過程如蛋白質的折疊 (Xu, Y., J. Shen, et al. (2005). "Conformational transition of amyloid b-peptide." Proc Natl Acad Sci U S A 102(15): 5403-5407)、配體一受體相互作用、以及離子通道的開合 時間范圍往往是從毫秒(1 ms=l(T3 s)到秒(Adcock, S. A. and J. A. McCammon (2006). "Molecular dynamics: survey of methods for simulating the activity of proteins." Chem Rev 106(5): 1589-1615 )。為了研究分子的各種運動及其引起的構象變化與生物功能之間的 關系，除了上述的常規(guī)分子動力學模擬外，拉伸分子動力學(Steered Molecular Dynamics, SMD) (Grubmiiller, H.， B. Heymann, et al. (1996). "Ligand binding: molecular mechanics calculation of the streptavidin-biotin rupture force." Science 271(5251): 997-999; Isralewitz, B., S. Izrailev, etal. (1997). "Binding pathway of retinal to bacterio-opsin: a prediction by molecular dynamics simulations." Biophvs J 73(6): 2972-2979 )、把向分 子動力學(Targeted Molecular Dynamics， TMD) ( Engels， M., E. Jacoby， et al. (1992). "The T<~>R structural transition of insulin; pathways suggested by targeted energy minimization." Protein Eng 5(7》669-677; Schlitter， J.， M. Engels, et al. (1994). "Targeted molecular dynamics: a new approach for searching pathways of conformational transitions." J Mol Graph 12(2): 84-89 )、主成分分子動力學(Essential Dynamics) ( Amadei， A., A. B. Linssen， et al. (1993). "Essential dynamics of proteins." Proteins 17(4): 412-425 )和構象漲躍方法(Conformational Flooding Method) (Grubmiiller, H. (1995). "Predicting slow structural transitions in macromolecular systems: conformational flooding." Phys Rev E 52(3): 2893-2906)等模擬方法中通過引入外力和減少蛋白質構象空間自由度 等策略，使模擬的時間范圍比常規(guī)分子動力學寬，因而可以模擬配體 和受體結合/解離(Shen, L.， J. Shen, et al. (2003). "Steered molecular dynamics simulation on the binding of NNRTI to HIV國l RT." Biophvs J 84(6): 3547-3563; Xu， Y.， J. Shen, et al. (2003). "How does huperzine A enter and leave the binding gorge of acetylcholinesterase Steered molecular dynamics simulations." J Am Chem Soc 125(37): 11340-11349; Li, W.， H. Liu， et al. (2005). "Possible pathway(s) of testosterone egress from the active site of cytochrome P450 2B1: a steered molecular dynamics simulation." Drug Metab Dispos 33(7): 910-919; Niu， C., Y. Xu， et al. (2005). "Dynamic mechanism of E2020 binding to acetylcholinesterase: a steered molecular dynamics simulation." J Phvs Chem B 109(49): 23730-23738; Li, W.， H, Liu， et al. (2007). "Possible pathway(s) ofmetyrapone egress from the active site of cytochrome P450 3A4: a molecular dynamics simulation." Drug Metab Dispos 35(4): 689-696 )、蛋白質的折疊和去折疊等生物現象。以上方法均是針對不考慮生物大分子發(fā)生轉動運動進行的模擬， 而這些方法具有研究對象和研究運動的局限性。由于近來對離子通道 家族(ligand-gated ion channel superfamily)的研究發(fā)展，尤其是隨著 對煙堿樣乙酰膽堿受體(nicotinic acetylcholine receptor, nAChR)三維 結構的解沖斤 (Unwin, N. (2005). "Refined structure of the nicotinic acetylcholine receptor at 4A resolution." J Mol Biol 346(4): 967-989)，為 深入研究這類分子的門控機制提供了基礎 (Purohit, P., A. Mitra, et al. (2007). "A stepwise mechanism for acetylcholine receptor channel gating." Nature 446(7138): 930-933; Szarecka, A.， Y. Xu， et al. (2007). "Dynamics of heteropentameric nicotinic acetylcholine receptor: implications of the gating mechanism." Proteins 68(4): 948-960 ), —系列基于分子結構才莫態(tài) 分析的研究表明乙酰膽堿受體的結構特征表明它有進行轉動運動的可 負fe (Taly， A.， M. Delarue， et al. (2005). "Normal mode analysis suggests a quaternary twist model for the nicotinic receptor gating mechanism." Biophys J 88(6): 3954-3965; Cheng, X.， B. Lu, et al. (2006). "Channel opening motion of a7 nicotinic acetylcholine receptor as suggested by normal mode analysis." J Mol Biol 355(2): 310-324; Szarecka, A" Y. Xu， et al. (2007). "Dynamics of heteropentameric nicotinic acetylcholine receptor: implications of the gating mechanism." Proteins 68(4): 948-960), 這就需要對這類分子的轉動加以研究和模擬，以達到研究結構運動(繞 軸轉動) 一功能(離子通道開合)之間相互作用的效果。然而，目前尚未有針對模擬生物大分子的區(qū)域或者全體的轉動運動的模擬方法， 尤其是還沒有提出一種運用角動量守恒的方法來進行才莫擬的方法。發(fā)明內容為實現生物大分子的模擬，從而達到研究分子運動—觀察構象變 化—實現生理功能的效果，本發(fā)明的目的主要是提供一種模擬生物大 分子的區(qū)域或者整體的轉動的方法，尤其是一種運用角動量守恒的方 法進行模擬的方法，可以在常規(guī)分子動力學基礎上模擬轉動運動，進 而達到研究分子構象變化一功能的目的。為實現該目的，主要是采用角動量守恒的原理，根據分子的構象 得到瞬間轉動慣量，通過實時調整分子的角速度來達到模擬轉動運動。具體采用以下步驟(a )對生物大分子的指定區(qū)域整體沿主軸方向施加一個初始角速度；(b) 根據角動量守恒原理，通過計算該區(qū)域的瞬間轉動慣量，調 整瞬間角速度的大小，使角動量始終保持不變；(c) 通過模擬該區(qū)域的轉動，引起整個生物大分子的構象變化， 包括離子通道的開合變化，進而研究結構構象-功能之間的關系。本發(fā)明的 一個優(yōu)選實施方式主要包括以下步驟l)選擇要模擬旋轉的區(qū)域，得到這些區(qū)域的原子的數目(w)、類 型和位置信息(xj,，z,)。對于不同類型的原子有不同的質量(w,)， 而質量是影響轉動慣量的因素之一。同時，對于這些原子根據需要， 可以采用不同的權重系數(G)來強調或者弱化其在形成整體轉動慣量的作用，優(yōu)選地。=1。2) 設計旋轉軸的位置和方向。根據原子的坐標位置，可以得到系享c,w,x,2 ;c,w 鄉(xiāng)克合勺專爭動'貫量(/) ， /= Sc;w,.乂2 ^c,附,yA冬c，,'z, 本c,w, gc 2由此轉動慣量，通過求解特征值問題，就可以得到三個正交的特征向量，即為該體系的三個相互垂直的主軸(ZJ，Z),對應地，每個特征向量都對應有一個特征值，即為各個主軸方向的轉動慣量(/》/y,/z)。特別地，對于離子通道蛋白而言，其中一^f艮主軸方向和通道軸(pore axis)是重合的。3) 根據研究對象的特征和可支配計算資源的大小，合理設計轉動 的初始角速度3，優(yōu)選地27txl()4rad.ps"，再由公式￡ = /3，計算得出角 動量Z。以繞主軸Z軸方向轉動為例，設計初始角速度為&,才艮據公式 ￡z =/z& ,得到在沖莫擬中一直保持不變的沿主軸Z軸方向的角動量&。4) 最后運動的粒子隨著時間變化(xj,.,z,.) — ，因此瞬g。附X2 本c",々；;c,w,々；間轉動'lf量<formula>formula see original document page 10</formula>￡z=/z&，按需要設計初始角速度&，然后通過實時調整角速度 來保持角動量4守恒。由于分子是時刻都在運動的，分子內原子的位置也發(fā)生變化(A乂A) — ,由公式(1)，可以得到^ —/'z。<formula>formula see original document page 10</formula>因此，角速度也要實時調整^ — &。通過上述方法，我們就可以對給定區(qū)域或者整個分子繞主軸的轉 動運動進行才莫擬了，并且通過調節(jié)初始角速度&,達到模擬各種轉速 的轉動。本發(fā)明中發(fā)明者基于角動量守恒原理，在常規(guī)分子動力學模擬的 方法上，發(fā)展了模擬生物大分子繞主軸轉動的方法。 一方面，我們加 速了在常規(guī)分子動力學模擬的時間尺度內難以觀測到的轉動運動，凸現了這種轉動對分子構象的影響；另一方面，通過保持沿轉動方向的 角動量守恒的方式，考慮實時的分子構象(即各原子坐標位置)，調 整整體角速度的大小，避免了破壞分子內部結構。我們的模擬方法運 用于離子通道家族的生物大分子，對研究通道的門控機制和位置具有 重要意義和良好的應用前景。


圖1為以平行于細胞膜的方向視角，在煙堿樣乙酰膽堿受體的膜 外區(qū)域(黑色)沿通道軸(Z軸)向按角動量守恒施加轉動；圖2為圖1的俯^L圖。
具體實施方式
基于傳統(tǒng)分子動力學方法的旋轉動力學方法的建立主要涉及四個 步驟1)選定模擬的區(qū)域對于不同的體系結構，按照研究的目的和需要，選取部分或者全部的原子(")，進而獲得原子的類型(m,.)和坐標(義,，乂,z,);再根據 關系，決定每個權重系數(c,.);最后，得到整個旋轉區(qū)域的轉動慣量 (/)。具體計算方法如公式(1)所示2c,m,x,2 S。m^x,),以圖1和圖2所示為例，我們選取了實施例煙堿樣乙酰膽堿受體的膜外 區(qū)域做為轉動的模擬區(qū)域(黑色部分)。2) 計算慣性主軸和沿主軸方向的轉動分量。對(l)中3x3矩陣求解特征值問題，解得3個特征矢量和3個特征 值，即為該轉動慣量的三個互相垂直的主軸(《F,Z)和對應的主軸方 向的轉動分量(/x/k/z)。再根據結構的特征和通道軸的取向，確定模 擬繞哪根主軸轉動。尤其對于離子通道家族的蛋白質而言，其分子是 由三元、.四元和五元亞基的聚合而成，繞通道的軸線對稱，因此，結 構的其中一個主軸(如圖1所示設為Z軸)即為通道軸。3) 設計初始角速度和計算角動量。根據計算資源和模擬體系的大小，參考實際生理過程的尺度，設 計沿Z軸方向的初始角速度&。舉例來說，如果繞Z軸轉動一圏的周期 是10納秒(lns=103ps)，那么角速度大小就是27^10-^&(1，-1。根據 角動量的定義可以得到(3.1)其中，巧為粒子測主軸Z的距離。這樣，第一項就是("2>,,2>2=/>2，、^ 」而由于結構關于z軸的對稱性，第二項可消去。因此，我們可以得到，Zz = /Z5Z (3.2)由步驟l)中轉動區(qū)域的原子信息和步驟2)中解特征值的方法，可以 得到/z;再由設計的初始角速度&，即可獲得繞主軸Z軸轉動的角動量4了。4)調整角速度&，使角動量4不變。在模擬轉動的過程中，需要實時調整角速度&，以保持角動量^ 不變。這是由于兩方面原因引起的 一方面，由于分子的構象一直是 處于動態(tài)的變化中，轉動區(qū)域的原子坐標也發(fā)生變化 0c,j^,)—(x',.，力z',),根據公式(1.1 )，繞主軸Z軸的轉動慣量也發(fā)生了 改變/z —/'z;另一方面，模擬的體系往往是耗散系統(tǒng)，角動量4會隨 著時間有衰減。因此，需要在模擬的過程中，根據4 =/z& (4.1 )實時調整角速度&,以保持在模擬轉動的過程中角動量&不變。然后在實施例中，釆用我們實驗室已發(fā)表的文獻(Xu, Y.， F.丄 Barrantes， et al. (2005). "Conformational dynamics of the nicotinic acetylcholine receptor channel: a 35-ns molecular dynamics simulation studv." J Am Chem Soc 127(4): 1291-1299)中的方法，對 在包脂的水溶液環(huán)境里;f莫擬了 20納秒的分子動力學的煙堿樣乙酰膽堿 受體結構，選取其膜外區(qū)域的結構以初始角速大小為5xl(T4 rad.ps"，方向為逆時針(俯視)繞離子通道軸向進行轉動了模擬(如圖1和圖2 所示)。本發(fā)明中我們在使用Gromacs程序包(Berendsen, H.丄C., D. van der Spoel' et al. (1995). "GROMACS: a message-passing parallel molecular dynamics implementation." Comput Phvs Commun 91(1-3): 43-56)進行常規(guī)分子動力學模擬基礎上,使用對選 定區(qū)域內的系統(tǒng)按角動量守恒的方式進行繞主軸轉動的模擬。以煙堿 樣乙酰膽堿受體為實施例，對膜外部分施加了逆時針方向(俯視)的 轉動后，離子通道的門控位置發(fā)生了明顯的打開(如下表所示)，得到 了很好結果。離子通道的門起始半徑(A)轉動后半徑增加的半徑擴大的面積控位置(A)(A)(%)L2511.702.500.8116V2552.052.750.781以上對本發(fā)明的實施例煙堿樣乙酰膽堿受體的膜外區(qū)域當前優(yōu)選 的實施方式進行了說明。本發(fā)明不僅僅局限于以上實施例，我們認為 考察本文提出的應用角動量守恒原理模擬生物大分子運動的方法后將 會對本領域技術人員提示其他的修改、變形和改變。因此，所有這些范圍內。
權利要求
1、一種應用角動量守恒原理模擬生物大分子運動的方法，主要包括以下步驟(a)對生物大分子的指定區(qū)域整體沿主軸方向施加一個初始角速度；(b)通過調整該方向角速度的大小來控制轉動；(c)通過模擬該區(qū)域的轉動，引起整個生物大分子的構象變化。
2、 根據權利要求1所述的方法，其特征在于在模擬體系轉動過 程中，保持角動量不變，不會破壞體系的平衡和扭曲體系的內在關聯。
3、 根據權利要求1所述的方法，其特征在于通過觀察離子通道 的開合變化，研究離子通道蛋白的門控機制。
4、 根據權利要求2所述的方法，其特征在于保持角動量不變， 要計算上述系統(tǒng)的轉動慣量/為<formula>formula see original document page 2</formula>
5、 根據權利要求2所述的方法，其特征在于運用保持角動量不 變，實時調整瞬間角速度6 —以保持轉動。
6、 根據權利要求l或4所述的方法，其特征是對于持續(xù)運動的 粒子團，它們的轉動慣量隨著它們位置的變化也發(fā)生著變化/—/',本發(fā)明通過轉動慣量/這個物理量，來反映和衡量體系內所有粒子的整體 性變化。
7、 一種應用角動量守恒原理模擬生物大分子運動的方法，主要包 括以下步驟l)對于不同的體系結構，按照研究的目的和需要，選取部分或者 全部的原子"，進而獲得原子的類型叫.和坐標(<formula>formula see original document page 3</formula>)再根據關系， 決定每個權重系數ci; 最后，通過如下公式(1)得到整個旋轉區(qū)域的 轉動慣量/:<formula>formula see original document page 3</formula>2) 計算慣性主軸和沿主軸方向的轉動分量，對(1)中3x3矩陣 求解特征值，解得3個特征矢量和3個特征值，即為該轉動慣量的三個 互相垂直的主軸(《KZ)和對應的主軸方向的轉動分量(/乂/y，/z),再 根據結構的特征和通道軸的取向，確定模擬繞哪根主軸轉動；3) 設計初始角速度并計算角動量，根據研究對象的特征和計算資 源的大小，設計初始角速度3，再根據Z-W的關系計算得出角動量Z;4) 通過實時調整角速度來保持角動量Z守恒，當模擬轉動繞主軸Z 方向轉動時，設計初始角速度為&，根據公式得到￡z =，由于分 子是時刻都在運動的，分子內原子的位置也由(;c,j^ )變化至(), 由公式(1)，可以得到轉動慣量/z變化至/'z，根據角動量守恒 Zz=/^z=/;&,實時調整角速度&變化至g。
8、根據權利要求7所述的應用角動量守恒原理模擬生物大分子運動的方法，其特征在于，所述步驟U中權重系數c尸l;所述步驟3) 中初始角速度6 =2ti x l (T4 rad.ps"。
全文摘要
本發(fā)明提供一種模擬生物大分子運動的方法，適用于研究離子通道開合和其他因此引起的構象的變化。首先，對生物大分子的指定區(qū)域整體沿主軸方向施加一個初始角速度；然后，通過調整該主軸方向的角速度的大小來控制轉動，優(yōu)選地根據角動量守恒原理，通過計算該區(qū)域的瞬間轉動慣量，調整瞬間角速度的大小，使角動量始終保持不變；最終，通過模擬該區(qū)域的轉動，引起整個生物大分子的構象變化，優(yōu)先地分析離子通道的開合變化，進而研究結構構象-功能之間的關系。
文檔編號G06F19/00GK101276382SQ20081003602
公開日2008年10月1日 申請日期2008年4月14日 優(yōu)先權日2008年4月14日
發(fā)明者劉信力, 王希誠, 蔣華良 申請人:中國科學院上海藥物研究所;大連理工大學