新型多苯并咪唑环芳识别氨基酸甲酯的分子动力学和密度泛函研究

娄永锋，肖 蓉，向明礼，杨明理

(1.四川大学 a.华西医院生物治疗国家重点实验室;b.华西医院纳米生物医学技术与膜生物学研究所，成都 610041;2.四川大学 化工学院，成都 610065)

环芳是一类含有芳香环的环状有机化合物［1］。环芳由于可能存在某种隙间共振效应(interstitial resonance effect)引起的芳香性改变以及在不同环中分别引入吸、供电子基团后将形成分子内配合物［2］等诱人特性而备受关注。自从1949年被首次发现［3］以来，环芳被广泛应用于手性识别［4－6］、超分子［7］、分子催化［8］、抗菌［9］等领域。咪唑环芳更为特殊:1)咪唑既能作为质子受体又能作为质子供体，具有特殊亲核、亲电性［6］，是生物体系的重要基团［10］;2)含咪唑环的一些化合物具有抗炎活性［11］，并对细胞色素P450酶有抑制活性;3)咪唑衍生物的合成颇具挑战性。由于上述诸多因素［12］，使得咪唑环芳的研究更为瞩目。环芳的研究工作主要集中在合成品种不同、功能各异的环芳以及发展新的合成方法［6，13］方面，对于环芳识别氨基酸甲酯的分子作用机制的认识［6］，还有待深入。本文采用分子动力学与密度泛函理论相结合的研究方法，对笔者最新合成的一个多苯并咪唑环芳受体(A)［6］与一系列氨基酸甲酯配体(L)间的相互作用进行了探讨。采用以密度泛函理论计算结果为核心参数的组合策略，来表达该多苯并咪唑环芳对一系列氨基酸甲酯的识别，获得了与实验数据一致的研究结果。

图1 多苯并咪唑环芳结构

1 计算方法

如图1所示，该多苯并咪唑环芳(MBICP)受体由2条醚链桥接2 个2，2’-二苯并咪唑构成。

MBICP所识别的5个氨基酸甲酯(AAME)配体分别为甘氨酸甲酯(LG)、丙氨酸甲酯(LA)、脯氨酸甲酯(LP)、亮氨酸甲酯(LL)和缬氨酸甲酯(LV)。因尚无该MBICP受体及其与AAME形成的复合物的晶体结构，故而对于该受体和这几个配体的单体三维结构，笔者采用ChemBio3D构建并用Gaussian 03［14］在B3LYP/6-31G(d)理论水平上进行优化得到。

在选择性的分子络合中，受、配体结合位点间的互补性和复合物中配体的空间位阻是影响选择性络合的关键因素，而静电作用、范德法力、氢键等也是形成分子复合物的推动力［15］。因此，在没有晶体结构的情况下，如何获得考虑了上述因素的复合物构象，对于用量子化学研究MBICP受体与AAME配体间的作用和识别机制，是很关键的。

本文采用基于模拟退火的分子动力学程序CDOCKER［16］，将AAME配体与MBICP受体对接。对接时，选用CHARMM力场，把受体的几何中心默认为结合位点，球半径参数设置为5 A。。每个AAME配体与受体对接，都产生1 000个随机构象。对所产生的随机构象进行逐渐升温和模拟退火动力学模拟。对接和动力学模拟过程中的其余参数采用默认值。

要成功预测弱相互作用体系的结构、能量性质需要精确的计算。鉴于此，笔者把动力学模拟后得到的每种复合物的1 000构象，按能量高低顺序排序，仅从中选出10个能量值为负值的最低能量构象。利用Turbomole［17］程序在BLYP/def-TZVP［18］水平进行优化。再将优化收敛后得到的最低能量构象在相同理论和基组下作NBO分析。对于之前结构优化后得到的MBICP受体和AAME配体，也在该条件下进行单点能计算，以便获得受体与配体间的相互作用能。

2 结果与讨论

2.1 多苯并咪唑环芳对氨基酸甲酯识别的复杂性

分子识别过程中所形成的复合物的最高占据轨道能级EHOMO、最低空轨道能级EHOMO、形成复合物过程中受、配体间所转移的电荷量CT都是衡量复合物稳定性的重要指标。研究表明，结构相似的分子，最高占有轨道与最低空轨道的能量差越大，则分子越稳定;受体与配体间电荷转移量越大，则配合物越稳定。

