王長生， 李曉蕾， 高希嬋

(遼寧師范大學(xué) 化學(xué)化工學(xué)院，遼寧 大連 116029)

應(yīng)用可極化偶極-偶極作用模型快速計算分子間非鍵作用強度

王長生， 李曉蕾， 高希嬋

(遼寧師范大學(xué) 化學(xué)化工學(xué)院，遼寧 大連 116029)

非鍵作用廣泛存在于生物體系中，對生物分子的結(jié)構(gòu)和性能有重要影響，快速準確預(yù)測這些非鍵作用強度具有重要科學(xué)意義.為了實現(xiàn)快速計算分子間非鍵作用強度，提出了可極化偶極-偶極作用模型.為了進一步測試可極化偶極-偶極作用模型的適用范圍，將該模型應(yīng)用于更多的含有有機和生物小分子的非鍵復(fù)合物體系.計算結(jié)果表明：可極化偶極-偶極作用模型的預(yù)測結(jié)果很好地重復(fù)了高精度從頭計算計算結(jié)果，其預(yù)測能量時的計算精度至少與M06-2X和M06-2X-D3方法計算精度相當(dāng)，優(yōu)于AMBER99和AMOEBA方法.說明可極化偶極-偶極作用模型可望在生物分子體系的分子模擬領(lǐng)域得到應(yīng)用.

鍵偶極；偶極-偶極作用；極化作用；非鍵作用；分子模擬

非鍵作用廣泛存在于核酸和蛋白質(zhì)等生物大分子體系中，對生物大分子的結(jié)構(gòu)、性能以及分子識別過程有重要影響[1-4].要準確模擬這些復(fù)雜的生命過程、正確描述其中的反應(yīng)機制，首先必須準確預(yù)測生物分子體系中這些非鍵作用強度.耦合簇方法因其計算結(jié)果的準確性備受關(guān)注，其中CCSD(T)/CBS方法更被譽為計算化學(xué)的金字招牌[5].雖然耦合簇方法在計算非鍵作用強度時可以獲得高精度的計算結(jié)果，然而由于其計算過程需要消耗大量的磁盤空間和時間，一般只能適用于含有幾十個原子的分子體系.近年來，密度泛函理論方法快速發(fā)展，其中，M06-2X以及M06-2X-D3方法在計算非鍵作用時可以獲得與CCSD(T)方法符合較好的計算結(jié)果[6-7].雖然密度泛函方法相較于耦合簇方法在計算效率上要快很多，但是現(xiàn)階段只能適用于包含幾百個原子的體系，想要模擬含成千上萬個原子的生物大分子仍是不切實際的.分子力學(xué)方法能夠用于生物大分子體系的模擬，如AMBER99是一種常見的分子力場，其在計算靜電作用時采用固定點電荷模型，這使得它們在進行分子模擬時，忽略了環(huán)境因素對模擬結(jié)果的影響[8].為了克服這一限制，極化力場應(yīng)運而生并迅速發(fā)展.譬如AMOEBA是一種常用的極化力場，該力場是使用點偶極模型來估算極化作用[9].多年來，分子力場已經(jīng)有了很大進步，但是其計算精度仍有很大提升空間.因此發(fā)展一種快速準確預(yù)測生物大分子體系中非鍵作用強度的理論方法十分必要.近年來，本課題組建立了可極化偶極-偶極作用模型[10-12]，該模型可以快速準確預(yù)測生物分子體系中氫鍵、堆積、T型和X-H…π作用強度.本文將可極化偶極-偶極作用模型應(yīng)用于更多的包含有機和生物小分子的非鍵復(fù)合物體系中，進一步測試可極化偶極-偶極作用模型的適用范圍.

1 理論模型

可極化偶極-偶極作用模型[10-12]包含4項作用：靜電作用(Ees)、極化作用(Epol)、范德華作用(Evdw)以及軌道作用(Eorb).

1.1 靜電作用

對于2個電中性分子形成的非鍵復(fù)合物，本文使用偶極-偶極作用來估算復(fù)合物分子間的靜電相互作用，公式中，μi是孤立分子I中的第i個化學(xué)鍵偶極，μj是孤立分子J中的第j個化學(xué)鍵偶極，r，θ，θ′和ζ是描述2個化學(xué)鍵偶極μi和μj間相對位置的參數(shù).

1.2 極化作用

當(dāng)2個分子彼此靠近時，由于分子間的彼此極化，化學(xué)鍵偶極大小隨之變化，此時分子中的化學(xué)鍵偶極與孤立分子中的化學(xué)鍵偶極的差值稱為誘導(dǎo)偶極δμ.本文通過式

δμ=(q-q0)·d

