陳 強(qiáng)，王登文，鐵治欣，2，洪 亮

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

0 引 言

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

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

QSAR預(yù)測模型主要分為基于統(tǒng)計(jì)分析方法的預(yù)測模型和基于機(jī)器學(xué)習(xí)算法的預(yù)測模型。例如：El Ghalia Hadaji以多元線性回歸構(gòu)建QSAR預(yù)測模型；Afaf Zekri以多元線性逐步回歸構(gòu)建QSAR預(yù)測模型；Lu Yang基于遺傳算法的多元線性回歸構(gòu)建QSAR預(yù)測模型；Svetnik Vladimir以隨機(jī)森林算法構(gòu)建QSAR預(yù)測模型；代志軍以支持向量機(jī)回歸構(gòu)建QSAR預(yù)測模型；楊杰元以BP神經(jīng)網(wǎng)絡(luò)算法構(gòu)建QSAR預(yù)測模型；Li Jingshan以梯度下降樹決策樹（GBDT）構(gòu)建QSAR預(yù)測模型。雖然或多或少實(shí)現(xiàn)了預(yù)測，但是基于統(tǒng)計(jì)分析的方法隨著變量急劇增多也變得束手無策。為了提高基于機(jī)器學(xué)習(xí)算法的預(yù)測精度，本文提出了基于改進(jìn)的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ò)訓(xùn)練算法，使其具有廣泛的實(shí)用性，能夠?qū)崿F(xiàn)輸入到輸出的非線性映射。BPNN是單向傳播的多層前向神經(jīng)網(wǎng)絡(luò)（結(jié)構(gòu)如圖1所示），由輸入層（個(gè)節(jié)點(diǎn)）、輸出層（個(gè)節(jié)點(diǎn)）和多個(gè)隱含層組成。

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

1.2 SVR預(yù)測算法

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

1.3 XGBoost預(yù)測算法

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

2 QSAR模型

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

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

2.1 基于改進(jìn)PCA的特征提取

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

圖2 改進(jìn)的PCA算法流程Fig.2 Improved PCA algorithm

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

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

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

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

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

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

由文獻(xiàn)［20］可知，當(dāng)（）＝時(shí)，效果最好。

2.1.2 隨機(jī)森林（RF）特征選擇法

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

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

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

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

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

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

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

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

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

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

由式（6）計(jì)算得到20個(gè)主要特征變量見表1。

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

2.1.4 確定模型輸入變量

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

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

（3）計(jì)算累計(jì)貢獻(xiàn)率并確定主成分個(gè)數(shù)。

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

表2 新特征累計(jì)貢獻(xiàn)率Tab.2 Cumulative variance contribution rate of new features

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

2.2 ISSA－BPNN

傳統(tǒng)的BPNN對權(quán)值和閾值較敏感，存在收斂速度慢和極易陷入局部最優(yōu)的問題。因此，本文通過改進(jìn)的麻雀搜索算法（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)者是為了獲得食物的同時(shí)，為所有加入者提供覓食的方向；偵察者選擇安全第一為目標(biāo)，在發(fā)現(xiàn)危險(xiǎn)的情況下，提醒種群放棄食物。

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

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

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

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

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

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

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

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

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

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

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

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

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

3 實(shí)驗(yàn)結(jié)果與分析

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

3.1 預(yù)測模型的評價(jià)指標(biāo)

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

3.2 結(jié)果分析

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

圖4 4種模型的pIC50預(yù)測值與真實(shí)值對比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ù)測化合物的生物活性，本文提出了一種基于改進(jìn)的PCA和改進(jìn)的麻雀搜索算法優(yōu)化BP神經(jīng)網(wǎng)絡(luò)（ISSA－BPNN）預(yù)測模型，其具有良好的尋優(yōu)能力。

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