有機(jī)化合物生物富集因子的計(jì)算機(jī)預(yù)測(cè)研究

孫露，陳英杰，吳曾睿，李衛(wèi)華，劉桂霞，Philip W. Lee，唐赟

華東理工大學(xué)藥學(xué)院 上海市新藥設(shè)計(jì)重點(diǎn)實(shí)驗(yàn)室，上海 200237

有機(jī)化合物生物富集因子的計(jì)算機(jī)預(yù)測(cè)研究

孫露，陳英杰，吳曾睿，李衛(wèi)華，劉桂霞，Philip W. Lee，唐赟*

華東理工大學(xué)藥學(xué)院 上海市新藥設(shè)計(jì)重點(diǎn)實(shí)驗(yàn)室，上海 200237

有機(jī)化合物在生物體內(nèi)的富集，通常用生物富集因子(bioconcentration factor，簡(jiǎn)稱BCF)來(lái)表達(dá)，這是化合物生態(tài)環(huán)境毒性評(píng)估的重要指標(biāo)。為合理預(yù)測(cè)有機(jī)化合物是否易于生物富集，首先從美國(guó)環(huán)保局網(wǎng)站收集了624個(gè)具有不同BCF值的化合物，然后采用7種分子指紋結(jié)合5種機(jī)器學(xué)習(xí)方法(包括支持向量機(jī)、C4.5決策樹(shù)、k最近鄰法、隨機(jī)森林法和樸素貝葉斯法)，構(gòu)建了化合物BCF的分類預(yù)測(cè)模型，所有模型均采用獨(dú)立外部驗(yàn)證集進(jìn)行驗(yàn)證。其中，使用ChemoTyper分子指紋結(jié)合支持向量機(jī)方法得到的二分類模型，整體預(yù)測(cè)準(zhǔn)確度最好，達(dá)到了85.4%。通過(guò)采用信息增益、頻率分析等方法，進(jìn)一步確定了化合物中易于引起生物富集的關(guān)鍵子結(jié)構(gòu)，包括芳基氯、二芳基醚、氯代烷烴等。研究中所用到的方法為有毒化學(xué)品的生態(tài)風(fēng)險(xiǎn)評(píng)價(jià)提供了良好可靠的預(yù)測(cè)工具。

生物富集因子；計(jì)算機(jī)預(yù)測(cè)；二分類模型；警示子結(jié)構(gòu)；環(huán)境毒理學(xué)

生物從周?chē)h(huán)境中吸收并累積某種元素或難分解的化合物，從而導(dǎo)致生物體內(nèi)該物質(zhì)的濃度超過(guò)環(huán)境中濃度的現(xiàn)象叫做生物富集。生物富集常用生物富集因子[1](bioconcentration factor，簡(jiǎn)稱BCF)來(lái)表達(dá)。生物富集因子是化學(xué)品在水生生物和水體之間的平衡分配過(guò)程，其在生物體內(nèi)(平衡)濃度與其水體中(平衡)濃度的比值稱為生物富集因子，它反映水生生物對(duì)水體中有機(jī)物的吸收儲(chǔ)存能力，是評(píng)價(jià)有機(jī)污染物生物累積性的重要指標(biāo)。無(wú)論是確定持久性有機(jī)污染物(POPs)，還是確定持久生物累積性有毒(PBT)污染物的清單，BCF都是一個(gè)不可或缺的參數(shù)。

BCF具有標(biāo)準(zhǔn)的測(cè)試方法[2]，但實(shí)驗(yàn)測(cè)定BCF成本高、周期長(zhǎng)，以及人力、經(jīng)費(fèi)、時(shí)間的限制，不可能對(duì)進(jìn)入環(huán)境中的每種化學(xué)品都進(jìn)行生物富集實(shí)驗(yàn)測(cè)定。因此，對(duì)化學(xué)品的生物富集因子進(jìn)行計(jì)算機(jī)預(yù)測(cè)研究便顯得尤為重要。迄今只有有限的有機(jī)化合物具有實(shí)測(cè)的BCF值。在環(huán)境科學(xué)技術(shù)中，定量結(jié)構(gòu)-活性關(guān)系(QSAR)[3]指關(guān)聯(lián)有機(jī)污染物的分子結(jié)構(gòu)與其理化性質(zhì)、環(huán)境行為和毒理學(xué)參數(shù)(統(tǒng)稱為活性)的定量預(yù)測(cè)模型。QSAR可以彌補(bǔ)基礎(chǔ)數(shù)據(jù)的缺失，降低昂貴的測(cè)試費(fèi)用，減少動(dòng)物實(shí)驗(yàn)。

在生態(tài)的角度上BCF的實(shí)驗(yàn)數(shù)據(jù)是很重要的，從監(jiān)管的角度上來(lái)說(shuō)也是很重要的。歐盟法規(guī)《化學(xué)品的注冊(cè)、評(píng)估、授權(quán)和限制》(Registration, Evaluation, Authorization and Restriction of Chemicals, REACH)[4]要求每個(gè)化合物都要有與之對(duì)應(yīng)的BCF值，從而推動(dòng)了有機(jī)化合物BCF預(yù)測(cè)模型的建立。目前文獻(xiàn)中已經(jīng)報(bào)道了很多BCF預(yù)測(cè)模型，從最簡(jiǎn)單的log BCF與log KOW的線性模型[5-7]開(kāi)始，有基于2D分子描述符的QSAR模型[8]、基于基線的BCF模型[9]、基于代謝的動(dòng)力學(xué)模型[10]、基于拓?fù)渲笖?shù)的模型[11]、基于分子電性距離矢量(Molecular Electronegativity Distance Vector, MEDV)的模型[12]以及混合模型[13-14](幾組模型的組合)等。同時(shí)也有很多預(yù)測(cè)模型整合到商業(yè)或者免費(fèi)的軟件中，例如ACD Labs[15]、VEGA[16]、OECD QSAR ToolBox[17]以及EPI suit[18]等。

