周悅，李瀟嵐，程薇，馮艷，王艷,2*

（1.上海交通大學醫(yī)學院公共衛(wèi)生學院，上海 200025）（2.上海交通大學醫(yī)學院附屬第九人民醫(yī)院，上海 200011）

神經(jīng)毒性是外源性化合物引起的神經(jīng)系統(tǒng)結構或功能的損害，而經(jīng)口攝入是眾多可導致神經(jīng)毒性中毒癥狀的化學品重要的暴露方式之一[1,2]。如今隨著食品工業(yè)的發(fā)展，食品種類增多，食品中添加劑使用的增加，許多食品化學污染物如重金屬、食品添加劑、農(nóng)藥獸藥殘留物等，均可誘導發(fā)生氧化應激、神經(jīng)細胞發(fā)育受損、神經(jīng)遞質含量變化等各種神經(jīng)系統(tǒng)損傷，嚴重威脅人體健康[3-5]。因此，神經(jīng)毒性評估成為了食品污染物暴露毒性評價必不可少的環(huán)節(jié)。傳統(tǒng)的神經(jīng)毒性評估主要通過動物實驗測試進行。但動物實驗周期長且成本高昂，在高通量篩選和人類中樞神經(jīng)系統(tǒng)的外推性方面存在一定的局限性[6]。因此，迫切需要開發(fā)替代方法和工具，來評估化學品的神經(jīng)毒性。

定量構效關系（Quantitative StructureActivity Relationship，QSAR）是目前常用的毒性預測模型之一，利用化合物的物理化學、幾何等結構信息，預測其可能的生物活性[7]。分子描述符分為2D和3D描述符。雖然3D描述符在模擬分子受體配體對接，處理空間化學方面可能有更大的優(yōu)勢，但2D描述符具有無需估算生物活性構象，計算成本低，數(shù)據(jù)處理快的優(yōu)點[8]，同時預測結果依然優(yōu)異[9,10]。QSAR已被成功運用于經(jīng)口毒性判斷和內(nèi)分泌干擾物篩選[11]，但針對神經(jīng)毒性，相關研究多聚焦于同源化合物，納入同屬數(shù)據(jù)集建立，如多氯聯(lián)苯類（Polychlorinated Biphenyls，PCBS）[12]、芳基氨基脲類化合物等[13]，通過同源化合物間的相似結構構建線性模型，分析特定取代基或分子結構與其神經(jīng)毒性間關聯(lián)[14]。因食品中的化學性污染物來源廣泛、種類繁多，化合物的結構特征多樣，神經(jīng)毒性癥狀各異且毒性作用靶點各不相同，化合物可通過影響神經(jīng)分化成熟、神經(jīng)元遷移/空間定向、突觸形成/可塑性、神經(jīng)遞質、神經(jīng)營養(yǎng)功能、神經(jīng)免疫系統(tǒng)、神經(jīng)內(nèi)分泌、離子穩(wěn)態(tài)干擾等功能誘發(fā)神經(jīng)毒性[15]。因此，傳統(tǒng)QSAR線性模型不適用于食品化學污染物分子結構與神經(jīng)毒性關聯(lián)的建立。機器學習算法中的監(jiān)督學習，在學習框架的基礎上半自動化創(chuàng)建分析模型，并通過交叉驗證和迭代改進以細化和優(yōu)化預測模型的準確性，具有“多指標、多關聯(lián)、多交叉”的特點，尤其適用于大規(guī)模樣本的歸納研究[16]。

考慮到該模型的適用范圍和高通量篩選的需要，本研究將選則2D描述符用于神經(jīng)毒性預測模型的建立。通過支持向量機、隨機森林、神經(jīng)網(wǎng)絡等機器學習算法，構建食品化學污染物神經(jīng)毒性分類預測模型。據(jù)此闡明影響化合物神經(jīng)毒性的主要分子結構因素。進一步了解化合物性質，為后續(xù)化合物的毒性預警及風險評估開展提供參考依據(jù)。

1 材料與方法

1.1 數(shù)據(jù)庫的建立與神經(jīng)毒性資料獲取