來估算誘導(dǎo)偶極，其中，q0是孤立分子的原子分數(shù)電荷；q是復(fù)合物中對應(yīng)的原子分數(shù)電荷；d是偶極內(nèi)正負電荷間的距離.單體和復(fù)合物的原子分數(shù)電荷使用AM1方法計算得到.

1.3 范德華作用

使用經(jīng)典的Lennard-Jones 12-6函數(shù)計算2個分子之間的范德華作用.公式中

1.4 軌道作用

非鍵復(fù)合物的相互作用能大多包含電荷轉(zhuǎn)移的貢獻，即從質(zhì)子受體的非鍵軌道至質(zhì)子供體的反鍵軌道之間的電荷轉(zhuǎn)移對相互作用能的貢獻，這種電荷轉(zhuǎn)移的強度依賴于非鍵軌道與反鍵軌道間重疊積分的大小.基于這樣的理解，本文采用包含重疊積分的函數(shù)來近似處理該項作用.軌道間重疊積分的計算借助于類氫離子波函數(shù)來完成，同時采用雙曲正切函數(shù)形式來描述軌道能量與重疊積分的關(guān)系.公式中，Dm、am、m、pm和qm是依賴于非鍵作用類型的參數(shù).

通過模型分子，確定了可極化偶極-偶極作用模型所需參數(shù)，結(jié)果列于表1.

表1 本文確定的參數(shù)

2 應(yīng) 用

將本文模型和表1中的參數(shù)應(yīng)用于計算多個非鍵復(fù)合物的分子間平衡距離和相互作用能，并與高精度從頭算方法計算結(jié)果進行比較.同時，以文獻中報道的CCSD(T)/CBS方法計算結(jié)果為標準，進一步比較了本文模型、密度泛函理論M06-2X和M06-2X-D3方法、非極化分子力場AMBER99方法、極化分子力場AMOEBA方法的計算結(jié)果.通過對比以檢驗本文模型的合理性以及參數(shù)的可靠性.圖1給出了4個體系的計算結(jié)果.

圖1 非鍵復(fù)合物的結(jié)構(gòu)和相互作用能Fig.1 The structures and interaction energies of nonbonded complexes結(jié)構(gòu)：括號內(nèi)數(shù)值，MP2方法；括號外數(shù)值，本文模型；單位：nm.作用能：IECCSD(T)/CBS，CCSD(T)/CBS；IEM06-2X，M06-2X/aug-cc-pVDZ；IEM06-2X-D3，M06-2X-D3/aug-cc-pVDZ；IEAMBER99，AMBER99力場方法；IEAMOEBA，AMOEBA力場方法；IEthis work，本文模型；CP，包含基組重疊誤差校正(圖1中IE單位：kJ·mol-1)

2.1 預(yù)測分子間平衡距離

2.2 預(yù)測相互作用能

由圖1可知，本文將CCSD(T)/CBS方法計算復(fù)合物D1的相互作用能-86.32 kJ·mol-1視為標準值.本文模型計算的相互作用能為-87.53 kJ·mol-1，與標準值相差1.21 kJ·mol-1；不包含基組重疊誤差校正的M06-2X/aug-cc-pVDZ方法計算的相互作用能為-84.65 kJ·mol-1，與標準值相差1.67 kJ·mol-1；包含基組重疊誤差校正的M06-2X/aug-cc-pVDZ方法計算的相互作用能為-80.38 kJ·mol-1,與標準值相差5.94 kJ·mol-1;不包含基組重疊誤差校正的M06-2X-D3/aug-cc-pVDZ方法計算的相互作用能為-86.28 kJ·mol-1,與標準值相差0.04 kJ·mol-1；包含基組重疊誤差校正的M06-2X-D3/aug-cc-pVDZ方法計算的相互作用能為-82.01 kJ·mol-1，與標準值相差4.31 kJ·mol-1；AMBER99力場計算的相互作用能為-68.68 kJ·mol-1,與標準值相差17.64 kJ·mol-1；AMOEBA力場計算的相互作用能為-75.20 kJ·mol-1,與標準值相差11.12 kJ·mol-1.再如，使用CCSD(T)/CBS方法計算復(fù)合物D4的標準相互作用能為-6.94 kJ·mol-1.本文模型計算的相互作用能同樣為-6.94 kJ·mol-1；使用M06-2X/aug-cc-pVDZ方法(包含/不包含基組重疊誤差校正)計算的相互作用能分別為-6.10和-8.15 kJ·mol-1,與標準值分別相差0.84和1.21 kJ·mol-1；使用M06-2X-D3/aug-cc-pVDZ方法(包含/不包含基組重疊誤差校正)計算的相互作用能分別為-7.40和-9.49 kJ·mol-1，與標準值分別相差0.46和2.55 kJ·mol-1；AMBER99力場計算的相互作用能為-6.56 kJ·mol-1，與標準值相差0.38 kJ·mol-1；AMOEBA力場計算的相互作用能為-6.23 kJ·mol-1，與標準值相差0.71 kJ·mol-1.以上計算結(jié)果表明:本文模型可以準確計算氫鍵、堆積、T型和X—H…π復(fù)合物的相互作用能，計算精度至少與M06-2X和M06-2X-D3方法計算精度相當(dāng)，優(yōu)于AMBER99和AMOEBA方法.