本研究首先收集具有實(shí)驗(yàn)測(cè)定BCF值的化合物，然后使用7種不同的分子指紋結(jié)合5種機(jī)器學(xué)習(xí)方法來(lái)構(gòu)建化合物BCF的定性分類預(yù)測(cè)模型，并使用多種方法來(lái)識(shí)別導(dǎo)致化合物發(fā)生生物富集的關(guān)鍵子結(jié)構(gòu)作為預(yù)警。本研究所用方法也可用于環(huán)境毒理學(xué)中的其它毒性端點(diǎn)預(yù)測(cè)評(píng)價(jià)。

1 材料與方法(Materials and methods)

1.1 數(shù)據(jù)集準(zhǔn)備

以美國(guó)環(huán)保局Estimation Program Interface (EPI) Suite[18]的BCFBAF程序中的Non-Ionic Training set[19]作為訓(xùn)練集，validation set作為外部驗(yàn)證集。其中訓(xùn)練集中有466個(gè)化合物，外部驗(yàn)證集中有158個(gè)化合物。依照生物富集的分類標(biāo)準(zhǔn)[20]，將化合物分為低生物富集、中度生物富集以及高度生物富集。數(shù)據(jù)集的統(tǒng)計(jì)結(jié)果如表1所示。由于高生物富集化合物較少，因此在建模時(shí)將高生物富集和中度生物富集歸為一類“易于生物富集”，將低生物富集歸為一類“不易生物富集”，構(gòu)建二分類定性預(yù)測(cè)模型。

表1 數(shù)據(jù)集分布Table 1 Data sets and chemical toxicity categories

1.2 分子指紋計(jì)算

本研究使用7種分子指紋來(lái)表達(dá)化合物的分子結(jié)構(gòu)，其中6種分子指紋分別是Fingerprint (FP，1024位)、EState fingerprint (EStFP，79位)、Extende fingerprint (ExtFP，1024位)、MACCS keys (MACCS，166位)、PubChem fingerprint (PubFP，881位)和Substructure fingerprint (FP4，307位)，使用PaDEL-Descriptor[21]軟件計(jì)算得到。第7種分子指紋使用ChemoTyper (CT，729位)軟件[22]計(jì)算得到。在計(jì)算分子指紋之前數(shù)據(jù)集中所有化合物的SMILES先通過(guò)ChemAxon Standardizer[23]處理，其中的設(shè)置參數(shù)如下[24]：add explicit hydrogens, aromatize, clean 2D, remove fragment。

1.3 建模方法

本研究使用5種機(jī)器學(xué)習(xí)方法構(gòu)建二分類模型。這5種機(jī)器學(xué)習(xí)方法分別是支持向量機(jī)(Support Vector Machine, SVM)，最近鄰居法(k-Nearest Neighbors, k-NN)，樸素貝葉斯分類器(Na?ve Bayes, NB)，隨機(jī)森林(Random Forest, RT)，決策樹(shù)(Decision Tree, C4.5)。這些方法均在Orange[25]軟件中實(shí)現(xiàn)，Orange是一個(gè)基于Python腳本的數(shù)據(jù)挖掘和機(jī)器學(xué)習(xí)軟件套裝。為了測(cè)試模型的有效性和準(zhǔn)確性，首先對(duì)訓(xùn)練集采用10倍交叉驗(yàn)證，檢驗(yàn)?zāi)Ｐ偷聂敯粜?；然后進(jìn)行外部驗(yàn)證集驗(yàn)證，檢驗(yàn)?zāi)Ｐ偷念A(yù)測(cè)準(zhǔn)確性。

支持向量機(jī)(SVM)：本算法[26]是1995年由Vapnik和Cortes提出的一種統(tǒng)計(jì)學(xué)習(xí)算法。支持向量機(jī)通過(guò)對(duì)輸入樣本空間進(jìn)行非線性映射轉(zhuǎn)換，將輸入空間變換為一種高維空間，從而利用線性分類平面來(lái)描述非線性的分類邊界。這種非線性變換是通過(guò)核函數(shù)(Kernel Functions)實(shí)現(xiàn)的，在本研究中我們選取的核函數(shù)是RBF(Radial Basis Function)核函數(shù)，RBF核函數(shù)的優(yōu)點(diǎn)包括可以將樣本非線性的映射到高維空間，從而處理非線性問(wèn)題，另外其只含有一個(gè)參數(shù)，形式簡(jiǎn)單。在Orange軟件中，c值設(shè)定為10，g值設(shè)定為0.00212，且不勾選“Normalize data”選項(xiàng)。

