农药对蚯蚓急性毒性的分类建模研究

楚留意，丛建业，刘 洋，魏宠芝，任月英

（兰州交通大学环境与市政工程学院，兰州 730070）

蚯蚓被誉为“生态系统的工程师”［1］，是土壤生态系统的重要组成部分［2］。蚯蚓在土壤改良、有机物质分解、土壤养分提高和促进农作物生长等方面发挥着至关重要的作用［3］。蚯蚓的数量可以作为评价土壤健康状况和衡量农业可持续发展的一个非常有效的指标［4-5］。因此，欧洲联盟和经济合作与发展组织选择蚯蚓作为化学品环境危险的非目标生物［6］。

中国作为全球农药生产和使用的大国，农药在预防生物虫害和农业生产安全保障方面起着非常重要的作用［2］，但是农药长期大规模的使用也对我们的生态环境产生了一定的负面影响［7］。中国农药管理相关部门一直关注着农药的使用对生态环境产生的影响，但是针对农药对鸟类、鱼类影响的研究较多，农药对蚯蚓等土壤动物的影响的相关研究还未受到重视［8］。近些年来，国内已经有一些科研工作者开展了关于农药对蚯蚓毒性的相关研究，主要针对的是杀虫剂和杀菌剂［9-10］。蚯蚓在土壤中能通过皮肤接触或者食用被污染的食物而直接接触到农药［11-12］，因此导致蚯蚓死亡或者影响其生长［13］。土壤中残留的农药会通过不同的途径和毒性机制产生毒性作用，在蚯蚓毒性测试实验中最大的问题就是确定哪种化合物会对蚯蚓产生致死效果。但是，要想对市面上的数10万种化学品逐一进行毒性测试是非常耗时的。利用定量构效关系（quantitative structure-activity relationship，QSAR）方法建模来预测化合物的毒性则具有一定的优势，但该方面对于蚯蚓毒性的研究则报道较少。目前可见文献报道的两例相关研究均来自于印度贾达普大学的乔伊教授等人，针对化合物对蚯蚓的毒性，分别建立了一个分类模型［14］和回归模型［15］。分类模型中共有163种物质，用随机分组方法将117种物质作为训练集，建立了包含8 个结构参数的线性判别方法（linear discriminant analysis，LDA）模型。该模型可以正确划分83.96%的训练集物质，对测试集数据的划分正确率为71.74%，但是对于测试集数据的灵敏度较低（仅为58.82%）。

在本研究中，使用DRAGON 描述符［16］来表征化合物的结构信息。考虑到生物毒性／活性与结构参数之间的关系比较复杂，用非线性方法一般更能准确地表述这种关系。因此，在本文中除采用LDA选择结构参数重新建模外，我们还尝试采用SVM来建立非线性模型用于农药对蚯蚓毒性的分类研究，并将二者结果进行对比。本研究的目的是建立一个准确可靠的对蚯蚓毒性物质分类的QSAR 模型，同时分析影响此类物质对蚯蚓毒性作用的重要结构信息特征。此外，采用主成分分析（principle component analysis，PCA）方法进行数据分组，使得训练集样本具有更广泛的代表性。与前人的工作相比，本研究在数据集分组以及建模方法上均有改进。

1 材料与方法

1.1 实验数据

本研究中数据来源于文献［14］，共包含163种物质（化合物编号6因结构信息不足摈弃不用），其种类涵盖了烷烃、烯烃、醇、醛、酮、酰胺、酚以及酯等。毒性数据用LC50（mg／kg）表示；根据OECD关于“蚯蚓急性毒性实验”［17］，把实验浓度为1 000 mg／kg作为极限实验浓度。故本文中物质毒性类别与文献［14］相同，将163种物质划分为两类。其中“1”代表有毒化合物（LC50＜1 000 mg／kg）共106个，“0”代表无毒化合物（LC50＞1 000 mg／kg），共57个。

在QSAR 的研究中最关键的一步是训练集和测试集的划分，因为所建立模型的泛化能力依赖于数据集的划分［18］。主成分分析（PCA）是一种基于降维思想产生的高维数据处理方法［19］。本文应用意大利Talete公司研发的Dragon软件［20］计算得到的分子结构描述符对全部化合物作了主成分分析。取本征值大于1时，前3个主成分的总方差达到了51.78%，各自的方差贡献率分别为38.3%、8.316%和5.164%。根据主成分分析的结果，我们挑选了131个化合物作为训练集，用来调节模型的结构和训练参数；剩余的32个化合物作为测试集，用于测试模型的泛化能力。从图1中可以看出样本分布比较均匀，训练集中的样本能够代表测试集中样本的结构特点和信息，说明训练集和测试集的数据划分是合理的。

图1 主成分分析载荷图Fig.1 Principal component analysis load diagram

1.2 方法

线性判别分析是一种模式判别方法，LDA 理论可参考相关文献［21－22］。

支持向量机由数学家Vapnik于1995年首先提出来［23-24］，是一类可用于回归和分类的有监督的机器学习模型。它在解决非线性、小样本及高维模式识别等问题中具有独特的优势［25］，如今SVM 已被广泛使用于函数拟合等其他机器学习问题中［26］。关于SVM 算法理论可参考相关文献［26-27］，这里不再赘述。本研究中SVM 算法由R 语言实现。

2 结果与讨论

2.1 LDA模型结果