3 結(jié) 論

可極化偶極-偶極作用模型可以較好地描述有機和生物小分子體系中的氫鍵作用、堆積作用、T型作用和X—H…π作用，模型的計算精度可以達到M06-2X和M06-2X-D3方法計算精度或者更好，明顯優(yōu)于AMBER99和AMOEBA方法的計算精度.

[1] GUCKIAN K M,SCHWEITZER B A,REN R X F,et al.Experimental measurement of aromatic stacking affinities in the context of duplex DNA[J].J Am Chem Soc,1996,118(34):8182-8183.

[2] DILL K A,MACCALLUM J L.The protein-folding problem,50 years on[J].Science,2012,338(6110):1042-1046.

[3] 王長生,李蕾.精氨酸側(cè)鏈與生物小分子間離子氫鍵強度的快速計算[J].遼寧師范大學(xué)學(xué)報(自然科學(xué)版),2017,40(1):47-51.

[4] 黃翠英,王長生.準確快速計算含多肽酰胺和核酸堿基的氫鍵復(fù)合物的氫鍵強度和氫鍵作用勢能曲線[J].化學(xué)進展,2012,24(6)1214-1226.

[5] HALDAR S,GNANASEKARAN R,HOBZA P.A comparison of ab initio quantum-mechanical and experimentalD0binding energies of eleven H-bonded and eleven dispersion-bound complexes[J].Phys Chem Chem Phys,2015,17(40):26645-26652.

[6] ZHAO Y,TRUHLAR D G.Density functional with broad applicability in chemistry[J].Acc Chem Res,2008，41(2):157-167.

[7] GOERIGK L,KRUSE H,GRIMME S.Benchmarking density functional methods against the S66 and S66×8 datasets for non-covalent interactions[J].ChemPhysChem,2011,12(17):3421-3433.

[8] PATON R S，GOODMAN J M.Hydrogen bonding andπ-stacking:how reliable are force fields.a critical evaluation of force field descriptions of nonbonded interactions[J].J Chem Inf Model，2009，49(4):944-955.

[9] SHI Y,XIA Z,ZHANG J,et al.Polarizable atomic multipole-based AMOEBA force field for proteins[J].J Chem Theory Comput,2013,9(9):4046-4063.

[10] LI S S,HUANG C Y,HAO J J,et al.A polarizable dipole-dipole interaction model for evaluation of the interaction energies for N—H…O=C and C—H…O=C hydrogen-bonded complexes[J].J Comput Chem,2014，35(6):415-426.

[11] 王長生,黃翠英.核酸堿基間氫鍵距離和作用能的快速預(yù)測[J].遼寧師范大學(xué)(自然科學(xué)版),2016,39(3):368-372.

[12] GAO X C,HAO Q,WANG C S.Improved polarizable dipole-dipole interaction model for hydrogen bonding，stacking，T-shaped and X—H…πinteractions[J].J Chem Theory Comput,2017,13(6):2730-2741.

Rapidcalculationoftheintermolecularnonbondedinteractionstrengthbyusingthepolarizabledipole-dipoleinteractionmodel

WANGChangsheng，LIXiaolei，GAOXichan

(School of Chemistry and Chemical Engineering， Liaoning Normal University， Dalian 116029， China)

Nonbonded interactions are known to play key roles in determining the structures，properties，and functions of biomolecules.The ability to accurately and efficiently simulate these nonbonded interactions is critical for the computer simulation of biological processes.The polarizable dipole-dipole interaction model has been formulated in our laboratory to rapidly and accurately estimate the nonbonded interaction strength in biosystems.In this paper，the polarizable dipole-dipole interaction model is applied to systems containing organic and small biomolecules to further check the model’s utility.The results prove that the polarizable dipole-dipole interaction model reproduced the equilibrium intermolecular distances and interaction energies in good agreement with those yielded by the high qualityabinitiomethods.The calculation precision for interaction energies is at least equal to M06-2X and M06-2X-D3 methods，and outperforms the AMBER99 and AMOEBA methods，demonstrating the model may serve as a new tool for biological process simulations.

bond dipole；dipole-dipole interaction；polarization interaction；nonbonded interaction；molecular modeling

O641.121

A

2017-08-15

國家自然科學(xué)基金資助項目(21573098)

王長生(1963- )，男，遼寧大連人，遼寧師范大學(xué)教授，博士，博士生導(dǎo)師.

1000-1735(2017)04-0479-05

10.11679/lsxblk2017040479