CN102507630A - 一种基于分子结构和环境温度预测化学物质与臭氧氧化反应速率常数的方法 - Google Patents

Abstract

本发明公开了一种基于分子结构和环境温度预测化学品与臭氧反应速率常数(k_O3)的方法，属面向化学品环境监管的定量结构-活性关系(QSAR)技术领域。k_O3影响化学品在大气中滞留时间(半衰期)，是评估化学品的环境持久性的重要参数。本发明采用了量子化学描述符、Dragon计算的描述符、分子结构碎片表征化合物的分子结构特征，构建了预测k_O3的QSAR模型。所构建的不同环境温度下k_O3的QSAR模型，经内部验证和外部验证证实了其有良好的稳健性和预测能力，并基于Williams图法表征了QSAR模型的应用域，明确了模型的应用范围。

Description

一种基于分子结构和环境温度预测化学物质与臭氧氧化反应速率常数的方法

技术领域

本发明涉及不同环境温度下臭氧与有机化学物质氧化反应速率常数(k_O3)的快速预测方法，属于面向环境风险评价的定量结构与活性关系(QSAR)技术领域。 

背景介绍 

大多数的工业化学品具有一定的生化稳定性，由光所引发的直接降解以及自由基反应是其在环境中转化的主要途径。在大气层中有机污染物的持久性(半衰期)可以通过与羟基自由基、硝基自由基和臭氧进行化学反应的速率常数等信息进行计算。大气对流层中的有机污染物，日间主要是通过羟基自由基和臭氧降解来去除，夜间主要是通过硝基自由基反应降解。其中，化合物与臭氧反应的降解速率常数(k_O3)是表征化学品在环境中持久性的关键参数。 

目前对于化学品环境持久性评价多基于实验测试，如光解活性和活性氧物种(ROS)氧化活性测试等，但完全依靠实验测试来获得化学品的环境持久性数据，具有巨大的财政压力。同时，新合成的日常使用有机化学品正以每年500～1000种的速度增加，对这些化学品一一进行实验测试，无法满足环境管理的要求，亟需发展成本低且快速的化学品环境持久性评估方法。基于定量结构与活性关系(QSAR)，可以从化学品的分子结构出发，预测化学品与臭氧反应的速率常数，弥补基础数据的缺失，用于评估化学品的光稳定性。 

已报导关于k_O3的QSAR模型，但这些模型在实用性、预测能力与稳健性、应用域等方面仍存在一些问题。Fatemi等(Fatemi，M.H.Prediction of ozone tropospheric degradation rate constant of organic compounds by using artificial neural networks.Analytica Chimica Acta.2006，556：355-363)建立了基于人工神经网络(ANN)的137种化合物与臭氧反应的QSAR模型，该模型采用非线性的方法，主要用于预测298K条件下有机化合物与臭氧反应的k_O3，模型的透明度不强，不能预测不同温度下的k_O3。Pompe和Veber(Pompe，M.，Veber，M.，Prediction of rate constants for the reaction of O₃ with different organic compounds.Atmospheric Environment.2001，35(22)：3781-3788)建立了基于多元线性回归(MLR)的117种不同种类化合物与臭氧反应的QSAR模型，不能预测多个温度下的k_O3。Jiang等(Jiang，J.L.，Yue，X.A.，Chen，Q.F.Determination of ozonization reaction rate constants of aromatic pollutants and QSAR study.Bulletin of Environmental Contamination and Toxicology.2010，85：568-72)建立了基于密度泛函理论(DFT)的39种芳香族化合物与臭氧反应的QSAR模型，该模型选择的化合物种类单一，应用域偏小，需扩大应用域的范围。 

构建实用性强且算法简单透明的k_O3温度依附性预测模型，并对其预测能力和稳健性进 行评价，在此基础上对应用域范围进行表征，对有机化学品的生态风险评价与监管工作具有重要意义。 

发明内容

本发明要解决的技术问题是提供一种不同环境温度下有机化合物与臭氧反应速率常数的简便、快速、低成本、应用范围广的预测方法。 

为解决以上技术问题本发明采用的技术方案如下： 

