基于多數據源與機器學習的藥物副作用預測

杜 瑤

（上海理工大學管理學院，上海 200093）

0 引言

藥物不良反應（Adverse Drug Reaction，ADR）通常定義為在正常藥物劑量下，在預防、診斷或治療中產生的有害和非預期的影響［1］。這種影響會給患者身心造成危害，嚴重的甚至會威脅到患者生命安全。近年來，因藥物副作用所引發的藥物安全問題備受關注。在美國，每年有200 萬人受到藥物副作用的影響，其中有10 萬人因嚴重的藥物副作用而死亡［2］。嚴重的藥物副作用也是新藥開發失敗和已上市藥物撤回的一個重要因素［3］。由于藥物研發周期過長以及藥物研發所需的人力物力和經濟成本非常高昂，所以提前預測藥物副作用可起到指導藥物開發的重要作用，針對已上市藥物的副作用進行預測還可有效指導醫生進行臨床治療。因此，及時、準確地預測藥物副作用已成為國內外研究的熱點問題［4］。

近年來，眾多研究人員也提出了一些預測藥物不良反應的方法。例如，Liang 等［5］在2019 年利用多視圖學習重要特征進行藥物副作用預測；Hu 等［6］在2018 年利用堆疊深度異構信息網絡嵌入方法進行藥物副作用預測；Emir等［7］在2017 年將生物醫學知識以圖的形式生成不同特征集以預測藥物副作用；Zheng 等［8］在2018 年利用藥物特征的逆相似性構建可靠的負樣本以預測藥物副作用。以上都是使用計算方法對藥物副作用進行預測。因此，利用數據挖掘和機器學習等計算方法分析海量的生物醫學數據，實現對藥物副作用的精準預測已成為了該領域的熱門研究方法［9］。很多將藥物自身特征作為研究的主要對象，這些特征包括藥物的化學結構、靶點蛋白、分子通路等。綜上所述，整合藥物副作用相關信息，利用當下流行且高效的計算方法開發一套藥物副作用預測工具，既可以幫助藥物研發人員進行新藥研發，也可以幫助患者規避臨床治療中產生副作用的風險，具有重要的科學意義及廣闊的應用前景［4］。

1 預測方法

1.1 原理

通常認為在藥物的作用下，如果基因表達量與患病情況下的基因表達量變化情況相反，則認為藥物起到了積極作用，如果基因表達量與患病情況下的基因表達量變化情況相同，則認為藥物可能沒有起到積極作用，或者說藥物更可能產生副作用，因為其導致了基因的非正常表達，這也是本文的核心思想。此前，也有不少學者使用基因相關信息作為研究對象建立預測模型。例如，Xiang 等［1］利用基因與藥物副作用的關聯網絡進行藥物副作用預測；Uner等［10］將基因表達作為特征，使用深度學習框架進行藥物副作用預測；Fukuzaki 等［11］利用統計相似性模型，基于藥物通路和基因表達譜信息預測藥物副作用；Wang 等［12］在2016 年利用LINCS L1000 經典基因表達數據集結合化學結構等信息進行藥物副作用預測。但以上預測方法都未將基因表達量變化情況與疾病聯系起來。在疾病治療過程中，用藥后即會引起體內基因的表達量發生變化，一旦出現副作用，體內某些基因就會出現非正常表達，也即是說副作用與基因表達始終是聯系在一起的。因此，基于基因表達的副作用預測方法就是將藥物副作用與基因表達聯系起來，找到兩者之間的關聯性。本文利用患病前后和用藥前后重要基因的表達量變化情況作為解釋變量，因變量副作用標簽來自權威的SIDER 數據庫，再利用當下流行且高效的機器學習算法建立模型進行藥物副作用預測。

1.2 工作內容

