基于范數指數構效關系預測單端孢霉烯真菌毒素毒性

劉 英，尹京晨，計惠芬，賈青竹

(1. 天津科技大學海洋與環境學院，天津 300457；2. 浙江東洋環境工程有限公司，湖州 313000)

單端孢霉烯族毒素(trichothecenes)是一類化學性質相近的真菌毒素，由鐮刀菌屬、單端孢屬和頭孢霉屬等真菌產生[1-2]，其結構屬于倍半枯環氧化物，具有1個三環核心，在C12和C13位含有1個環氧化物，這是一個毒性必需基團(核心結構見圖 1)．該類毒素不僅污染小麥和玉米等禾谷類作物，也會通過食物鏈危害肉、奶、蛋等畜產品，WTO已將其列為天然存在的最危險的食品污染物之一[3]．而且，由于該類毒素對真核細胞具有多重抑制作用，通過對蛋白質、DNA 和 RNA 合成的抑制，對線粒體功能、細胞分裂和膜功能的抑制以及促進凋亡基因的表達，對人畜健康產生免疫抑制，導致其具有典型的“致癌、致畸、致突變”作用，從而嚴重威脅人畜健康[4]．

圖1 單端孢霉烯核結構Fig. 1 Structure of the trichothecene nucleus

定量構效關系(quantitative structure-activity rela-tionship，QSAR)是研究化合物結構與生物活性相關性的典型方法，廣泛應用于新藥物設計以及對于毒性機制的理解探究領域[5-8]．盡管已有學者對單端孢霉烯族毒素的結構與活性關系進行了考察，但是至今關于這一系列毒素的 QSAR 研究依然很少[5-6].Steinmetz等[6]基于分子比較力場分析(comparative molecular field analysis，CoMFA)方法和比較分子相似因子分析(comparative molecular similarity indices analysis，CoMSIA)方法對32個大環內酯物單端孢霉烯族毒素和 15個單端孢霉烯族毒素的毒性進行了3D(3 dimension)QSAR研究，結果表明單端孢霉烯族毒素類物質的取代基以及分子構象對其生物活性都有一定影響．

本課題組之前提出系列范數指數描述符已經成功地預測了有機物多種活性，包括麻醉性污染物的水生毒性、雜環化合物的藥理學和毒理學活性等[7-10].先前的研究工作表明該系列范數描述符也許可以從本質上體現分子結構，從而反映化合物的多個物化性質.

在本課題組之前工作的基礎上，本工作又提出新的描述符，并據此建立單端孢霉烯真菌毒素的毒性預測模型，對 35個單端孢霉烯真菌毒素的毒性進行了計算，采用留一法交叉驗證(LOO)、外部驗證和 Y隨機驗證手段對模型進行了驗證，并對模型的應用域(AD)進行了評價．本工作深入研究此類化合物分子結構對其毒性的影響，將有助于人們采取更科學的毒素防控及清除機制．

1 材料與方法

1.1 樣本集

Jarvis等[11-12]測定了大環內酯物真菌毒素在延長接種p-388白血病小鼠生命中的能力，這些數據可以從根本上反映大環內酯物真菌毒素對于小鼠白血病細胞的影響，并將測定活性定量為R，即R＝100×測試集動物存活天數/對照組動物存活天數．在此工作中，包含大環內酯物真菌毒素毒性數據的化合物分子(代表性分子結構見圖2)從文獻中[6,13-14]獲得．

圖2 兩種大環內酯物真菌毒素結構Fig. 2 Structure of two macrolide trichothecene mycotoxins

1.2 分子結構優化

利用軟件HyperChem 7.0構建分子的3D結構，對有機物分子結構進行優化，具體采用量子化學中的從頭算法ab initio在ST0-3G進行能量最低優化，分子結構信息中包括分子中原子的種類、數量、分支度，原子之間的鍵連接形式以及穩定結構下的原子電荷．

1.3 范數指數描述符和模型

在分子結構優化基礎上，為了區分原子的特點和種類，本工作提出了性質矩陣．同時，歐氏空間距離矩陣是一種能夠全面、有效的描述分子的方法，因此，本工作提出了不同的矩陣．

下面列出上述具體矩陣．

