高 闖，唐 冕，趙 亮，2*

（1.廣西大學計算機與電子信息學院，南寧 530004；2.湖北醫(yī)藥學院太和醫(yī)院，湖北十堰 442000）

（?通信作者電子郵箱S080011@e.ntu.edu.sg）

0 引言

B 細胞表位是抗原表面上可以被抗體的免疫反應識別的特定區(qū)域［1］。識別抗原與抗體相結(jié)合的表位對現(xiàn)代藥物和疫苗的開發(fā)起到至關(guān)重要的作用［2-3］。通過實驗方法（如X 射線晶體學［4］）檢測表位雖然準確，但需要消耗大量的時間和資源，因此，探索有效并且可靠的計算方法進行抗原上的B細胞表位識別，具有重要的現(xiàn)實意義；同時，蛋白質(zhì)數(shù)據(jù)庫（Protein Data Bank，PDB）［5］中提供的抗原-抗體結(jié)構(gòu)數(shù)據(jù)也為預測算法的研究提供了幫助。

目前，已經(jīng)開發(fā)出許多可用于預測表位的計算方法。Kringelum 等［6］提出的方法使用了來自376 種抗原-抗體復合物的107 個抗原數(shù)據(jù)的結(jié)構(gòu)和幾何特征；Sun 等［7］通過分析161 種免疫球蛋白復合物的理化性質(zhì)和結(jié)構(gòu)特征，識別表位和非表位之間的差異。隨著機器學習算法的發(fā)展，這些技術(shù)也被用于設計表位預測算法。例如，Zhang等［8］結(jié)合常規(guī)特征和從3D 結(jié)構(gòu)中提取的特殊特征，將隨機森林算法用作分類器，設計B 細胞構(gòu)象表位的預測模型。此外基于圖的方法也被應用于表位預測，如CEP（Conformational Epitope Prediction server）［9］。但是，這些方法大部分只關(guān)注單個或多個分離的表位，對抗原中存在的重疊表位的預測效果不理想。與免疫反應相關(guān)的蛋白質(zhì)不止存在單個或多個分離的結(jié)合位點情況，還存在多個重疊的結(jié)合位點的情況［10］。Zhao 等［11］提出的模型已經(jīng)證明重疊的子圖挖掘方法可以提高對抗原表位的識別能力，特別是對于多個重疊表位；但是，該模型在子圖擴展階段需要設置閾值，這會降低模型的泛化能力。

針對現(xiàn)有表位預測方法對重疊表位預測能力不佳的問題，本文提出了一種基于局部度量（Local Metric，L-Metric）［12］的重疊子圖發(fā)現(xiàn)算法用于表位預測的模型SGLMEP（overlapping SubGraph mining based on L-Metric for Epitope Prediction）。該模型分為三個主要步驟：1）利用抗原-抗體復合物中的表面原子構(gòu)建氨基酸殘基圖；2）利用馬爾可夫聚類算法（Markov CLustering algorithm，MCL）［13］將氨基酸殘基圖劃分為互不重疊的種子子圖，并利用重疊子圖發(fā)現(xiàn)算法對種子子圖擴展以得到重疊子圖；3）利用圖卷積神經(jīng)網(wǎng)絡（Graph Convolutional neural Network，GCN）［14］和全連接網(wǎng)絡（Fully Connected Network，F(xiàn)CN）［15］構(gòu)建的分類器對子圖進行分類。

1 氨基酸殘基圖的構(gòu)建

