以網(wǎng)絡(luò)拓?fù)渚嚯x為特征的有效藥物組合預(yù)測

任 彪

（上海理工大學(xué) 管理學(xué)院，上海 200093）

0 引言

將兩種或以上藥物組合在一起治療某種疾病，能夠產(chǎn)生協(xié)同作用，治療效果遠(yuǎn)大于單個(gè)藥物作用的效果，即為所謂的組合用藥［1］。組合用藥作用于多個(gè)靶點(diǎn)，避免了可能引起的反饋機(jī)制，減少了藥物毒副作用［2］?；谶@些優(yōu)勢，可明確組合藥物對于復(fù)雜疾病，如癌癥、心血管疾病等具有更好的治療效果，因此組合用藥具有重要的研究價(jià)值［3］。傳統(tǒng)藥物組合的發(fā)現(xiàn)是基于高通量實(shí)驗(yàn)篩選，通過分析高維實(shí)驗(yàn)數(shù)據(jù)，判斷藥物之間是否會產(chǎn)生協(xié)同作用，進(jìn)而確定藥物組合［4］。但如今藥物總體數(shù)量龐大，藥物配對的數(shù)量超出了高通量實(shí)驗(yàn)?zāi)軌虺惺艿姆秶?。組學(xué)技術(shù)的出現(xiàn)可幫助研究者解決這一難題，其結(jié)合目前積累的大量數(shù)據(jù)和經(jīng)驗(yàn)，有效利用生物信息學(xué)和計(jì)算工具，挖掘出可靠的信息及各生物信息之間的潛在關(guān)系［5］。通過該方法既能發(fā)現(xiàn)藥物的作用機(jī)制，又能降低研究成本、縮短研究周期。

近年來，各種針對藥物組合的預(yù)測方法發(fā)展迅速，其中一些利用組合特征之間的相似性建立機(jī)器學(xué)習(xí)模型以預(yù)測有效的藥物組合。如Xu 等［6］整合生物學(xué)（靶蛋白、疾病通路）、化學(xué)（二維結(jié)構(gòu)、化合物相互作用的可能性）及藥理學(xué)信息（藥物療效），運(yùn)用隨機(jī)梯度提升算法預(yù)測藥物組合。此外，對于非線性的動(dòng)態(tài)作用數(shù)據(jù)，基于數(shù)學(xué)模型的方法具有不錯(cuò)的效果。如針對特定分子的常微分方程（ODE）模型，以及針對分子運(yùn)作過程的Petri 網(wǎng)絡(luò)模型等［7-9］。但這些方法都有不足之處，例如基于數(shù)學(xué)模型的方法需要豐富的動(dòng)力學(xué)參數(shù)，否則會降低預(yù)測效果［1］。

除上述方法外，基于網(wǎng)絡(luò)的方法也得到了廣泛應(yīng)用，利用網(wǎng)絡(luò)探索目標(biāo)在網(wǎng)絡(luò)上的關(guān)聯(lián)情況，預(yù)測潛在的協(xié)同藥物組合。如Zou 等［10］描述藥物靶點(diǎn)與其鄰居在PPI 網(wǎng)絡(luò)中的相互作用，區(qū)分協(xié)同藥物組合；Chen 等［11］構(gòu)建一個(gè)通路—通路的相互作用（WWI）網(wǎng)絡(luò)，描述復(fù)雜的通路—通路之間的關(guān)系，從基于通路的角度探討藥物協(xié)同作用［12］。這些方法考慮兩個(gè)藥物的靶點(diǎn)在網(wǎng)絡(luò)中的拓?fù)潢P(guān)系，但并沒有考慮其與疾病單元之間的關(guān)系。故本文引入藥物—藥物—疾病單元三者在蛋白質(zhì)相互作用網(wǎng)絡(luò)中的拓?fù)渚嚯x作為特征，結(jié)合生物學(xué)和藥理學(xué)相似性特征，利用支持向量機(jī)、邏輯回歸和隨機(jī)森林3 種機(jī)器學(xué)習(xí)算法構(gòu)建藥物組合預(yù)測模型。若藥物能夠作用于更多疾病單元中的基因，即兩個(gè)藥物需要同時(shí)作用于疾病單元，且兩者作用于疾病單元的部分盡可能不重合［13］，則能獲得更好的治療效果。將這三者在蛋白質(zhì)相互作用網(wǎng)絡(luò)中的拓?fù)潢P(guān)系轉(zhuǎn)化為兩兩之間的距離，以此作為特征進(jìn)行預(yù)測，最終提高預(yù)測效果。