為了討論上述矩陣的模，本工作使用了兩類范數指數 norm(TP，1)和 norm(TP，2)，其定義式如下：norm(TP，1)：Fro函數

在式(8)中，操作符(.×)的含義是標量相乘．通過范數計算得到48個范數描述符，為避免描述符過多而造成模型過度擬合，對建模所用描述符進行優化篩選，最后優選了6個范數描述符，建立預測模型如下：

模型中系數b見表1．

表1 模型系數Tab. 1 Parameters of this model

1.4 模型評價

通過下列參數評價本工作建立模型的穩定性：相關系數的平方(R2)、平均相對誤差值．為了進一步驗證模型的準確性，還使用了留一交叉驗證法(Q2)、外部驗證法、Y隨機驗證法進行驗證．

同時，本工作選擇了應用域驗證本模型的應用范圍[15]，其定 義為式(15)：

在上述等式中，h是杠桿值，x是變量值．

通過與閾值(h*)比較，如果h＜h*，表明在這個區域內模型的預測性是穩定的；反之，這個模型的預測能力是較差的[16-17]．閾值的定義為式(16)

式中：p′是自變量數量加1；n是訓練集數量．

2 結果與討論

2.1 大環內酯物真菌毒素毒性的QSAR模型

利用新建模型式(9)對 35個大環內酯物真菌毒素的毒性進行了預測，圖3是大環內酯物真菌毒素毒性的實驗值log10R與預測值對比散點圖．

圖3 log10R預測值和實驗值相關性Fig. 3 Correlation between log10R predicted by this model and the experimental data

由圖3可知：所有毒素的毒性預測點和實驗點均位于對角線上及附近，表明本模型計算結果與實驗值有很好的一致性．本預測模型相關統計數據、和 F值分別為 0.925,4、0.996,6和 78.95，說明本模型計算結果的精確性較好．

圖 4是殘差與實驗值的對比圖．從圖 4可以看出所有毒素的毒性預測殘差都分布在-0.3～0.3，這表明本工作的系統誤差小，模型預測能力穩定；同時，可以說明本工作提出的變量能有效地反映分子結構與毒性之間的關系．

圖4 殘差與實驗值對比圖Fig. 4 Plot of the residual vs．experimental data from this model

本工作與Steinmetz等[6]對比結果見表2．由表 2可以看出文獻中 CoMSIA方法和本工作得到的相關系數 R2都大于 0.900，R2均在同一數量級，但是本工作交叉驗證的相關系數遠大于 CoMSIA方法，說明本工作的模型更穩定．對于F值而言，本工作的F值同樣遠大于文獻中兩種方法得出的 F值，更進一步證明本模型的準確性高、穩定性好．

表2 本方法和文獻方法預測結果對比Tab. 2 Statistics comparison of reference methods and this work

本工作采用多元線性回歸分析法，建立關于大環內酯物真菌毒素的結構與毒性關系模型．Steinmetz等[6]利用CoMFA方法研究結果表明，空間位阻和靜電因素對毒素的活性有重大影響，其中靜電因素的貢獻更大；同時，相比之下包含氫鍵供體/受體的CoMSIA方法計算結果更為精確．本工作在計算描述符的優化確定中，首先利用歐氏距離以及相關距離矩陣考慮了大環內酯物真菌毒素分子中各個組成原子的空間分布；同時，為了實現大環內酯物真菌毒素分子中不同原子對其毒性的貢獻影響，提出了原子屬性矩陣即增廣矩陣，其中包括電負性、原子質量、原子支化度和原子電荷．由此，根據以上相關矩陣建立的范數描述符能從分子水平上反映真菌毒素分子結構特點以及真菌毒素與生物體之間的相互作用力．綜上，與前人研究相比，本工作優點在于：能給出確定的模型形式，方便其他研究者驗證或使用；模型預測準確性和穩定性強．

在此以第13個化合物為例，計算該化合物的毒性.化合物Baccharin B8的結構如圖5所示．

首先，按照式(1—3)的定義式，根據化學圖計算出該化合物的矩陣 T1、T2和性質矩陣 PE；然后依據式(6—8)，計算出增廣矩陣 TP1,a,b、TP2,a,b和 TP3,a,b,.表 3是上述所有矩陣的范數描述符的值．應用式(9)以及表2和表3中的參數，第13個化合物的毒性預測值計算過程如下：