本文首要進行的工作是準備好高質量的數據，本文數據全部來自官方數據庫。從官方數據庫中下載所需的原始數據，包括患病基因表達數據集、用藥基因表達數據集、人類重要基因集以及藥物副作用數據集。本文以乳腺癌細胞作為研究對象，因此從數據集中選擇與乳腺癌相關數據，在數據集預處理過程中，按照需要對原始數據進行合理篩選。經過層層篩選后得到所需數據，最終選定了401種藥物和142 個基因用于模型構建，標簽則是一個表示副作用有無的401 維向量。

在完成數據預處理后，還要篩選適合本文數據類型的機器學習算法建立模型。這里根據數據類型及數據量大小選擇了隨機森林算法（Random Frost，RF），并且選擇K 近鄰（k-Nearest Neighbor，KNN）、決策樹（Decision Tree）以及樸素貝葉斯（Naive Bayesian，NB）3 個算法作為對比。算法確定之后，再利用處理好的特征數據進行模型訓練與測試，然后對比每種機器學習模型所得結果，最終發現隨機森林在幾種算法中表現最好。

2 特征工程

2.1 數據收集及預處理

GEO 和cmap 兩個基因表達數據庫是本文的核心數據庫，其中基因受到疾病影響所產生的基因表達量變化情況來自GEO 數據庫（https：//www.ncbi.nlm.nih.gov/），這是美國國立衛生研究院1988 年創立的美國國立生物技術信息中心，目的是給分子生物學家提供一個信息存儲與處理的系統。在數據庫中找到所需的樣本組，這里選擇的是12 個乳腺細胞樣本，分為6 個正常對照組和6 個患病實驗組。Value 值即為樣本某個探針所對應基因的表達量，如表1所示。基因的差異表達分析是利用R 語言的Limma 包實現的，這是一個專門用于基因差異分析的包［13］。

Table 1 Breast cell sample data表1 乳腺細胞樣本數據

本文選擇直接在數據庫官網在線對兩組樣本的基因表達量差異進行分析，由于這里的差異分析同樣是使用R語言的Limma 包得出的結果，所以該結果真實可靠，一共得到2 萬多個基因的差異表達。以下呈現了前5 個基因的差異分析結果，如表2 所示。

Table 2 Results of difference analysis表2 差異分析結果

但這里得到的差異分析結果并不能直接使用，還要對其進行篩選。在統計學上，adj.P.Val 的值越小，表示結果越好。通常情況下，本文設置adj.P.Val 的閾值為0.05，取小于等于0.05 的值［14］。logFC 表示差異量，絕對值越大，表示兩組樣本的基因表達量差異越大，大于零表示后者的表達量大于前者，小于零表示后者的表達量小于前者。其值由如下公式求出：

式中，x為實驗組樣本的基因表達量，y 為對照組的基因表達量。這里log（FC）的閾值設置為1，篩選出實驗組基因表達量為對照組基因表達量至少兩倍的基因作為差異表達基因。log（FC）大于零的部分基因作為在疾病作用下上調的基因，小于零的部分基因作為在疾病作用下下調的基因。在對數據進行篩選的過程中，還需要用到一個對重要基因進行篩選的數據庫（http：//tubic.tju.edu.cn/deg/blast.php？type=single&db=e）。這些基因的表達量變化會對人體產生較大影響，因此本文重點關注這些基因的表達量變化情況。將重要基因數據集與之前預處理的基因集合取交集，便得到患病之后重要基因的表達量變化情況。

接下來是對用藥之后可體現基因表達量變化的相關數據進行處理，這里選擇的是cmap 數據庫（https：//portals.broadinstitute.org/cmap/）中藥物對疾病細胞的處理結果。cmap 數據庫是一個基因表達譜數據庫，利用藥物小分子對疾病細胞進行處理，得到一個前后基因表達量變化譜。為了與前面的疾病相對應，同樣選擇乳腺癌細胞的處理結果。本文將ratio 矩陣中作用于乳腺癌細胞的部分篩選出來，進行相應處理后，形成一個藥物與基因表達量的矩陣，每個藥物與每個基因作一次映射。ratio 矩陣中大于1 的值表示用藥后基因表達量上調，反之下調。最后同樣將基因集與此前處理的疾病情況下上下調基因結果取交集，得到最終的基因集合。

