基于改進的PCA和ISSA－BPNN的定量構(gòu)效關(guān)系預(yù)測模型

陳 強，王登文，鐵治欣，2，洪 亮

（1浙江理工大學 信息學院，杭州 310018；2浙江理工大學 科技藝術(shù)學院，浙江 紹興 312369；3浙江傳媒學院 媒體工程學院，杭州 310018）

0 引 言

據(jù)近幾年全球癌癥統(tǒng)計數(shù)據(jù)表明，乳腺癌癥發(fā)病率和死亡病例逐漸增加，其防治須引起人們高度重視。臨床、流行病學和生物學證據(jù)表明，雌激素參與了乳腺癌的發(fā)生和發(fā)展。雌激素化合物的大多數(shù)生理功能，在基因調(diào)控水平上主要由雌激素受體（ER）調(diào)節(jié)，這些蛋白質(zhì)在細胞核中發(fā)揮作用，控制著各種器官系統(tǒng)的關(guān)鍵生理功能，并通過與相關(guān)的DNA調(diào)控序列相結(jié)合，來調(diào)節(jié)特定靶基因的轉(zhuǎn)錄。雌激素受體α亞型（Estrogen receptors alpha，ERα）在乳腺癌病中起著至關(guān)重要的作用，但在正常乳腺上皮細胞中極少被表達。通過使用選擇性雌激素受體調(diào)節(jié)劑（SERM）和雌激素受體降解劑（SERD），可用來降低ERα的穩(wěn)定性。

目前，在藥物研發(fā)中，為了節(jié)約時間和成本，通常采用建立化合物活性預(yù)測模型的方法，來篩選潛在活性化合物。這種定量構(gòu)效關(guān)系（Quantitative Structure Activity Relationship，QSAR）方法是一種預(yù)選工具，旨在減少化合物的數(shù)量，并增加選擇候選藥物的可能性。其以一系列分子結(jié)構(gòu)描述符作為自變量，化合物的生物活性作為因變量建立模型，根據(jù)可測量的物理、化學參數(shù)，精確預(yù)測化合物的生物活性，或者對已有活性化合物的結(jié)構(gòu)進行優(yōu)化，QSAR本質(zhì)上是數(shù)據(jù)驅(qū)動模型。近年來，人工智能、機器學習、大數(shù)據(jù)等技術(shù)的發(fā)展，為QSAR帶來了挑戰(zhàn)和機遇，通過成千上萬的化學結(jié)構(gòu)數(shù)據(jù)集，為藥物的生物活性和安全性進行更精確的回歸和分類預(yù)測帶來了可能，對推動中國化學品的管理有著重要的意義。

QSAR預(yù)測模型主要分為基于統(tǒng)計分析方法的預(yù)測模型和基于機器學習算法的預(yù)測模型。例如：El Ghalia Hadaji以多元線性回歸構(gòu)建QSAR預(yù)測模型；Afaf Zekri以多元線性逐步回歸構(gòu)建QSAR預(yù)測模型；Lu Yang基于遺傳算法的多元線性回歸構(gòu)建QSAR預(yù)測模型；Svetnik Vladimir以隨機森林算法構(gòu)建QSAR預(yù)測模型；代志軍以支持向量機回歸構(gòu)建QSAR預(yù)測模型；楊杰元以BP神經(jīng)網(wǎng)絡(luò)算法構(gòu)建QSAR預(yù)測模型；Li Jingshan以梯度下降樹決策樹（GBDT）構(gòu)建QSAR預(yù)測模型。雖然或多或少實現(xiàn)了預(yù)測，但是基于統(tǒng)計分析的方法隨著變量急劇增多也變得束手無策。為了提高基于機器學習算法的預(yù)測精度，本文提出了基于改進的PCA和ISSA－BPNN的預(yù)測模型。

1 相關(guān)預(yù)測方法

1.1 BP神經(jīng)網(wǎng)絡(luò)預(yù)測算法

BP神經(jīng)網(wǎng)絡(luò)（BPNN）結(jié)構(gòu)簡單，使用方便，非循環(huán)多級網(wǎng)絡(luò)訓練算法，使其具有廣泛的實用性，能夠?qū)崿F(xiàn)輸入到輸出的非線性映射。BPNN是單向傳播的多層前向神經(jīng)網(wǎng)絡(luò)（結(jié)構(gòu)如圖1所示），由輸入層（個節(jié)點）、輸出層（個節(jié)點）和多個隱含層組成。

