硝基芳烃对梨形四膜虫急性毒性的CoMFA模型

唐自强，冯长君

1. 徐州工程学院材料与化学工程学院，徐州 221018 2. 徐州技师学院，徐州 221151

分子结构是决定物质理化性质及其在环境中迁移转化行为和生态毒理学效应的内因，因此，结构与性质之间具有内在联系，以数学方程形式予以表达，就形成定量结构-活性相关关系(quantitative structure-activity relationship, QSAR)。经统计诊断为质量良好的QSAR模型，可为有机污染物环境毒性评估提供参考信息，以降低昂贵的测试费用和减少动物实验[1-8]。

硝基芳烃类有机物结构稳定，种类繁多且复杂，难以降解，属高毒污染物。因此，这类化合物一直备受关注，有些硝基芳烃已被美国环境保护局和中国生态环境保护部列入优先控制污染物黑名单。为预测硝基芳烃类化合物生物毒性，闫秀芬等[9]和顾云兰等[10]采用量化方法研究硝基芳烃对梨形四膜虫(Tetrahymenapyriformis)的急性毒性(以半数生长抑制浓度的负对数值(pIC50)表征)；何琴等[11]和堵锡华等[12]基于人工神经网络方法给出45种硝基芳烃对梨形四膜虫急性毒性的良好QSAR模型；金飙等[13]以量化参数给出了36种硝基芳烃对梨形四膜虫急性毒性的人工神经网络模型；Liao等[14]以三维全息原子相互作用矢量场的(three-dimensional holographic vector of atomic interaction field, 3D-HoVAIF)方法建立了上述45种硝基芳烃对梨形四膜虫急性毒性的QSAR模型。这些研究多为QSAR的二维方法，且所用三维-定量构效关系(3D-QSAR)方法[15-17]中也没有涉及使用最为广泛的比较分子力场分析法(comparative molecular force field analysis, CoMFA)。因此，笔者以CoMFA方法建立了45种硝基芳烃对梨形四膜虫急性毒性pIC50值[14]的3D-QSAR模型，以揭示硝基芳烃分子周围立体场和静电场分布对其急性毒性的影响，探讨其生物毒性作用的分子机理。

1 数据来源(Data source)

45种硝基芳烃对梨形四膜虫的急性毒性，以IC50表示，文献[14]给出其pIC50。硝基芳烃分子的具体结构及pIC50值如表1所示。

表1 硝基芳烃的分子结构与急性毒性pIC50值Table 1 The molecular structures and pIC50 value of acute toxicity of nitroaromatic compounds

2 建模方法(Modeling method)

使用Tripos公司最新版本Sybylx2.1.1分子模拟软件完成45种硝基芳烃[14]对梨形四膜虫的急性毒性的3D-QSAR分析及建立CoMFA模型。

2.1 化合物低能构象的确定及分子叠合

一般采用分子的最低能量构象代替其生物活性构象。首先使用Sybylx2.1.1软件中的Sketch Molecule模块构建45个硝基芳烃分子[14]的初始三维空间结构，然后通过Minimize模块，选取Tripos力场，加Gasteiger-Huckel电荷。以此进行所有分子的分子力学优化，获取相应最低能量构象，并以此构象作为分子叠合的生物活性构象。

将45个硝基芳烃分子按照从低至高排序，遵循5∶1原则选取分子2、7、12、17、22、27、32、37、42和45为测试集(test set)，其中，45号为模板分子；余下36个分子作为训练集(training set)，亦含45号模板分子。在CoMFA研究中，所谓模板分子是指总体生物活性最强的分子。常用的是基于公共骨架的原子契合叠合方法，即保证所有分子取向的一致性，使分子间相互重叠时的均方根偏差最小。以对梨形四膜虫的急性毒性pIC50值最大的2,3,4,5-四氯硝基苯分子中的公共骨架，运用Align Database模块分别对训练集和测试集中分子予以叠合。训练集的叠合图如图1所示(测试集的叠合图与图1相同，故予以省略)，由图1可知，公共骨架的原子契合非常密切，呈现良好的叠合效果。

2.2 CoMFA模型的建立

采用Tripos标准力场，对训练集和测试集的叠合分子周围每个网格点上的立体场(steric, St)和静电场(electrostatic, El)的场能值进行计算，采用偏最小二乘法(partial least squares, PLS)方法进行建模[18-19]。第一步，用逐一剔除法(leave-one-out, LOO)予以交叉验证，即以pIC50为因变量，场能值为自变量建模，得交叉验证系数(Q2)和最佳主成分数(N)。以Q2衡量模型的稳健性与预测能力，Q2>0.5的模型，才具有95%的可信度[20]，所含机会相关的可能性<5%。模型中的化合物数(m)与变量数(即最佳主成分数)之比称为样变比(Sv)。只有满足Sv≥5的模型，才可能具有统计学意义且或然性较低[21-22]。第二步，用非交叉验证的回归分析，得到非交叉验证判定系数(R2)、估计标准误差(SD)和统计方差比(F)。一般认为，判定系数(R2)>0.80的方程是合理的，具有良好的拟合能力[23]。最后采用View CoMFA模块，以三维等势图直观反映出立体场和静电场对化合物pIC50的贡献。