逐步线性判别分析得到了包含12个描述符的最佳模型（见表1）。其统计参数如下：n＝131，AWilk’s-λ＝0.393，BEigenvalue＝1.546，OCanonicalCorr＝0.779，SChisquare＝114.932，p＜0.000 1。从表1中最右一列可知，所有描述符的方差膨胀因子（variance inflation factor，VIF）值都小于5，说明它们之间是相互独立的。这些描述符都能反映了不同的分子结构特征，其具体含义见表1。

本研究中涉及到的分子结构呈多样化，所以化合物对蚯蚓的毒性作用与分子结构之间的关系就显得较为复杂。在所选择的12个描述符中，Sp MAD_L（spectral mean absolute deviation from Laplace matrix）属于二维矩阵描述符，NssS（number of atoms of type ssS）是原子型电子状态指数描述符，表述分子中含-S-片段的个数［20］。F03［O-O］（frequency of O-O at topological distance 3）、B04［CO］（presence／absence of C-O at topological distance 4）和F09［F-Cl］（frequency of F-Cl at topological distance 9）是二维原子对描述符，表示分子结构中特定键或结构片段对毒性的影响［20］。F03［OO］与F09［F-Cl］分别指在拓扑距离3对应的O-O结构片段的频率和在拓扑距离9对应的F-Cl结构片段的频率，这些结构片段的存在意味着会增强电负性从而物质表现的毒性更强，这一点从其在模型中系数为正可以解释。GGI5（topological charge index of order 5）、GATS8i（Geary autocorrelation of lag 8 weighted by ionization potential）、JGI2（mean topological charge index of order 2）和MATS8i（Moran autocorrelation of lag 8 weighted by ionization potential）是二维自相关描述符，由原子物理化学性质加权化的含氢拓扑分子结构计算而得，均和分子内部的电子信息和分子几何形状有关［20］。对于GATS8i而言，值越小表明存在较强的空间自相关，所以取值在0～1之间意味着是正相关，而大于1时意味着负相关，值为零则不存在相关。而对MATS8i则情况相反，该参数值为正意味着存在正相关，取值为负则为负相关。GGI5等于分子中5级最大拓扑距离值的电荷指数，而JGI2则由GGI5除以分子中对应级数总电荷指数计算得到。MEcc（molecular eccentricity）属于几何描述符，表述分子形状信息；Mor14s（signal 14／weighted by I-state）和Mor13p（signal 13／weighted by polarizability）是3D-MoRSE 描述符（3D-molecule representation of structures based on electron diffraction），这些描述符都是剔除了氢原子，考虑其他原子对分子几何结构的影响，包括极性和离子势的作用［20］。综上所述，这些描述符主要包含了分子几何形状、电荷分布特征以及极性等多方面的信息，其中对于分子毒性影响最大的是GATS8i和GGI5。

该模型对于训练集的准确率可达91.60%，无毒物质有81个被正确区分，有毒物质有39个被正确区分。模型对于测试集准确率为75.00%，无毒物质有19个被正确区分，有毒物质有5个被正确区分；统计结果见表2。

表2 LDA和SVM 模型的统计验证参数Tab.2 Statistical validation parameters of LDA and SVM models

2.2 SVM 结果

SVM 模型的结果取决于容量因子C，核函数K以及γ。但是要如何选择最优的参数，目前还没有明确的理论指导。本研究中采用留一法交互检验的结果确定最优模型的参数为C＝334，Ggamma＝0.036。该模型对训练集的准确率可达96.18%，其中83个无毒物质被正确区分，43个有毒物质被正确区分；对测试集的准确率达到了81.25%，其中20个无毒物质被正确区分，6个有毒物质被正确区分。

2.3 LDA和SVM 结果比较

上述可知，SVM 模型的预测结果要明显优于LDA 模型。比较两模型的ROC 曲线（见图2）可以看出，LDA 模型训练集和测试集的曲线下的面积（AUC）分别为0.964 和0.705；而SVM 模型则分别为0.991和0.732，均优于LDA 模型。这意味着有机物对蚯蚓毒性作用的因素是复杂的，并非都与毒性作用成线性相关；而应用同样的描述符，非线性模型更能准确地识别出毒性作用与结构因素之间的关系。

图2 LDA和SVM 模型的ROC曲线图Fig.2 ROC curve of LDA and SVM models

此外我们还对两个模型对有毒物质和无毒物质的预测结果做了比较。SVM 模型对无毒物质和有毒物质预测的准确率分别是92.79% 和94.23%，均高于LDA模型的88.50%和88.00%；SVM 和LDA 模型的结果均优于文献中的80.19% 和80.70%。

在评价物质毒性危害的时候，错误的否定要比错误的肯定更难令人接受，因前者会直接导致危害，而后者起到了警示作用［22］。SVM 模型对有毒物质预测的准确率达到了94.23%，高于LDA 的88.00%和文献中的80.70%。在这种意义上说，我们在本文中建立的模型有利于风险评估。

两个模型的敏感度、特异性、灵敏度等见表2。

3 结论

本文基于分子结构描述符，应用线性判别分析（LDA）和支持向量机（SVM）对163种有机物做了蚯蚓急性毒性的分类研究，所得两模型的预测结果都是令人满意的。与前人工作相比较，本文所建立的LDA 模型结果有所提高，而SVM 模型的结果更优于LDA 模型结果。这说明化合物对蚯蚓的毒性作用是复杂的，非线性方法更能准确地描述化合物的结构因素与化合物毒性之间的关系。