采用定量结构-活性关系模型预测有机化合物的土壤或沉积物吸附系数

技术领域

本发明属于生态风险评价测试策略领域，涉及一种采用定量结构-活性关系模型预测有机 化合物的土壤/沉积物吸附系数的方法。

背景技术

土壤/沉积物吸附系数是表征有机化合物在土壤/沉积物和水相之间的吸附、分配行为及其 环境归趋的重要参数，用分配系数K_p表示，K_p=C_S/C_W，其中C_S和C_W分别表示有机化合物 在土壤/沉积物和水中达到分配平衡时的浓度。由于土壤/沉积物的类型多样，成分复杂，为了 更加准确地表征分配行为，引入了标准化的分配系数K_oc，K_oc=K_p/X_oc，X_oc表示土壤/沉积物 中有机碳的含量(kg·L^-1)。有机化合物在土壤/沉积物和水相之间的分配程度，不仅影响有机化 合物的迁移行为，还影响有机化合物的水解、光解及生物降解等过程。因此，获得有机化合 物的K_oc值，对于评价该化合物的生态风险性具有重要意义。

通常，可以通过实验的方法如高效液相色谱法或批量平衡法测得化合物的K_oc值。截至 目前，已开展的研究工作中仅有800余种有机化合物具有实验测定的K_oc值。据估计，目前 人类日常使用的有机化学品达14万种以上，并且新的日常使用有机化学品以每年500～1000 种的速度增长。如果仅靠实验方法对这些有机化学品的K_oc值进行一一测定，不仅耗时、耗 力，而且，实验数据还具有时间滞后性，不能满足有毒有害化学品风险管理的“预先防范原则”。 因此，有必要发展一种非实验的方法快速预测有机化合物的K_oc值。

预测有机化合物K_oc值的模型主要有如下两类：第一类是根据有机化合物的其他理化参 数如正辛醇/水分配系数(K_ow)、水溶解度(S_w)等建立模型估算K_oc值，但是，该类模型仅仅对 logK_ow在0.5-7.5范围内的化合物预测较为准确，推广应用受到限制；第二类是根据有机化合 物的分子结构与K_oc之间的内在联系建立定量结构-活性关系(QSAR)模型估算K_oc值，该类模 型不受logK_ow值的限制，能够高效快捷地对有机化合物的K_oc值进行预测。

由于QSAR技术有助于实现有毒有害化学品污染管理的“预先防范原则”，能够减少或替 代相关实验，弥补实验数据的缺失、降低实验费用，在世界各国有毒有害化学品生态风险性 评价和管理方面得到了广泛的开发。2004年经济合作与发展组织(OECD)正式确定了QSAR 模型发展和使用的导则，具体如下：(1)具有明确定义的环境指标；(2)具有明确的算法； (3)定义了模型的应用域；(4)有适当的拟合度，稳定性和预测能力；(5)最好能够进行机 理解释。

截止目前，已有很多研究人员运用QSAR技术成功建立了有机化合物K_oc值的预测模型。 如文献“J.Environ.Sci.Health.,Part A,2012,47:441–449.”针对209种多氯反式-偶氮苯类化 合物，采用遗传算法-人工神经网络建立QSPR模型对该类化合物的logK_oc值进行预测，模型 具有较好的拟合能力和预测能力，但是该算法不透明，不利于机理解释。文献“Environ.Sci. Technol.,2006,40(22):7005–7011.”收集了571种有机化合物的logK_oc值，采用多元线性回归 (MLR)分析建立模型，虽然算法透明，但采用了29个参数，不便于使用。文献“QSAR Comb Sci, 2009,28(5):561–567.”收集了同样的571种有机化合物的logK_oc值，采用最小二乘支持向量 机的方法建立模型，模型的拟合能力和预测能力较MLR方法建立的模型有很大的提高，但 是该算法不够透明，不便于应用。文献“Chemosphere,2012,86:634–640.”仅使用由分子结构 计算的logP值作为变量，针对logP在0.5–7.5范围的化合物建立线性模型预测logK_oc值， 对于logP<0.5或logP>7.5的化合物建立非线性模型进行预测，大大减少了描述符的个数， 方便使用，但是，模型的拟合能力不如以往的模型。综上所述，目前现有模型存在算法不透 明或模型参数多的缺点，并且不能对溴代阻燃剂、真菌毒素和植物雌激素的logK_oc值进行准 确的预测。因此，有必要构建一个数据集涵盖溴代阻燃剂、真菌毒素和植物雌激素、具有明 确算法、便于应用推广的QSAR模型，并且，依照OECD导则对模型进行机理解释和表征。

发明内容

本发明目的是发展一种高效、快速、简洁且应用范围广的有机化合物K_oc值的预测方法。 该方法可直接根据有机化合物的分子结构预测K_oc值，进而了解目标化合物的迁移和归趋， 为化合物的生态风险性评价和管理提供必要的基础数据。