1 特征工程

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

本文采用的樣本數(shù)據(jù)主要來源于Drug Combination Database（DCDB，version 2.0），DCDB 包含1 363 種藥物組合，這里選取用于治療高血壓的86 對藥物組合作為研究對象。Drugbank 數(shù)據(jù)庫是一個(gè)提供生物信息和化學(xué)信息的數(shù)據(jù)庫，其具有詳細(xì)的藥物數(shù)據(jù)和全面的藥物目標(biāo)信息［14］。從drugbank 數(shù)據(jù)庫中獲取藥物的靶標(biāo)信息、二維分子結(jié)構(gòu)信息與ATC 編碼信息，其中靶標(biāo)的蛋白質(zhì)序列信息來自UniProt 數(shù)據(jù)庫（http：//www.uniprot.org/），高血壓的致病基因來自O(shè)MIM 數(shù)據(jù)庫［15］，包括目前所有已知的遺傳病以及超過15 000 個(gè)基因信息。蛋白質(zhì)相互作用信息源自Menche 等［16］提供的蛋白質(zhì)互作文件和HuRI-人類蛋白互作組數(shù)據(jù)庫（http：//www.interactome-atlas.org/），最新公布的HuRI-人類蛋白互作組數(shù)據(jù)庫中包含了9 094 個(gè)獨(dú)立蛋白構(gòu)成的64 006 對相互作用關(guān)系［17］。

1.2 正負(fù)樣本構(gòu)建

將每個(gè)藥物的化學(xué)結(jié)構(gòu)信息、ATC 編碼信息、靶標(biāo)信息與蛋白質(zhì)序列信息進(jìn)行整合，并與86 對藥物進(jìn)行組合匹配，若藥物組合缺失以上任意一個(gè)信息，該組合將被刪除。通過整合刪選，最終留下81 對治療高血壓的藥物組合，以此作為正樣本。這81 對藥物組合包含69 種獨(dú)立藥物，將其兩兩組合，共得到2 346 種組合。排除DCDB 中已出現(xiàn)的組合，從剩下的組合中挑選81 對組合作為負(fù)樣本。

1.3 蛋白質(zhì)相互作用網(wǎng)絡(luò)構(gòu)建

將參考文件和HuRI-人類蛋白互作組數(shù)據(jù)庫中的蛋白質(zhì)相互作用數(shù)據(jù)取并集之后，共獲得由15 063 個(gè)獨(dú)立蛋白質(zhì)所構(gòu)成的187 371 對相互作用關(guān)系。利用uinprot 在線mapping 工具，將這些蛋白質(zhì)名稱統(tǒng)一轉(zhuǎn)化成Entrez_ID，以便于構(gòu)建蛋白質(zhì)網(wǎng)絡(luò)，并保證靶標(biāo)及致病基因的ID 相同。將15 063 個(gè)蛋白看成15 063 個(gè)節(jié)點(diǎn)，每一對相互作用視為一條邊，即可構(gòu)成一張蛋白質(zhì)相互關(guān)系的網(wǎng)絡(luò)圖。

1.4 特征構(gòu)建

1.4.1 化學(xué)結(jié)構(gòu)相似性

若兩個(gè)藥物化學(xué)結(jié)構(gòu)的相似性分?jǐn)?shù)越高，則其更有可能產(chǎn)生相似作用，這里調(diào)用第三方包RDkit 計(jì)算兩個(gè)藥物化學(xué)結(jié)構(gòu)的相似性。首先利用獲得的每個(gè)藥物分子的二維結(jié)構(gòu)信息（smiles）調(diào)用RDkit 計(jì)算對應(yīng)的二維MACCS 分子指紋［18］，其是一種基于SMARTS、長度為167 的分子指紋，然后利用Tanimoto 系數(shù)計(jì)算二維MACCS 分子指紋的相似度。兩個(gè)藥物的Tanimoto 相似度定義如下：