從現(xiàn)有文獻[17,18]、政府通告的違法食品添加物或新檢出物、生活常見藥物和食品安全污染物中選擇了107種化合物組建數(shù)據(jù)庫?；衔锝M成包括20種殺蟲劑（18.69%）、14種阻燃劑（13.08%）、16種多環(huán)芳烴（14.95%）、19種食品包裝化工原材料（17.76%）、6種食品添加劑或污染物（5.61%）及32種藥物（29.91%）。選擇現(xiàn)有人類、非人靈長類動物、哺乳動物實驗以及體外細胞實驗中，有兩篇及以上出版物表明其具有神經(jīng)毒性的化合物做陽性判定。選擇現(xiàn)有文獻中無神經(jīng)毒性相關實驗，實驗結果矛盾或僅單細胞系體外實驗陽性的化合物做陰性判定[19]。根據(jù)以上判定原則，最終獲得包含各類神經(jīng)毒性機制，如阻礙神經(jīng)元發(fā)育成熟、干擾神經(jīng)遞質傳遞、影響突觸形成等化合物57種，及明確無神經(jīng)毒性化合物50種。數(shù)據(jù)庫于2020年5月完成建立，于2021年5月對化合物文獻檢索結果進行更新。

1.2 數(shù)據(jù)集的建立

按8:2比例隨機劃分訓練集及測試集，訓練集用于模型構建，測試集用于模型驗證[20,21]。

1.3 描述符的獲取與特征選擇

分子描述符是化合物物理化學性質的量化表征，其將分子中的結構化學信息，轉換成特定的指標或某些標準化實驗的結果，其選擇一定程度上決定了模型的優(yōu)劣[22,23]。2D分子描述符一般包括：（1）組成描述符（如分子量、鹵素類型）；（2）拓撲描述符（如Wiener指數(shù)、Randic指數(shù)）；（3）電子描述符（如極性指數(shù)、Hammett常數(shù)）；（4）幾何參數(shù)描述符（如分子體積、溶劑接觸面積）等[24,25]。本研究從pubchem網(wǎng)站獲取化合物的2D分子結構信息，導入MOLD2（2.0）軟件后計算得到777個分子描述符。考慮到疏水性在經(jīng)口毒性預測、藥物篩選方面的重要作用，使用Episuite（version 4.1 US EPA 2012）及VEGA（core version1.2.1 IRFMN）軟件，計算常見于相關毒性預測模型中的正辛醇-水分配系數(shù)（Octanol-Water Partition Coefficient，Kow）納入研究[26,27]，并針對以上778個計算變量開展特征選擇。

特征選擇主要由以下五個步驟組成。（1）刪去恒定、半恒定（相同取值比例＞0.85）描述符后剩余472個計算變量。（2）二元logistics回歸完成描述符的初步篩選。（3）結合已有文獻經(jīng)驗，補充納入現(xiàn)有研究中廣泛運用于預測藥物吸收、代謝、毒性的相關描述符，包括：分子量、柔性可旋轉鍵數(shù)、Kow、拓撲極性表面積等[24,28]。（4）將472個描述符和訓練集89種化合物帶入五種機器學習模型并分析計算變量重要性，納入重要性排名前五位者。（5）刪去重復后將以上合計52種計算變量納入隨機森林模型進行擬合，采用遞歸特征消除法進行重復構建和計算變量篩選，每次構建刪除其中最不重要的變量[29]。最終納入如表1所示5種描述符，主要覆蓋化合物可旋轉鍵數(shù)、電拓撲、分子極性相關性質，包括1個原子數(shù)和鍵數(shù)描述符、2個理化性質描述符、2個拓撲描述符。根據(jù)“Rule-of-Thumb”規(guī)則，2D-QSAR模型中訓練集（N）與描述符（M）的比值（N/M）應不小于5，本研究所得最優(yōu)模型的N/M為17.80＞5，滿足此項規(guī)則，表明描述符的選擇數(shù)量適當[30]。

表1 構建QSAR模型所使用描述符Table 1 The descriptors used in QSAR model

1.4 模型的構建

共采用5種機器學習算法[25,28,31]：