最近鄰居法(k-NN)：又稱k-近鄰法，是著名的模式識(shí)別和統(tǒng)計(jì)學(xué)習(xí)方法之一[27]，被廣泛應(yīng)用于文本分類、模式識(shí)別、圖形圖像以及空間分布等領(lǐng)域。該方法主要根據(jù)特征空間(描述空間)中最接近的樣本進(jìn)行分類，其基本思想是首先在多維向量中尋找與待分類樣本最接近的k個(gè)鄰居，然后根據(jù)這k個(gè)臨近點(diǎn)的類別決定待分類樣本所屬的種類。本研究中選取的k值為9，距離公式選擇歐幾里得距離，同時(shí)選取了距離權(quán)重這個(gè)選項(xiàng)。

樸素貝葉斯分類器(NB)：本方法[28]同樣在化學(xué)分類模型中有著廣泛的應(yīng)用。它是貝葉斯分類器中簡(jiǎn)單而有效的一種，通過(guò)計(jì)算樣本屬于不同種類的概率，具有最大概率的類便是該樣本所屬的類。樸素貝葉斯分類器的優(yōu)點(diǎn)是分類過(guò)程中占用的計(jì)算資源很少，分類結(jié)果很穩(wěn)定，魯棒性很好。在Orange軟件中使用默認(rèn)參數(shù)設(shè)置。

決策樹(shù)(C4.5)：本方法在眾多的模式識(shí)別方法中是最經(jīng)典和最古老的方法之一。C4.5是一種以實(shí)例為基礎(chǔ)的歸納學(xué)習(xí)算法，從無(wú)序的訓(xùn)練樣本中，歸納出分類的標(biāo)準(zhǔn)，其基本構(gòu)成包括決策結(jié)點(diǎn)，分支和葉結(jié)點(diǎn)。決策樹(shù)的目標(biāo)是根據(jù)簡(jiǎn)單的幾個(gè)變量(描述符)輸入建立一種簡(jiǎn)單的規(guī)則預(yù)測(cè)一個(gè)目標(biāo)值。樸素貝葉斯算法在Orange軟件中使用的是默認(rèn)參數(shù)設(shè)置。

隨機(jī)森林(RF)：本方法是一個(gè)包含多個(gè)決策樹(shù)的分類器，并且其輸出的類別是由個(gè)別樹(shù)輸出的類別的眾數(shù)而定。在Orange軟件中，number of trees in forest參數(shù)設(shè)置為100。

1.4 預(yù)測(cè)模型評(píng)價(jià)

所有的模型均采用10倍交叉驗(yàn)證(10-fold cross validation)和外部驗(yàn)證集驗(yàn)證。通過(guò)計(jì)算如下變量的化合物數(shù)目：真陽(yáng)性(True Positives, TP)、真陰性(True Negatives, TN)、假陽(yáng)性(False Positives, FP)和假陰性(False Negatives, FN)。然后計(jì)算每個(gè)模型的敏感性、特異性和預(yù)測(cè)準(zhǔn)確率。敏感性(Sensitivity, SE)表示的是對(duì)中/高生物富集化合物預(yù)測(cè)正確的百分比；特異性(Specificity, SP)表示的是對(duì)低生物富集化合物預(yù)測(cè)正確的百分比；整體預(yù)測(cè)正確率(Accuracy, Q)表示的是中/高生物富集化合物和低生物富集化合物都預(yù)測(cè)正確的百分比。SE、SP和Q值均在0到1之間，值越高，表明模型的預(yù)測(cè)能力越好。馬修斯相關(guān)系數(shù)(MCC)代表模型的整體預(yù)測(cè)能力，MCC的范圍處于-1與+1之間。當(dāng)某個(gè)預(yù)測(cè)模型的MCC值等于+1時(shí)，這就代表著這個(gè)模型的預(yù)測(cè)能力達(dá)到了最優(yōu)水平；相反，如果某個(gè)預(yù)測(cè)模型的MCC值等于-1，那么它得到的則是最差的預(yù)測(cè)結(jié)果。具體的計(jì)算方程如下[29]：

另外，受試工作特征曲線(ROC曲線)也被用來(lái)表征模型的預(yù)測(cè)能力[30]。ROC曲線能形象地呈現(xiàn)出模型的預(yù)測(cè)能力，一般曲線越往上偏，表明模型預(yù)測(cè)的精準(zhǔn)度越高。工作特性曲線下面積(AUC)值越大，表示模型預(yù)測(cè)能力越高。最理想的模型AUC=1，隨機(jī)產(chǎn)生的模型AUC=0.5。

1.5 警示子結(jié)構(gòu)識(shí)別

本研究使用信息增益(Information gain，簡(jiǎn)稱IG)[31]和子結(jié)構(gòu)片段頻率分析[32]來(lái)獲得警示子結(jié)構(gòu)。如果一個(gè)子結(jié)構(gòu)在高/中生物富集化合物中出現(xiàn)的頻率比低生物富集化合物中出現(xiàn)的頻率更高，這個(gè)子結(jié)構(gòu)就可以被認(rèn)為是易于產(chǎn)生生物富集的特征子結(jié)構(gòu)，需要警惕[33]。警示子結(jié)構(gòu)是由毒性機(jī)理總結(jié)而來(lái)，是一種很重要的預(yù)測(cè)工具[34]。生物富集化合物中，片段出現(xiàn)的頻率定義為如下公式計(jì)算：