式中，a 代表藥物A 分子指紋中的分子位數(shù)，b 代表藥物B 分子指紋中的分子位數(shù)，c 代表A 和B 中相同的位數(shù)。

1.4.2 藥物療效相似性

藥物療效相似性是指藥物在治療疾病時(shí)功能相似，其在預(yù)測協(xié)同藥物組合中具有重要作用。ATC 編碼是解剖治療學(xué)及化學(xué)分類系統(tǒng)編碼，可表示藥物療效信息。編碼由7 位組成，共分為5 級：第1 級為一位字母，第2 級為兩位數(shù)字，第3 級為一位字母，第4 級為一位字母，第5 級為兩位數(shù)字，分別表示藥物在解剖學(xué)、治療學(xué)、藥理學(xué)、化學(xué)、化合物上的分類。對于有多個(gè)ATC 編碼的藥物，需要計(jì)算每一對ATC 編碼的相似性Sk(A,B)，最后選擇最大的Sk(A,B)值作為藥物療效相似性的數(shù)值。藥物療效相似性定義如下：

式中，ATCK代表藥物ATC 編碼。

1.4.3 靶標(biāo)序列相似性

利用萊溫斯坦距離［19］計(jì)算兩個(gè)藥物靶標(biāo)序列的相似性，其可用來計(jì)算兩個(gè)不等長序列的相似性，原理是計(jì)算出將序列A 變換成序列B 需要變換的次數(shù)。

式中，a、b 為兩個(gè)數(shù)組（字符串），i、j為數(shù)組下標(biāo)，這里可直接調(diào)用Python 中的固定模塊進(jìn)行計(jì)算。

1.4.4 靶標(biāo)單元與疾病單元距離計(jì)算

首先需要將藥物靶標(biāo)信息和高血壓致病基因信息轉(zhuǎn)化成Entrez_ID，并將其映射到相互作用網(wǎng)絡(luò)中，刪除蛋白質(zhì)網(wǎng)絡(luò)中的獨(dú)立節(jié)點(diǎn)和不在網(wǎng)絡(luò)中的靶標(biāo)。藥物靶標(biāo)或致病基因在網(wǎng)絡(luò)中往往位于相鄰位置，所以認(rèn)為藥物—藥物—疾病單元三者之間的位置關(guān)系滿足互補(bǔ)性的拓?fù)潢P(guān)系時(shí)，對于預(yù)測協(xié)同藥物具有很大幫助［20］。利用d（X，Y）表示藥物X 與疾病單元Y 之間距離，具體公式如下：

式中，d（X，Y）代表藥物靶蛋白x 與疾病蛋白y 在網(wǎng)絡(luò)中的最短距離，可采用Python 計(jì)算得到。

1.4.5 兩個(gè)靶標(biāo)單元之間距離計(jì)算

計(jì)算兩靶標(biāo)之間距離與計(jì)算靶標(biāo)單元到疾病單元的距離相似。由于兩個(gè)藥物靶標(biāo)單元中含有的蛋白數(shù)量為同一數(shù)量級，故可通過比較靶標(biāo)單元之間最短距離和兩個(gè)靶標(biāo)單元在網(wǎng)絡(luò)中半徑的均值來計(jì)算靶標(biāo)單元之間距離。具體公式如下：

2 機(jī)器學(xué)習(xí)模型

監(jiān)督學(xué)習(xí)分為學(xué)習(xí)和預(yù)測兩個(gè)過程，其任務(wù)是學(xué)習(xí)一個(gè)模型，使模型對于任意給定的輸入會預(yù)測最好的結(jié)果作為輸出［21］。模型是通過對訓(xùn)練數(shù)據(jù)集的學(xué)習(xí)獲得的，之后再對測試樣本集進(jìn)行預(yù)測。監(jiān)督學(xué)習(xí)利用的數(shù)據(jù)集都由特征及其對應(yīng)標(biāo)簽組成，將測試集的特征向量作為模型輸入，并將得到的結(jié)果與特征向量原有標(biāo)簽進(jìn)行對比，從而獲得模型預(yù)測的準(zhǔn)確度。例如，每對藥物組合具有5 個(gè)特征，其對應(yīng)1 個(gè)標(biāo)簽，即如果該組合有效則為1，否則為0。機(jī)器學(xué)習(xí)算法會對訓(xùn)練集中的特征和標(biāo)簽進(jìn)行學(xué)習(xí)，分析一組特征是有效還是無效，并形成自己的判別標(biāo)準(zhǔn)，從而判斷測試的數(shù)據(jù)是否有效。