相互作用能ΔEint则更是衡量复合物稳定性的重要参数。ΔEint是指把单个受体(A)和单个配体(L)结合在一起，形成复合物的能量［19］，其定义如下:

MBICP-AAME复合物结合常数的实验测定值［6］和相关参数计算值如表1所示，各项指标EHOMO、ELUMO、CT和ΔEint与实验测定的复合物结合常数之间没有简单的对应关系。甚至在EHOMO、ELUMO、CT与ΔEint之间，也没观察到任何关联。这清楚地表明，多苯并咪唑环芳对氨基酸甲酯的识别情况比较复杂，用单一的量子化学计算参数，很难表达和认识其间的作用机制。

2.2 组合策略及其对结合常数的计算

量子化学精确计算是深入理解主客体识别形成复合物作用机制的重要途径。但不能忽视的是:实验数据［6］是在氯仿溶剂环境下获得的。另外，预组织、疏水相互作用以及形成复合物的单体间的空间位阻和结合位点间的互补置放等，都是影响分子识别的重要因素［15］。目前量子化学尚未能对溶剂效应、疏水相互作用等进行精确的计算。为此，笔者提出了以密度泛函理论计算的相互作用能ΔEint(见表1)为核心参数的组合策略，以期能对实验测定的结合常数进行准确计算。

表1 MBICP-AAME复合物结合常数K［6］的实验值、计算值及其相关参数

显然，MBICP受体与5个AAME配体所形成的复合物的稳定性的差异源于各配体之间的差异。因为预组织和疏水相互作用等因素对分子间选择性识别有重要贡献，所以，如何体现各配体对预组织和疏水相互作用贡献的差异就显得很关键。为此，用ChemBio3D计算了优化后的AAME配体的脂水分配系数P和溶剂可及表面积SASA(计算数据见表1)。P可以作为表达各配体对疏水相互作用贡献的描述符。由预组织原理［15］可知，在形成复合物前，受体和配体间组织得愈好，溶剂化程度愈低，则复合物愈稳定。而根据非静电相互作用的连续模型，溶剂化自由能与SASA间有很好的线性关系。因此，SASA可以用来表达AAME配体的溶剂化作用。

由于影响MBICP受体与AAME配体间的相互识别的因素很多，个别描述符难以恰当地表述这种受、配体间的相互作用并与实验测定的结合常数定量关联，但也不能为了建立定量关系模型而无限制地加入描述符。数据的相关性分析表明，AAME配体的SASA与复合物的结合常数K间的相关度较高，而其P与结合常数间的相关性则差一些。这说明，对该MBICP受体而言，复合物的稳定性，更多地来自于AAME配体对预组织的贡献。

因此，选定ΔEint作为核心参数，组合AAME配体的SASA，构建了MBICP识别AAME形成复合物结合常数的定量计算模型为:

其中，α、β和γ为待定参数。用表1中的数据，对(2)式进行多项式回归最小二乘拟合求解。则得:－lgKcalc=－2.464 ×10－3× ΔEint+2.356 ×10－3×SASA– 3.664 ，计算值(见表1)与实验值［6］之间有良好的相关性，相关系数达0.961。

2.3 识别机制的认识

2.3.1 预组织效应等因素的贡献

从表1可知，配体LG的SASA最小。根据非静电相互作用的连续模型和预组织原理［15］，该配体的溶剂化自由能最小，溶剂化程度最低，对复合物的稳定性贡献最大。而配体LP和LL的SASA较大，溶剂化自由能较大，溶剂化程度较高，因而复合物的稳定性较差。这与实验结果是一致的。

用Turbomole在BLYP/def-TZVP水平计算的各MBICP-AAME复合物的最低能量构象，如图2所示。各配体与MBICP受体间的相互作用模式不尽相同。配体LV、LP、LL、LA，其氮端上的H与受体醚链上的O形成了氢键相互作用。在5个配体中，LG的体积最小。也只有LG在与MBICP形成复合物时，配体骨架与受体的苯并咪唑环平面平行，配体的酯甲基、羰基部分伸向溶液，而将α-氨基和α-碳嵌入苯并咪唑环平面之间。这是一种很恰当的空间互补作用模式。配体受体间有最大的范德法作用接触面，而与溶剂作用的接触面最小，溶剂化作用最小，因而A-LG复合物最稳定。

图2 多苯并咪唑环芳受体(A)与氨基酸甲酯配体(L)间相互作用的作用模式

2.3.2 电荷转移的影响

