qmmm組合理論及其發(fā)展

1凝聚態(tài)體系的化學(xué)表征從帶頭算的方法是化學(xué)和化學(xué)生物學(xué)研究的重要武器。它能準(zhǔn)確、詳細(xì)地提供信息，如結(jié)構(gòu)、能量變化和反應(yīng)活性，對測試結(jié)果起到補(bǔ)充和驗(yàn)證的作用，有助于識別現(xiàn)象的本質(zhì)。經(jīng)過幾十年的努力，隨著計(jì)算機(jī)科學(xué)技術(shù)的快速發(fā)展和計(jì)算算法的不斷優(yōu)化，使用頭上算方法可以有效地處理包括100個(gè)非氫原子系統(tǒng)。然而，由于頭算方法的長期復(fù)雜性通常存在于氮n7（n是原子數(shù)）的數(shù)量。從現(xiàn)在到將來，它仍然無法使用頂部算術(shù)方法來處理包含1000多個(gè)非氫原子的系統(tǒng)。從頭算方法在體系大小方面受到的限制嚴(yán)重地阻礙了人們對于凝聚態(tài)體系以及生命化學(xué)體系的理論研究,而目前有關(guān)凝聚態(tài)體系(如納米材料)以及生命化學(xué)體系(如酶與分子機(jī)器)的研究正逐漸成為化學(xué)領(lǐng)域中的重點(diǎn).為此,科研人員往往不得不采用較為粗糙的、純屬經(jīng)驗(yàn)性質(zhì)的分子力學(xué)方法來處理大型的化學(xué)體系.這些分子力學(xué)方法簡單地將體系用彈簧力加范德華力束縛起來,的確能夠提供一些有關(guān)分子構(gòu)象以及靜電勢分布的信息.但是,除了計(jì)算精度上的問題以外,分子力學(xué)方法還有一個(gè)根本的缺陷,就是不能有效地描述化學(xué)鍵的形成和斷裂過程.因此,當(dāng)我們關(guān)心的是體系的“化學(xué)”問題(也即體系中的化學(xué)鍵形成與斷裂問題)時(shí),傳統(tǒng)的分子力學(xué)方法是無能為力的.一個(gè)巧妙地解決上述一系列矛盾的方法是將量子力學(xué)從頭算方法與分子力學(xué)方法有機(jī)地結(jié)合起來,對于體系中重要的、涉及到化學(xué)鍵變化的部分使用量子力學(xué)從頭算方法來描述,而對于體系中次重要的、不直接涉及到化學(xué)鍵變化的部分使用分子力學(xué)方法來描述.這種方法就是量子力學(xué)和分子力學(xué)聯(lián)用方法,簡稱QM/MM方法.該方法在科學(xué)計(jì)算理論中有著簡單而又正確的邏輯,也就是對于主要矛盾部分進(jìn)行精算,而對于次要矛盾部分進(jìn)行粗算.該方法在計(jì)算的可行性方面也邁出了一大步.原則上講,在目前的條件下使用QM/MM方法處理含有上萬個(gè)非氫原子的體系是完全可行的.最早有關(guān)QM/MM方法的研究由Warshel等人在上個(gè)世紀(jì)70年代開展的.經(jīng)過隨后近30年的改進(jìn),QM/MM方法已經(jīng)從理論探討逐漸走上實(shí)踐應(yīng)用.直到今天,QM/MM方法仍然尚未發(fā)展成為一種成熟可靠的、規(guī)范化和商業(yè)化的理論工具,我們應(yīng)該看到這一方法在未來化學(xué)研究中的巨大潛力和長遠(yuǎn)價(jià)值.因此,學(xué)習(xí)、掌握并發(fā)展QM/MM方法,在目前有著重要的意義,機(jī)不可失.本文扼要地闡述了QM/MM方法的基本理論、實(shí)現(xiàn)方法以及使用范圍,希望能夠?yàn)檫@方面的深入研究提供有益的信息.2哈密頓算符及能量使用QM/MM方法,第一步是將整個(gè)體系劃分為兩個(gè)部分.其中一個(gè)部分我們稱之為活性中心,它是我們將要考察的關(guān)鍵部位,由那些直接參與化學(xué)鍵形成或斷裂的原子構(gòu)成.對于活性中心我們采用了量子力學(xué)從頭算的方法來描述,因此該部分又稱之為QM區(qū)域.活性中心之外的所有原子構(gòu)成了外圍環(huán)境.對于外圍環(huán)境我們采用了分子力學(xué)方法來描述,因此該部分又稱之為MM區(qū)域.描述用復(fù)合勢場處理的體系時(shí),使用有效哈密頓Heff算符.在這個(gè)公式中,可以通過求解薛定諤方程,得到給定分子結(jié)構(gòu)的體系能量Etotal,?effΨ(ri,Rα,RM)=Etotal(Rα,RM)Ψ(ri,Rα,RM)(1)其中,ri和Rα為電子和核子在QM區(qū)的系數(shù);RM為原子在MM區(qū)的系數(shù).整個(gè)體系的有效混合哈密頓算符可表示為:?Heff=?HQΜ(QΜ)+?HΜΜ(ΜΜ)+?HQΜ/ΜΜ(QΜ/ΜΜ)(2)第一項(xiàng)為QM區(qū)的哈密頓算符,可以用標(biāo)準(zhǔn)的量子力學(xué)形式表示.例如用abinitio處理時(shí),?HQΜ=-12∑i?i2+∑ij1rij-∑iαΖαriα+∑αβΖαΖβrαβ(3)其中,Zα為α原子核的電荷;rij為粒子i和j之間的距離.下標(biāo)(i,j)和(α,β)分別為電子和原子核.第二項(xiàng)為MM區(qū)哈密頓算符HMM,而且只在MM原子間有影響,體系的能量只跟原子坐標(biāo)有關(guān),可以表示為鍵的伸縮Ebond、鍵彎曲Eangle、鍵扭曲Edihedral以及非鍵相互作用,如靜電作用和Lennard-Jones相互作用.?HΜΜ=Ebond+Eangle+Edihedral+ELennard-Jones+EClulomb=∑iΚb(bi-b0)2+∑iΚθ(θi-θ0)2+∑iΚφcos[n(?i-?0)]+∑i>jQiQjrij+∑i>jεij[(rij0rij)12-2(rij0rij)6](4)第三項(xiàng)為QM區(qū)與MM區(qū)之間相互作用的哈密頓算符,它構(gòu)成了QM/MM方法的核心,但是它依賴于QM和MM的形式.對此我們將作詳細(xì)的介紹.當(dāng)我們在劃分QM/MM區(qū)域時(shí),必然涉及邊緣邊界問題.在模擬液相體系時(shí),一般選取數(shù)百個(gè)溶劑分子組成的正方體或長方體箱子作為MM區(qū),溶質(zhì)分子則自然成為QM體系,而邊界區(qū)域則通過溶劑箱的周期平移得到,把邊界區(qū)域看成連續(xù)介質(zhì)模型.這樣通過式(1)～(4),整個(gè)體系的能量(Etotal)由下式給出:Etotal=EQΜ(QΜ)+EΜΜ(ΜΜ)+EQΜ/ΜΜ(QΜ/ΜΜ)(5)式中,EQM(QM)為QM區(qū)域的能量,取決于所需的精度,這一部分能量可以由某一種量子力學(xué)方法(如半經(jīng)驗(yàn)方法、Hartree-Fock方法、M?ller-Plesset微擾法以及DFT方法等)給出;EMM(MM)為MM區(qū)域的能量,通常可以使用某一力場(如AMBER和CHARMM力場)計(jì)算得到;EQM/MM(QM/MM)為QM區(qū)域與MM區(qū)域之間的相互作用能.一旦總能量知道了,原子之間的力參數(shù)很容易通過能量對原子坐標(biāo)的微分得到:Fα=-?Etotal?Rα(6)FΜ=-?Etotal?RΜ(7)得到總的體系能量以及相關(guān)的力參數(shù)后,就可以處理分子幾何優(yōu)化、分子動力學(xué)研究以及MonteCarlo模擬等問題了.對EQM/MM(QM/MM)的計(jì)算則取決于研究的體系,對不同的體系要使用不同的計(jì)算方法.①Q(mào)M區(qū)域與MM區(qū)域之間沒有共價(jià)相連.這種情況多見于溶液體系,其中一個(gè)溶質(zhì)分子(QM區(qū)域)被大量溶劑分子(MM區(qū)域)包圍,而溶質(zhì)分子與溶劑分子之間并沒有共價(jià)相連.在這種體系中,QM區(qū)域與MM區(qū)域之間的相互作用能主要來自于它們之間的靜電吸引(或排斥)和范德華作用,也即:EQΜ/ΜΜ(QΜ/ΜΜ)=Eelectrostatic(QΜ/ΜΜ)+EvdW(QΜ/ΜΜ)(8)其中范德華作用能(EvdW(QM/MM))可以用經(jīng)驗(yàn)的公式(如12～6次Lennard-Jones勢場)求算,其表達(dá)式為:EvdW=S∑s=1Μ∑m=14εms[(σmsRms)12-(σmsRms)6](9)σms、εms為Lennard-Jones參數(shù),且σms=(σmσs)1/2,εms=(εmεs)1/2,σm、σs可以從MM力場中獲得,εm、εs為QM原子的經(jīng)驗(yàn)參數(shù).而靜電相互作用能(Eelectrostatic(QM/MM))可以首先使用分子力場方法求出QM區(qū)域和MM區(qū)域中每個(gè)原子所帶的點(diǎn)電荷,然后使用庫侖定律得到總的靜電相互作用,其表達(dá)式為:Eelectrostatic(QΜ/ΜΜ)=-S∑s=1Ν∑i=1qSris+S∑s=1Μ∑m=1qSΖmRms(10)式中,qS為MM原子上的點(diǎn)電荷;Zm為QM原子的核電荷;ris、Rms分別為MM點(diǎn)電荷與QM電子、原子核間的距離.另一種更為精確的求算Eelectrostatic(QM/MM)的方法是將MM區(qū)域中每個(gè)原子所帶的經(jīng)驗(yàn)點(diǎn)電荷直接寫進(jìn)QM區(qū)域的Schr?dinger方程里.由此我們得到有效Hamilton方程:Ηeff=-12Νeff∑i?i2+Νeff∑i≠j1rij-Νeff∑i∑α∈QΜΖαrαi-Νeff∑i∑β∈ΜΜqβriβ+∑α1≠α2∈QΜΖα1Ζα2rα1α2+∑α∈QΜ?β∈ΜΜΖαqβrαβ(11)式中,Neff為QM區(qū)域總的電子數(shù)目;rij為電子i到電子j的距離;Zα為QM區(qū)域中原子α的核電荷;rαi為電子i到原子核α的距離;qβ為MM區(qū)域中原子β所帶的經(jīng)驗(yàn)點(diǎn)電荷;rβi為電子i到原子β的距離;rα1α2為原子核α1到原子核α2的距離;rαβ為原子核α到原子β的距離.顯而易見,使用上式求解QM區(qū)域的Schr?dinger方程時(shí),Eelectrostatic(QM/MM)已經(jīng)被之間包括在體系總的能量中,即:EQΜ(QΜ)+Eelectrostatic(QΜ/ΜΜ)=?ψ|Ηeff|ψ?(12)②QM區(qū)域與MM區(qū)域之間有共價(jià)相連.這種情況較前一種情況更為常見,如在酶化學(xué)中反應(yīng)中心與整個(gè)蛋白質(zhì)分子并非是兩個(gè)彼此獨(dú)立的單元,而是通過共價(jià)鍵相互連接在一起.在這種情況下我們不能將QM區(qū)域與MM區(qū)域之間的化學(xué)鍵簡單切斷,否則造成末端為自由基的分子片斷,從而導(dǎo)致荒謬的計(jì)算結(jié)果.為了飽和QM區(qū)域中與MM區(qū)域相鍵連的原子(以下稱之為QM末端原子)的化學(xué)價(jià),人們提出了以下幾種方法.第一,使用氫原子來飽和QM末端原子,也就是由Singh、Field等人發(fā)展的連接原子法,以及在此基礎(chǔ)上由Morokuma及其合作者發(fā)展的IMOMM方法和INIOM方法.這些方法簡明易行,但也有先天的缺陷,而且它對連接原子的選擇有嚴(yán)重的依賴性.產(chǎn)生缺陷的主要原因來自于氫原子獨(dú)特的成鍵能力(即使用1s軌道成鍵)及其較小的電負(fù)性.這使得QM末端原子不僅在成鍵后電子云分布與真實(shí)情況有著較大的差異,而且其電荷轉(zhuǎn)移情況(即由QM區(qū)域向MM區(qū)域轉(zhuǎn)移電子,或者QM區(qū)域從MM區(qū)域得到電子)往往失真.第二,使用有效原子對QM末端原子進(jìn)行“蓋帽”.該原子帶有一定數(shù)目的電子、一定量的核電荷以及一定的有效電勢場.使用一定數(shù)目的電子以滿足兩條基本要求,一是該有效原子只能與QM末端原子(通常是sp3雜化的碳原子)形成一根σ單鍵,二是該有效原子的外層電子云基本上保持球?qū)ΨQ.使用一定量的核電荷則是為了使得該有效原子的總電荷平衡.最后,使用一定的有效電勢場可以調(diào)節(jié)QM區(qū)域與MM區(qū)域之間的電荷轉(zhuǎn)移.值得指出的是,在目前主要有兩種方法來實(shí)現(xiàn)有效電勢場.一種方法是固定軌道法,即在有效原子上人為地添加一些系數(shù)一定的已經(jīng)填充好電子的原子軌道,比較有代表性的是Rivail及其合作者的LSCF方法,以及在此基礎(chǔ)上由Gao發(fā)展的廣義雜化軌道理論(GHO).由于在組成分子軌道時(shí)這些原子軌道的系數(shù)并不被優(yōu)化,這相當(dāng)于人為地固定了一些具有一定空間分布的負(fù)電荷勢場在有效原子周圍.另一種方法是有效核勢場法.在這種方法中,人為地強(qiáng)加一些描述電勢場的函數(shù)在有效原子的核上.通過調(diào)節(jié)函數(shù)的性質(zhì)和參數(shù)可以很好地模擬再現(xiàn)QM區(qū)域與MM區(qū)域之間的電荷轉(zhuǎn)移,這種方法比較有代表性就是Yang及其合作者的pseudobond法.顯而易見,使用有效原子的方法的根本思路就是最大程度地再現(xiàn)QM區(qū)域與MM區(qū)域之間的共價(jià)鍵特征.盡管我們可以使用多種方法來添加有效原子,但實(shí)質(zhì)上各種方法都是通過某種參數(shù)調(diào)節(jié)來實(shí)現(xiàn)最優(yōu)的擬合.因此,好的有效原子方法不僅要擬合程度高,而且要參數(shù)化過程簡單.3u3000改進(jìn)的yps-mbr算法根據(jù)以上的基本原理和簡便參數(shù)化的思想,我們在這里具體地給出了一種QM/MM的實(shí)現(xiàn)方法.該方法的核心是解決QM區(qū)域與MM區(qū)域之間的共價(jià)相連問題.而針對QM區(qū)域與MM區(qū)域之間沒有共價(jià)相連這種簡單情況,我們假設(shè)下述的方法也可以直接應(yīng)用.首先,考慮到計(jì)算精度與計(jì)算速度的關(guān)系,我們決定采用DFT方法(B3LYP)來描述QM區(qū)域.這種方法與Hartree-Fock方法有著近似的時(shí)間復(fù)雜度,卻能夠較好描述體系的電子相關(guān)能.其次,由于我們的QM/MM方法主要針對于酶化學(xué),我們決定采用CHARMM力場.這種力場在處理多肽與蛋白質(zhì)體系時(shí)有著較好的能力.最后,為了簡化參數(shù)化的過程,我們將“蓋帽”原子定義為擬鹵素原子.該原子有7個(gè)核電荷,7個(gè)價(jià)電子,一套核勢場函數(shù).顯然,該原子能且僅能與sp3雜化的碳原子形成一根σ單鍵.形成單鍵后,其8(7+1)個(gè)價(jià)電子保證了球?qū)ΨQ性,并最大程度地模擬了飽和成鍵后的碳原子.我們先將分子體系做如圖1的劃分.圖1中,Active區(qū)表示反應(yīng)的活性中心,用QM來處理,Environment區(qū)表示反應(yīng)的活性中心以外的部分,用MM來處理.在pseudobond的方法中,Y原子(MM)將被單自由價(jià)的邊界原子Yps所取代,而Yps原子就是前面提到的擬鹵素原子,形成的X-Yps將會擬合原來的X-Y鍵.這個(gè)時(shí)候?qū)╕ps原子在內(nèi)的Active區(qū)域表示為QM體系,而將剩下的部分表示為MM體系.則系統(tǒng)能量式中EQM/MM(QM/MM)可以由下式給出:EQΜ/ΜΜ(QΜ/ΜΜ)=Eelectrostatics(QΜ/ΜΜ)+EvdW(QΜ/ΜΜ)+Emm-bonded(QΜ/ΜΜ)(13)其中第三項(xiàng)表示QM部分與MM部分成健的作用能,用的是分子力場的方法.所以,在計(jì)算的過程中EvdW(QM/MM)、Emm-bonded(QM/MM)和Emm(MM)是使用分子力場的方法來計(jì)算的,而用量子化學(xué)方法計(jì)算的只有EQM(QM)和Eelectrostatics(QM/MM)可以用式(12)來得到,其中有效Hamiltonian算符描述的是一組有效電子Neff,它包括Active部分的所有電子以及邊界原子Yps的價(jià)電子,Heff可以由下式給出:Ηeff=-12Νeff∑i?i2+Νeff∑i≠j1rij-Νeff∑i∑α∈QΜΖαrαi-Νeff∑i∑β∈ΜΜqβrβi+Νeff∑i∑γ∈YpsVeffYps(riγ)+∑α1≠α2∈QΜΖα1Ζα2rα1α2+∑α∈QΜ?β∈ΜΜΖαqβrαβ(14)其中,Zα為QM部分的原子核電荷;qβ為MM部分原子的點(diǎn)電荷;VeffYps為邊界原子的有效核電勢(effectivepseudopotential).所以,要想得到Heff,其中有兩項(xiàng)是未知的,就是qβ和VYpseff,求解這兩項(xiàng)的方法如下:①要得到qβ,需要注意的兩個(gè)問題是:MM體系的部分電荷的總和加上QM體系的電荷應(yīng)該等于整個(gè)體系的電荷;QM計(jì)算中怎么樣降低由于存在非常靠近QM體系的MM原子而產(chǎn)生大的靜電作用.解決的方法就是將直接和QM體系成鍵的MM原子以及和QM體系通過兩個(gè)鍵相連的原子的點(diǎn)電荷視為0,這種方法的原理是利用了分子力場中的電荷形式.通過這種定義就解決了如上提到的兩個(gè)問題.②要得到有效核電勢VeffYps.Yps原子有效核電勢的作用就是Yps-X鍵可以模擬原來的sp3雜化的X-Yσ單鍵,也就是說Yps-X鍵應(yīng)該和原來的X-Y單鍵有相似的鍵長、鍵的強(qiáng)度以及對QM部分的影響,其表達(dá)式如下:{Veff(r)=VeffL(r)+L-1∑l=0[Veffλ(r)-VeffL(r)]×∑m|λm??λm|VeffL(r)=1raLe-bLr2Veffλ(r)-VeffL(r)=aλe-bλr2(λ=0,1,?,L-1)(15)其中L為基組中的最大λ值;ai和bi分別為適當(dāng)?shù)膮?shù).有效核電勢的參數(shù)化過程就是讓Yps-X鍵的參數(shù)和原來的X-Y鍵的參數(shù)擬合得最好,這樣在QM的計(jì)算過程中使得Yps和MM部分的作用是相當(dāng)?shù)?當(dāng)?shù)玫竭@些參數(shù)后再對其他的化學(xué)體系進(jìn)行檢驗(yàn).而且值得注意的是有效核電勢的參數(shù)化過程不是在分子力場的計(jì)算中得到的,而是在純QM計(jì)算過程中得到的.這種算法的和其他方法相比,明顯的優(yōu)勢就是參數(shù)的通用性很強(qiáng),而且對不同體系計(jì)算得到的結(jié)果與用純QM計(jì)算得到的結(jié)果是相近的.4目標(biāo)函數(shù)的確定PseudobondQM/MM計(jì)算的流程如圖2.具體的步驟如下:①從蛋白質(zhì)、核酸等數(shù)據(jù)庫下載相應(yīng)的分子晶體結(jié)構(gòu)文件,然后通過一系列的能量最小化、平衡步驟,從中選擇幾個(gè)比較合理的分子結(jié)構(gòu),對每一個(gè)分子進(jìn)行如下的步驟.②先用Tinker程序計(jì)算MM部分的能量、電荷以及能量梯度等數(shù)據(jù),它們將在下面的計(jì)算過程中起到非常重要的作用.③用Gaussian98程序?qū)M部分進(jìn)行優(yōu)化,在這個(gè)優(yōu)化過程中MM部分是固定不動的.并且這個(gè)過程中Yps原子用F原子來表示,它的贗勢場目前只有3-21G和6-31G*兩個(gè)基組下的參數(shù),這些參數(shù)在文獻(xiàn)中有詳細(xì)的說明.④用tinker中的Newton程序?qū)M部分進(jìn)行構(gòu)型的優(yōu)化,同樣在這個(gè)過程中QM部分是固定不動的.循環(huán)進(jìn)行上述三個(gè)過程直到系統(tǒng)總的能量以及RMSD(均方根偏差)均收斂為止,其中RMSD的計(jì)算表達(dá)式如下:RΜSD=√1Νatoms∑i=1(ri-ri0)2(16)它表明不同構(gòu)象的兩個(gè)分子間的構(gòu)象差異.至此,對此構(gòu)象分子的優(yōu)化過程結(jié)束.⑤要研究一個(gè)反應(yīng)的機(jī)理,也就是反應(yīng)路徑的搜索步驟,必須選取合適的反應(yīng)坐標(biāo),反應(yīng)坐標(biāo)可以是兩個(gè)原子間的距離、兩組原子的距離差、三個(gè)原子間的夾角以及四個(gè)原子間的二面角等參數(shù).一旦反應(yīng)坐標(biāo)確定后,就可以對反應(yīng)路徑進(jìn)行一系列的研究.在這個(gè)過程中對反應(yīng)坐標(biāo)上的每一個(gè)構(gòu)型都要進(jìn)行上述②～④步驟的優(yōu)化.在優(yōu)化的過程中為了讓反應(yīng)坐標(biāo)穩(wěn)定在我們預(yù)先設(shè)定的某個(gè)值上,必須在體系的能量函數(shù)中引入約束能這一能量項(xiàng),它的表達(dá)式如下:ERestrain=k(R-s)2(17)其中,k為力常數(shù);R為反應(yīng)坐標(biāo);s為反應(yīng)坐標(biāo)變化的步長.這樣做的目的就是迫使分子體系向預(yù)定反應(yīng)路徑的方向進(jìn)行優(yōu)化.通過如上的步驟就可以得到整個(gè)反應(yīng)路徑的能量-反應(yīng)坐標(biāo)曲線.⑥得到了反應(yīng)路徑,就相當(dāng)于掌握了反應(yīng)的整個(gè)過程,還可以對這個(gè)反應(yīng)路徑進(jìn)行一系列的計(jì)算,比如可以計(jì)算整個(gè)反應(yīng)過程中的自由能壘,反應(yīng)中的自由能變化,這個(gè)可以用自由能微擾來得到,自由能微擾就是將整個(gè)反應(yīng)路徑分成很多的點(diǎn),然后計(jì)算相鄰兩點(diǎn)間的自由能,這樣反應(yīng)路徑中任何兩點(diǎn)間的自由能變化等于它們之間所有相鄰兩點(diǎn)間的自由能變化的疊加,對反應(yīng)路徑劃分的越細(xì),得到的自由能也越準(zhǔn)確.在步驟⑤中,我們剛好得到了反應(yīng)路徑-反應(yīng)坐標(biāo)的曲線,通過對這個(gè)曲線進(jìn)行自由能微擾的計(jì)算,就可以得到整個(gè)反應(yīng)過程的自由能變化和自由能壘等重要參數(shù),而這些參數(shù)是有實(shí)驗(yàn)數(shù)據(jù)可以考證的,所以這樣反過來又可以檢驗(yàn)得到的這個(gè)反應(yīng)路徑是否合理.在進(jìn)行自由能微擾的計(jì)算過程中值得注意的一點(diǎn)就是對反應(yīng)路徑上的每一個(gè)分子假定QM部分對整個(gè)分子體系自由能的貢獻(xiàn)是固定的,它們有前面得到的固定的分子幾何構(gòu)型和電荷分布.5生物溶液體系的模擬在此,我們重點(diǎn)討論一下利用復(fù)合QM/MM方法,進(jìn)行的分子模擬研究.目前,最主要的應(yīng)用是溶液體系中的化學(xué)現(xiàn)象模擬和生物分子的研究;QM/MM方法的最大優(yōu)勢就是能夠處理龐大而又復(fù)雜的體系.這主要包括兩個(gè)方面,即溶液體系和生物大分子的處理[5,22,23,24,25,26,27,28,29,30,31].5.1不同兩因素下ch2cl/cl-形成混合物時(shí)能量變化利用QM/MM方法,可以很好地分析整個(gè)反應(yīng)過程中參與反應(yīng)原子的電荷分布,進(jìn)一步可以說明溶劑分子對反應(yīng)的影響.而且在模擬溶液體系中的化學(xué)行為時(shí),QM/MM邊界劃分顯的比較明顯,參與反應(yīng)的分子用QM的方法處理,而反應(yīng)體系周圍的溶劑分子將用MM的方法處理;有些情況下,溶劑分子也參與反應(yīng)過程,特別是水作為溶劑的反應(yīng)過程,QM/MM的邊界劃分將顯得比較復(fù)雜.1976年,Warshel和Levitt介紹了QM/MM,并且用它研究了溶菌酶的作用機(jī)理.Singh和Kollman發(fā)展了用abinitio分子軌道理論和AMBERMM力場相組合的QM/MM方法,即abinitio/AMBER的QM/MM方法,并且研究了氯甲烷/氯負(fù)離子的交換反應(yīng),溶液中以及氣相中的聚醚的質(zhì)子化作用(protonationofpolyethers).氣態(tài)條件下,CH3Cl+Cl-形成復(fù)合物時(shí),能量有36kJ/mol的降低;但是動力學(xué)顯示,穩(wěn)定的對稱復(fù)合物之間有50kJ/mol的能壘.當(dāng)考慮溶劑化效應(yīng)時(shí),情況變得很簡單,但是沒有找到穩(wěn)定極小,只有一個(gè)簡單能壘為113kJ/mol.Singh和Kollman在QM部分采用了abinitio的方法來優(yōu)化反應(yīng)途徑中各種結(jié)構(gòu),發(fā)現(xiàn)了一個(gè)更高的反應(yīng)能壘,以及兩個(gè)沒有氣態(tài)條件下穩(wěn)定的化合物.但是這種方法的缺點(diǎn)是abinitio的QM方法太耗時(shí),并且準(zhǔn)確的同時(shí)也排除了分子動力學(xué)和MonteCarlo模擬中勢能的使用,并且僅限于靜態(tài)勢能曲線的研究.為了解決這個(gè)問題,Karplus及合作者發(fā)展了AM1/MNDO的半經(jīng)驗(yàn)QM方法和CHARMM的MM力場的組合方法.并且也進(jìn)行了氯甲烷/氯負(fù)離子交換反應(yīng)的研究.后來由Warshel等人提出的EVB(EmpiricalValenceBond)方法,在處理溶液體系中得到了廣泛的應(yīng)用,并且研究了許多關(guān)于溶液和酶的反應(yīng).Voth及其合作者對EVB方法做了進(jìn)一步的擴(kuò)展,并且應(yīng)用于水溶液體系中的質(zhì)子遷移反應(yīng)的模擬.由Gao等人在EVB的基礎(chǔ)上發(fā)展的MOVB(MolecularOrbitalandValenceBond)方法.對反應(yīng)體系PES(PotentialEnergySurfaces)的處理有很好效果.除了對溶劑化能進(jìn)行了詳細(xì)的研究外,Gao通過QM/MM方法對溶液中溶劑與溶質(zhì)間的相互作用的研究為溶劑極性不同對溶質(zhì)光譜性質(zhì)的影響提供了有力的理論依據(jù).Thompson及其合作者也對溶劑對分子光譜特性的影響進(jìn)行了研究,主要包括羰基化合物中v→n*光譜位移.這項(xiàng)工作盡管采用了AM1半經(jīng)驗(yàn)的方法,但是它包含了MM原子的極化作用,并且發(fā)展了一套完整的方法來計(jì)算MM原子的誘導(dǎo)偶極與QM原子上電荷的相互作用.利用這種方法,他們研究了Phodopseudomonasviridis光合成反應(yīng)中心的細(xì)菌葉綠素(bacteriochlorophyll)激發(fā)態(tài)的光譜特性.他們還研究了金屬K+和二甲醚(dimethylether)復(fù)合物的光譜特性,但是沒有包括MM原子的極化作用.Bakowies和Thiel也發(fā)展了考慮MM原子極化作用的QM/MM集合的方法,研究了一些氣態(tài)條件下的反應(yīng).由Morokuma及其合作者提出的ONION模型在處理大分子Diels-Alder反應(yīng)時(shí),給出了和試驗(yàn)數(shù)據(jù)較接近的結(jié)果.早期的QM/MM方法一般使用半經(jīng)驗(yàn)的方法,但是,隨著計(jì)算機(jī)的發(fā)展以及算法的優(yōu)化,使得對QM部分一般使用DFT方法或者abinitio方法.5.2qm/mm邊緣QM/MM復(fù)合的方法為我們模擬酶催化反應(yīng)提供了可能.而且,可以這么說,QM/MM方法將來最大的應(yīng)用是來解決酶催化反應(yīng)的模擬.所以近年來QM/MM復(fù)合方法在酶催化反應(yīng)中得到了廣泛的應(yīng)用.在這個(gè)體系中,酶催化的活性中心一般用QM處理,而酶分子的其他部分以及溶劑分子一般用MM的方法處理.與溶液體系中的反應(yīng)相比,酶催化反應(yīng)的模擬有兩點(diǎn)很不一樣.其一,QM/MM邊緣的區(qū)分將會變得很不明顯,我們要盡量多的用QM的方法來處理活性中心;其二,對于蛋白質(zhì)分子,由于QM/MM間有化學(xué)鍵的斷裂,而在溶液體系中這種情況是比較少的.所以在酶催化反應(yīng)體系中,QM/MM的處理方法和溶液體系有很大的不同.Warshel和Levitt在這方面作出了開創(chuàng)性的貢獻(xiàn),早期他們研究了溶菌酶lysozyme催化的糖苷鍵的glycosidic裂解過程中,碳正陽離子中間體的穩(wěn)定性,結(jié)果表明靜電作用的穩(wěn)定,對酶活性的保持起著重要的作用.在他們對酶作用機(jī)理的研究中,得出了影響酶反應(yīng)的各種因素,包括反應(yīng)中化學(xué)鍵的重排、電荷的重新分布、酶與酶作用物的空間及靜電相互作用、酶以及周圍水分子由于極化作用而導(dǎo)致的介電效應(yīng).并且對酶催化反應(yīng)中QM/MM區(qū)域的劃分提出了不少建議.隨后Warshel提出了經(jīng)驗(yàn)價(jià)鍵理論(EVB)方法,利用此方法,他們研究了一系列的酶催化反應(yīng),包括絲氨酸蛋白酶和胰蛋白酶在溶液中的催化反應(yīng).Antonczak及其合作者通過Rivail提出的LSCF方法對嗜熱菌蛋白酶的作用機(jī)理進(jìn)行了研究,并且總結(jié)了酶所在的環(huán)境對酶反應(yīng)過程的影響因素.6使用示例6.1催化問題的研究4-草酰丁烯酸酯互變異構(gòu)酶(4-OT)可以催化非軛的α-酮酸轉(zhuǎn)變?yōu)楣曹椀耐之悩?gòu)體,例如可以將2-羰基-4-烯己二酸酯通過2-羥基-2,4-二烯己二酸酯中間產(chǎn)物轉(zhuǎn)變?yōu)?-羰基-3-烯己二酸酯.反應(yīng)過程如圖3.4-OT由6個(gè)亞單元組成,其中每一個(gè)亞單元由62個(gè)氨基酸殘留基團(tuán)構(gòu)成.4-OT含有6個(gè)活性中心,他們分布在每一個(gè)亞單元中Pro-1周圍的疏水區(qū)域,同時(shí)Arg-11′、Arg-39″、Arg-61′也參與同分異構(gòu)反應(yīng)的催化過程.由于4-OT尺寸很小,而且它可以在沒有任何金屬離子的作用下進(jìn)行質(zhì)子遷移反應(yīng)的催化,所以對它的催化機(jī)理的研究也就顯得非常的重要.實(shí)驗(yàn)結(jié)果顯示反應(yīng)機(jī)理是一個(gè)廣義酸堿反應(yīng),而且其中的Pro-1氨基酸殘留基團(tuán)起到廣義堿的作用,但是廣義酸的角色一直受到爭議.反應(yīng)中心如圖4.理論計(jì)算得到Arg39″的pKa在9.8～10.3,可以起到廣義酸的作用,但是實(shí)驗(yàn)證明上面3個(gè)Arg氨基酸殘留基團(tuán)沒有一個(gè)pKa在9附近,所以尋找廣義酸基團(tuán)便成了這個(gè)工作的重點(diǎn).之前有人提出在這個(gè)反應(yīng)過程中,如圖4所示,Arg-39″和一個(gè)水分子通過和過渡態(tài)中的氧負(fù)離子形成兩個(gè)氫鍵來穩(wěn)定過渡態(tài).在這個(gè)工作中,體系一共加入2407個(gè)水分子,總共由13197個(gè)原子,在分子動力學(xué)平衡的過程中,發(fā)現(xiàn)有兩個(gè)定序的水分子與底物和Arg-61′形成氫鍵.為了尋找反應(yīng)過程中的廣義酸,在這個(gè)工作中對反應(yīng)體系進(jìn)行四組不同QM/MM界線的劃分,然后分別對這四組體系進(jìn)行QM/MM的計(jì)算,在這個(gè)過程中先在HF/3-21G的基礎(chǔ)上對體系進(jìn)行構(gòu)象優(yōu)化,然后在B3LYP/6-31G*的基礎(chǔ)上對體系進(jìn)行單點(diǎn)能的計(jì)算,然后進(jìn)行自由能微擾的計(jì)算,得出整個(gè)反應(yīng)的自由能能壘為67.76kJ/mol,比實(shí)驗(yàn)值約54.35kJ/mol低14.21kJ/mol左右.而這個(gè)實(shí)驗(yàn)值也是一個(gè)近似的結(jié)果.并且通過對反應(yīng)能壘的分解,可以得到單個(gè)氨基酸殘留基團(tuán)對過渡態(tài)穩(wěn)定的貢獻(xiàn).通過如上的對比研究可以得知在這個(gè)反應(yīng)中沒有氨基酸殘留基團(tuán)或者水分子起到廣義酸的角色,而是通過兩個(gè)定序的水分子和Arg-39″通過靜電相互作用來穩(wěn)定底物分子的負(fù)電荷,而且在這個(gè)反應(yīng)中Arg-11′在穩(wěn)定過渡態(tài)的過程中起到了非常重要的作用,同時(shí)這些結(jié)果也使以往在尋找廣義酸過程中解釋不了的現(xiàn)象得到了比較滿意的解釋.6.2glu334在氫鍵上的作用機(jī)理乙酰膽堿酯酶(簡稱AChE)對于催化神經(jīng)遞質(zhì)乙酰膽堿(ACh)的水解有著很高的效率,反應(yīng)生成膽堿ROH和乙酸.反應(yīng)路線分為兩步,?；磻?yīng)和脫?；磻?yīng).AchE的反應(yīng)核心是一個(gè)催化三聯(lián)體(Ser203,His447,Glu334)和氧負(fù)離子孔(oxyanionhole).其中Ser203氨基酸殘基起著親核攻擊的作用,而His447起著廣義酸堿催化的作用,但是Glu334氨基酸殘基的作用機(jī)理還不是太清楚,雖然氧負(fù)離子孔能通過氫鍵來穩(wěn)定高能中間體和過渡態(tài)分子,但是在整個(gè)反應(yīng)過程中這些氫鍵是否始終存在也是這個(gè)工作所要解決的問題.首先為了研究Glu334的作用機(jī)理,在劃分體系的過程中有兩種方案,如圖5所示.一種是將Glu334氨基酸殘基規(guī)劃到QM部分,一種則是將Glu334氨基酸殘基規(guī)劃到MM部分,其他部分都是一樣的.在QM/MM的計(jì)算過程中,以Ser203氧原子為中心20?為半徑的球外的原子在分子優(yōu)化的過程中是固定不動的,這樣大大