圖1 BP神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)圖Fig.1 BP neural network structure

1.2 SVR預(yù)測算法

支持向量機回歸（SVR）是將支持向量機分類（SVM）算法應(yīng)用于回歸預(yù)測中，兩者不同的是：SVM將間隔之內(nèi)的空間樣本算入損失函數(shù)中，以達到分類的目的；而SVR則是將間隔之外的空間樣本算入損失函數(shù)中，以達到回歸的目的。對于非線性SVR模型，使用核函數(shù)將數(shù)據(jù)映射到高維空間，而后進行回歸預(yù)測。由于徑向基核函數(shù)（RBF）應(yīng)用廣泛且具有較好的回歸效果，因此本文選擇RBF作為SVM分析的核函數(shù)。

1.3 XGBoost預(yù)測算法

XGBoost（Extreme Gradient Boosting）是在Boosting算法基礎(chǔ)上進行改良的，在預(yù)測精度以及訓練速度方面有較大的突破，屬于GBDT的范疇，并且也是一種前向特征的算法，本質(zhì)上是由許多回歸和分類的決策樹組成。XGBoost相較于GBDT而言：前者加入正則項防止過擬合，對目標損失函數(shù)進行二階泰勒展開，從而增加了精度，根據(jù)最佳切分點進行葉子節(jié)點分裂優(yōu)化計算，從而優(yōu)化結(jié)果。

2 QSAR模型

本文實驗數(shù)據(jù)集源自乳腺癌治療靶標ERα時，得到的1 974個化合物作為ERα生物活性數(shù)據(jù)樣本。其中包括729個分子描述符信息和（實際QSAR建模中，一般采用來表示生物活性值，即因變量），值越大表明生物活性越高。

由于變量的數(shù)量比較多，本文首先提出基于改進的PCA特征選擇算法，對模型的輸入變量進行篩選，然后提出ISSA－BPNN算法對BPNN算法進行改進。

2.1 基于改進PCA的特征提取

改進的PCA算法流程如圖2所示。首先對數(shù)據(jù)進行標準化，然后在729個分子描述符信息中，用基于Pearson、MIC和RF的加權(quán)得分算法得到前20個特征變量，最后基于PCA算法提取4個新特征代替原特征，作為模型的主要輸入變量。

圖2 改進的PCA算法流程Fig.2 Improved PCA algorithm

2.1.1 最大互信息系數(shù)法（MIC）

MIC是一種通過繪制變量散點圖計算兩個變量的互信息，來衡量變量間關(guān)聯(lián)程度的算法。其實現(xiàn)步驟如下：

（1）散點圖網(wǎng)格化，計算互信息值。給定個有序?qū)?shù)據(jù)集（，），將數(shù)據(jù)集劃分為的網(wǎng)格，方向和方向的網(wǎng)格數(shù)分別為、?；バ畔⒅档挠嬎闳缡剑?）：

式中，（，）為與之間的聯(lián)合概率密度，（）和（）分別為和的邊緣概率密度。

（2）互信息值歸一化，如式（2）：

（3）變換網(wǎng)格劃分情況，選擇不同尺度下互信息的最大值作為值，如式（3）：

由文獻［20］可知，當（）＝時，效果最好。

2.1.2 隨機森林（RF）特征選擇法

隨機森林（Random Forest，RF）實質(zhì)是一個包含多個決策樹的組合分類器。其通過特征隨機置換前后的誤差分析，計算每個特征重要度得分，分值越高，特征越重要，從而進一步確定特征排序。隨機森林結(jié)合Bootstrap重采樣技術(shù)和決策樹，構(gòu)建一個包含多個基本分類器的樹型分類器集合，采用簡單多數(shù)投票的方法得到結(jié)果。

RF中決策樹數(shù)目為N，原始數(shù)據(jù)集有個特征，單特征X（1，2，…，）基于誤差分析的特征重要性度量，按以下步驟計算：

（1）計算第棵決策樹相應(yīng)的袋外數(shù)據(jù)OOB的袋外錯誤樣本數(shù)；

（3）重復(fù)步驟（1）、（2）得到：

（4）由式（4）計算特征X的重要性得分。