电荷转移的情况比较复杂。为使问题简化，仅选取5个复合物的最低能量构象进行讨论。通常，电荷转移是有利于分子识别的，而且电荷转移量与复合物的稳定性是正相关的。但表1的电荷转移数据显示，在形成复合物的过程中有少量电子从受体MBICP向配体LV、LP、LL转移;而对于配体LA和LG，则是少量电子向受体MBICP转移。

仔细探究各配体在形成复合物前后分子内部的电荷变化情况时，发现了一个更加趣的现象。在LG形成复合物的过程中，其氮端氨基上带负电的氮原子，负电荷增加了 0.017 2，电荷更负;带正电的氢原子，正电荷也增加了0.008 3，电荷更正(见表2)。而在碳端，带负电的羰基氧和带负电的甲基碳，都在形成复合物后减少了负电荷。从BLYP/def-TZVP水平优化收敛得到的最低能量构象看(图2 A-LG)，氨基是与受体最为靠近的基团。在这个基团上，电荷分离加强，极化增大，有利于配体与受体间的静电相互作用，从而对体系的稳定做出贡献。同时，在伸向溶剂的碳端，电荷分离被弱化，极性减小，改善了疏水性，更有利于复合物的稳定性。因此，A-LG复合物在5个复合物中最为稳定。

表2 不同配体分子中特定原子上NBO电荷的变化a

同样，在LA中也观察到氮端氨基的电荷分离加强的现象。带负电的氨基上的氮原子，负电荷增加了0.015 3，电荷更负;氨基上带正电的氢原子，正电荷也增加了0.029 1，电荷更正。配体与受体间的氢键作用得到加强，有利于复合物的稳定。伸向溶液的碳端羰基氧和甲基碳上的负电荷略有增加，碳端极性增加，疏水效应减弱。但因碳端电荷转移的量较之于氮端小了一个数量级，因而在总体上表现为电荷转移有利于复合物的稳定性。A-LA的电荷转移量小于A-LG，因而稳定性不如A-LG。这与实验结果是一致的。

总体而言，在其他3个配体中电荷从MBICP受体向配体转移，使得这3个配体中伸向溶液的碳端上带负电的羰基氧和带负电的甲基碳，负电荷显著增加。其结果是碳端极性显著增大，导致疏水作用减弱，从而对复合物的稳定性产生负贡献。虽然氮端氨基上的氮原子所带负电荷和氢原子所带正电荷也有增加的趋势，但增加的量远不如在LA和LG中那样(见表2)。故而它们对复合物稳定性的贡献相对较小。因此，这3个配体与MBICP受体形成的复合物稳定性都要差一些。这也同样与实验相符。

就个别而言，在这3个配体中，配体LV的氮端电荷增加最多，且碳端极化最小。所以，A-LV也是这3个相关复合物中最稳定的。

3 结语

为研究新型多苯并咪唑环芳(MBICP)受体对一系列氨基酸甲酯 (AAME)配体的识别作用，用ChemBio3D构建了受体和配体分子单体，并用Gaussian 03在B3LYP/6-31G(d)水平上优化结构，得到了受体和各配体的三维结构。用基于模拟退火的分子动力学程序CDOCKER，将结构优化后的AAME配体与MBICP受体对接。对于每一个复合物，从随机产生的并经过升温和退火动力学模拟的1 000个构象中，筛选出10个能量最低的复合物构象，作为复合物的初始构象。用Turbomole在BLYP/def-TZVP水平对每个复合物的10个合理初始构象进行优化。将收敛后得到的最低能量构象作为复合物的结构，并在相同理论和基组下作NBO分析。对此前得到的MBICP受体和AAME配体的优化结构，在该条件下进行单点能计算，获得了受体与配体间的相互作用能。

常用于表达和解释分子识别能力的量子化学计算参量如相互作用能、受体与配体间的电荷转移等，在本体系中都与实验测定的MBICP-AAME复合物的结合常数间没有简单关联。各量子化学计算参数间也没有关联性。这表明MBICP和AAME间的相互识别与作用机制是相当复杂的。

采用组合策略，以密度泛函理论计算所获得的受、配体间相互作用能ΔEint作为核心参数，引入能表达预组织效应的SASA参数，建立了表达MBICP-AAME复合物结合常数的定量关系模型。基于该模型计算的复合物结合常数，能较好地与实验测定值吻合。表明了预组织效应对该MBICP受体与AAME配体间的相互作用，具有重要贡献。