2.1 支持向量機(jī)

支持向量機(jī)（Support Vector Machines，SVM）于1964 年被提出后得到了迅速發(fā)展，并在多個(gè)領(lǐng)域獲得青睞，目前已應(yīng)用于人像識別、文本分類等多個(gè)模式識別問題中［22］。其是監(jiān)督學(xué)習(xí)中的二分類模型，目的是找到一個(gè)超平面分割正負(fù)樣本，并遵循間隔最大化原則。這也是其與感知機(jī)的區(qū)別所在，最終轉(zhuǎn)化為一個(gè)求解凸二次規(guī)劃問題。超平面，通俗的講就像二維平面中的直線或三維空間中的一個(gè)平面，以此類推。在對超平面沒有限制時(shí)，現(xiàn)實(shí)中可以找到無數(shù)個(gè)這樣分割正負(fù)樣本的直線或平面，但當(dāng)限定間隔最大時(shí)，則只能找到唯一滿足條件的超平面。通過線性函數(shù)能對樣本分類，則稱這些數(shù)據(jù)樣本線性可分，可利用支持向量線性分類器對其進(jìn)行分類。然而，對于那些非線性樣本，可通過核函數(shù)和軟間隔最大化形成非線性支持向量機(jī)，本質(zhì)上還是將非線性問題轉(zhuǎn)化成線性問題進(jìn)行處理。本文綜合考慮后采用了線性支持向量機(jī)。

2.2 邏輯斯蒂回歸

邏輯斯蒂回歸（Logistic Regression）屬于統(tǒng)計(jì)學(xué)習(xí)算法中的經(jīng)典分類方法，其是一種線性模型［23］。邏輯斯蒂回歸在線性回歸模型基礎(chǔ)上，使用Sigmod 函數(shù)估計(jì)預(yù)測結(jié)果概率P（y | x）的大小。邏輯斯蒂回歸并沒有對數(shù)據(jù)分布進(jìn)行建模，其不知道數(shù)據(jù)具體分布，而是直接求解分類超平面。在本文的應(yīng)用中利用邏輯斯蒂回歸模型直接計(jì)算藥物組合是有效組合的概率，范圍在0～1 之間，大于0.5 時(shí)則被認(rèn)為有效。

2.3 隨機(jī)森林算法

隨機(jī)森林算法（Random Frost）是基于Bagging 集成學(xué)習(xí)理論的代表算法，由Leo 于2001 年提出［24］。作為一種新興的、高度靈活的機(jī)器學(xué)習(xí)算法，具有廣闊的應(yīng)用前景，在準(zhǔn)確率方面相當(dāng)具有優(yōu)勢。隨機(jī)森林是指通過集成學(xué)習(xí)思想將多棵樹集成的一種算法，其基本單元是決策樹，而其本質(zhì)屬于機(jī)器學(xué)習(xí)中的集成學(xué)習(xí)（Ensemble Learning）方法。通俗的講，在分類問題中，每棵決策樹都是一個(gè)分類器，N 棵樹就會有N 種結(jié)果，隨機(jī)森林集成了所有投票結(jié)果，將投票次數(shù)最多的類別視為最終結(jié)果。隨機(jī)森林能夠有效運(yùn)行在大數(shù)據(jù)集上，但當(dāng)數(shù)據(jù)量較小時(shí)會出現(xiàn)過擬合現(xiàn)象，其還能夠處理具有高維特征的輸入樣本，且不需要降維，并能評估各個(gè)特征在分類問題上的重要性。

3 模型構(gòu)建與結(jié)果分析