第 13個化合物的 log10R預測值是 2.363，實驗值是2.367，相對偏差是0.17%,．

圖5 Baccharin B8的結構Fig. 5 Structure of Baccharin B8

表3 第13個化合物的相關變量值Tab. 3 Variable values of the thirteenth chemical

2.2 留一交叉驗證

留一交叉驗證法就是從所有的數據中每次選擇一個分子作為一次測試集，用剩余的分子建立模型預測測試集分子的毒性．如果留一交叉驗證的相關系數大于0.5，則認為模型是合理的[18]．圖6表明：留一交叉驗證法的毒性預測值與實驗值有較好吻合度，表明本工作建立的模型的內部預測能力較高．

圖6 留一交叉驗證預測值和實驗值相關性Fig. 6 Correlation between the data predicted by LOOCV and the experimental data

圖 7是兩種模型預測結果的標準誤差(SE)分布情況對比圖．由圖 7可知：交叉驗證的誤差趨勢和本方法建立的模型的誤差趨勢一致，說明本工作建立的模型穩定、預測能力好．同時，留一交叉驗證法的結果具有較高 Q2值(0.878,9)，以上所有驗證結果均表明本工作提出的方法可以穩定地預測大環內酯物真菌毒素的毒性．

圖7 本模型和LOO模型毒性預測標準誤差分布Fig. 7 Standard error distributions of the toxicity predicted by this model and LOO-CV

2.3 外部驗證

模型的內部驗證只能說明該方法對于類似物的活性預測能力好，如果用該方法預測新的不同種類的化合物的活性，則結果有可能是不可靠的．因此，需要進行外部驗證．外部驗證就是將全部的分子分為訓練集和測試集，訓練集用于建立模型，測試集用來驗證模型．有效的外部驗證能確保本工作提出的方法對測試集分子預測的可預見性和適用性．

在本工作中，訓練集和測試集的分子與文獻[6]相同，訓練集和測試集的分子毒性預測值見圖 8，其結果是0.925,4、是 0.996,6．這表明預測值和實驗值有很好的一致性，進而說明本工作建立的模型不僅可以很好地預測同類化合物毒性，而且對多種類化合物具有較高的預測性．

圖8 訓練集和測試集的實驗值和預測值的散點圖Fig. 8 Scatter plot predicted and experimental data of the training set and test set

2.4 Y隨機驗證

Y隨機驗證的原理是將實驗值隨機擾亂、變量 不變，再建立模型，得到的結果與最初建立模型的結果進行對比．

在本工作中，進行了 5次隨機驗證，新模型預測R2和 Q2的結果與原模型的對比結果見表 4．由表 4可以得出，5次Y隨機驗證得出的新模型的R2和Q2都很低，甚至為0，說明新模型的預測能力低．由此證明，本工作建立的原始模型是穩定的，并非偶然建立.

表4 Y隨機驗證結果Tab. 4 Results of Y randomization test of this model

2.5 應用域驗證

對于構效關系模型統計參數的驗證來說，模型的應用域評價依然是最困難的一個評價標準．應用域評價是模型能否準確預測新化合物性質或活性的重要指標．目前，已經有多種關于模型應用域評價的方法[19-21]．本工作利用杠桿方法檢測計算模型的應用域，結果如圖9所示．從圖 9可以看出，對于35個有機物，只有 1個有機物的杠桿值小于-3，其余所有分子的杠桿值均在規定的閾值內，說明本方法建立的模型應用域較廣泛．

圖9 本模型的應用域Fig. 9 Application domain of this model

3 結 語

本工作基于有機毒素化學結構構建的距離矩陣和性質矩陣，提出了一個新范數描述符，并據此建立了預測單端孢霉烯真菌毒素毒性的 QSAR模型．研究結果表明：本工作提出的范數描述符能夠有效描述化合物毒性與其結構之間的關系；本模型毒性預測值與實驗值有很好的一致性，其中模型的相關性系數R2均大于 0.9，F值為 78.95，留一交叉驗證測試(Q2為 0.878,9)、Y隨機驗證以及與文獻的對比結果均證明本模型穩定性更好、預測能力更高、準確性更好.同時，應用域驗證表明本模型計算結果準確、穩定、可靠、適用范圍廣．因此，本工作提出的方法將可能推廣用于化學品的環境風險評價、生態健康評價等領域.