受体、配体间的电荷转移，似乎无法解释MBICP对AAME的识别能力。但用Turbomole程序在BLYP/def-TZVP水平分析形成复合物前后AAME氮端原子和碳端原子的NBO电荷变化，发现形成复合物前后有电荷分离现象。利用该分析结果，能为MBICP对AAME的识别能力的强弱，提供合理的解释。MBICP受体与AAME配体间的作用模式分析表明，配体体积的大小也是影响复合物稳定性的重要因素。

［1］CRAM D J，CRAM J M.Cyclophane chemistry:bent and battered benzene rings［J］.Acc Chem Res，1971，4(6):204 －213.

［2］CRAM D J，STEINBERG H.Macro Rings.I.Preparation and Spectra of the Paracyclophanes［J］.J Am Chem Soc，1951，73(12):5691 －5704.

［3］BROWN C J，FARTHING A C.Preparation and Structure of Di-p-Xylylene［J］.Nature，1949，164:915 － 916.

［4］PU L.Enantioselective Fluorescent Sensors:A Tale of BINOL［J］.Acc Chem Res，2011，:DOI:10.1021/ar200048d.

［5］YANG L，QIN S，SU X，et al.1，1’-Binaphthyl-based imidazolium chemosensors for highly selective recognition of tryptophan in aqueous solutions［J］.Org Biomol Chem，2010，8(2):339 －348.

［6］XU X，WANG X，WU A，et al.Synthesis and molecular recognition of novel multiimidazole cyclophanes［J］.J Heterocyclic Chem，2009，46(6):1137－1141.

［7］CAMACHO D H，SALO E V，ZILLER J W，et al.Cyclophane-based highly active late-transition-metal catalysts for ethylene polymerization［J］.Angew Chem Int Ed，2004，43:1821 －1825.

［8］POYATOS M，MATA J A，PERIS E.Complexes with Poly(N-heterocyclic carbene)Ligands:Structural Features and Catalytic Applications［J］.Chem Rev，2009，109(8):3677 －3707.

［9］MELAIYE A，SUN Z，HINDI K，et al.Silver(I)-midazole Cyclophane gem-Diol Complexes Encapsulated by Electrospun Tecophilic Nanofibers:Formation of Nanosilver Particles and Antimicrobial Activity［J］.J Am Chem Soc，2005，127(7):2285 － 2291.

［10］YUAN Y，GAO G，JIANG Z，et al.Synthesis and selective anion recognition of imidazolium cyclophanes［J］.Tetrahedron，2002，58(44):8993 －8999.

［11］HERN N P，CABRERA M，LAVAGGI M L，et al.Discovery of new orally effective analgesic and anti-inflammatory hybrid furoxanyl N-acylhydrazone derivatives［J］.Bioorg Med Chem，2012，20(6):2158 －2171.

［12］GARRISON J C，SIMONS R S，TALLEY J M，et al.Synthesis and structural characterization of an imidazolium-linked cyclophane and the silver complex of an N-heterocyclic carbene-linked cyclophane［J］.Organometallics，2001，20(7):1276 － 1278.

［13］GABUTTI S，SCHAFFNER S，NEUBURGER M，et al.Planar chiral asymmetric naphthalenediimide cyclophanes:synthesis，characterization and tunable FRET properties［J］.Org Biomol Chem，2009(7):3222 －3229.

［14］FRISCH M J，TRUCKS G W，SCHLEGE H B，et al.Gaussian 03，Revision C.02［Z］.Wallingford CT:Gaussian Inc，2004.

［15］CRAM D J.The design of molecular hosts，guests，and their complexes［J］.J Incl Phenom，1988，6(4):397 －413.

［16］WU G，ROBERTSON D，BROOKS C，et al.Detailed analysis of grid-based molecular docking:A case study of CDOCKER-A CHARMm-based MD docking algorithm［J］.J Comput Chem，2003，24(13):1549 －1562.

［17］TREUTLER O，AHLRICHS R.Efficient molecular numerical integration schemes［J］.J Chem Phys，1995，102:346 －354.

［18］EICHKORN K，TREUTLER O，OEHM H，et al.Auxiliary basis sets to approximate Coulomb potentials［J］.Chem Phys Lett，1995，240(4):283 －289.

［19］MO Y，SUBRAMANIAN G，GAO J，et al.Cation-π Interactions:An Energy Decomposition Analysis and Its Implication in δ-Opioid Receptor-Ligand Binding［J］.J Am Chem Soc，2002，124:4832 － 4837.