其中，N子結(jié)構(gòu)類別表示的是在中/高生物富集化合物中包含有此片段的總數(shù)目，N總數(shù)表示數(shù)據(jù)集中所有化合物數(shù)目，N子結(jié)構(gòu)總數(shù)是表示包含該子結(jié)構(gòu)片段的所有化合物數(shù)目，N類別是表示在中/高生物富集化合物的數(shù)目。

我們也使用了免費(fèi)軟件KNIME (Konstanz Information Miner)[35]中的MoSS模塊和ChemoTyper軟件搜尋一系列分子數(shù)據(jù)集中頻繁出現(xiàn)的分子結(jié)構(gòu)片段。KNIME的MoSS模塊中，“minimum fragment size”是一個(gè)重要的參數(shù)，經(jīng)過(guò)測(cè)試比較，我們把這個(gè)值設(shè)定為8；“maximum fragment size”設(shè)定為100，“minimum focus support in %”設(shè)定為8，“maximum complement support in %”設(shè)定為20，其他參數(shù)使用默認(rèn)設(shè)置。ChemoTyper軟件可以根據(jù)化合物的結(jié)構(gòu)信息搜索警示子結(jié)構(gòu)。

2 結(jié)果(Results)

2.1 數(shù)據(jù)集分析

影響模型質(zhì)量好壞的關(guān)鍵因素是數(shù)據(jù)集的質(zhì)量，本研究中使用的是BCFBAF軟件中的數(shù)據(jù)集。通過(guò)化合物的分子量(molecule weight)和ALogP來(lái)定義數(shù)據(jù)集(包含訓(xùn)練集以及外部驗(yàn)證集)學(xué)空間分布，如圖1所示。從圖中可以看出，外部驗(yàn)證集化合物的化學(xué)空間與訓(xùn)練集分子的化學(xué)空間類似，分布在同一區(qū)域內(nèi)。為了進(jìn)一步研究數(shù)據(jù)集的應(yīng)用域，我們另使用5種物理描述符(AlogP、分子質(zhì)量、溶解度、氫鍵受體數(shù)目和氫鍵供體數(shù)目)雷達(dá)圖[36]來(lái)定義(圖2)。如圖2中所示，分子質(zhì)量的最小值為68.074，最大值為959.168；AlogP的最小值為-1.661，最大值為到10.874；溶解度的最小值為-15.489，最大值為1.094；氫鍵受體數(shù)目為0到9；氫鍵供體數(shù)目為0到3，這些數(shù)據(jù)說(shuō)明我們的數(shù)據(jù)集具有較大的應(yīng)用域。

2.2 模型構(gòu)建

本研究使用7種分子指紋結(jié)合5種機(jī)器學(xué)習(xí)方法構(gòu)建了35個(gè)二分類模型。通過(guò)分析訓(xùn)練集10倍交叉驗(yàn)證的結(jié)果(表2)發(fā)現(xiàn)，5種機(jī)器學(xué)習(xí)方法在預(yù)測(cè)生物富集方面存在一定的差異性，使用不同的分子指紋的建模結(jié)果也略有差異。其中使用SVM和k-NN這兩種機(jī)器學(xué)習(xí)方法時(shí)，模型的魯棒性較好。例如，在MACCS-kNN和MACCS-SVM模型中，Q值分別為0.856和0.850，AUC值分別為0.869和0.884，MCC值分別為0.586和0.546，均高于使用其他機(jī)器學(xué)習(xí)方法建立模型的結(jié)果。

2.3 外部驗(yàn)證集驗(yàn)證

通過(guò)分析訓(xùn)練集的10倍交叉驗(yàn)證結(jié)果發(fā)現(xiàn)使用SVM和k-NN機(jī)器學(xué)習(xí)方法時(shí)模型的魯棒性較好，因此對(duì)這兩種機(jī)器學(xué)習(xí)方法構(gòu)建的模型進(jìn)行外部驗(yàn)證集驗(yàn)證(見(jiàn)表3)。在這7種分子指紋中，使用CT、ExtFP、PubFP和MACCS這四種分子指紋描述分子特征時(shí)模型的預(yù)測(cè)結(jié)果較好，即模型CT-SVM、ExtFP-SVM和PubFP-SVM具有最好的預(yù)測(cè)準(zhǔn)確度。在這四個(gè)模型中，Q值分別為0.854、0.842、0.835和0.810，SE值分別為0.9、0.918、0.918和0.882，SP值分別為0.75、0.667、0.646和0.646；AUC值分別為0.910、0.910、0.911和0.890(圖3)，MCC值分別為0.654、0.614、0.597和0.541。通過(guò)對(duì)比分析Q、SE、SP、AUC以及MCC值，發(fā)現(xiàn)CT-SVM模型為最優(yōu)模型。

2.4 警示子結(jié)構(gòu)識(shí)別

通過(guò)信息增益分析，我們找到了10個(gè)潛在的具有致生物富集的警示子結(jié)構(gòu)，分別是芳基氯化物(arylchloride)、二芳基醚(diarylether)、芳基溴化物(arylbromide)、鹵素縮醛類似物(halogen acetal like)、氯代烷烴(chloroalkene)、稠環(huán)(annelated rings)、羧酸酯(carboxylic ester)、叔碳(tertiary carbon)、橋環(huán)(bridged rings)、仲碳(secondary carbon)。這些子結(jié)構(gòu)的命名均基于FP4命名規(guī)則，具體結(jié)構(gòu)及其信息增益值見(jiàn)表4。