(1)为了保证用于建立模型的数据准确性，对从文献收集的实验测定值进行评估和分析，删除离群点。首先，对同一化合物的多个实验测定值进行统计学评估，删除偏离平均值的数据点；其次，对同一化合物在不同温度下的logk_O3对1/T作图分析，删除偏离线性的数据点。最终得到129种有机化合物共264个logk_O3数据，温度范围从178K到364K。根据化合物与臭氧反应的机理分析，计算26个量子化学描述符、1481个Dragon描述符和12个分子结构碎片，加入1/T温度描述符。化合物包括单烯烃、二烯烃及多烯烃、卤代烯烃、环状化合物、芳香族化合物、含氧化合物、含氮化合物等。以4∶1的比例随机划分训练集和验证集。 

(2)联合应用逐步回归和PLS回归筛选训练集的描述符组合，构建预测k_O3的QSAR模型，具体步骤如下： 

第一步，采用逐步回归方法进行初筛变量，在多元线性回归(MLR)模型中，膨胀因子(VIF)小于10； 

第二步，基于PLS回归，筛选变量的过程如下：先建立包含所有初筛得到的变量的PLS模型，选择模型相关系数的平方(R²)和累计交叉验证系数(Q² _CUM)为参考值，R²和Q² _CUM是由Simca6.0软件计算得出；分别去除每个描述符一次，以R²和Q² _CUM数值升高最大的模型为基准模型，重复上述步骤直至R²和Q² _CUM不再升高为止。若有几个方程的R²和Q² _CUM数值相近，则比较经自由度调整后的决定系数(R² _adj)，选择R² _adj最大的模型为最优模型。 

最优模型的表达式： 

logk_O3＝-12.542-493.3×(1/T)+0.41722×E_HOMO+0.4443×electrophility+0.66971×n_C＝C

-0.26128×qC_max+0.74783×BELm2+4.8412×Mor32v+0.35198×H3u+0.38372×n_＝CHR

-1.7438×n_NH2+0.4576×n_＝CR2-1.1235×n_BM+0.28542×n_CIRCLE    (1) 

其中，1/T表示温度的倒数；E_HOMO表示最高占据分子轨道能；electrophility表示亲电性；n_C＝C表示碳碳双键的个数；qC_max表示碳原子最大正电荷；BELm2表示BCUT描述符；Mor32v表示3D-Morse描述符；H3u表示Getaway描述符；n_＝CHR表示＝CHR的个数；n_NH2表示氨基的个数；n_＝CR2表示＝CR₂的个数；n_BM表示苯环上甲基的个数；n_CIRCLE表示脂肪环的个数。 

经内部验证和外部验证证实k_O3预测模型的稳健性和预测能力。模型的拟合能力由R² _adj 和均方根误差(RMSE)来表征，R² _adj＝0.849，RMSE＝0.562，表明模型具有较好的拟合能力。稳健性用内部验证的交叉验证系数评价，交叉验证系数Q² _CUM＝0.838，表明该QSAR模型具有较好的稳健性。模型的预测能力采用未参与建模的50个外部数据评价，由外部预测相关系数的平方(Q² _ext)来表征，Q² _ext＝0.878，表明模型具有较好的预测能力。其中，R² _adj，RMSE，Q² _ext由下列公式得出： 

$R_{\mathrm{adj}}^2=1-\frac{{\mathrm{\Σ}}_{i=1}^n{(y_i-{\hat{y}}_i)}^2/(n-p-1)}{{\mathrm{\Σ}}_{i=1}^n{({\hat{y}}_i-\stackrel{\‾}{y})}^2/(n-1)}$

$\mathrm{RMSE}=\sqrt{\frac{{\mathrm{\Σ}}_{i=1}^n{(y_i-{\hat{y}}_i)}^2}{n}}$

$Q_{\mathrm{ext}}^2=1-\frac{{\mathrm{\Σ}}_{i=1}^{n_{\mathrm{EXT}}}{({\hat{y}}_i-y_i)}^2}{{\mathrm{\Σ}}_{i=1}^{n_{\mathrm{EXT}}}{(y_i-{\stackrel{\‾}{y}}_{\mathrm{EXT}})}^2}$

其中，y_i和 

分别为第i个化合物的实验值和预测值， 

为训练集所有数据点的实验值的平均值， 

为验证集所有数据点的实验值的平均值，n为训练集数据点个数，p为描述符个数，n_EXT为验证集的数据点个数。 