3 结果与讨论(Results and discussion)

3.1 CoMFA的3D-QSAR模型

训练集的CoMFA模型：交叉验证部分中Q2=0.738>0.5，N=5，Sv=m/N=7.2>5，表明模型具有很好的预测能力和稳健性；非交叉验证部分中R2=0.932>0.8，显示急性毒性与场能之间具有良好相关性；在95%显著性水平下，训练集的CoMFA模型的F临界值为F0.05(5,30)=2.53。而模型的F值为82.40，是其临界值2.53的32倍之多，表明所建模型密切相关。该模型的SD为0.201，仅是最大差值(ΔS=pIC50 max-pIC50 min=2.74-0.14=2.60)的7.7%，<10%，说明不存在离异点[24]。利用训练集的模型对测试集中分子的pIC50进行预测，以检验其预测能力。测试集的pIC50预测值如表1所示，与相应实验值基本吻合，平均绝对误差仅为0.28，表明模型确实具有较好的预测能力。

3.2 CoMFA等势图

图2(a)为训练集中以pIC50最高的分子45号为模板的CoMFA模型的三维立体作用分布图，其中，绿色区域表示增大取代基体积有利于提高化合物的急性毒性，黄色区域为减小立体位阻有利于提高急性毒性。由立体场的空间分布可知，在硝基的2、3和6位的3个位置上引入体积较大基团，相应化合物的pIC50升高。例如，在硝基的2、3和6位的3个位置上，溴取代的化合物急性毒性pIC50值高于氯取代的。分子17号强于12号，26号强于16号，29号强于23号。尚需指出，图2(a)中未出现黄色区域，这可能因为在苯环上的所有取代基的体积都比梨形四膜虫体内靶标提供的空穴要小，不存在空间位阻。

图1 训练集的三维叠合图Fig. 1 The three dimensional view of all the aligned molecules in training set

图2 CoMFA模型的立体场(a)与静电场(b)等势图注：CoMFA表示比较分子力场分析。Fig. 2 CoMFA contour maps for (a) steric field and (b) electrostatic fieldNote: CoMFA stands for comparative molecular force field analysis.

图2(b)为45号分子周围的静电场分布图，蓝色区域为增大基团正电性有利于化合物活性的提高，红色区域则表示引入带负电荷的基团对活性有利。由静电场的空间分布可知，仅有一块蓝色区域在硝基的2、3位，且覆盖苯环2、3位上的碳原子；红色区域较多，表现出硝基的2、3位之间以及较远离3位的2个红色区域对急性毒性影响较大。如在硝基的2、3位的2个位置上，氟(负电性最强的取代基)取代的化合物急性毒性pIC50值强于氯(负电性较弱)取代的。分子34号强于21号，37号强于32号。

训练集模型的立体场和静电场对pIC50值的贡献率依次为39.0%和61.0%，表明静电场强于立体场。其中，静电场主要为电荷库伦力、氢键和配位作用。苯环上连有负电基团(—F、—Cl)等，使苯环上正电荷量增加，易进攻梨形四膜虫体内酶、蛋白质等生物大分子的重要功能基团(如带有负电荷的巯基、羟基和氨基等),从而影响梨形四膜虫的正常生理功能，导致生物中毒。另外，苯环上的负电基团可与靶标内金属离子形成配位键。对于立体场，一般包含空间位阻和疏水作用。其中，取代基体积越大，相应分子疏水性越强，脂溶性越好，越易被细胞膜富集而产生生物毒性。顾云兰等[10]以量化参数建立了化合物的三元数学模型，模型中的变量包括苯环上净电荷增量(ΔQR)、分子体积(V)和最低空轨道能。显然，本研究结果与之基本一致。

综上所述，本研究结果表明：

(1)基于CoMFA方法研究硝基芳烃对梨形四膜虫急性毒性pIC50值的三维定量构效关系：N=5，Sv=7.2，Q2=0.738，R2=0.932，F=82.40，显示良好的相关性、稳定性和预测能力。

(2)所建CoMFA模型的三维等势图较好地解释了硝基芳烃对梨形四膜虫急性毒性pIC50值的影响，主要表现在苯环的2、3-位上键合体积较大的负电基团(如—F、—Cl等)利于急性毒性的增强。