圖1 訓(xùn)練集化合物和外部驗(yàn)證集化合物的空間分布Fig. 1 Diversity analysis of chemicals in the training set and validation set

圖2 5個(gè)物理描述符的雷達(dá)分析圖 Fig. 2 The radar chart of five physicochemical descriptors (AlogP, Molecular Weight, Solubility, H-Acceptors and H-Donors)

圖3 模型CT-SVM、ExtFP-SVM、MACCS-SVM和 PubFP-SVM的受試者工作特性曲線圖Fig. 3 Representation of receiver operating characteristics (ROC) curve for the validation set in model CT-SVM, ExtFP-SVM, MACCS-SVM and PubFP-SVM

表2 訓(xùn)練集的10倍交叉驗(yàn)證結(jié)果Table 2 The performance of 10-fold cross validation in training set

表3 SVM和k-NN構(gòu)建模型的外部驗(yàn)證集結(jié)果Table 3 The performance of models using SVM and k-NN methods for validation set

表4 警示子結(jié)構(gòu)的頻率分析和信息增益結(jié)果Table 4 The common substructure alerts identified using IG methods and frequency analysis

通過(guò)KNIME軟件的MoSS模塊檢索生物富集化合物結(jié)構(gòu)，得到了10個(gè)子結(jié)構(gòu)，分別為1-氯-3-甲基苯、1,3-二氯苯、異丙苯、1,2,4-三氯苯、1,2,3-三氯苯、1,2-二氯苯、硝基苯、苯甲醚、乙苯和1,4-二氯苯。上述子結(jié)構(gòu)及包含對(duì)應(yīng)結(jié)構(gòu)的代表性化合物分子見(jiàn)表5。從表5中可以看出，這10個(gè)子結(jié)構(gòu)在生物富集化合物中出現(xiàn)的頻率最大，同樣給了我們警示的作用。

使用ChemoTyper軟件根據(jù)數(shù)據(jù)集的結(jié)構(gòu)信息搜尋數(shù)據(jù)集中的警示子結(jié)構(gòu)，找到10個(gè)警示子結(jié)構(gòu)，見(jiàn)表6。表6中“X”和“？”均表示鹵素原子，當(dāng)“X”和“？”為氯原子時(shí)，這10個(gè)警示子結(jié)構(gòu)為1,4-二氯苯、1,2,4-三氯苯、1,2,3-三氯苯、3-氯酚、二苯基甲烷、1,3,5-三氯苯、1-氯-2-(3-氯苯氧基)苯、1,2,3-三氯、1,2-二氯乙烯和2-氯-2-甲基丙烷。這些子結(jié)構(gòu)中，有7個(gè)為芳香族結(jié)構(gòu)片段，且在苯環(huán)上均連接有鹵素原子；3個(gè)為烷烴類結(jié)構(gòu)，均為鹵代烷烴。

表5 MoSS分析結(jié)果以及包含其警示結(jié)構(gòu)的代表性化合物Table 5 The results of MoSS searched for structural alerts and representative structures

3 討論(Discussion)

3.1 模型結(jié)果分析

在我們的研究中，使用了5種不同的機(jī)器學(xué)習(xí)方法(SVM、C4.5、RF、k-NN和NB)。通過(guò)分析表2中的10倍交叉驗(yàn)證結(jié)果(Q、SE、SP、AUC以及MCC值)，可以看出整體趨勢(shì)上，在使用同一分子指紋描述分子特性時(shí)，SVM和k-NN兩種機(jī)器學(xué)習(xí)方法顯示了良好的預(yù)測(cè)精度。MCC值代表模型的整體預(yù)測(cè)能力，在使用SVM和k-NN建模時(shí)，模型的MCC值要明顯高于其他機(jī)器學(xué)習(xí)方法。例如模型CT-kNN和CT-SVM的MCC值分別為0.512和0.493，明顯優(yōu)于其他機(jī)器學(xué)習(xí)方法。

SVM具有很強(qiáng)的擬合非線性關(guān)系的能力，并在一定程度上成為預(yù)測(cè)準(zhǔn)確度的“黃金標(biāo)準(zhǔn)”。k-NN算法之所以在10倍交叉驗(yàn)證預(yù)測(cè)準(zhǔn)確率較高是由其算法的特殊性和生物富集這個(gè)特殊的毒性端點(diǎn)所決定的。一個(gè)化合物被預(yù)測(cè)為易于富集化合物還是不易富集化合物，主要根據(jù)它附近的鄰居化合物的富集與否所決定。被分為相同類的化合物之間結(jié)構(gòu)具有相似性。由于這個(gè)因素，如果數(shù)據(jù)庫(kù)包含的化合物數(shù)量足夠大和化合物結(jié)構(gòu)足夠多樣，那么以k-NN作為建模方法建立的模型去預(yù)測(cè)化合物的生物富集因子，就能夠獲得很高的預(yù)測(cè)準(zhǔn)確度。

表6 ChemoTyper分析結(jié)果Table 6 The results of ChemoTyper analysis searched for structural alerts

注：“X”和“？”均表示鹵素。

Note: “X” and “?” means halogen.