2.2 特征構建

在最終數據集的處理中，本文將患病與用藥之后同時上調的基因標記為1，同理，將患病與用藥之后同時下調的基因也標記為1，反向變化的則標記為0。最終的解釋變量形式如表3 所示。一個藥物樣本的形式為一串0、1 數字的組合，表示在藥物作用下特定基因的表達量變化情況。這里考慮同時上調和下調都屬于患病與用藥之后的同向變化，因此把同時上調和下調的基因組合到一個矩陣里，將藥物作用下同向變化的不管是同時上調還是同時下調的基因都標記為1，而反向變化的不管是先下調再上調，還是先上調再下調，都標記為0。假設此為矩陣T，如果=1，表示基因j在疾病和藥物i的作用下，表達量發生了同向變化。

因變量也即本文的數據標簽，來自權威的藥物及相關副作用資源庫SIDER 數據庫，里面收集整理了1 430 個藥物、5 880 個副作用信息以及99 423 個藥物—副作用關系［4］。這里將此1 430 種藥物與之前cmap 數據庫中作用于乳腺癌細胞的1 241 種藥物取交集，得到401 種藥物集合。這401 種藥物便是解釋變量的樣本數量，標簽則是這401 種藥物所對應的有無某副作用的情況，數據形式是一個401 維的向量，也用0 和1 表示，1 表示藥物有該副作用，0 則表示藥物沒有該副作用。

Table 3 Results of data processing表3 數據處理結果

3 分類器選擇

3.1 選擇依據

構建好模型特征之后，選擇一種合適且高效的算法也是本文工作的重要內容，目前最流行的是使用機器學習算法建立模型進行分類預測等任務。如今，機器學習、深度學習等技術已得到了廣泛應用，相關算法也因不斷優化使得代碼更加通俗易懂、易于實現，可根據數據類型、數據量大小以及任務內容進行算法選擇。

本文數據特征是一個藥物樣本的142 個基因在患病前后和用藥前后的基因表達量變化情況，其是由0 和1 的布爾類型數據組成的，并非具體某個數值。通過分析，可判斷這種數據類型不適合類似線性回歸的相關算法。由于復雜的算法往往需要大量數據才能發揮效果，本文數據量不多，而且是一個二分類問題，因此需要盡量選擇簡單、高效的分類算法。本文選擇的目標分類器是隨機森林，其是一種集成學習算法，可整合多棵決策樹，從中得到最優結果，另外還選擇了決策樹算法與隨機森林算法進行對比，驗證隨機森林融合多棵決策樹取最優結果的特點，最后使用具有概率特點的樸素貝葉斯作為本文的分類器。最終，本文選擇了隨機森林、K 近鄰、決策樹、樸素貝葉斯4種機器學習算法。

3.2 分類器介紹

3.2.1 隨機森林

隨機森林算法（Random Frost）是基于Bagging 集成學習理論的代表算法，由Breiman［15］于2001 年提出。其是利用bootsrap 重抽樣方法從原始樣本中抽取多個樣本，對每個bootsrap 樣本進行決策樹建模，然后組合多棵決策樹的預測，通過投票得出最終預測結果。一棵決策樹的分類能力可能很弱，但在隨機產生大量決策樹并組成隨機森林后，每個樣品都逐一通過一棵樹分類決策，最后組合的結果將更接近于正確分類［16］。大量理論與實證研究都證明了RF 具有很高的預測準確率，對異常值和噪聲具有很好的容忍度，且不容易出現過擬合。隨機森林通過在每個節點處隨機選擇特征進行分支，每棵分類樹之間的相關性得到最小化，故對多元共線性不敏感，從而提高了分類精確性及抗噪聲能力［17］。隨機森林是一種自然的非線性建模工具，也是目前數據挖掘領域的熱門研究方向之一。