本发明的技术方案如下：

首先，搜集得到有机化合物的logK_oc值，如果一种化合物的logK_oc值有多个来源，则取 其平均值，最终得到813种有机化合物的logK_oc值，其中有机化合物包括多溴联苯醚、多氯 联苯、酞酸酯、多环芳烃及其取代物、杂环化合物及其衍生物、卤代烷烃、卤代烯烃、有机 酸、酯、醚、酮、醇、苯酚、苯胺等。对上述的813种有机化合物的logK_oc值划分为训练集 和验证集，训练集包括609种有机化合物，验证集包括204种有机化合物。训练集中的有机 化合物用于构建模型，验证集中的有机化合物用于模型构建后的外部验证。采用去一法对上 述构建的模型进行内部验证。

上述模型所使用的描述符均为Dragon描述符，用3545个Dragon描述符和训练集中的有 机化合物的logK_oc值进行逐步回归分析，得到模型的线性关系式如下：

logK_oc=0.063×MLOGP2+0.332×WiA_Dt+0.260×nHM-0.002×H_D/Dt+0.338×O-061- 1.037×HATS4v-0.803×P-117+1.011×nR=CRX-0.123×F05N-O+1.185×B08Br-Br-1.868 ×R3e+-0.537×B03N-S-0.227×CATS2D_05_NL+0.220×F02S-S+0.627×nRCN+0.546(1) 其中，MLOGP2表示Moriguchi辛醇-水分配系数的平方；WiA_Dt表示由迂回矩阵得到的 Wiener-like指数；H_D/Dt表示由迂回矩阵得到的Harary-like指数；nHM表示重原子个数； O-061表示硝基上的氧原子碎片数；HATS4v和R3e+是GETAWAY描述符，HATS4v与分子的 范德华体积有关，R3e+与分子尺寸、电负性相关；P-117表示分子中X3-P=X结构存在与否， 存在取1，不存在取0；nR=CRX表示分子中R=CRX结构个数；F05N-O表示分子中N-O在 拓扑距离5处出现的频率；B08Br-Br表示分子中Br-Br结构存在与否，存在取1，不存在取0； B03N-S表示分子中N-S结构存在与否，存在取1，不存在取0；CATS2D_05_NL是CATS2D 描述符，与分子的亲脂性有关；F02S-S表示分子中S-S在拓扑距离2处出现的频率；nRCN 表示亚硝基个数。

最终筛选出15个Dragon描述符用于模型构建，每个描述符的变量膨胀因子(VIF)均小于 10，自变量与因变量组成的矩阵M_YX以及自变量矩阵M_X的K相关指数K_YX和K_X分别为0.316 和0.284，ΔK=K_YX-K_X=0.032>0，表明模型不存在多重相关性；模型的相关系数的平方(R²) 为0.853，均方根误差(RMSE)为0.473，说明模型具有较好的拟合能力；模型的内部交叉验证 系数(Q²_Loo)为0.847，说明模型的稳健性较好；外部验证决定系数(Q²_ext)为0.801，验证集的均 方根误差(RMSE_ext)为0.550，表明该模型具有良好的外部预测能力。相对现有的预测模型，该 模型的应用域更广，涵盖了溴代阻燃剂、真菌毒素和植物雌激素。该模型采用了15个描述符， 便于模型的应用推广。

采用Williams图对模型的应用域进行表征。当化合物的Hat(h_i)值大于警戒值(h^*)时，说 明化合物不在模型应用域内。h_i和h^*由如下公式计算：

h_i=x_i^T(X^TX)^-1x_i        (2)

h^*=3(k+1)/n           (3)

其中x_i是第i个化合物的描述符矩阵；x_i^T是x_i的转置矩阵；X是所有化合物的描述符矩阵； X^T是X的转置矩阵；(X^TX)^-1是矩阵X^TX的逆；k是模型中变量的个数。模型的h^*为0.079， 因此，该模型适用于对h_i小于0.079的化合物logK_oc值的预测。

本发明的有益效果是采用本发明方法可以通过分子结构特征快速预测有机化合物的 logK_oc值。该方法不仅简单快捷、成本低廉，而且节省了实验测定所需的人力、物力和财力。 本发明涉及的logK_oc预测方法的建立和验证严格依据OECD规定的QSAR模型发展和使用导 则，因此，使用该发明专利的logK_oc的预测结果，可以为化学品监管提供数据支持，对化学 品的生态风险性评价具有重要意义。

本发明提供的方法具有如下特点：

1.建模过程采用了透明的算法-MLR，此外，15个描述符用于预测模型的构建，模型简洁、 易于解释，便于应用推广；

2.模型的应用域广，涵盖溴代阻燃剂、真菌毒素和植物雌激素等多种有机化合物，可用于预 测不同种类有机化合物的logK_oc值，为化学品的生态风险性评价和监管提供数据支持；