當(dāng)使用相同的建模方法，不同的分子指紋作為屬性變量的建模結(jié)果略有差異。如表2中所示，當(dāng)使用k-NN建模時(shí)，模型MACCS-kNN、ExtFP-kNN、EStFP-kNN、FP-kNN、FP4-kNN、PubFP-kNN和CT-kNN的Q值分別為0.856、0.803、0.794、0.824、0.779、0.845和0.824。使用MACCS、ExtFP、FP、PubFP和CT這五個(gè)分子指紋的建模結(jié)果要優(yōu)于EStFP、FP4。EstFP分子指紋的長(zhǎng)度為79，在這79個(gè)分子片段中僅有35個(gè)片段用來(lái)描述模型的結(jié)構(gòu)信息；FP4分子指紋中含有307個(gè)片段，但用來(lái)描述模型結(jié)構(gòu)信息的片段僅有93位。大量的信息丟失是導(dǎo)致預(yù)測(cè)準(zhǔn)確度低的重要因素，可能是其不具備優(yōu)秀的特征來(lái)表征多樣數(shù)據(jù)集中的分子結(jié)構(gòu)特性。

3.2 與EPI中的預(yù)測(cè)軟件進(jìn)行比較

BCFBAF是由EPA研究開(kāi)發(fā)預(yù)測(cè)BCF的一款軟件，現(xiàn)已整合到EPI Suite中。本研究建模使用的數(shù)據(jù)為BCFBAF軟件中的non-ionic training set，外部驗(yàn)證集也為BCFBAF軟件中的validation set。BCFBAF軟件對(duì)外部驗(yàn)證集的預(yù)測(cè)結(jié)果為Q值為0.854、SE值為0.9以及SP值為0.75。通過(guò)與BCFBAF軟件的預(yù)測(cè)準(zhǔn)確度的比較顯示，我們構(gòu)建的模型CT-SVM與BCFBAF的預(yù)測(cè)準(zhǔn)確度與其相當(dāng)，從而進(jìn)一步證明了本研究所用建模方法的可行性。

3.3 警示子結(jié)構(gòu)分析

影響生物富集的因素有很多，例如生物物種特性、化合物的性質(zhì)、化合物的濃度和作用時(shí)間，以及環(huán)境因素等都是影響生物富集的主要因素。本研究采用信息增益方法、KNIME中的MoSS模塊以及ChemoTyper軟件找到了一系列警示子結(jié)構(gòu)，期望從化合物結(jié)構(gòu)上尋找易于引起生物富集的原因。化合物的穩(wěn)定性和脂溶性是引起生物富集的重要條件。例如DDT化學(xué)穩(wěn)定性強(qiáng)，為脂溶性物質(zhì)，易被吸收和積累在脂肪中。類似的化合物有機(jī)氯農(nóng)藥、多氯聯(lián)苯、甲基汞等化合物。我們使用信息增益的方法找到的10個(gè)警示子結(jié)構(gòu)中，有3個(gè)含有鹵素原子(見(jiàn)表4)，分別為芳基氯、芳基溴和氯代烷烴；使用KNIME中的MoSS模塊找到的10個(gè)警示子結(jié)構(gòu)中，有6個(gè)片段含有鹵素原子(見(jiàn)表5)；使用ChemoTyper軟件找到的10個(gè)警示子結(jié)構(gòu)中，9個(gè)片段含有鹵素原子(見(jiàn)表6)。以上警示子結(jié)構(gòu)均和文獻(xiàn)報(bào)道的已知生物富集化合物的結(jié)構(gòu)特征相吻合，證明了我們方法的可靠性。

從結(jié)果多樣性上分析以上三種方法找到的警示結(jié)構(gòu)發(fā)現(xiàn)，使用KNIME中的MoSS模塊找到的10個(gè)警示子結(jié)構(gòu)均為芳香族化合物，使用ChemoTyper軟件找到的10個(gè)警示子結(jié)構(gòu)中有7個(gè)子結(jié)構(gòu)為芳香族化合物，而使用信息增益的方法找到的10個(gè)警示子結(jié)構(gòu)屬于不同類的化合物。因此使用信息增益的方法尋找到的警示子結(jié)構(gòu)從結(jié)構(gòu)多樣化的角度，優(yōu)于KNIME中的MoSS模塊和ChemoTyper軟件。

本研究中我們使用7種不同的分子指紋結(jié)合5種機(jī)器學(xué)習(xí)方法構(gòu)建了具有高預(yù)測(cè)準(zhǔn)確度的二分類生物富集因子預(yù)測(cè)模型，使用10倍交叉驗(yàn)證的方法驗(yàn)證模型的魯棒性。其中四個(gè)模型(CT-SVM、ExtFP-SVM、PubChem-SVM和MACCS-SVM)對(duì)中/高生物富集和低生物富集化合物都具有很高的預(yù)測(cè)準(zhǔn)確度，從而保證了模型具有一定的實(shí)用性。另外，與BCFBAF軟件的預(yù)測(cè)結(jié)果相對(duì)比，發(fā)現(xiàn)我們構(gòu)建的模型CT-SVM與BCFBAF軟件的預(yù)測(cè)結(jié)果相當(dāng)，從而進(jìn)一步證明了使用分子指紋描述分子特征構(gòu)建分類模型方法的可行性。在本研究的最后，采用信息增益子結(jié)構(gòu)碎片分析、KNIME軟件的MoSS模塊分析和ChemoTyper軟件分析了生物富集的特權(quán)子結(jié)構(gòu)碎片和警示結(jié)構(gòu)，對(duì)生態(tài)系統(tǒng)安全評(píng)估具有一定的指導(dǎo)意義。生物富集過(guò)程并非一個(gè)簡(jiǎn)單、機(jī)械的分配過(guò)程，它受到很多因素的制約和影響，例如生物物種的特性、污染物的性質(zhì)、污染物的濃度及其作用時(shí)間，以及環(huán)境等都是影響生物富集的因素。因此只有應(yīng)用多參數(shù)分析的方法，在大量的實(shí)驗(yàn)數(shù)據(jù)的基礎(chǔ)上，才能尋找出更為合理的估算方法。