3.2.2 K 近鄰

K 近鄰算法（KNN）是基于實例的代表算法，于1967 年由Cover 等［18］提出，是一種用于分類與回歸的統計方法，主要思想是通過測量不同特征值之間的距離進行分類。如果一個樣本在特征空間中的k 個最相似或在特征空間中最鄰近樣本中的大多數屬于某一類別，則該樣本也屬于該類別。在KNN 算法中，所選擇的鄰居都是已正確分類的對象。該方法在確定樣本類別時依據最鄰近的一個或幾個樣本類別決定待分樣本所屬類別［17］。在本文工作中，觀察離它最近的一個或幾個藥物的基因表達特征標簽屬于哪一類，則目標樣本也屬于此類。

3.2.3 決策樹

決策樹，顧名思義，其類似于一棵樹，利用樹的結構對數據記錄進行分類，樹的一個葉結點即代表某個條件下的一個記錄集，根據記錄字段的不同取值建立樹的分支，在每個分支子集中重復建立下層結點和分支，便可生成一棵決策樹［19］。本文的數據特征為離散變量，非常適合使用決策樹對數據作分類預測，同理也非常適用于隨機森林算法。

3.2.4 樸素貝葉斯

1960 年，Maron 等［20］首先提出樸素貝葉斯分類方法，其是一種基于概率模型的分類方法。概率模型即為貝葉斯概率公式，即：

其中，P（C|X）為條件X 下C 的后驗概率，P（C）為C的先驗概率，P（X|C）為條件C 下X 的后驗概率，P（X）表示X 的先驗概率［21］。本文中C 為類別變量，X 為樣本變量，則貝葉斯概率公式就是求解樣本屬于某一類別的概率，然后選擇概率值最大的類別作為模型分類結果。

3.3 模型實現

本文模型是基于目前流行的Python 機器學習庫sklearn 實現的。由于樣本數量少，因此本文隨機選擇90%的樣本作為訓練集，剩下的10% 作為測試集。在機器學習任務中，需要盡量保證正負樣本的均衡，才能使模型具有較好的泛化性能。因此，從整理好的SIDER 數據庫中選擇正負樣本相對均衡的5 種副作用應用于本文模型，得到的測試集準確率結果如圖1 所示（彩圖掃OSID 碼可見）。

Fig.1 Results of model prediction圖1 模型預測結果

從圖1 中可以很直觀地看出，在副作用預測任務中，隨機森林算法得到的準確率明顯優于K 近鄰算法，并且在5 種副作用測試結果中都是最優的，這也說明本文數據在隨機森林中的表現相比其它3 種算法更好。在副作用皮疹的測試結果中，隨機森林的準確率達到90.24%，而K 近鄰只有58.54%。在副作用口腔炎的測試結果中，隨機森林的準確率為78.05%，但也依然優于其它3 種算法。這一結果驗證了本文的研究思路，證明了隨機森林分類器優越的分類能力，也充分說明基因表達量變化情況與藥物副作用之間具有很強的相關性。

4 結語

本文以患病前后與用藥前后的基因表達量變化情況為研究對象，分析了基因表達量變化情況與藥物副作用的關聯關系，提出一種藥物副作用預測方法。由于數據集中某些數據的缺失，以及對數據進行的層層篩選都導致了本文數據集的大幅減少，這對模型的構建會產生影響，因此數據集的擴展是本文需要突破的一個方面。由于生物信息數據存在重復樣本的特點，若在數據預處理過程中遇到重復樣本時采用平均處理方法，以及兩次基因表達量測試是在不同平臺上進行的，這都會導致數據出現誤差。另外，本文只使用了基因表達這一單一特征以及簡單的機器學習算法，因此對特征的融合以及算法選擇也是本文下一步要突破的方向。