本文采用兩個步驟構(gòu)建氨基酸殘基圖：第一，利用抗原鏈中氨基酸殘基的表面原子構(gòu)建原子圖；第二，將表面原子圖升級為氨基酸殘基圖。表面原子的定義是指可及表面積（Accessible Surface Area，ASA）不小于10 ?2的原子［16］。為了得到抗原主鏈和側(cè)鏈上的表面原子，本文使用工具NACCESS（Non-commercial atomic solvent ACCESSible area calculation）［17］計算每個原子的ASA，探針大小采用默認設置。根據(jù)表面原子的坐標構(gòu)建KD（K-Dimensional）樹［18］（利用空間劃分在k維空間中存儲點的一種數(shù)據(jù)結(jié)構(gòu)），計算由表面原子構(gòu)成的不大于閾值（設置為8 ?）的邊。通過此步驟可以得到抗原的原子圖表示。若兩個氨基酸殘基中的原子至少有一組邊，則定義兩個殘基存在邊。按照此規(guī)則，將原子圖表示升級為氨基酸殘基圖。當兩個氨基酸殘基的側(cè)鏈最外層原子之間存在邊時，則兩個殘基也一定存在邊，因此，為了減少計算量，在構(gòu)建原子圖表示時僅利用主鏈的表面原子與側(cè)鏈的最外層原子。

2 基于L-Metric的重疊子圖發(fā)現(xiàn)算法

對重疊社區(qū)的檢測已被廣泛用于發(fā)現(xiàn)社交網(wǎng)絡中的社區(qū)挖掘［19-20］。本文基于子圖劃分和社區(qū)挖掘思想提出了基于L-Metric的重疊子圖發(fā)現(xiàn)算法，包括子圖發(fā)現(xiàn)和子圖擴展兩個階段。

2.1 子圖發(fā)現(xiàn)

MCL 是一種子圖劃分的算法，它基于消息流傳遞發(fā)現(xiàn)圖中成員的不重疊分組［13］，因此，首先使用該算法將氨基酸殘基圖劃分為互不重疊的種子子圖。MCL的機制是模擬帶權(quán)重的無向圖中的消息流傳遞過程，需要輸入由邊權(quán)重組成的轉(zhuǎn)移矩陣。本文基于Zhao 等［16］提出的邊界邊思想，設計邊權(quán)重的計算公式。

本文計算氨基酸殘基（組成生命體中蛋白質(zhì)的20 種氨基酸）r1和r2構(gòu)成的邊的權(quán)重W(r1，r2)，由兩部分組成，即由r1和r2構(gòu)成的邊在表位和非表位中的頻率的χ2檢驗和對數(shù)函數(shù)生成，具體計算式為：

其中：α為超參數(shù)，通過網(wǎng)格搜索法優(yōu)化后為0.4；Wχ2(r1，r2)和的計算式如式（2）、（3）。

其中：c∈{c1，c2}，c1為表位殘基，c2為非表位殘基為在殘基圖x中標簽是c的殘基r1和r2構(gòu)成的邊的出現(xiàn)頻率；為標簽是c的殘基r1和r2構(gòu)成的邊出現(xiàn)頻率的期望。的計算式為：

其中：P為由氨基酸殘基圖所組成的數(shù)據(jù)集；的計算如式（5）。

其中：γ和φ為超參數(shù)，通過網(wǎng)格搜索法優(yōu)化后分別為4和1。

在殘基圖中，邊界邊（由表位殘基和非表位殘基構(gòu)成的邊）的存在會影響對表位與非表位殘基的鑒別。利用上述公式，計算在c1為邊界殘基和c2為表位殘基、c1為邊界殘基和c2為非表位殘基的兩種情況下邊的權(quán)重，這兩個權(quán)重中較大的為邊的新權(quán)重。當權(quán)重不小于θ（超參數(shù)，通過網(wǎng)格搜索法優(yōu)化后為0.5）時，則這條邊被當作邊界邊從殘基圖中刪除。

2.2 子圖擴展

L-Metric可以度量局部社區(qū)網(wǎng)絡中的成員關(guān)系緊密度，它的優(yōu)點是不需要全局網(wǎng)絡信息，且不需要設置任何參數(shù)。本文基于L-Metric對種子子圖進行擴展，以檢測重疊子圖。

圖1 描述了在氨基酸殘基圖G中的一個子圖的定義。D是G的一個子圖內(nèi)的點構(gòu)成的集合，它被分為兩部分：邊界點集合B（與D外的點存在連接關(guān)系的點構(gòu)成的集合）和核心點集合C（僅與D中的點存在連接關(guān)系的點構(gòu)成的集合）。鄰居點集合N是D外與邊界點存在連接關(guān)系的點構(gòu)成的集合。