參考文獻：

［1］ Zhou T，He J W，Gong J. Microbial transformation of trichothecene mycotoxins[J]. World Mycotoxin Journal，2008，1(1)：23-30.

［2］ 鄒忠義，賀稚非，李洪軍，等. 單端孢霉烯族毒素轉化降解研究進展[J]. 食品科學，2010，31(19)：443-448.

［3］ 周闖，何成華，司慧民，等. 2012年國內飼料及原料霉菌毒素污染調查分析[J]. 畜牧與獸醫，2014，46(1)：81-84.

［4］ 許偉，耿芳芳，范夢雪，等. 脫氧雪腐鐮刀菌烯醇毒性的研究進展[J]. 生物學雜志，2016，33(1)：78-81.

［5］ Betina V. Structure-activity relationships among mycotoxins[J]. Chemico-Biological Interactions，1989，71(2/3)：105-146.

［6］ Steinmetz W E，Rodarte C B，Lin A. 3D QSAR study of the toxicity of trichothecene mycotoxins[J]. European Journal of Medicinal Chemistry，2009，44(11)：4485-4489.

［7］ Zhu Z，Wang Q，Jia Q，et al. Structure-property relationship for the pharmacological and toxicological activity of heterocyclic compounds[J]. Acta Physico-chimica Sinica，2014，30(6)：1086-1090.

［8］ Qian H，Kanwal S，Jia Q，et al. Norm index-based quantitative structure-activity relationship to predict βcyclodextrin complex binding constants[J]. Acta Physico-chimica Sinica，2015，31(5)：893-898.

［9］ Wang Q，Jia Q，Yan L，et al. Quantitative structuretoxicity relationship of the aquatic toxicity for various narcotic pollutants using the norm indexes[J]. Chemosphere，2014，108：383-387.

［10］ Jia Q，Cui X，Li L，et al. A Quantitative structure-activity relationship for high affinity 5-HT1A receptor ligands based on norm indexes[J]. Journal of Physical Chemistry B，2015，119(51)：15561-15567.

［11］ Jarvis B B，Stahly G P，Pavanasasivam G. Isolation and characterization of the trichoverroids and new roridins and verrucarins[J]. Journal of Medicinal Chemistry，1982，47(6)：1117-1124.

［12］ Jarvis B B，Vrudhula V M，Midiwo J O，et al. New trichoverroids from myrothecium verrucaria：Verrol and 12，13-deoxytrichodermadiene[J]. Journal of Medicinal Chemistry，1983，48(15)：2576-2578.

［13］ Cramer R D，Patterson D E，Bunce J D. Comparative molecular field analysis(CoMFA). 1. Effect of shape on binding of steroids to carrier proteins[J]. Journal of the American Chemical Society，1988，110(18)：5959-5967.

［14］ Klebe G，Abraham U，Mietzner T. Molecular similarity indices in a comparative analysis(CoMSIA)of drug molecules to correlate and predict their biological activity[J]. Journal of Medicinal Chemistry，1994，37(24)：4130-4146.

［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]. American Theological Library Association，2005，33(2)：155-173.

［16］ Gramatica P. Principles of QSAR models validation：Internal and external[J]. QSAR & Combinatorial Science，2007，26(5)：694-701.

［17］ Gramatica P，Giani E，Papa E. Statistical external validation and consensus modeling：A QSPR case study for Kocprediction[J]. Journal of Molecular Graphics and Modelling，2007，25(6)：755-766.

［18］ Mihali? Z，Trinajsti? N. A graph-theoretical approach to structure-property relationships[J]. Journal of Chemical Education，1992，69(9)：701.

［19］ Tropsha A. Best practices for QSAR model development，validation，and exploitation[J]. Molecular Informatics，2010，29(6/7)：476-488.

［20］ Tropsha A，Golbraikh A. Predictive QSAR modeling workflow，model applicability domains，and virtual screening[J]. Current Pharmaceutical Design，2007，13(34)：3494-3504.

［21］ Jaworska J，Nikolova-Jeliazkova N，Aldenberg T. QSAR applicabilty domain estimation by projection of the training set descriptor space：A review[J]. American Theological Library Association，2005，33(5)：445-459.