[1] Arnot J A, Gobas F A. A review of bioconcentration factor (BCF) and bioaccumulation factor (BAF) assessments for organic chemicals in aquatic organisms [J]. Environmental Reviews, 2006, 14(4): 257-297

[2] OECD. OECD Guidelines for Testing of Chemicals. TG 305: Bioaccumulation in Fish: Aqueous and Dietary Exposure [R]. OECD, 2012

[3] Cherkasov A, Muratov E N, Fourches D, et al. QSAR modeling: Where have you been? Where are you going to? [J]. Journal of Medical Chemistry, 2014, 57(12): 4977-5010

[4] European Chemicals Agency. Regulation (EC) No 1907/2006 of the European Parliament and of the Council of 18 December 2006 concerning the Registration, Evaluation, Authorisation and Restriction of Chemicals (REACH), establishing a European Chemicals Agency, amending Directive 1999/45/EC and repealing Council Regulation (EEC) No 793/93 and Commission Regulation (EC) No 1488/94 as well as Council Directive 76/769/EEC and Commission Directives 91/155/EEC, 93/67/EEC, 93/105/EC and 2000/21/EC [S]. European Chemicals Agency, 2007

[5] Weisbrod A V, Burkhard L P, Arnot J, et al. Workgroup report: Review of fish bioaccumulation databases used to identify persistent, bioaccumulative, toxic substances [J]. Environmental Health Perspectives, 2007, 115(2): 255-261

[6] Sa?an M T, Erdem S S, Ozpinar G A, et al. QSPR study on the bioconcentration factors of nonionic organic compounds in fish by characteristic root index and semiempirical molecular descriptors [J]. Journal of Chemical Information and Modeling, 2004, 44(3): 985-992

[7] Neely W B, Branson D R, Blau G E. Partition coefficient to measure bioconcentration potential of organic chemicals in fish [J]. Environmental Science & Technology, 1974, 8(13): 1113-1115

[8] Roy K, Sanyal I, Roy P P. QSPR of the bioconcentration factors of non-ionic organic compounds in fish using extended topochemical atom (ETA) indices [J]. SAR and QSAR in Environmental Research, 2006, 17(6): 563-582

[9] Dimitrov S, Dimitrova N, Parkerton T, et al. Base-line model for identifying the bioaccumulation potential of chemicals [J]. SAR and QSAR in Environmental Research, 2005, 16(6): 531-554

[10] Stadnicka J, Schirmer K, Ashauer R. Predicting concentrations of organic chemicals in fish by using toxicokinetic models [J]. Environmental Science & Technology, 2012, 46(6): 3273-3280

[11] Khadikar P V, Singh S, Mandloi D, et al. QSAR study on bioconcentration factor (BCF) of polyhalogented biphenyls using the PI index [J]. Bioorganic & Medicinal Chemistry, 2003, 11(23): 5045-5050

[12] Cui S H, Yang J, Liu S S, et al. Predicting bioconcentration factor values of organic pollutants based on MEDV descriptors derived QSARs [J]. Science in China Series B: Chemistry, 2007, 50(5): 587-592

[13] Zhao C, Boriani E, Chana A, et al. A new hybrid system of QSAR models for predicting bioconcentration factors (BCF) [J]. Chemosphere, 2008, 73(11): 1701-1707

[14] Gissi A, Nicolotti O, Carotti A, et al. Integration of QSAR models for bioconcentration suitable for REACH [J]. The Science of the Total Environment, 2013, 456-457: 325-332

[15] ACD Labs. ACD Labs homepage [OL]. [2014-12-04]. http://www.acdlabs.com/home/

[16] VEGA. VEG homepage [OL]. [2014-12-04]. http://www.vega-qsar.eu/

[17] OECD. QSAR ToolBox [OL]. [2014-12-04]. http://www.qsartoolbox.org/

[18] Estimation Program Interface (EPI) Suite [OL]. [2014-12-04]. http://www.epa.gov/opptintr/exposure/pubs/episuite.htm

[19] US Environmental Protection Agency. EPI Suite Data [OL]. [2014-12-04]. http://esc.syrres.com/interkow/EpiSuiteData.htm

[20] Costanza J, Lynch D G, Boethling R S, et al. Use of the bioaccumulation factor to screen chemicals for bioaccumulation potential [J]. Environmental Toxicology and Chemistry, 2012, 31(10): 2261-2268

[21] Yap C W. PaDEL-descriptor: An open source software to calculate molecular descriptors and fingerprints [J]. Journal of Computational Chemistry, 2011, 32(7): 1466-1474

[22] ChemoTyper Community. ChemoTyper website [OL]. [2014-12-04]. https://chemotyper.org/

[23] ChemAxon. ChemAxon website [OL]. [2014-12-04]. http://www.chemaxon.com

