量子化学计算解析手性共价有机框架材料6色谱固定相的手性拆分机理

胡园园, 张忠杰, 黄 露

(闽江学院海洋学院化学工程系, 福建 福州 350108)

共价有机框架材料(COFs)是通过共价键和有机结构单元构建的新型多孔聚合物材料。目前,COFs已被广泛应用于能量存储[1,2]、光电器件[3,4]、传感器[5,6]、催化[7,8]、样品前处理[9,10]等领域。进入21世纪以来,手性药物发展迅速,手性研究已经成为化学、材料、生物等学科的研究热点,手性COFs也因此应运而生[11]。作为一种新型的手性分离材料,手性COFs已被成功用于液相色谱[12]、气相色谱[13]、毛细管电色谱[14]拆分手性对映体。然而,由于手性COFs材料出现时间较短,尚未有关于其手性识别机理的报道。

自色谱法被用于拆分手性对映体以来,科研工作者就一直在尝试用各种方法解释手性拆分机理。除了通过实验进行手性拆分机理研究[15,16],还可以通过化学计算从分子水平上研究分子之间的弱相互作用,从而揭示手性固定相识别对映体的机理。以环糊精为例,人们采用分子力学[17]、分子动力学[18]和蒙特卡洛[19]等多种分子模拟方法研究其与客体的包合作用及手性识别机理。与上述分子模拟方法相比,量子化学计算是一种精确度更高的计算方法,已被成功用于研究杯状芳烃[20]、手性金纳米通道[21]、手性离子液体[22]、手性CdSe纳米片[23]、手性CdTe量子点[24]等与手性对映体之间的相互作用。目前,最常用的量子化学计算程序是Guassian,但由于分子间弱相互作用体系通常都比较大,Guassian计算非常耗时,而ORCA程序使用密度拟合近似(resolution-of-the-identity, RI)进行加速后,其计算速度明显优于Guassian,非常适合计算大体系弱相互作用。

手性共价有机框架材料6(chiral covalent organic frameworks 6, CCOF6)是一种新型手性COFs,已被成功用于液相色谱拆分手性对映体[12]。本文使用AutoDock程序进行分子对接模拟,获得CCOF6与4对对映体(即1-苯基-2-丙醇、1-苯基-1-戊醇、1-苯基-1-丙醇、1-(4-溴苯基)乙醇)相互作用的初始构型。采用ORCA计算相互作用体系的结合自由能,阐明对映体的结合自由能差与其选择性因子之间的关系。采用Multiwfn程序中的独立梯度模型(independent gradient model, IGM)对ORCA计算获得的最终构型进行分析,使用视觉分子动力学(visual molecular dynamics, VMD)程序将体系弱相互作用图形化,并对手性拆分机理进行深度探讨。

图 1 CCOF6与4种外消旋化合物的化学结构Fig. 1 Chemical structures of CCOF6 and four racematesThe chemical structure of CCOF6 was drawn according to Ref. [12].

1 原理和方法

1.1 分子模型的建立

采用ChemOffice 2016程序构建CCOF6及4种外消旋化合物(化学结构式见图1)的三维分子模型。为了简化计算,选取CCOF6基本结构单元中易发生相互作用的手性结构部分(图1虚线框内结构)作为CCOF6的简化分子结构。

1.2 分子对接

将CCOF6作为受体分子,另外4对手性对映体作为配体分子,采用AutoDock 4.2.6[25]分别对上述受体分子及配体分子进行对接。采用拉马克遗传算法结合基于网格的能量评估方法,根据受体、配体分子中所有原子类型的原子间相互作用势,计算格点盒子(盒子大小为1.5 nm×2.25 nm×1.5 nm,格点间距为0.037 5 nm)。共进行100次独立对接,其他参数使用默认值。

1.3 结合自由能计算

使用ORCA程序[26]对CCOF6及4对手性对映体进行结构优化以及能量计算。采用B3LYP泛函[27],带DFT-D3[28]校正,轨道基组为def2-TZVP[29], def2/J[30]作为RI-J的辅助基组,RIJCOSX用来加速计算。利用该方法分别计算CCOF6、手性对映体(CE)及结合物的自由能,可以得到CCOF与对映体的结合自由能(ΔG),即:

ΔG=GCCOF-CE-(GCCOF+GCE)

(1)

进一步求得对映体之间的结合自由能差(|ΔΔG|),即:

|ΔΔG|=|ΔGR-ΔGS|

(2)

另外,根据半经验结合自由能评价函数公式:

ΔG=-RTlnK

(3)

以及选择性因子与分配系数的关系(假定KR>KS):

(4)

可知对映体之间的结合自由能差(|ΔΔG|)与α的关系为:

|ΔΔG|=RTlnα

(5)

1.4 独立梯度模型分析

2017年,Lefebvre等[31]提出了一种叫作IGM的方法,重在展现弱相互作用区域及其特征。用Multiwfn进行IGM分析[32],函数表达式如下:

(6)

(7)

δg(r)=gIGM(r)-g(r)

(8)

其中δg是三维实空间函数,i是原子序号,▽ρ是梯度矢量,abs(▽ρ)表示▽ρ矢量的每个分量都取绝对值,||表示对矢量取模。

2 结果与讨论

2.1 分子对接

首先,按照1.1节建立CCOF6和对映体的分子模型,使用ORCA在B3LYP-D3/def2-TZVP计算水平上对上述分子结构进行优化。其次,按照1.2节进行分子对接。如图1所示,4种手性化合物均含有一个羟基基团,此处为手性化合物与CCOF6发生相互作用的关键部位;而CCOF6的手性结构单元含有4个亚氨基、4个羰基、两个羟基和两个醚键,即其与手性化合物可能发生相互作用的部位有4处。结合分子对接结果,选取CCOF6在上述4处与对映体发生对接的构型中能量最低的构型用于下一步的结合自由能计算。

表 1 CCOF6与4对对映体的结合自由能差、选择性因子及分离度

图 2 CCOF6与4种外消旋化合物的相互作用构型Fig. 2 Interaction configurations of CCOF6 with four racemates

按照1.3节,计算2.1节中获得的相互作用构型的结合自由能,比较结合自由能的大小,最终确定CCOF6与4对对映体的相互作用构型并用Avogadro程序进行查看(见图2)。从图2中可以看出,除S-1-苯基-1-丙醇是以羟基和CCOF6的醚键发生相互作用,其他对映体皆为羟基与CCOF6的羰基发生相互作用。表1列出了CCOF6与4对对映体的结合自由能及结合自由能差,计算方式有3种:AutoDock; ORCA(考虑溶剂化效应); ORCA(不考虑溶剂化效应)。根据文献[12],对映体拆分实验是在正己烷/异丙醇(99∶1, v/v)的流动相中进行的,因此在计算时加入正己烷作为溶剂。

根据文献[12],对映体的选择性因子由大到小依次为:1-苯基-1-丙醇>1-苯基-1-戊醇>1-苯基-2-丙醇>1-(4-溴苯基)乙醇。AutoDock和不考虑溶剂化效应ORCA的结合自由能差大小次序与选择性因子大小次序不相符;而考虑溶剂化效应ORCA的结合自由能差大小次序与选择性因子大小次序完全一致。由公式(5)可知,结合自由能差应与选择性因子的ln值呈线性关系。将考虑溶剂化效应ORCA的结合自由能差对lnα作图,获得的线性方程为|ΔΔG|=158.32lnα-30.82,其相关系数(r2)为0.94。这进一步证明了考虑溶剂效应的ORCA计算方法的精确性。由于计算时只选取了CCOF6的手性结构单元作为受体分子,因此拟合曲线的斜率和截距与真实值之间存在较大差距。另外,由于分离度Rs除了与选择性因子有关,还与柱效和保留因子有关,因此结合自由能差与Rs不一定呈现比例关系。

图 3 IGM等值面图及散点图Fig. 3 Images of independent gradient model (IGM) isosurfaces and scatter plots

2.3 弱相互作用图形化

图 4 CCOF6与对映体相互作用的IGM图Fig. 4 IGM graphs showing interactions of CCOF6 with four racemates

