杜 瑤

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

0 引言

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

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

1 預(yù)測方法

1.1 原理

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

1.2 工作內(nèi)容

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

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

2 特征工程

2.1 數(shù)據(jù)收集及預(yù)處理

GEO 和cmap 兩個基因表達(dá)數(shù)據(jù)庫是本文的核心數(shù)據(jù)庫，其中基因受到疾病影響所產(chǎn)生的基因表達(dá)量變化情況來自GEO 數(shù)據(jù)庫（https：//www.ncbi.nlm.nih.gov/），這是美國國立衛(wèi)生研究院1988 年創(chuàng)立的美國國立生物技術(shù)信息中心，目的是給分子生物學(xué)家提供一個信息存儲與處理的系統(tǒng)。在數(shù)據(jù)庫中找到所需的樣本組，這里選擇的是12 個乳腺細(xì)胞樣本，分為6 個正常對照組和6 個患病實驗組。Value 值即為樣本某個探針?biāo)鶎?yīng)基因的表達(dá)量，如表1所示?；虻牟町惐磉_(dá)分析是利用R 語言的Limma 包實現(xiàn)的，這是一個專門用于基因差異分析的包［13］。

Table 1 Breast cell sample data表1 乳腺細(xì)胞樣本數(shù)據(jù)

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

Table 2 Results of difference analysis表2 差異分析結(jié)果

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

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

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

2.2 特征構(gòu)建

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

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

Table 3 Results of data processing表3 數(shù)據(jù)處理結(jié)果

3 分類器選擇

3.1 選擇依據(jù)

構(gòu)建好模型特征之后，選擇一種合適且高效的算法也是本文工作的重要內(nèi)容，目前最流行的是使用機器學(xué)習(xí)算法建立模型進(jìn)行分類預(yù)測等任務(wù)。如今，機器學(xué)習(xí)、深度學(xué)習(xí)等技術(shù)已得到了廣泛應(yīng)用，相關(guān)算法也因不斷優(yōu)化使得代碼更加通俗易懂、易于實現(xiàn)，可根據(jù)數(shù)據(jù)類型、數(shù)據(jù)量大小以及任務(wù)內(nèi)容進(jìn)行算法選擇。

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

3.2 分類器介紹

3.2.1 隨機森林

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

3.2.2 K 近鄰

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

3.2.3 決策樹

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

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 為樣本變量，則貝葉斯概率公式就是求解樣本屬于某一類別的概率，然后選擇概率值最大的類別作為模型分類結(jié)果。

3.3 模型實現(xiàn)

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

Fig.1 Results of model prediction圖1 模型預(yù)測結(jié)果

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

4 結(jié)語

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