3.建模过程严格依照OECD关于QSAR模型的构建和使用导则，所建模型具有良好的拟合 能力、稳健性和预测能力。

附图说明

图1为训练集logK_oc的实测值与预测值的拟合图，训练集化合物为609种。

图2为验证集logK_oc的实测值与预测值的拟合图，验证集化合物为204种。

图3为模型的Williams图，●表示训练集化合物，Δ表示验证集化合物，警戒值h^*为0.079。

具体实施方式

实施例1

给定一个化合物尼古丁，要预测其logK_oc值。首先根据尼古丁的结构信息，使用Gaussian 09软件包对其进行结构优化，基于高斯优化的结构，采用Draogon6.0软件计算出nHM， WiA_Dt，H_D/Dt，HATS4v，R3e+，nRCN，nR=CRX，O-061，P-117，CATS2D_05_NL， B03N-S，B08Br-Br，F02S-S，F05N-O和MLOGP2的值分别为0，5.955，156.433，0.159，

0.069，0，0，0，0，0，0，0，0，0和1.619。根据公式(2)计算的h值为0.005(<0.079)， 所以该化合物在模型应用域内，将以上描述符的值代入公式(1)，得到logK_oc的预测值为2.02， 其实验测定的logK_oc值为2.01，预测值和实验值的数据非常相符。

实施例2

给定一个化合物邻苯二甲酸二乙酯，要预测其logK_oc值。首先根据邻苯二甲酸二乙酯的 结构信息，使用Gaussian09软件包对其进行结构优化，基于高斯优化的结构，采用Draogon 6.0软件计算出nHM，WiA_Dt，H_D/Dt，HATS4v，R3e+，nRCN，nR=CRX，O-061，P-117， CATS2D_05_NL，B03N-S，B08Br-Br，F02S-S，F05N-O和MLOGP2的值分别为0，5.85， 203.902，0.102，0.066，0，0，0，0，0，0，1，0，0和6.641。根据公式(2)计算的h值为 0.003(<0.079)，所以该化合物在模型应用域内，将以上描述符的值代入公式(1)，得到logK_oc的预测值为2.27，其实验测定的logK_oc值为1.84，预测值和实验值的数据非常相符。

实施例3

给定一个化合物2-氯联苯醚，要预测其logK_oc值。首先根据2-氯联苯醚的结构信息，使 用Gaussian09软件包对其进行结构优化，基于高斯优化的结构，采用Draogon6.0软件计算 出nHM，WiA_Dt，H_D/Dt，HATS4v，R3e+，nRCN，nR=CRX，O-061，P-117，CATS2D_05_NL， B03N-S，B08Br-Br，F02S-S，F05N-O和MLOGP2的值分别为1，6.462，192.8，0.24，0.111， 0，0，0，0，0，0，0，0，0和19.499。根据公式(2)计算的h值为0.007(<0.079)，所以该化 合物在模型应用域内，将以上描述符的值代入公式(1)，得到logK_oc的预测值为3.34，其实验 测定的logK_oc值为3.47，预测值和实验值的数据非常相符。

实施例4

给定一个化合物甲基对硫磷，要预测其logK_oc值。首先根据甲基对硫磷的结构信息，使 用Gaussian09软件包对其进行结构优化，基于高斯优化的结构，采用Draogon6.0软件计算 出nHM，WiA_Dt，H_D/Dt，HATS4v，R3e+，nRCN，nR=CRX，O-061，P-117，CATS2D_05_NL， B03N-S，B08Br-Br，F02S-S，F05N-O和MLOGP2的值分别为2，5.242，184.6，0.258, 0.062，0，0，2，0，0，0，0，0，1和2.865。根据公式(2)计算的h值为0.020(<0.079)，所 以该化合物在模型应用域内，将以上描述符的值代入公式(1)，得到logK_oc的预测值为2.79， 其实验测定的logK_oc值为3，预测值和实验值的数据非常相符。

实施例5

给定一个化合物2,4-二氯苯氧基乙酸，要预测其logK_oc值。首先根据2,4-二氯苯氧基乙酸 的结构信息，使用Gaussian09软件包对其进行结构优化，基于高斯优化的结构，采用Draogon 6.0软件计算出nHM，WiA_Dt，H_D/Dt，HATS4v，R3e+，nRCN，nR=CRX，O-061，P-117， CATS2D_05_NL，B03N-S，B08Br-Br，F02S-S，F05N-O和MLOGP2的值分别为2，5.231， 139.967，0.169，0.139，0，0，0，0，3，0，0，0，0和5.543。根据公式(2)计算的h值为 0.111(>0.079)，所以该化合物在模型应用域外，将以上描述符的值代入公式(1)，得到logK_oc的预测值为1.76，其实验测定的logK_oc值为1.66，预测值和实验值的数据仍然非常相符，说 明对于应用域外的化合物，本模型的预测结果也具有一定的参考价值。