運(yùn)用上述3 類模型與計(jì)算出的特征進(jìn)行有效藥物組合預(yù)測，選擇一倍于正樣本的負(fù)樣本構(gòu)建模型。根據(jù)簡單交叉驗(yàn)證原理，從數(shù)據(jù)集中選取75%的樣本作為訓(xùn)練集，25%的樣本作為測試集。本文利用AUC（ROC 曲下面積）、Ac?curacy、F-measure、Recall、precision 作為模型評價(jià)標(biāo)準(zhǔn)，將未加入3 個(gè)距離特征與加入距離特征后的模型預(yù)測結(jié)果進(jìn)行比較。

3.1 ROC 曲線

ROC（Receiver Operating Characteristic）在機(jī)器學(xué)習(xí)領(lǐng)域用來評判分類、檢測結(jié)果好壞，是一種重要與常見的統(tǒng)計(jì)分析方法［25］?；煜仃囀荝OC 曲線繪制的基礎(chǔ)，主要涵蓋4 個(gè)指標(biāo)：TP（將正例正確預(yù)測為正例）、FN（將正例錯(cuò)誤預(yù)測為負(fù)例）、PF（將負(fù)例錯(cuò)誤預(yù)測為正例）、TN（將負(fù)例正確預(yù)測為負(fù)例）。在ROC 曲線圖中，每個(gè)點(diǎn)以對應(yīng)的FPR值為橫坐標(biāo)，TPR 值為縱坐標(biāo)，其中FPR 為假陽性率，TPR為真陽性率。具體公式如下：

利用4 個(gè)指標(biāo)還可得到精確率（Precision）、召回率（Re?call）和準(zhǔn)確度（Accuracy），進(jìn)一步計(jì)算得到F1 評分，其是精確率與召回率的加權(quán)平均。具體公式如下：

3.2 ROC 曲線下面積

AUC（曲線下面積）可用來評價(jià)分類器區(qū)分正負(fù)類的性能，AUC 的值在0～1 之間，該值越接近1，說明分類效果越好。當(dāng)AUC 的值等于0.5 時(shí)，模型能夠正確分類的概率只有一半，類似于隨機(jī)拋硬幣的概率，此時(shí)分類效果很差。利用TPR 值和FPR 值，可在對應(yīng)坐標(biāo)平面中確定一個(gè)點(diǎn)，由各個(gè)點(diǎn)確定ROC 曲線。例如，給定一個(gè)閾值為0.5，意味著若分類模型對每個(gè)實(shí)例判斷為正類的概率大于等于0.5，則判斷其為正類，小于0.5 則判斷為負(fù)類，在這種情況下計(jì)算出TPR 值和FPR 值，則可確定一個(gè)點(diǎn)（FPR，TPR）；同樣，將閾值設(shè)置為0.6 時(shí)，又可確定一個(gè)（FPR，TPR）點(diǎn)，這樣眾多點(diǎn)即可構(gòu)成一條ROC 曲線。理想情況是曲線經(jīng)過（0，1）點(diǎn)，即FPR 為0，TPR 為1，但在現(xiàn)實(shí)中很難達(dá)到，所以只要曲線越接近這個(gè)點(diǎn)越好。通過計(jì)算曲線下方面積，面積越大則曲線越靠近（0，1）點(diǎn)，以此判斷預(yù)測結(jié)果的好壞。

3.3 預(yù)測過程與結(jié)果

由于藥物化學(xué)結(jié)構(gòu)相似性、藥物療效相似性與藥物靶標(biāo)序列相似性已在藥物組合預(yù)測方面表現(xiàn)突出，故本文首先選擇這3 個(gè)特征作為輸入，分別利用當(dāng)下比較流行的機(jī)器學(xué)習(xí)算法進(jìn)行預(yù)測。利用不同模型進(jìn)行預(yù)測，并對比上述各項(xiàng)評價(jià)模型性能指標(biāo)，最終支持向量機(jī)、邏輯斯蒂回歸、隨機(jī)森林算法展示出不錯(cuò)的預(yù)測效果。具體指標(biāo)如表1 所示，ROC 曲線如圖1 所示（彩圖掃OSID 碼可見，下同）。之后為體現(xiàn)藥物—藥物—疾病單元三者在蛋白質(zhì)相互作用網(wǎng)絡(luò)中的距離在預(yù)測有效藥物組合中的作用，本文將3個(gè)計(jì)算出的距離作為特征加入，保持各項(xiàng)參數(shù)不變，再進(jìn)行上述同樣的操作，得到的各項(xiàng)指標(biāo)如表2 所示，ROC 曲線如圖2 所示。