圖1 子圖描述Fig.1 Description of subgraph

本文中度量L的計算方法為：

式中，Lin和Lex分別為內(nèi)部和外部的度量，具體計算式為：

式中：Win和Wex為子圖內(nèi)部和外部權(quán)重；W(ri，rj)為由式（1）計算得到的邊權(quán)重。

式中：Ein為由D中的點構(gòu)成的邊集；Eex為B中的點與N中的點構(gòu)成的邊集。

因此，當新加入的點滿足第1）種和第2）種情況時，該點被加入后形成新的D。當沒有其他點被添加進D時，則得到擴展后的重疊子圖。擴展后的重疊子圖包含表位殘基和非表位殘基。本文將表位殘基數(shù)不少于3 個的重疊子圖標記為表位子圖，其他標記為非表位子圖。

3 基于深度學習的子圖分類模型

3.1 模型總體架構(gòu)

抗原表位預測的難點之一是特征選擇問題［21］。深度學習是一種無特征學習算法，它能夠自動提取特征，降低因?qū)I(yè)知識不足而導致忽略某些重要特征的可能性。本文構(gòu)建了基于深度學習的子圖分類模型，總體架構(gòu)如圖2 所示，主要由基于圖卷積神經(jīng)網(wǎng)絡（GCN）構(gòu)建的特征提取器和基于全連接網(wǎng)絡（FCN）構(gòu)建的分類器兩部分構(gòu)成。

圖2 分類模型總體架構(gòu)Fig.2 Overall architecture of classification model

圖卷積神經(jīng)網(wǎng)絡（GCN）是卷積神經(jīng)網(wǎng)絡（Convolutional Neural Network，CNN）在圖數(shù)據(jù)上的擴展，它可以有效地利用圖數(shù)據(jù)的全局信息表示圖節(jié)點的特征。本文構(gòu)建的特征提取器包含兩層GCN。假設GCN 的輸入子圖G=(A，Z)，A∈Rn×n表示子圖帶權(quán)重的鄰接矩陣，Z∈Rn×d表示節(jié)點的特征矩陣，其中，n表示圖中節(jié)點的數(shù)量，d表示節(jié)點特征向量的維度。第i∈{1，2}層的GCN通過式（13）得到第i層節(jié)點特征矩陣：

由兩層FCN 構(gòu)成的分類器，以特征提取器得到的節(jié)點特征矩陣為輸入。第j∈{1，2}層FCN的映射公式為：

分類器得到子圖節(jié)點的類別屬性（表位殘基或者非表位殘基），然后利用子圖類別定義（子圖中含有不少于3 個表位殘基的子圖為表位子圖）對子圖進行分類。

3.2 基于代價敏感的損失函數(shù)

由于表位與非表位數(shù)據(jù)的數(shù)量極不平衡，表位子圖的數(shù)量遠遠少于非表位子圖的數(shù)量。在本文使用的訓練數(shù)據(jù)中，表位子圖的數(shù)量與非表位子圖的數(shù)量之比為1∶6。為了避免與此問題相關(guān)的性能損失，本文在訓練過程使用焦點損失函數(shù)（Focal Loss function，F(xiàn)L）［23］作為網(wǎng)絡結(jié)構(gòu)的損失函數(shù)，計算式為：

式中，λ為超參數(shù)，通過網(wǎng)格搜索優(yōu)化后為2。pt和αt的計算式如下：

式中：p∈[0，1]為正類的概率；α∈[0，1]為正類在訓練樣本中所占比例；y∈{-1，+1}為子圖標簽類型（其中，-1 表示負類即非表位子圖，+1表示正類即表位子圖）。

4 實驗和結(jié)果分析

4.1 實驗數(shù)據(jù)