(4)要加強地質(zhì)人才培養(yǎng)，引進人才激勵機制。地質(zhì)人才是理論創(chuàng)新、技術(shù)變革的重要基礎(chǔ)，也是進行深部找礦的關(guān)鍵，為此我們應(yīng)該加強地質(zhì)人才的培養(yǎng)，建立健全人才激勵機制，鼓勵更多的地質(zhì)人才去為深部找礦發(fā)展做貢獻。一方面要采取相應(yīng)的鼓勵機制，鼓勵更多的人報考地質(zhì)領(lǐng)域的相關(guān)專業(yè)。另一方面要保障地質(zhì)工作者的各方面待遇，在改善他們生活條件的同時也要加強隊伍培養(yǎng)和人才業(yè)務(wù)能力建設(shè)。鼓勵更多專業(yè)素養(yǎng)過硬、有吃苦耐勞精神的年輕人加入，同時也要完善相關(guān)的人才激勵機制，加強技能培訓力度，更好地激發(fā)他們的積極性和工作熱情，不斷提高他們的實踐能力。只有人才得到保證才能使我國的地質(zhì)事業(yè)更好地可持續(xù)地發(fā)展下去。

2.1.3 基于Pearson、MIC和RF的加權(quán)得分算法

由于各變量的數(shù)值量綱之間存在較大差異，為了消除量綱的影響，需要對數(shù)據(jù)進行標準化處理。本文采用Z－score標準化方法，對變量進行歸一化處理，如式（5）：

Pearson和MIC反映了自變量與因變量之間的線性和非線性關(guān)系，而RF是以特征重要度計算值來表示自變量與因變量的相關(guān)性。加權(quán)得分由式（6）計算得到：

其中，g r a d e表示第（1，2，3，…，729）個分子描述符的加權(quán)分；P表示第個自變量與因變量的Pearson系數(shù)絕對值；MIC表示第個自變量與因變量的最大互信息系數(shù)絕對值；RF表示第個自變量與因變量的特征重要度計算值，和均應(yīng)在0和1之間（本文取025）。

由式（6）計算得到20個主要特征變量見表1。

表1 加權(quán)得分分子描述符顯著性排序Tab.1 Significance ranking of weighted score molecular descriptors

2.1.4 確定模型輸入變量

PCA算法的原理是以原始特征的線性組合方式，得到新特征來代替原特征，從而達到降維的效果。根據(jù)方差越大新特征越重要的原則，對個主成分按照貢獻率進行排序，再從中提取個主成分來代表全部數(shù)據(jù)，最后將新特征作為QSAR模型的輸入值。算法流程如下：

（1）計算數(shù)據(jù)的協(xié)方差矩陣。假設(shè)原始數(shù)據(jù)集為，其協(xié)方差矩陣記為；

（3）計算累計貢獻率并確定主成分個數(shù)。

將表1中的20個特征變量由PCA算法特征提取后，得到新特征的貢獻率見表2。

表2 新特征累計貢獻率Tab.2 Cumulative variance contribution rate of new features

由表2可知，前4個新特征已包含原始特征95%以上的信息。故本文取前4個新特征代替原特征作為模型的主要輸入變量。

2.2 ISSA－BPNN

傳統(tǒng)的BPNN對權(quán)值和閾值較敏感，存在收斂速度慢和極易陷入局部最優(yōu)的問題。因此，本文通過改進的麻雀搜索算法（improved sparrow search algorithm，ISSA）來優(yōu)化BPNN的權(quán)值和閾值。

SSA是根據(jù)麻雀覓食并逃避捕食者的行為而提出的群智能優(yōu)化算法，其模擬了麻雀群覓食的過程。在SSA中有3種狀態(tài)，分別是發(fā)現(xiàn)者、加入者、偵察者。其中，適應(yīng)度值較好的發(fā)現(xiàn)者是為了獲得食物的同時，為所有加入者提供覓食的方向；偵察者選擇安全第一為目標，在發(fā)現(xiàn)危險的情況下，提醒種群放棄食物。

由于SSA容易陷入局部最優(yōu)，且全局搜索能力較弱，可將SSA中發(fā)現(xiàn)者和加入者位置更新公式分別改為式（7）、（8）。加入者以一定概率向發(fā)現(xiàn)者靠攏，保證了全局收斂。同時，后加入的麻雀要盡快飛到其他區(qū)域覓食。