Table 1 Indicators of prediction results of three models（without distance feature）表1 3 個(gè)模型（未加入距離特征）預(yù)測結(jié)果各項(xiàng)指標(biāo)

Table 2 Indicators of prediction results of three models（adding distance feature）表2 3 個(gè)模型（加入距離特征）預(yù)測結(jié)果各項(xiàng)指標(biāo)

Fig.1 Comparison of ROC curves without distance feature圖1 未加入距離特征的ROC 曲線對比

Fig.2 ROC curve comparison after adding distance feature圖2 加入距離特征后的ROC 曲線對比

3.4 結(jié)果分析

從表1、表2 對比中可以看出，在未加入距離特征進(jìn)行預(yù)測時(shí)，準(zhǔn)確率在0.7 左右，且AUC 值相對較低。在加入距離特征后，準(zhǔn)確率達(dá)到0.8 左右，AUC 值也有所提高。由于樣本數(shù)量較少，采用隨機(jī)森林算法預(yù)測時(shí)可能出現(xiàn)過擬合現(xiàn)象，并在實(shí)驗(yàn)中表現(xiàn)為結(jié)果不穩(wěn)定，所以本文以其他兩個(gè)模型的準(zhǔn)確率作為最終比較依據(jù)。支持向量機(jī)和邏輯斯蒂回歸的預(yù)測結(jié)果分別提升了12%與9%。

4 結(jié)語

組合藥物的使用無疑可幫助人們治療復(fù)雜疾病，同時(shí)利用計(jì)算機(jī)技術(shù)、組學(xué)和網(wǎng)絡(luò)技術(shù)幫助人們發(fā)現(xiàn)新的藥物組合，因此是一種行之有效的手段。該方法極大地縮小了搜索范圍，在小范圍內(nèi)再進(jìn)行實(shí)驗(yàn)檢驗(yàn)，更為安全、可靠。發(fā)現(xiàn)新的可靠特征也是準(zhǔn)確預(yù)測的關(guān)鍵之一，本文引入的網(wǎng)絡(luò)距離特征對于提高預(yù)測準(zhǔn)確率起到了重要作用，加入網(wǎng)絡(luò)距離特征后，預(yù)測準(zhǔn)確率提高了10%以上。本文使用的化學(xué)結(jié)構(gòu)特征、藥物療效特征、藥物靶標(biāo)序列特征基本都是比較經(jīng)典的，前人都曾經(jīng)采用過，故可保證其可信度。

此外，網(wǎng)絡(luò)距離特征原理上是希望藥物—藥物—疾病單元三者之間能夠滿足一定的特殊條件，即兩種藥物不僅能同時(shí)作用于疾病單元，且相互之間可保持足夠的距離。通俗一點(diǎn)講，即希望在作用于更多致病基因的同時(shí)，產(chǎn)生的不良反應(yīng)更小。本文直接運(yùn)用網(wǎng)絡(luò)中三者的距離，并未完全展現(xiàn)出這種互補(bǔ)性關(guān)系，所以下一步希望探索出能夠準(zhǔn)確描述這種拓?fù)潢P(guān)系的特征。本文利用的特征相對較少，其他文獻(xiàn)中還運(yùn)用了化合物相互作用、蛋白質(zhì)相互作用等特征，故本文之后也會結(jié)合生物學(xué)、藥理學(xué)、化學(xué)等多種類的數(shù)據(jù)特征，融合靜態(tài)數(shù)據(jù)與動(dòng)態(tài)數(shù)據(jù)進(jìn)行藥物組合預(yù)測。本文所用模型都是監(jiān)督模型，而負(fù)樣本是隨機(jī)構(gòu)建的，不排除當(dāng)中會有有效的藥物組合，若利用半監(jiān)督模型或其他更先進(jìn)的算法，相信在預(yù)測精度上會得到進(jìn)一步提升。