由于數(shù)據(jù)集的各個特征取值范圍存在較大差異，對描述符取值進行標準化處理：

式中：

X——化合物描述符參數(shù)數(shù)值；

μ——參數(shù)平均值；

σ——參數(shù)標準差。

1.4.1 貝葉斯網(wǎng)絡（Bayesian Network，BN）

貝葉斯算法適用于基于定量研究的不確定性推理，利用概率統(tǒng)計對未知樣品進行分類。變量間的關系強度由兩個概率組成，一為條件概率，可被父節(jié)點個數(shù)和局部概率影響，二為節(jié)點固定先驗概率。貝葉斯網(wǎng)絡要求各個變量在除父節(jié)點以外節(jié)點間條件獨立，現(xiàn)有研究中可見于方劑配伍規(guī)律等藥性研究[32]。使用Rstudio（version 4.0.2下同）“e1071”軟件包中的“naiveBayes”函數(shù)實現(xiàn)。

1.4.2 隨機森林（Random Forests，RF）

隨機森林算法因其易于理解、執(zhí)行力高等優(yōu)勢常見于藥物安全性評價等藥物分析領域[32]。研究中使用隨機森林算法的基分類器為cart算法，通過重采樣方法生成隨機自助訓練集，未被納入訓練集樣本稱為袋外樣本。構建多顆樹后選擇最佳模型輸出。使用Rstudio軟件“randomforest”軟件包進行模型構建和優(yōu)化，“pROC”軟件包進行圖形繪制，“varImpPlot”函數(shù)進行模型參數(shù)重要性可視化。最終構建模型節(jié)點中二叉樹變量個數(shù)（mtry）取值為3。構建基分類器模型數(shù)量（ntree）取值為500，最小子節(jié)點（nodesize）取值為1，由于特征參數(shù)較少，不對最大子節(jié)點數(shù)加以限制。

1.4.3 類神經(jīng)網(wǎng)絡（Artificial Neural Network，ANN）

類神經(jīng)網(wǎng)絡具有強大的擬合能力和自我學習能力，用于非線性映射回歸、多分類預測類研究[33]。從信息處理角度對人腦神經(jīng)元網(wǎng)絡進行抽象，包括接受外部信號數(shù)據(jù)的輸入單元、完成系統(tǒng)結果對外傳遞的輸出單元，和不能由系統(tǒng)外部觀察得到的隱藏單元。本研究采用ANN中的多層感知機反向傳播算法（BP）的三層人工神經(jīng)網(wǎng)絡進行擬合。使用Rstudio軟件“nnet”軟件包進行模型構建。同時為進一步增加模型穩(wěn)定性，避免過度擬合，使用SPSS MODELER（version18.0 IBM）軟件中的bagging算法模塊對原ANN模型進行優(yōu)化，合并基分類器規(guī)則設定為投票。bagging組件模型數(shù)為10。

1.4.4 支持向量機（Support Vector Machine，SVM）

支持向量機相較于類神經(jīng)網(wǎng)絡，更適用于高緯度、小樣本回歸或二分類研究[34]，在現(xiàn)有研究中常根據(jù)待測化合物特征加以選擇，構建毒性預測模型。作為一種常見的，對數(shù)據(jù)進行二元分類的廣義線性算法，其決策邊界是對學習樣本求解的最大邊距超平面。當線性不可分時，將松弛變量引入超平面，并在模型中產(chǎn)生另一個懲罰系數(shù)。使用Rstudio軟件“e1071”軟件包中的“svm”函數(shù)進行模型的構建，核函數(shù)為分類非線性SVM中使用最多的高斯徑向基核（RBF），運用“tune”函數(shù)試錯法將γ范圍設定為10^[-6,-1]，懲罰因子C范圍設定為10^[1,4]進行模型優(yōu)化，最終選擇最佳參數(shù)γ取值0.1，C取值1000。

1.4.5 K最近鄰值算法（K-NearestNeighbor，KNN）