其中，代表當前迭代次數(shù)；（0，1）和是服從標準正態(tài)分布的隨機數(shù)；是1的矩陣，代表維度；x是第個麻雀在第維的位置；∈［0，1］代表預(yù)警值；∈［05，1］代表安全值。

當≥時，表示發(fā)現(xiàn)者已經(jīng)發(fā)現(xiàn)捕食者，此時種群內(nèi)其它麻雀盡可能飛到其它安全地方進行覓食；當時，發(fā)現(xiàn)者可以廣泛搜索。是種群規(guī)模，x是當前全局最差的位置，x是當前發(fā)現(xiàn)者的位置，∈［1，1］表示加入者跟隨生產(chǎn)者尋找食物的概率。當2時，表示適應(yīng)度值較差的第個加入者處于挨餓狀態(tài)，需要盡快飛到其它區(qū)域繼續(xù)尋找食物來獲得能量。

偵察者的位置更新如式（9）：

式中，是［1，1］范圍內(nèi)的一個隨機數(shù)；是步長控制參數(shù)，其服從標準正態(tài)分布的隨機數(shù)；x表示當前的全局最佳位置；f、f和f分別代表當前麻雀的全局最差、全局最優(yōu)和個體適應(yīng)度。分母加上一個常數(shù)量，是為了防止分母出現(xiàn)0的情況。

本文提出的ISSA－BPNN流程如圖3所示，其實現(xiàn)步驟為：

圖3 ISSA－BPNN流程Fig.3 ISSA－BPNN flow chart

（1）初始化麻雀搜索算法；

（2）計算麻雀種群個體適應(yīng)度，并得到最佳位置、最差位置和最佳適應(yīng)度值、最差適應(yīng)度值；

（3）根據(jù)式（7）～（9）分別更新發(fā)現(xiàn)者、加入者和偵查者的位置信息，并更新適應(yīng)度值；

（4）若算法達到最大迭代次數(shù)或達到最初設(shè)定的收斂精度，則執(zhí)行步驟（5），否則返回步驟（2）；

（5）將得到的最優(yōu)值賦給BPNN的權(quán)值和閾值；

（6）使用BPNN進行學習，不斷調(diào)整直至達到訓練終止條件，最終實現(xiàn)預(yù)測輸出。

3 實驗結(jié)果與分析

依據(jù)上述方法對數(shù)據(jù)進行新特征選取后，將1 974個樣本按照7：3的比例劃分訓練集和測試集。訓練集用來擬合模型，測試集用來對模型的性能進行評價。驗證本文所提出模型的有效性，分別利用SVR、XGBoost、BPNN和ISSA－BPNN模型對前述數(shù)據(jù)集合進行預(yù)測。

3.1 預(yù)測模型的評價指標

本文采用平均絕對誤差（）、平均絕對百分比誤差（）和均方根誤差（）評價模型的預(yù)測精度。其計算公式分別為式（10）（12）：

3.2 結(jié)果分析

4種模型的預(yù)測值與真實值曲線對比如圖4所示，預(yù)測精度對比結(jié)果見表3。

圖4 4種模型的pIC50預(yù)測值與真實值對比Fig.4 Comparison of predicted pIC50 values and true values of four models

表3 4種模型預(yù)測精度對比Tab.3 Comparison of prediction accuracy of four models

由表3可知，ISSA－BPNN模型的、、均是最低的，表明ISSA－BPNN預(yù)測誤差值最小、穩(wěn)定性最高、效果最佳。其中，ISSA－BPNN模型的值較SVR模型提高了13.10%，較XGBoost模型提高了10.53%，較BPNN模型提高了20.22%。

4 結(jié)束語

為了更精確地預(yù)測化合物的生物活性，本文提出了一種基于改進的PCA和改進的麻雀搜索算法優(yōu)化BP神經(jīng)網(wǎng)絡(luò)（ISSA－BPNN）預(yù)測模型，其具有良好的尋優(yōu)能力。

算法中，利用改進的PCA算法提取模型的主要變量，再利用ISSA優(yōu)化BPNN的權(quán)值和閾值，改善了BPNN易陷入局部極值的缺點。通過實驗對比結(jié)果表明，基于ISSA－BPNN預(yù)測模型的預(yù)測精度最高，并具有較強的擬合能力和泛化能力。但是，由于訓練的數(shù)據(jù)量較少，導致模型的預(yù)測精度不是太高，后期研究可增加訓練數(shù)據(jù)來提高模型的預(yù)測精度。