化学品在给人们生产生活带来便利的同时，也给人体健康带来极大威胁。大量的挥发性有机化学品可通过直接排放或从土壤和水相挥发进入大气层。大气中有机化学品可通过物理过程移除，如干、湿沉降；也可通过化学过程降解，如直接光解或与大气氧化剂(OH、NO3和臭氧)反应。在大气对流层中有机物与臭氧的反应是其转化的重要途径。表征有机化学品与O3自由基反应的速率常数(KO3，cm3·mol-1·s-1)是反映有机污染物在大气中持久性能力的重要参数，是进行有机污染物生态风险评价的基础指标[1-2]。然而，目前KO3的实验数据较少，且实验耗时, 费力, 成本高，不能满足有机化学品生态风险评价的需求。因此，有必要开发出快速有效预测KO3的方法[3]。目前，定量结构-活性关系(QSAR)模型是用于获取KO3的重要方法。为指导各国构建满足化学品风险管理需求的QSAR模型，经济合作与发展组织于2007年发布了QSAR模型构建与验证的导则(简称导则)[4]，导则全面阐述了符合管理要求的QSAR模型应满足的标准。
目前，已有研究开发了关于臭氧反应速率常数的QSAR模型[5-6]。2007年，REN等[7]报道了116种有机化合物基于KO3的QSAR模型，该模型采用DUPLEX分类算法划分模型训练集和测试集，采用启发式方法(Heuristic Algorigthm)筛选最优描述符，并采用多元线性回归、支持向量机和投影寻踪回归方法构建预测模型。GRAMATICA等[8]基于遗传算法筛选最优描述符，采用多元线性回归算法构建了125种有机化合物基于KO3的QSAR模型，结果显示留一法交叉验证系数(QLOO2)达到82%~88%，外部验证决定系数(QEXT2)达到90%，均方根误差(RMSE)达到0.73。由美国国家环境保护局有毒物质污染防治办公室和Syracuse Research Corporation (SRC)公司共同开发的EPI(estimation programs interface)Suite软件，采用基团贡献法构建了112个烯烃和炔烃不饱和有机物基于KO3的QSAR模型，相关系数达到0.94，绝对平均残差达到0.35[9]。但上述模型并不满足导则要求，如缺少模型稳健性和预测能力表征，或未定义模型应用域[6-9]，不利于模型使用者评估需预测的有机化合物是否处于模型应用域内。因此，根据导则要求采用简单透明的遗传算法-多元线性回归(GA-MLR)算法构建基于KO3的新QSAR模型，并对模型进行拟合优度、稳健性、预测能力、应用域表征和机制解释。所构建的模型对实现环境行为参数预测软件化具有重要意义。
1 材料与方法 1.1 数据来源与处理烷烃、烯烃、芳香烃、含氧挥发性有机物和酚类等152种有机化学品的KO3数据来源于文献[10]。选择-lg KO3作为模型指标。为避免样本分布不均匀，采用KENNARD等[11]分组方法将数据集划分为训练集和验证集，将结构差异较大的样本选入训练集，其他与之相近的样本选入验证集，从而使代表性样本全部进入训练集。训练集有107种化学品，验证集有45种化学品。
1.2 分子结构描述符的计算分子结构描述符是用于反映分子结构信息的参数，根据分子结构按照一定理论或规则计算得到。笔者构建的模型采用的分子结构描述符为Dragon描述符。采用Hyperchem 7.0软件中MM+和AM1方法对选取的152种有机化学品结构进行优化[12]，采用Dragon 5.4软件计算优化后结构的描述符[13]，并对得到的1 664个描述符进行初步筛选，去掉常数项、近似常数项和高度相关的分子结构描述符，最终得到488个分子结构描述符。
1.3 模型的建立采用MobyDigs软件中遗传算法选择与-lg KO3高度相关的描述符[14]。由遗传算法变量筛选得到最优描述符，并采用多元线性回归(MLR)方法构建预测模型。遗传算法相关参数设置为种群大小(population size)为100，初始模型允许的最大变量数(maximum allowed variables)为7，变异均衡值(mutation trade-off，T)为0.5，交叉(crossover)和变异(mutation)概率均基于T值。当增加变量数目对结果影响不大时，得到8个最优描述符。
1.4 模型的表征与评价根据导则要求，对构建的QSAR模型进行内部验证(训练集的拟合优度和稳健性评估)和外部验证(验证集的预测能力评估)。采用实验值与预测值之间校正后的决定系数(Radj2)和均方根误差(RMSE，ERMS)表征模型拟合优度，采用留一法交叉验证系数(QLOO2)表征模型稳定性，采用外部检验参数(QEXT2)、验证集相关系数(REXT2)和验证集均方根误差(ERMS，EXT)等外部验证决定系数表征模型预测能力，基于杠杆值(leverage，hi)的Williams图定义模型应用域[15]。外部验证决定系数计算公式为
| $ {R_{{\rm{adj}}}}^2 = 1 - \frac{{\sum\limits_{i = 1}^n {{{({y_i} - {{\hat y}_i})}^2}/\left( {n - m - 1} \right)} }}{{\sum\limits_{i = 1}^n {{{({y_i} - \bar y)}^2}/\left( {n - 1} \right)} }}, $ | (1) |
| $ {E_{{\rm{RMS}}}} = \sqrt {\sum\limits_{i = 1}^n {{{({y_i} - {{\hat y}_i})}^2}/n} }, $ | (2) |
| $ {Q_{{\rm{LOO}}}}^2 = 1 - \sum\limits_{i = 1}^n {{{({y_i} - {{\hat y}_i})}^2}/\sum\limits_{i = 1}^n {{{({y_i} - {{\bar y}_{{\rm{TR}}}})}^2}} } , $ | (3) |
| $ {Q_{{\rm{EXT}}}}^2 = 1 - \sum\limits_{i = 1}^{{n_{{\rm{EXT}}}}} {{{({y_i} - {{\hat y}_i})}^2}\sum\limits_{i = 1}^{{n_{{\rm{EXT}}}}} {{{({y_i} - {{\bar y}_{{\rm{EXT}}}})}^2}} } 。$ | (4) |
式(1)~(4)中，n为化合物总种数；m为预测变量个数；yi和
Williams图是标准残差(δ)和hi值定义的模型应用域，其计算公式为
| $ \delta = {y_i} - {{\hat y}_i}/\sqrt {\sum\limits_{i = 1}^n {{{({y_i} - {{\hat y}_i})}^2}/\left( {n - m - 1} \right)} }, $ | (5) |
| $ {h_i} = {x_i}^{\rm{T}}{\left( {{{\mathit{\pmb{X}}}^{\rm{T}}}{\mathit{\pmb{X}}}} \right)^{ - 1}}{x_i}。$ | (6) |
式(5)~(6)中，xi为第i种化合物分子结构描述符的行向量；X为n×m的矩阵，构成训练集化合物的描述符空间。
当训练集中化合物hi值大于警戒值(h*)时，说明在数据集中该物质的子结构出现较少，会对模型预测结果有显著影响。h*值计算公式为
| $ {h^*} = 3\left( {m + 1} \right)/n。$ | (7) |
模型描述符意义、回归系数、回归系数偏差和标准回归系数见表 1。构建的GA-MLR回归方程为Y=17.898-0.371X1+0.334X2+0.215X3+0.193X4+0.426X5-0.453X6-0.260X7-0.308X8，n训练集=107，Radj，训练集2=0.784，QLOO2=0.744，ERMS，训练集=1.127，P < 0.000 1，nEXT=45，REXT2=0.664，QEXT2=0.761，ERMS，EXT=1.039。
|
表 1 臭氧自由基反应速率常数模型描述符物理化学意义和相应系数 Table 1 Physical chemistry meaning and corresponding coefficients of the descriptors used in the -lg KO3 MLR model |
GOLBRAIKH等[20]研究认为，QSAR模型可接受标准为Q2>0.50和R2>0.60。由图 1可知，笔者构建的模型拟合优度和稳健性较好，预测能力也较好。
|
图 1 臭氧自由基反应速率常数MLR模型的预测值与实验值的拟合图 Fig. 1 Plot of experimental versus predictive of -lg KO3 values of MLR model |
利用杠杆方法制作Williams图分析和评价模型应用范围，可以图形方式量化模型应用范围。模型对应用域内物质预测性能较好，而对应用域外物质预测性能差。采用Williams图表征的QSAR模型应用域见图 2。
|
虚线为警戒值(h*=0.252)。 图 2 臭氧反应速率常数模型的应用域Williams图 Fig. 2 Williams plot of the application domain of the -lg KO3 MLR model |
由图 2可知，数据集152种化合物中只有肼的h值≥h*(h*=0.252)，位于应用域范围外，为X离群点。所有化合物标准残差在-3~3范围内，即无Y离群点。因此，构建的QSAR模型可用于预测应用域内其他化合物-lg KO3值。
2.3 有机化学品与臭氧反应机制解释通过解释线性化合物臭氧反应速率的QSAR模型中所选描述符的物理化学意义，可以获得决定化合物臭氧反应速率的结构信息。描述符的相对重要程度由模型中每个描述符的标准回归系数决定。标准回归系数绝对值的大小表示对应描述符对臭氧反应速率影响程度的强弱，正负号表示对应描述符与臭氧反应速率相关性的正负。在模型的8个描述符中，nR=Cs和nR=Ct均为官能团数目描述符，且标准回归系数均为负值(表 1)，这表明nR=Cs和nR=Ct与臭氧的反应速率常数呈正相关。HATS2e为GETAWAY(geometry，topology and atom-weights assembly)类描述符，与-lg KO3呈较大负相关。HOMA为几何描述符，PW3为拓扑描述符，可通过计算分子中每个原子的path数目与walk数目的比值，再将这些比值求和后除以分子中的原子数目得到。由于path/walk独立于分子大小，所以PW3可以较好地表征分子形状。RDF035u为径向分布函数描述符，表示在一个半径为R的球形体内发现特定类似原子的概率。G1s为WHIM描述符，在模型中表征分子静电拓扑状态。H-050为以原子为中心的碎片描述符，表征与杂原子相连的氢原子个数。
2.4 模型验证和应用案例近年来，计算毒理学技术在欧美、日本和OECD得到大力发展。美国国家环境保护局研发了化学品理化性质/环境行为指标参数与预测模型软件EPI Suite，其中的AOPWIN模块采用基团贡献法预测有机化学品臭氧自由基反应常数。OECD允许使用QSAR方法弥补数据缺失，并于2008年发布了第1版QSAR Toolbox工具包。其中的臭氧反应速率数据主要来源于EPI Suite软件。与发达国家相比，我国在计算毒理学技术研发和应用方面具有较大差距。近年来我国已经启动化学品环境安全信息预测技术研究，在一定程度上填补了我国化学品固有属性预测技术的空白。其中，生态环境部南京环境科学研究所基于简化分子线性输入规范(SMILES)解析碎片拆分技术，开发了具有我国自主知识产权的化学品定量构效预测软件[17]。
将笔者研究数据集之外的20种有机化学品-lg KO3实验数据[18-19]与该模型和EPI Suite软件中AOPWIN模块的预测结果进行比较发现，20种有机化学品-lg KO3实验值与笔者模型预测值的决定系数(R2)达到0.794，与EPI Suite预测值的R2为0.695(表 2)。其中该模型15种化学品预测结果优于EPI Suite，EPI Suite软件5种化学品预测结果较好。SMILES是化学物质1维结构的线性表达，而2维和3维结构描述符可更全面地表达化学物质立体结构的空间形态。由于EPI Suite软件基于SMILES码碎片拆分，选取的结构碎片也许不能完全表达分子结构信息，同时也未给出模型应用域，因此笔者构建的模型弥补了EPI Suite软件的不足[20-21]。
|
|
表 2 20种有机化学品-lg KO3实验和预测数据比较 Table 2 Comparison of predicted results with experimental results for -lg KO3 of 20 organic chemicals |
分别以丙烯酸甲酯〔SMILES:OC(OC)CC〕和四氟乙烯〔SMILES:FC(F)C(F)F〕为例，预测有机化学品-lg KO3值。使用Dragon软件计算丙烯酸甲酯PW3、HOMA、RDF035u、G1s、HATS2e、nR=Cs、nR=Ct和H-050 8种描述符所得结果分别为0.281、-0.022、1.479、0.218、0.905、1、0和0，hat值为0.106，小于h*(h*=0.252)，在模型应用域内。模型预测丙烯酸甲酯-lg KO3值为18.078，实验值为17.978，残差为-0.100。使用Dragon软件计算四氟乙烯PW3、HOMA、RDF035u、G1s、HATS2e、nR=Cs、nR=Ct和H-050 8种描述符所得结果分别为0.267、0、1.379、1、1.096、0、0和0，hat值为0.113，小于h*，在模型应用域内。模型预测丙四氟乙烯-lg KO3值为18.755，实验值为19.036，残差为0.281。
3 结论该研究建立了包括烷烃、烯烃、芳香烃、含氧挥发性有机物和酚类152种有机化合物与臭氧反应速率常数预测模型。根据经济合作与发展组织关于QSAR模型构建与验证的导则要求，构建的有机化学品与臭氧反应速率预测模型拟合能力、稳健性和预测能力均较好，Williams图定义模型应用域(AD)结果也表明该模型应用域较广。模型机理研究结果表明分子芳香性、电负性和仲碳原子数目是影响有机化学品与臭氧自由基反应速率(KO3)的关键因素。综上所述，构建的有机化合物与臭氧自由基反应速率常数QSAR模型可以用于预测应用域内难以测定或未知有机化合物与臭氧自由基反应速率常数，评估其持久性，进而对有机污染物进行生态风险评价。
| [1] |
COOPER O R, PARRISH D D, ZIEMKE J, et al. Global Distribution and Trends of Tropospheric Ozone:An Observation-Based Review[J]. Elementa:Science of the Anthropocene, 2014, 2: 000029. (  0) |
| [2] |
ATKINSON R. Gas-Phase Tropospheric Chemistry of Volatile Organic Compounds:1.Alkanes and Alkenes[J]. Journal of Physical and Chemical Reference Data, 1997, 26(2): 215-290. DOI:10.1063/1.556012 ( 0) |
| [3] |
MEYLAN W M, HOWARD P H. A Review of Quantitative Structure-Activity Relationship Methods for the Prediction of Atmospheric Oxidation of Organic Chemicals[J]. Environmental Toxicology and Chemistry, 2003, 22(8): 1724. ( 0) |
| [4] |
Organisation for Economic Co-operation and Development (OECD).Guidance Document on the Validation of (Quantitative) Structure-Activity Relationship[(Q)SAR]Models[R].Paris, France: OECD, 2007.
( 0) |
| [5] |
MEYLAN W M, HOWARD P H. Computer Estimation of the Atmospheric Gas-Phase Reaction Rate of Organic Compounds With Hydroxyl Radicals and Ozone[J]. Chemosphere, 1993, 26(12): 2293-2299. DOI:10.1016/0045-6535(93)90355-9 ( 0) |
| [6] |
GRAMATICA P, CONSONNI V, TODESCHINI R. QSAR Study on the Tropospheric Degradation of Organic Compounds[J]. Chemosphere, 1999, 38(6): 1371-1378. ( 0) |
| [7] |
REN Y Y, LIU H X, YAO X J, et al. Prediction of Ozone Tropospheric Degradation Rate Constants by Projection Pursuit Regression[J]. Analytica Chimica Acta, 2007, 589(1): 150-158. ( 0) |
| [8] |
GRAMATICA P, PILUTTI P, PAPA E. QSAR Prediction of Ozone Tropospheric Degradation[J]. QSAR & Combinatorial Science, 2003, 22(3): 364-373. ( 0) |
| [9] |
U. S. Environmental Protection Agency.On-line User Documentation Available of SuiteTM-Estimation Program Interface With the AOPWIN[R].Washington DC, USA: Atmospheric Oxidation Program, 2000.
( 0) |
| [10] |
YU X L, YI B, WANG X Y, et al. Predicting Reaction Rate Constants of Ozone With Organic Compounds From Radical Structures[J]. Atmospheric Environment, 2012, 51: 124-130. ( 0) |
| [11] |
KENNARD R W, STONE L A. Computer Aided Design of Experiments[J]. Technometrics, 1969, 11(1): 137-148. ( 0) |
| [12] |
SHAO Y H, LIU J N, WANG M X, et al. Integrated QSPR Models to Predict the Soil Sorption Coefficient for a Large Diverse Set of Compounds by Using Different Modeling Methods[J]. Atmospheric Environment, 2014, 88: 212-218. ( 0) |
| [13] |
MAURI A, CONSONNI V, PAVAN M, et al. DRAGON Software:An Easy Approach to Molecular Descriptor Calculations[J]. Match Communication in Mathematical and in Computer Chemistry, 2006, 56(2): 237-248. ( 0) |
| [14] |
LIU H X, GRAMATICA P. QSAR Study of Selective Ligands for the Thyroid Hormone Receptor Β[J]. Bioorganic & Medicinal Chemistry, 2007, 15(15): 5251-5261. ( 0) |
| [15] |
NETZEVA T I, WORTH A P, ALDENBERG T, et al. Current Status of Methods for Defining the Applicability Domain of (Quantitative) Structure-Activity Relationships[J]. Alternatives to Laboratory Animals, 2005, 33(2): 155-173. DOI:10.1177/026119290503300209 ( 0) |
| [16] |
GOLBRAIKH A, SHEN M, XIAO Z Y, et al. Rational Selection of Training and Test Sets for the Development of Validated QSAR Models[J]. Journal of Computer-Aided Molecular Design, 2003, 17(2/3/4): 241-253. ( 0) |
| [17] |
生态环境部南京环境科学研究所.化学品定量构效预测软件: 2015SR018599[CP].2014-10-13.
( 0) |
| [18] |
MCGILLEN M R, CAREY T J, ARCHIBALD A T, et al. Structure-Activity Relationship (SAR) for the Gas-phase Ozonolysis of Aliphatic Alkenes and Dialkenes[J]. Physical Chemistry Chemical Physics, 2008, 10(13): 1757-1768. DOI:10.1039/b715394e ( 0) |
| [19] |
GROSJEAN D, GROSJEAN E, WILLIAMS E L. Rate Constants for the Gas-Phase Reactions of Ozone With Unsaturated Alcohols, Esters, and Carbonyls[J]. International Journal of Chemical Kinetics, 1993, 25(9): 783-794. DOI:10.1002/kin.550250909 ( 0) |
| [20] |
范德玲, 刘济宁, 王蕾, 等. 羟基自由基反应常数定量预测模型[J]. 环境化学, 2015, 34(10): 1924-1931. [ FAN De-ling, LIU Ji-ning, WANG Lei, et al. QSAR Model for Predicting Hydroxyl Radical Reaction Constant of Organic Chemicals[J]. Environmental Chemistry, 2015, 34(10): 1924-1931. DOI:10.7524/j.issn.0254-6108.2015.10.2015032605] ( 0) |
| [21] |
范德玲, 宋波, 刘济宁, 等. 化学品正辛醇空气分配系数定量预测模型研究[J]. 生态与农村环境学报, 2015, 31(2): 269-272. [ FAN De-ling, SONG Bo, LIU Ji-ning, et al. A Quantitative Structure-Activity Relationship Model for Prediction of Octyl Alcohol/Air Partition Coefficient of Organic Chemicals[J]. Journal of Ecology and Rural Environment, 2015, 31(2): 269-272.] ( 0) |