采用Williams图来表征模型的应用域。综合应用分子结构描述符的杆杠值(h_i)以及化合物分子结构描述符的预警值(h^*＝0.196)，以标准化的交叉验证残差为纵坐标，杆杠值为横坐标表征离群点。 

本发明的有益效果是采用本发明方法可以有效预测不同环境温度下，宽范围、多种类有机污染物与臭氧反应的速率常数，该方法成本低廉、简便而快速，能够节省大量的人力、物力、财力。该发明涉及的k_O3预测方法的建立和验证严格依照经济合作与发展组织(OECD)规定的QSAR模型发展和使用的导则，因此使用该发明专利的k_O3预测结果，可以为化学品监管工作提供重要的基础数据，并对生态风险评价具有重要的指导意义。 

本发明提供的k_O3预测模型有如下几个优点：(1)不同环境温度下(温度范围从178K到364K)有机物k_O3的预测模型，可用于计算大气对流层中不同温度区间内污染物被臭氧氧化降解的半衰期；(2)预测模型中的分子结构表述符可简单计算获得；(3)模型具有较高的稳定性和预测能力；(4)明确了模型最佳应用范围。 

附图说明

图1a是训练集logk_O3的实测值与预测值的对比图，训练集共110种有机化合物，214个logk_O3数据。 

图1b是验证集logk_O3的实测值与预测值的对比图，验证集共33种有机化合物，50个logk_O3数据。 

图2是表征logk_O3预测模型应用域的Williams图。 

具体实施方式

实施例1 

1-庚烯：采用Williams图法计算得到其杆杠值为0.0576＜h^*(预警值)＝0.196，标准残差(SE)＝0.2838＜3，说明此化合物在QSAR模型应用域内。采用MOPAC 2009的PM6算法，DRAGON软件(version 2.1)，以及考虑化合物的官能团及个数等特征最终计算出13个描述符。 

1-庚烯(296K)的logk_O3实验预测值为-16.76cm³molecule^-1s^-1。基于QSAR模型预测步骤如下： 

logk_O3＝-12.542-493.3×(0.003378)+0.41722×(-9.970)+0.4443×(-1.6584)+0.66971×1-

        0.26128×(-0.0611)+0.74783×1.684+4.8412×(-0.128)+0.35198×1.354+ 

        0.38372×1-1.7438×0+0.4576×0-1.1235×0+0.28542×0 

        ＝-16.92 

实施例2 

1，1-二氯乙烯：采用Williams图法计算得到其杆杠值为0.0616＜h^*(预警值)＝0.196，标准残差(SE)＝-3.12＜-3，说明此化合物在应用域外，使用该化合物的logk_O3的预测值要提高警惕。采用MOPAC 2009的PM6算法，DRAGON软件(version 2.1)，以及考虑化合物的官能团及个数等特征最终计算出了13个描述符。 

1，1-二氯乙烯(298K)的logk_O3实验预测值为-20.43cm³molecule^-1s^-1。基于QSAR模型预测步骤如下： 

logk_O3＝-12.542-493.3×(0.003356)+0.41722×(-10.225)+0.4443×(-2.5540)+0.66971×1 

        -0.26128×0.0855+0.74783×0.000+4.8412×0.058+0.35198×0.004+0.38372× 

        0-1.7438×0+0.4576×0-1.1235×0+0.28542×0 

        ＝-18.67 

实施例3 

莰烯：采用Williams图法计算得到其杆杠值为0.213＞h^*(预警值)＝0.196，标准残差(SE)＝-2.78＞-3，说明此化合物在应用域外，使用该化合物的logk_O3的预测值要提高警惕。采用MOPAC 2009的PM6算法，DRAGON软件(version 2.1)，以及考虑化合物的官能团及个数等特征最终计算出了13个描述符。 

莰烯(298K)的logk_O3实验预测值为-18.05cm³molecule^-1s^-1。基于QSAR模型预测步骤如下： 

logk_O3＝-12.542-493.3×(0.003356)+0.41722×(-9.663)+0.4443×(-1.5099)+0.66971×1-

        0.26128×0.1549+0.74783×1.705+4.8412×(-0.196)+0.35198×1.246+0.38372 

        ×0-1.7438×0+0.4576×1-1.1235×0+0.28542×2 

＝-16.48 

实施例4 