KNN算法依靠距離劃分未知樣本類別，通過尋找多個離測試樣本最鄰近的訓練樣本，以該樣本的多數(shù)類別作為最終類別，操作簡單結果容易解讀，可見于部分藥物分析模型[35]。但這一計算方法不適用于大樣本數(shù)據(jù)，且樣本不平衡時，模型預測偏差較大[36]。本研究通過SPSS.MODELER軟件KNN模塊進行內(nèi)部十折交叉驗證選擇最合適K值，距離計算方法為Euclidean法。最終模型選定K值為6時具有最低錯誤率。

1.5 模型的評價與驗證

本研究采用十折交叉檢驗進行數(shù)據(jù)集內(nèi)部魯棒性驗證，即將原始數(shù)據(jù)隨機平均分為10個子集，每個子集做測試集的同時，其余9個子集合并作為訓練集，進行10次訓練預測[37]。十次十折交叉驗證由Rstudio軟件完成。測試集化合物分類進行外部驗證，使用Rstudio軟件“confusionMatrix”函數(shù)輸出混淆矩陣，計算其中的真陽性（True Positive，TP）、真陰性（True Negative，TF）、假陽性（False Positive，F(xiàn)P）、假陰性（False Negative，F(xiàn)N），分別計算靈敏度（Sensitivity，SE）、特異度（Specificity，SP）、假陽性率（False Positive Rate，F(xiàn)PR）、假陰性率（False Negative Rate，F(xiàn)NR）、總準確性（Global Accuracy，GA）。通過“pROC”軟件包繪制受試者工作曲線（Receiver Operating Characteristic，ROC）并計算曲線下面積（The Area Under The ROC Curve，AUC）進行模型評價。

2 結果與討論

2.1 描述符共線性分析

5種描述符間相關性分析結果如表2所示，可見任意兩個描述符間相關系數(shù)＜0.75。描述符之間不存在明顯的共線性問題。

表2 5個分子描述符的相關矩陣分析Table 2 Correlation matrix analysis of five molecular descriptors

2.2 模型質量評價

在進行數(shù)據(jù)分割后，訓練集和測試集分別有89和18種化合物，且類別分布較為均勻無需進行數(shù)據(jù)不均衡處理。具體數(shù)據(jù)見圖1。

五種機器學習的表現(xiàn)各不相同，根據(jù)外部驗證的結果，五種模型中隨機森林模型綜合表現(xiàn)最優(yōu)。訓練集中靈敏度91.49%、特異度100.00%，GA和AUC均超過90%。在測試集中表現(xiàn)也較為穩(wěn)定，靈敏度、特異度分別為80.00%、87.50%，同時GA和AUC均超過80%。具體數(shù)據(jù)見表3。

表3 基于5個描述符對于神經(jīng)毒性預測訓練集及測試集的統(tǒng)計結果Table 3 Statistical results of neurotoxicity prediction training set and test set based on 5 descriptors

模型魯棒性檢驗結果顯示，在五種模型中，隨機森林模型表現(xiàn)最為穩(wěn)定，測試集平均準確率達到70.24%。其次為ANN和KNN模型，準確率分別為67.34%和64.50%。支持向量機和貝葉斯模型的魯棒性稍差，分別為61.28%和62.91%。

模型的預測準確性方面，隨機森林的模型表現(xiàn)也更優(yōu)，訓練集準確性為95.51%，測試集準確性為83.33%，該模型給出的訓練集、測試集ROC曲線如圖2所示，曲線下面積可達到0.99和0.85，是個較為理想的分類模型。

目前已有少量研究利用機器學習算法建立了的化合物毒性預測模型。如一項基于概率神經(jīng)網(wǎng)絡算法（PNN），納入47種有機溶劑作為研究對象，預測其神經(jīng)毒性的分類模型。訓練集、測試集中GA分別達0.91、0.93，曲線下面積分別達0.92、0.94[38]。一項基于RF算法，納入61種農(nóng)藥作為研究對象，預測其熊蜂急性接觸毒性的分類模型。訓練集、測試集中GA分別達0.83、0.87，曲線下面積分別達0.91、0.92[39]。相較于前兩種模型較為局限的化合物類別，本研究納入的食品污染物種屬各異，脫離了傳統(tǒng)的線性模型和化合物類型的限制。在擴大模型適用范圍的同時，保持了較高的總預測準確度和AUC值，訓練集、測試集AUC分別達0.99、0.85。而相較于現(xiàn)有的大樣本中藥的神經(jīng)保護作用預測模型（訓練集SE:0.67、SP:0.99；測試集SE:0.69、SP:0.98）[40]，極大的提高了模型靈敏度，訓練集、測試集靈敏度分別達0.91、0.8，適用于大樣本數(shù)據(jù)庫化合物神經(jīng)毒性的初步篩查。