本文從PDB 中檢索到808種抗體-抗原復合物，利用通用的數(shù)據(jù)選擇標準［10，16］從其中篩選出397 種復合物。然后，使用 工 具 CD-HIT（Cluster Database at High Identity with Tolerance）［24］去除序列相似性不小于0.9 的冗余數(shù)據(jù)以降低噪聲。最終，獲得了由254個復合物組成的實驗數(shù)據(jù)。

4.2 模型評估及分析

由于非表位子圖數(shù)據(jù)的數(shù)量遠多于表位子圖數(shù)據(jù)，為了實驗的準確性，本文采用F1值、召回率Re（Recall）和精確度Pre（Precision）作為模型表位子圖分類效果的評估指標，計算式為：

式中：TP（True Positive）表示實際為表位子圖被正確預測的數(shù)量；FN（False Negative）表示實際為表位子圖被錯誤預測的數(shù)量；FP（False Positive）表示實際為非表位子圖被錯誤預測的數(shù)量。

選擇當前主流的表位預測模型DiscoTope 2（Discontinuous epiTope prediction 2）［25］、ElliPro（Ellipsoid and Protrusion）［26］、EpiPred（Epitope Prediction server）［27］和Glep（overlapping Graph clustering-based B-cell epitope predictor）［10］與本文所提出的模型SGLMEP 進行比較。圖3 描述了各個表位預測模型的F1值分布。從圖3 中可以看出，SGLMEP 的F1值優(yōu)于其他模型，且擁有更高的性能下限。

圖3 表位預測模型的F1值比較Fig.3 F1-score comparison of epitope prediction models

表1 給出了各表位預測模型的平均F1值、召回率和精確度。與表位預測模型DiscoTope 2、ElliPro、EpiPred 和Glep 相比，本文所提出的模型SGLMEP 將平均F1值分別提高了267.3%、57.0%、65.4% 和3.5%。從表1 可以看出，模型SGLMEP 主要的優(yōu)勢在于擁有更高的召回率，與Glep 相比提高了18.3%。實驗結(jié)果說明，SGLMEP 能夠識別出更多的潛在表位殘基，對重疊表位的識別優(yōu)于Glep。雖然SGLMEP 的精確度不是十分突出，但仍然高于大多數(shù)模型。

表1 不同預測模型結(jié)果對比Tab.1 Result comparison of different prediction models

4.3 消融實驗結(jié)果及分析

本文所提出的模型SGLMEP 的核心是基于L-Metric 的重疊子圖發(fā)現(xiàn)算法。通過在相同實驗數(shù)據(jù)上進行消融實驗，驗證該算法對表位預測性能的影響。

通過圖4可以直觀地看出，基于L-Metric的重疊子圖發(fā)現(xiàn)算法的F1值高于未使用的情況。實驗結(jié)果表明，該重疊子圖發(fā)現(xiàn)算法對表位預測性能的提升是有效的。應用該算法的模型的優(yōu)勢主要是提高了召回率，即對未知表位殘基的識別能力。

表2 給出了消融實驗的詳細結(jié)果。通過表2 中數(shù)據(jù)可以看出，應用了基于L-Metric 的重疊子圖發(fā)現(xiàn)算法的模型相較于未使用子圖發(fā)現(xiàn)算法時的平均F1值、召回率分別提高了19.2%和38.9%。

表2 消融實驗的詳細結(jié)果Tab.2 Detailed results of ablation experiment

5 結(jié)語

為提高重疊表位預測性能，本文提出了基于L-Metric 的重疊子圖發(fā)現(xiàn)算法，自適應對種子子圖進行擴展。同時，利用深度學習算法自學習特征的特點設計分類模型，降低人工設計特征構(gòu)建分類器的難度。通過實驗驗證了本文所提出的表位預測模型SGLMEP 對抗原表位具有良好的預測性能。此外，通過消融實驗驗證了本文所提出的重疊子圖算法的有效性。在今后的研究工作中，將對重疊子圖發(fā)現(xiàn)算法進行優(yōu)化，旨在發(fā)現(xiàn)更多潛在表位的基礎(chǔ)上，進一步提高對表位識別的精確程度。