甲胺：采用Williams图法计算得到其杆杠值为1.0115＞h^*(预警值)＝0.19626，标准残差(SE)＝-0.54＞-3，在应用域外，模型有较好的范化能力。采用MOPAC 2009的PM6算法，DRAGON软件(version 2.1)，以及考虑化合物的官能团及个数等特征最终计算出了13个描述符。 

甲胺(296K)的logk_O3实验预测值为-19.67cm³molecule^-1s^-1。基于QSAR模型预测步骤如下： 

logk_O3＝-12.542-493.3×(0.003378)+0.41722×(-9.415)+0.4443×(-0.6806)+0.66971×0-

        0.26128×(-0.3390)+0.74783×0.750+4.8412×0.031+0.35198×0.083+0.38372 

        ×0-1.7438×1+0.4576×0-1.1235×0+0.28542×0 

        ＝-19.35 

实施例5 

亚硝酸乙酯：采用Williams图法计算得到其杆杠值为0.0658＜h^*(预警值)＝0.19626，标准残差(SE)＝0.2707＜3，说明此化合物在应用域内。采用MOPAC 2009的PM6算法，DRAGON软件(version 2.1)，以及考虑化合物的官能团及个数等特征最终计算出13个描述符。 

亚硝酸乙酯(310K)的logk_O3实验预测值为-18.80cm³molecule^-1s^-1。基于QSAR模型预测步骤如下： 

logk_O3＝-12.542-493.3×(0.003226)+0.41722×(-10.062)+0.4443×(-2.9546)+0.66971×0 

        -0.26128×0.0274+0.74783×0.909+4.8412×(-0.022)+0.35198×0.349+ 

        0.38372×0-1.7438×0+0.4576×0-1.1235×0+0.28542×0 

        ＝-18.96 。

Claims (2)

1.一种基于分子结构和环境温度预测化学物质与臭氧氧化反应速率常数的方法，其特征在于，

(1)首先，对同一化合物的多个实验测定值进行统计学评估，删除偏离平均值的数据点；其次，对同一化合物在不同温度下的logk_O3对1/T作图分析，删除偏离线性的数据点；最终得到129种有机化合物共264个logk_O3数据，温度范围从178K到364K；根据化合物与臭氧反应的机理分析，计算26个量子化学描述符、1481个Dragon描述符和12个分子结构碎片，加入1/T温度描述符；

(2)联合应用逐步回归和PLS回归筛选描述符组合，构建预测k_O3的QSAR模型；具体步骤如下：

第一步，采用逐步回归方法进行初筛变量，在多元线性回归(MLR)模型中，膨胀因子(VIF)小于10；

第二步，基于PLS回归，筛选变量的过程如下：先建立包含所有初筛得到的变量的PLS模型，选择模型相关系数的平方(R²)和累计交叉验证系数(Q² _CUM)为参考值，分别去除每个描述符一次，以R²和Q² _CUM数值升高最大的模型为基准模型，重复上述步骤直至R² _Y和Q² _CUM不再升高为止；若有几个方程的R²和Q² _CUM数值相近，则比较训练集的经自由度调整后的决定系数(R² _adj)，选择R² _adj最大的模型为最优模型；

最优模型的表达式：

logk_O3＝-12.542-493.3×(1/T)+0.41722×E_HOMO+0.4443×electrophility+0.66971×n_C＝C

-0.26128×qC_max+0.74783×BELm2+4.8412×Mor32v+0.35198×H3u+0.38372×n_＝CHR

-1.7438×n_NH2+0.4576×n_＝CR2-1.1235×n_BM+0.28542×n_CIRCLE    (1)

其中，1/T表示温度的倒数；E_HOMO表示最高占据分子轨道能；electrophility表示亲电性；n_C＝C表示碳碳双键的个数；qC_max表示碳原子最大正电荷；BELm2表示BCUT描述符；Mor32v表示3D-Morse描述符；H3u表示Getaway描述符；n_＝CHR表示＝CHR的个数；n_NH2表示氨基的个数；n_＝CR2表示＝CR₂的个数；n_BM表示苯环上甲基的个数；n_CIRCLE表示脂肪环的个数。

2.根据权利要求1所述的的方法，其特征在于，所述的化合物包括单烯烃、二烯烃及多烯烃、卤代烯烃、环状化合物、芳香族化合物、含氧化合物、含氮化合物。

Images