目前,最常用的弱相互作用图形化分析法是约化密度梯度(reduced density gradient, RDG),而IGM是受到RDG启发建立的。与RDG相比,IGM可以直接考察分子间相互作用而不被分子内相互作用干扰,且IGM给出的等值面图比RDG更为饱满。一般方式计算的初分子密度的梯度等于各个原子密度梯度的加和,IGM型密度梯度把每个原子的密度梯度先取绝对值再加和,二者求差,得到所谓的δg函数。δg可以划分为用于展现片段内相互作用的δg_intra和展现片段间相互作用的δg_inter。采用Multiwfn绘制δg_intra或δg_inter对sign(λ2)ρ的散点图来考察体系中的相互作用,其中sign(λ2)ρ为电子密度Hessian矩阵的第二大本征值λ2的符号。图3给出了(±)-1-苯基-1-丙醇的IGM等值面图和散点图,其中,右侧黑色散点图为分子内相互作用散点图,红色散点图为分子间相互作用散点图。从红色散点可以看出,CCOF6与S-1-苯基-1-丙醇的分子间作用力明显小于其与R-1-苯基-1-丙醇的分子间作用力。利用VMD绘制sign(λ2)ρ填色的δg_inter等值面图(图3左侧)。其中,氢键区域等值面为明显的蓝色,绿色扁片为电子密度很小的弱相互作用区域,一般为色散作用区域。其他3对对映体的IGM等值面图如图4所示。从这些图可以非常直观地看到CCOF6与对映体之间的氢键和色散力。另外,还可从等值面的厚薄大小初步判断相互作用力的大小。例如,CCOF6与S-1-苯基-1-丙醇的IGM等值面基本没有蓝色区域,这表明其相互作用力非常弱,从表1也可以看出,CCOF6与S-1-苯基-1-丙醇的结合自由能绝对值远低于CCOF6与其他对映体的结合自由能绝对值。CCOF6与R-1-苯基-2-丙醇的绿色区域明显少于CCOF6与R-1-苯基-2-丙醇的绿色区域,这预示着CCOF6与R-1-苯基-2-丙醇的相互作用力弱于其与S-1-苯基-2-丙醇的相互作用力。

2.4 手性拆分机理

结合2.2节和2.3节可知,除S-1-苯基-1-丙醇的羟基与CCOF6的醚键形成弱相互作用外,其他对映体的羟基则与CCOF6的羰基形成氢键(O-H…O),因此S-1-苯基-1-丙醇与CCOF6的结合力最弱(即结合自由能绝对值最小)。通常,固定相与对映体之间的结合力差别越大,对映体的选择性因子也越大。根据文献[12],(±)-1-苯基-1-戊醇、(±)-1-苯基-1-丙醇、(±)-1-苯基-2-丙醇、(±)-1-(4-溴苯基)乙醇的出峰时间分别为24和29 min、32和43 min、24和32 min、25和31 min。出峰时间顺序与结合自由能大小顺序没有相关性。这主要是因为影响出峰时间的因素除了物质与固定相的相互作用外,还有物质与流动相的相互作用。不考虑其他因素的影响,物质与固定相的结合力越弱,出峰时间越早,而与CCOF6结合力最弱的1-苯基-1-丙醇出峰时间最迟。这说明,正己烷/异丙醇流动相大大推迟了1-苯基-1-丙醇的出峰时间,即流动相对1-苯基-1-丙醇的洗脱能力最弱。而流动相对1-苯基-2-丙醇的洗脱能力次弱,因为无论是S型还是R型对映体,其与CCOF6的结合力均位于同型中的倒数第二位,而其出峰时间仅早于1-苯基-1-丙醇。

另外,由表1可知,除了(±)-1-苯基-2-丙醇外,其余3对对映体的S型与CCOF6的结合力均弱于R型与CCOF6的结合力。相同的实验条件下,同一对对映体的出峰时间主要取决于对映体与固定相的相互作用。对映体与固定相的结合力越强,出峰时间越迟;结合力越弱,出峰时间越早。因此尽管文献[12]并没有验证对映体的出峰顺序,根据对映体与CCOF6的结合自由能大小,可以预估:除了S-1-苯基-1-丙醇出峰时间迟于R-1-苯基-1-丙醇,其余对映体均为R型出峰时间迟于S型。

3 结论

在计算CCOF6色谱固定相与对映体相互作用的结合自由能方面,ORCA程序比AutoDock更为精确,但需要考虑溶剂化效应。色谱固定相与对映体的结合自由能差与对映体选择性因子存在相关性,结合自由能差越大,对映体的选择性因子也越大。然而,色谱固定相与对映体的结合自由能差与分离度不一定存在相关性。除了S-1-苯基-1-丙醇的羟基是与CCOF6的醚键形成弱相互作用,其他对映体的羟基均与CCOF6的羰基形成氢键相互作用,且S-1-苯基-1-丙醇与CCOF6的结合力最弱。结合对映体的出峰顺序及其结合自由能,可以推断:正己烷/异丙醇流动相对1-苯基-1-丙醇的洗脱能力最弱,对1-苯基-2-丙醇的洗脱能力次弱;R-1-苯基-1-丙醇出峰时间早于S-1-苯基-1-丙醇,其余R型对映体的出峰时间均迟于S型对映体。