应用可极化偶极-偶极作用模型快速计算分子间非键作用强度

王长生， 李晓蕾， 高希婵

(辽宁师范大学 化学化工学院，辽宁 大连 116029)

应用可极化偶极-偶极作用模型快速计算分子间非键作用强度

王长生， 李晓蕾， 高希婵

(辽宁师范大学 化学化工学院，辽宁 大连 116029)

非键作用广泛存在于生物体系中，对生物分子的结构和性能有重要影响，快速准确预测这些非键作用强度具有重要科学意义.为了实现快速计算分子间非键作用强度，提出了可极化偶极-偶极作用模型.为了进一步测试可极化偶极-偶极作用模型的适用范围，将该模型应用于更多的含有有机和生物小分子的非键复合物体系.计算结果表明：可极化偶极-偶极作用模型的预测结果很好地重复了高精度从头计算计算结果，其预测能量时的计算精度至少与M06-2X和M06-2X-D3方法计算精度相当，优于AMBER99和AMOEBA方法.说明可极化偶极-偶极作用模型可望在生物分子体系的分子模拟领域得到应用.

键偶极；偶极-偶极作用；极化作用；非键作用；分子模拟

非键作用广泛存在于核酸和蛋白质等生物大分子体系中，对生物大分子的结构、性能以及分子识别过程有重要影响[1-4].要准确模拟这些复杂的生命过程、正确描述其中的反应机制，首先必须准确预测生物分子体系中这些非键作用强度.耦合簇方法因其计算结果的准确性备受关注，其中CCSD(T)/CBS方法更被誉为计算化学的金字招牌[5].虽然耦合簇方法在计算非键作用强度时可以获得高精度的计算结果，然而由于其计算过程需要消耗大量的磁盘空间和时间，一般只能适用于含有几十个原子的分子体系.近年来，密度泛函理论方法快速发展，其中，M06-2X以及M06-2X-D3方法在计算非键作用时可以获得与CCSD(T)方法符合较好的计算结果[6-7].虽然密度泛函方法相较于耦合簇方法在计算效率上要快很多，但是现阶段只能适用于包含几百个原子的体系，想要模拟含成千上万个原子的生物大分子仍是不切实际的.分子力学方法能够用于生物大分子体系的模拟，如AMBER99是一种常见的分子力场，其在计算静电作用时采用固定点电荷模型，这使得它们在进行分子模拟时，忽略了环境因素对模拟结果的影响[8].为了克服这一限制，极化力场应运而生并迅速发展.譬如AMOEBA是一种常用的极化力场，该力场是使用点偶极模型来估算极化作用[9].多年来，分子力场已经有了很大进步，但是其计算精度仍有很大提升空间.因此发展一种快速准确预测生物大分子体系中非键作用强度的理论方法十分必要.近年来，本课题组建立了可极化偶极-偶极作用模型[10-12]，该模型可以快速准确预测生物分子体系中氢键、堆积、T型和X-H…π作用强度.本文将可极化偶极-偶极作用模型应用于更多的包含有机和生物小分子的非键复合物体系中，进一步测试可极化偶极-偶极作用模型的适用范围.

1 理论模型

可极化偶极-偶极作用模型[10-12]包含4项作用：静电作用(Ees)、极化作用(Epol)、范德华作用(Evdw)以及轨道作用(Eorb).

1.1 静电作用

对于2个电中性分子形成的非键复合物，本文使用偶极-偶极作用来估算复合物分子间的静电相互作用，公式中，μi是孤立分子I中的第i个化学键偶极，μj是孤立分子J中的第j个化学键偶极，r，θ，θ′和ζ是描述2个化学键偶极μi和μj间相对位置的参数.

1.2 极化作用

当2个分子彼此靠近时，由于分子间的彼此极化，化学键偶极大小随之变化，此时分子中的化学键偶极与孤立分子中的化学键偶极的差值称为诱导偶极δμ.本文通过式

δμ=(q-q0)·d

来估算诱导偶极，其中，q0是孤立分子的原子分数电荷；q是复合物中对应的原子分数电荷；d是偶极内正负电荷间的距离.单体和复合物的原子分数电荷使用AM1方法计算得到.

1.3 范德华作用

使用经典的Lennard-Jones 12-6函数计算2个分子之间的范德华作用.公式中

1.4 轨道作用

非键复合物的相互作用能大多包含电荷转移的贡献，即从质子受体的非键轨道至质子供体的反键轨道之间的电荷转移对相互作用能的贡献，这种电荷转移的强度依赖于非键轨道与反键轨道间重叠积分的大小.基于这样的理解，本文采用包含重叠积分的函数来近似处理该项作用.轨道间重叠积分的计算借助于类氢离子波函数来完成，同时采用双曲正切函数形式来描述轨道能量与重叠积分的关系.公式中，Dm、am、m、pm和qm是依赖于非键作用类型的参数.

通过模型分子，确定了可极化偶极-偶极作用模型所需参数，结果列于表1.

表1 本文确定的参数

2 应 用

将本文模型和表1中的参数应用于计算多个非键复合物的分子间平衡距离和相互作用能，并与高精度从头算方法计算结果进行比较.同时，以文献中报道的CCSD(T)/CBS方法计算结果为标准，进一步比较了本文模型、密度泛函理论M06-2X和M06-2X-D3方法、非极化分子力场AMBER99方法、极化分子力场AMOEBA方法的计算结果.通过对比以检验本文模型的合理性以及参数的可靠性.图1给出了4个体系的计算结果.

图1 非键复合物的结构和相互作用能Fig.1 The structures and interaction energies of nonbonded complexes结构：括号内数值，MP2方法；括号外数值，本文模型；单位：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 预测分子间平衡距离

2.2 预测相互作用能

由图1可知，本文将CCSD(T)/CBS方法计算复合物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方法计算复合物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.以上计算结果表明:本文模型可以准确计算氢键、堆积、T型和X—H…π复合物的相互作用能，计算精度至少与M06-2X和M06-2X-D3方法计算精度相当，优于AMBER99和AMOEBA方法.

3 结 论

可极化偶极-偶极作用模型可以较好地描述有机和生物小分子体系中的氢键作用、堆积作用、T型作用和X—H…π作用，模型的计算精度可以达到M06-2X和M06-2X-D3方法计算精度或者更好，明显优于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] 王长生,李蕾.精氨酸侧链与生物小分子间离子氢键强度的快速计算[J].辽宁师范大学学报(自然科学版),2017,40(1):47-51.

[4] 黄翠英,王长生.准确快速计算含多肽酰胺和核酸碱基的氢键复合物的氢键强度和氢键作用势能曲线[J].化学进展,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] 王长生,黄翠英.核酸碱基间氢键距离和作用能的快速预测[J].辽宁师范大学(自然科学版),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

国家自然科学基金资助项目(21573098)

王长生(1963- )，男，辽宁大连人，辽宁师范大学教授，博士，博士生导师.

1000-1735(2017)04-0479-05

10.11679/lsxblk2017040479