2.3 分類錯誤案例分析

如上文所述，隨機森林算法對數(shù)據(jù)庫中化合物具有較好的分類準確率，但預測模型中仍有一些化合物不能被正確的分類，分類錯誤的7個化合物的結構、錯分類型和描述符數(shù)值如表4、表5所示。

表4、表5數(shù)據(jù)表明，該隨機森林模型在特異度上表現(xiàn)良好，但預測化合物時存在一定的的假陰性現(xiàn)象。同時，描述符D564（質量加權Burden矩陣最大特征值）取值可能在一定程度影響模型對于化合物神經(jīng)毒性的判定，D564取值大的化合物更傾向于被判定為陽性，而D564取值小的化合物則相反。

表4 隨機森林算法中分類錯誤化合物Table 4 Misclassified compounds in random forest algorithm

化合物分類錯誤的原因可能由活性懸崖（化合物具有類似結構但活性不同）和骨架躍遷（化合物具有不同結構但活性相似）導致，如苯并(a)芘和苯并(b)芘。而化合物中質量加權Burden矩陣最大特征值在化合物神經(jīng)毒性推斷方面可能有比較重要的作用。為了進一步驗證該結論，我們選取了兩種已知神經(jīng)毒性，且D564描述符相差較大的數(shù)據(jù)庫外化合物作為參考帶入模型，具體化合物和描述符取值見表5。結果支持這一推論。本研究模型相較于現(xiàn)有研究中基于SVM算法的周圍神經(jīng)病變預測模型（訓練集、測試集GA均為0.89，未見AUC指標）[41]，具有更好的模型解讀能力，強調了質量加權Burden矩陣最大特征值在化合物神經(jīng)毒性預測中的作用。未來研究中我們將進一步增加數(shù)據(jù)庫中化合物數(shù)量，突破單一模型限制，構建集成模型以提高化合物神經(jīng)毒性預測準確性和普適性，有望用于指導高風險化合物的毒性測試篩選策略的制定。

表5 錯誤分類化合物及庫外參考化合物描述符取值Table 5 Misclassified compounds and values of off-library reference compound descriptors

3 結論

3.1 本研究收集整理了具有影響神經(jīng)分化成熟、神經(jīng)元遷移/空間定向、突觸形成/可塑性、神經(jīng)遞質、神經(jīng)營養(yǎng)功能、神經(jīng)免疫系統(tǒng)、神經(jīng)內(nèi)分泌、離子穩(wěn)態(tài)干擾等誘發(fā)神經(jīng)毒性機制的食品污染物57種，無神經(jīng)毒性污染物50種共107種食品污染物，運用5種機器學習算法，建立了食品化學污染物神經(jīng)毒性預測模型。

3.2 運用R和SPSS軟件對778個計算變量進行了篩選，最終確定將5個分子描述符納入模型，即D015可旋轉鍵數(shù)比例、D142 Balaban均方頂點距離指數(shù)、D152碳-sp3的平均原子極化率、D175 Wiener最大路徑指數(shù)和D564質量加權Burden矩陣最大特征。

3.3 5種機器學習算法中，基于隨機森林算法建立的分類模型具有最好的預測效果，外部驗證中訓練集靈敏度91.49%、特異度100.00%，總準確性95.51%，測試集靈敏度80.00%、特異度87.50%，總準確性83.33%，曲線下面積分別達0.99和0.85。內(nèi)部驗證中模型魯棒性達70.24%。模型預測結果中共有7種化合物不能被正確分類，造成這一結果的原因可能是活性懸崖或骨架躍遷現(xiàn)象的存在。在五種分子描述符中，質量加權Burden矩陣最大特征值對預警化合物高風險結構特征貢獻最佳。