[24] Fourches D, Muratov E, Tropsha A. Trust, but verify: On the importance of chemical structure curation in cheminformatics and QSAR modeling research [J]. Journal of Chemical Information and Modeling, 2010, 50: 1189-1204

[25] Orange website. Orange website [OL]. [2014-12-04]. http://orange.biolab.si/

[26] Cortes C, Vapnik V. Support-Vector Networks [J]. Machine Learning, 1995, 20(3): 273-379

[27] Itskowitz P, Tropsha A. kappa Nearest neighbors QSAR modeling as a variational problem: Theory and applications [J]. Journal of Chemical Information and Modeling, 2005, 45(3): 777-785

[28] Watson P. Na?ve Bayes classification using 2D pharmacophore feature triplet vectors [J]. Journal of Chemical Information and Modeling, 2008, 48(1): 166-178

[29] Cheng F X, Yu Y, Zhou Y D, et al. Insights into molecular basis of cytochrome p450 inhibitory promiscuity of compounds [J]. Journal of Chemical Information and Modeling, 2011, 51(10): 2482-2495

[30] Baldi P, Brunak S, Chauvin Y, et al. Assessing the accuracy of prediction algorithms for classification: An overview [J]. Bioinformatics, 2000, 16(5): 412-424

[31] Shen J, Cheng F X, Xu Y, et al. Estimation of ADME properties with substructure pattern recognition [J]. Journal of Chemical Information and Modeling, 2010, 50(6): 1034-1041

[32] Jensen B F, Vind C, Padkjaer S B, et al. In silico prediction of cytochrome P450 2D6 and 3A4 inhibition using Gaussian kernel weighted k-nearest neighbor and extended connectivity fingerprints, including structural fragment analysis of inhibitors versus noninhibitors [J]. Journal of Medicinal Chemistry, 2007, 50(3): 501-511

[33] Kruhlak N L, Contrera J F, Benz R D, et al. Progress in QSAR toxicity screening of pharmaceutical impurities and other FDA regulated products [J]. Advanced Drug Delivery Reviews, 2007, 59(1): 43-55

[34] Benigni R, Bossa C. Structure alerts for carcinogenicity, and the Salmonella assay system: A novel insight through the chemical relational databases technology [J]. Mutation Research, 2008, 659(3): 248-261

[35] KNIME. KNIME website [OL]. [2014-12-04]. http://www.knime.org/

[36] Cheng F X, Shen J, Xu Y, et al. In silico prediction of Tetrahymena pyriformis toxicity for diverse industrial chemicals with substructure pattern recognition and machine learning methods [J]. Chemosphere, 2011, 82(11): 1636-1643

◆

InSilicoPrediction of Chemical Bioconcentration Factor

Sun Lu, Chen Yingjie, Wu Zengrui, Li Weihua, Liu Guixia, Philip W. Lee, Tang Yun*

Shanghai Key Laboratory of New Drug Design, School of Pharmacy, East China University of Science and Technology, Shanghai 200237, China

5 December 2014 accepted 9 January 2015

Bioconcentration is an important endpoint in evaluation of chemical adverse effects on ecosystems. In this study, in silico methods were used to predict chemical bioconcentration factor (BCF). At first a data set containing 624 chemicals with BCF values was collected from the Estimation Program Interface Suite of the U. S. Environmental Protection Agency. Using seven fingerprints to represent the molecules, binary classification models were developed with five machine learning methods, including support vector machine (SVM), C4.5 decision tree (C4.5 DT), k-nearest neighbors (kNN), random forest (RF), and Na?ve Bayes (NB). Reliable predictive models were then obtained and validated by 10-fold cross validation and external validation set. Among them, the model built by SVM with ChemoTyper fingerprint performed best, with predictive accuracy up to 85.4%. Moreover, some substructures were identified to be key for bioconcentration via several methods, such as arylchloride, diarylether, chloroalkene, and so on. The approaches used in this study provide a useful tool for environmental risk assessment of chemicals.

bioconcentration factor; in silico prediction; binary classification models; substructural alerts; environmental toxicology

國(guó)家自然科學(xué)基金(No. 81373329)；學(xué)科創(chuàng)新引智計(jì)劃即111計(jì)劃(No. B07023)

孫露(1989-)，女，碩士，研究方向?yàn)橛?jì)算機(jī)輔助藥物設(shè)計(jì)、藥物信息學(xué)和計(jì)算毒理學(xué)，E-mail: sunlu900326@yeah.net；

*通訊作者(Corresponding author)， E-mail: ytang234@ecust.edu.cn

10.7524/AJE.1673-5897.20141205001

2014-12-05 錄用日期：2015-01-09

1673-5897(2015)2-173-10

X171.5

A

唐赟(1968-)，男，博士，教授，主要研究方向?yàn)橛?jì)算機(jī)輔助藥物設(shè)計(jì)、藥物信息學(xué)、計(jì)算生物學(xué)和計(jì)算毒理學(xué)，已發(fā)表學(xué)術(shù)論文100余篇。

孫露, 陳英杰, 吳曾睿, 等. 有機(jī)化合物生物富集因子的計(jì)算機(jī)預(yù)測(cè)研究[J]. 生態(tài)毒理學(xué)報(bào), 2015, 10(2): 173-182

Sun L, Chen Y J, Wu Z R, et al. In silico prediction of chemical bioconcentration factor [J]. Asian Journal of Ecotoxicology, 2015, 10(2): 173-182 (in Chinese)