基于快速多示例多標(biāo)記學(xué)習(xí)的G蛋白偶聯(lián)受體生物學(xué)功能預(yù)測

吳建盛 馮巧遇 袁京洲 胡海峰 周家特 高 昊

1(南京郵電大學(xué)地理與生物信息學(xué)院 南京 210023)2 (南京郵電大學(xué)通信與信息工程學(xué)院 南京 210003) (jansen@njupt.edu.cn)

G蛋白偶聯(lián)受體(G protein-coupled receptors，GPCRs)具有7個跨膜結(jié)構(gòu)域，是細(xì)胞信號傳導(dǎo)過程中的重要蛋白質(zhì)，它在人類視覺、嗅覺、味覺、神經(jīng)傳遞、新陳代謝以及免疫調(diào)節(jié)等各項正常生理活動和疾病過程中都發(fā)揮著重要的作用[1].另外，大約40%的現(xiàn)代藥物都以GPCRs作為靶點[2].準(zhǔn)確了解GPCRs的生物學(xué)功能是理解它們參與的生物學(xué)過程及其藥物作用機制的關(guān)鍵.

蛋白質(zhì)通常含有多個保守結(jié)構(gòu)域，同時具有多種生物學(xué)功能，每個結(jié)構(gòu)域可以獨立或者與相鄰結(jié)構(gòu)域相互作用執(zhí)行生物學(xué)功能.我們以前的研究發(fā)現(xiàn)，蛋白質(zhì)功能預(yù)測可抽象為多示例多標(biāo)記學(xué)習(xí)(multi-instance multi-label learning, MIML)問題，每個蛋白質(zhì)對應(yīng)于多示例多標(biāo)記學(xué)習(xí)中的一個樣本對象，每個結(jié)構(gòu)域?qū)?yīng)于一個示例，每個功能對應(yīng)于一個標(biāo)記[3].

Zhang & Zhou(2006)提出了多示例多標(biāo)記學(xué)習(xí)(MIML)框架[4].在該學(xué)習(xí)框架下，每個訓(xùn)練樣本由多個示例表示且同時隸屬于多個概念標(biāo)記.近年來，多示例多標(biāo)記學(xué)習(xí)模型由于其強大的表示能力得到了研究者的廣泛關(guān)注，各種算法被提了出來.其中的代表有基于退化框架下的算法MIMLSVM[4],MIMLBOOST[4],MIMLNN[5]和SISL-MIML[6]；基于正則化機制的算法D-MIMLSVM算法[5]和M3MIML算法[7]；將單示例數(shù)據(jù)恢復(fù)成MIML形式數(shù)據(jù)的方法INSDIF[8]；生成式方法DBA[9]；MIML距離度量學(xué)習(xí)[10]；MIML示例預(yù)測算法Rankloss-Sim[11].

在過去的研究中，有很多的計算學(xué)方法用來預(yù)測蛋白質(zhì)生物學(xué)功能[12].在這些預(yù)測方法中，用于表征蛋白質(zhì)的屬性很多，有氨基酸序列信息、進化信息、基因組分信息、蛋白質(zhì)-蛋白質(zhì)相互作用信息、蛋白質(zhì)結(jié)構(gòu)信息、基因芯片信息等[12].蛋白質(zhì)的生物學(xué)功能有多種描述方式，其中基因本體學(xué)(gene ontology， GO)使用最為廣泛[13].基因本體學(xué)從3個方面來描述蛋白質(zhì)的生物學(xué)功能，即分子功能(molecular function)、生物學(xué)過程(biological process)和細(xì)胞組分(cellular component)[13].GPCRs一般都位于細(xì)胞膜上，因此，本文中我們只考慮GPCRs的分子功能和生物學(xué)過程.

本文中，我們利用快速多示例多標(biāo)記學(xué)習(xí)方法MIMLfast[14]，基于新的混合特征，對GPCRs的GO分子功能和生物學(xué)過程進行了預(yù)測.實驗結(jié)果表明，本文的模型取得了很好的預(yù)測性能.

1 數(shù)據(jù)與方法

1.1 數(shù) 據(jù)

本文先從UniProt數(shù)據(jù)庫[15]中下載7tmrlist文件，該文件含有3 052個GPCRs的UniProt ID號，我們通過這些ID號得到所有GPCRs的氨基酸序列.然后，使用NCBI的Blastclust程序[16]對GPCRs氨基酸序列進行去冗余處理(相似性小于90%)， 通過保留聚類結(jié)果中每行的第1個蛋白質(zhì)得到非冗余的GPCRs樣本數(shù)據(jù)集.

我們從UniProt-GOA ftp站點[15]下載得到gene_association.goa_ref_uniprot文件，并通過上面非冗余GPCRs數(shù)據(jù)集的UniProt ID號得到其GO術(shù)語(剔除Evidence Code 為IEA的GO術(shù)語).然后，我們從基因本體學(xué)數(shù)據(jù)庫[13]下載go.obo文件，并運行數(shù)據(jù)庫提供的obo2csv.py程序，得到go.obo.F.is_a和go.obo.P.is_a文件.文件分別提供了每個分子功能和生物學(xué)過程GO術(shù)語的父節(jié)點，這樣可以得到GPCRs樣本含有的GO術(shù)語及其所有父節(jié)點GO術(shù)語，并得到GPCRs數(shù)據(jù)集的分子功能和生物學(xué)過程的GO術(shù)語標(biāo)記空間.

本文中，我們刪除樣本個數(shù)特別多(大于1 000)和特別少 (小于4)的GO術(shù)語；對分子功能，得到非冗余GPCRs樣本1 327個，GO術(shù)語206個；對生物學(xué)過程，得到非冗余GPCRs樣本1 331個，GO術(shù)語1 406個.

1.2 結(jié)構(gòu)域特征

將非冗余GPCRs氨基酸序列文件提交到NCBI的Batch CD-Search服務(wù)器[17]，得到GPCRs蛋白質(zhì)的保守結(jié)構(gòu)域.對于每個結(jié)構(gòu)域，本文考慮以下7個特征信息：三聯(lián)氨基酸信息、氨基酸關(guān)聯(lián)信息、進化信息、二級結(jié)構(gòu)關(guān)聯(lián)信息、信號肽信息、無序殘基信息和物化屬性.

1) 三聯(lián)氨基酸信息.本文中，依據(jù)氨基酸的偶極矩和側(cè)鏈體積，把20種氨基酸分為6類:①A類.有Ala,Gly,Val.②B類.有Ile,Leu,Phe,Pro.③C類.有Tyr,Met,Thr,Ser,Cys.④D類.有His,Asn,Gln,Tpr.⑤E類.有Arg,Lys.⑥F類.有Asp,Glu[18-19].對于每個結(jié)構(gòu)域，根據(jù)其氨基酸序列計算三聯(lián)體出現(xiàn)頻率(Triad)[19]，

(1)

其中，a,b,c∈{A,B,C,D,F};Na bc表示三聯(lián)氨基酸的個數(shù);n為氨基酸序列長度.對每個結(jié)構(gòu)域，Triad特征為216維.

2) 氨基酸關(guān)聯(lián)信息.依據(jù)上面的6類氨基酸A,B,C,D,E,F，氨基酸關(guān)聯(lián)信息(amino acid correla-tion, AAC)[20]為

(2)

其中，i,j∈{A,B,C,D,E,F}；Pi和Pj分別表示第i和j類氨基酸在結(jié)構(gòu)域中出現(xiàn)的頻率；Pi j(k)表示第i和j類氨基酸在結(jié)構(gòu)域中間隔k個殘基的聯(lián)合出現(xiàn)頻率，k∈{2,4,8,16}.對每個結(jié)構(gòu)域，其AAC特征的維數(shù)為144.

3) 進化信息.本文中，位置特異性得分矩陣(position-specific scoring matrics, PSSMs)被用來描述結(jié)構(gòu)域的進化信息.我們通過psiblast程序[16]比對NCBI非冗余蛋白質(zhì)數(shù)據(jù)集nr得到PSSMs矩陣(3輪迭代和期望值設(shè)為0.001).對每個結(jié)構(gòu)域，其產(chǎn)生的PSSMs矩陣包含42n個元素，其中n為氨基酸序列長度. 然后，我們通過標(biāo)準(zhǔn)的logistic函數(shù)將PSSMs中的元素a歸一化到0～1之間，表達式為

(3)

因為結(jié)構(gòu)域的氨基酸序列長度往往不一樣，其得到的PSSMs元素個數(shù)將不同.本文中，我們把每個氨基酸殘基對應(yīng)的42個PSSMs元素當(dāng)著一個示例，每個結(jié)構(gòu)域可表示為由n個示例組成的包，其中n為結(jié)構(gòu)域的氨基酸序列長度.我們通過多示例算法miFV[21]將示例包轉(zhuǎn)化為單個向量.對每個結(jié)構(gòu)域，最終得到的PSSMs特征維數(shù)為84.

4) 二級結(jié)構(gòu)關(guān)聯(lián)信息.本文中，通 過PSIPRED在線分析工具[22]完成結(jié)構(gòu)域的二級結(jié)構(gòu)預(yù)測.這里，二級結(jié)構(gòu)包括3類:螺旋(helix)、折疊(sheet)和轉(zhuǎn)角(coil).然后，我們計算了結(jié)構(gòu)域中二級結(jié)構(gòu)關(guān)聯(lián)信息(secondary structure element correlation, SSC)[20]，

(4)

其中，i,j∈{helix,sheet,coil}；Pi和Pj分別表示第i和j種二級結(jié)構(gòu)元素在結(jié)構(gòu)域中出現(xiàn)的頻率；Pi j(k)表示第i和j種二級結(jié)構(gòu)元素間隔k個殘基的聯(lián)合出現(xiàn)頻率,k∈{2,4,8,16}.對每個結(jié)構(gòu)域，其SSC特征的維數(shù)為72.

5) 信號肽信息.我們通過SignalP[23]程序從氨基酸序列中預(yù)測信號肽信息.對每個結(jié)構(gòu)域，其得到的信號肽信息包含5n個元素，其中n為氨基酸序列長度.因為每個結(jié)構(gòu)域的氨基酸序列長度往往不一樣，我們通過多示例算法miFV[21]將其轉(zhuǎn)化為單個向量.對每個結(jié)構(gòu)域，其信號肽信息SignalP特征維數(shù)為84.

7) 物化屬性.我們通過SciDBMaker程序[25]從氨基酸序列中預(yù)測物化屬性信息，并通過標(biāo)準(zhǔn)的logistic函數(shù)將每個元素歸一化到0～1之間.對每個結(jié)構(gòu)域，其物化屬性Pychem特征維數(shù)為59.

綜上所述，對每個結(jié)構(gòu)域，總的特征維數(shù)為743.

1.3 快速多示例多標(biāo)記學(xué)習(xí)

1.3.1 算法框架

本文使用的快速多示例多標(biāo)記學(xué)習(xí)MIMLfast[14]為一個2層的分類模型.在第1層中算法MIMLfast從GPCRs的原始特征空間學(xué)習(xí)出一個低維的子空間，該子空間是被所有GO標(biāo)記共享；在第2層中，算法MIMLfast在基于GO標(biāo)記共享的子空間為每個GO標(biāo)記學(xué)習(xí)一個分類模型.這2層模型通過交替優(yōu)化來擬合訓(xùn)練數(shù)據(jù).

首先，對于GPCRs樣本中的單一示例x，其在第l個GO標(biāo)記上的分類模型為

(5)

其中，W0是一個b×d的共享子空間矩陣，目的是將原始GPCRs樣本的特征維度由d維降到b維；而wl是對應(yīng)于第l個GO標(biāo)記的模型權(quán)重向量.這里W0對應(yīng)為模型的第1層，而wl對應(yīng)為模型的第2層.

根據(jù)國務(wù)院辦公廳《關(guān)于進一步改革完善藥品生產(chǎn)流通使用政策的若干意見》（2017［13］號）的要求，促進綜合醫(yī)改試點的省、區(qū)、市和公立醫(yī)院改革，試點城市要率先推行“兩票制”，鼓勵和帶動其他地區(qū)實行“兩票制”，爭取到2018年在全國推開。

其次，如果GPCRs樣本含有某個GO標(biāo)記，是由于其中的某個示例含有該GO標(biāo)記，即一個GPCRs的GO標(biāo)記實際上由它的示例中最可能是正例的那個決定.因此，我們定義GPCRs樣本在第l個GO標(biāo)記上的預(yù)測值如下：

(6)

其中，示例x稱為GPCRs樣本X在第l個GO標(biāo)記上的關(guān)鍵示例.對于GPCRs示例包X以及它的一個相關(guān)標(biāo)記yl，R(X,l)的定義如下：

(7)

定義GPCRs樣本X在GO標(biāo)記yl上的排序錯誤如下：

(8)

由此可知，當(dāng)GO標(biāo)記yl被排的越靠后，排序錯誤ε越大.那么，整個GPCRs樣本集上的排序錯誤為

(9)

(10)

其中，I[·]是指示函數(shù)，當(dāng)括號內(nèi)為真時返回1，否則返回0.由于排序錯誤ε(X,l)是非凸非連續(xù)的，其優(yōu)化問題是NP難問題，因此難于對式(10)進行直接優(yōu)化.為解決此問題，我們采用hinge loss來代替：

(11)

其中，當(dāng)q≥0時，|q|+=q；否則，|q|+=0.

1.3.2 算法求解

(12)

(13)

其中γt為SGD算法更新的步長.

1.4 評價指標(biāo)

本文采用10倍交叉驗證來評估模型，即把GPCRs樣本隨機分為10等分，將其中9等分用來訓(xùn)練模型，剩下1等分用來測試模型，重復(fù)10次，以保證每個GPCRs樣本都被測試過1次.本文采用多示例多標(biāo)記學(xué)習(xí)的評價指標(biāo)，即Hamming loss(HL),Ranking loss(RL),OneError(OE)和Average precision(AP)[5].HL用于考察多標(biāo)記分類器分類錯誤的程度，評估樣本在單個標(biāo)記上的真實標(biāo)記和預(yù)測標(biāo)記之間的誤差率；RL用于考察在樣本的預(yù)測標(biāo)記排序序列中存在的排序錯誤，即樣本的無關(guān)標(biāo)記位于相關(guān)標(biāo)記之前的次數(shù)；OE用于考察在樣本的真實標(biāo)記的排序序列中，序列最前端即排在第1的標(biāo)記不在樣本的預(yù)測標(biāo)記集合中的次數(shù)；AP用來考察在樣本的預(yù)測標(biāo)記排序序列中，位于樣本的某個相關(guān)標(biāo)記之前的標(biāo)記仍為相關(guān)標(biāo)記的比例，反映了預(yù)測標(biāo)記的平均精確度,AP越大，表示模型性能越好,其他評價指標(biāo)相反.

2 結(jié)果與討論

2.1 GPCR結(jié)構(gòu)域各種特征的貢獻

表1顯示了GPCRs結(jié)構(gòu)域各種特征對模型預(yù)測性能的貢獻.由表1可以看出，對GO分子功能， GPCRs結(jié)構(gòu)域特征對模型AP的貢獻由大到小依次為Conjoint Triad>PSSMs>AAC>DISOPRED>SSC>Pychem>SignalP；而對GO生物學(xué)過程，其特征的貢獻由大到小依次為Conjoint triad>AAC>PSSMs>SSC>DISOPRED>Pychem>SignalP.當(dāng)使用所有的特征時(對應(yīng)于ALL列)，我們的模型得到了最好的性能.例如，對GO分子功能，模型的AP=0.886 8；對GO生物學(xué)過程，模型的AP=0.881 8.

Table 1 Contributions of Domain Features on Predictive Performance表1 結(jié)構(gòu)域特征對模型預(yù)測性能的貢獻

Notes: “↑” indicates the larger the value, the better the performance; “↓” indicates the smaller the value, the better the performance; the best results on each evaluation criterion are highlighted in boldface.

2.2 與多示例多標(biāo)記方法的比較

我們與4種多示例多標(biāo)記學(xué)習(xí)方法進行了比較，分別是MIMLRBF[27]，ENMIMLNNBP[3]，MIMLKNN[5]，MIMLSVM[4].MIMLfast和4種對比方法均采用參考文獻中的默認(rèn)參數(shù).對MIMLfast方法，共享空間維度設(shè)為100；對MIMLRBF方法，縮放因子設(shè)為0.08；對MIMLKNN方法，簇的數(shù)量設(shè)為樣本數(shù)量的40%；對MIMLSVM方法，高斯核半徑r=0.2； 對ENMIMLNNBP方法，學(xué)習(xí)率設(shè)為0.4.表2顯示了我們的模型與多示例多標(biāo)記學(xué)習(xí)方法的比較結(jié)果(基于所有的特征).其中↑表示評價指標(biāo)值越大，性能越好，↓則相反；最優(yōu)的結(jié)果用粗體標(biāo)注.從表2可知，在GPCRs的GO分子功能和生物學(xué)過程預(yù)測上，我們的方法均優(yōu)于4種基于多示例多標(biāo)記學(xué)習(xí)的預(yù)測方法.例如，對GO分子功能，我們模型比次優(yōu)的ENMIMLNNBP方法高0.071 9(AP值)；對GO生物學(xué)過程，我們模型比次優(yōu)的MIMLRBF方法高0.136 7(AP值).

效率是多示例多標(biāo)記學(xué)習(xí)的一個瓶頸，因此對方法的效率進行檢驗非常關(guān)鍵.本文的實驗全部在一臺具有4×2.60 GHz CPU 和16 GB內(nèi)存的PC機上完成.表3給出了各個多示例多標(biāo)記學(xué)習(xí)算法在分子功能和生物學(xué)過程數(shù)據(jù)集上的時間開銷比較.從表3可以看出，我們的方法在2個數(shù)據(jù)集上都是效率最高的，遠(yuǎn)遠(yuǎn)少于比較方法的時間開銷.Huang等人(2014)論文中的多個結(jié)果顯示，本文基于的MIMLfast方法的時間復(fù)雜性與樣本數(shù)量呈現(xiàn)一個線性相關(guān)，在處理這么大規(guī)模數(shù)據(jù)具有很大的優(yōu)勢[14].

Table 2 Comparison with Multi-Instance Multi-label Learning Based Methods表2 與多示例多標(biāo)記學(xué)習(xí)方法的比較

Notes: “↑” indicates the larger the value, the better the performance; “↓” indicates the smaller the value, the better the performance; the best results on each evaluation criterion are highlighted in boldface.

Table 3 Runtime Comparison with Multi-Instance Multi-label Learning Based Methods表3 與多示例多標(biāo)記學(xué)習(xí)方法的時間效率比較 s

2.3 與多標(biāo)記學(xué)習(xí)方法的比較

本文與2種多標(biāo)記學(xué)習(xí)方法進行比較，分別是BPMLL[28]和ML-KNN[29]，它們均采用參考文獻中的默認(rèn)參數(shù).即對BPMLL方法，學(xué)習(xí)率設(shè)為0.05；對ML-KNN方法，近鄰樣本個數(shù)設(shè)為10.我們將多示例轉(zhuǎn)化為單示例采用了3種方法：miFV[21]，miVLAD[21]和平均(Means).平均的方法就是將每個樣本包的所有示例求平均，得到單示例向量.表4顯示了我們的模型與多標(biāo)記學(xué)習(xí)算法的比較結(jié)果.其中↑代表評價指標(biāo)越大，性能越好，↓則相反；最優(yōu)的結(jié)果用粗體進行標(biāo)注.表4顯示我們的算法優(yōu)于所有的多標(biāo)記學(xué)習(xí)方法，例如，我們模型比最優(yōu)的多標(biāo)記學(xué)習(xí)BPMLL(基于miFV單示例化方法)在GO分子功能上高0.288 1(AP值)和在GO生物學(xué)過程上高0.310 9(AP值).這或許表明將多示例轉(zhuǎn)換為單示例過程中會損失較多信息.

Table 4 Comparison with Multi-Label Learning Based Methods表4 與多標(biāo)記學(xué)習(xí)方法的比較

Notes: “↑” indicates the larger the value, the better the performance; “↓” indicates the smaller the value, the better the performance; the best results on each evaluation criterion are highlighted in boldface.

2.4 與CAFA蛋白質(zhì)功能預(yù)測方法的比較

Critical Assessment of Functional Annotation(CAFA)是國際上最權(quán)威針對蛋白質(zhì)功能進行注釋的比賽，到現(xiàn)在已經(jīng)舉行了3屆.本文中，我們與多個CAFA比賽中蛋白質(zhì)功能預(yù)測方法進行比較，包括FFPred[30],Argot[31],INGA[32]，我們通過這3個方法的在線預(yù)測平臺得到GPCRs的預(yù)測結(jié)果.為了對比的公平性，我們把平臺得到的預(yù)測結(jié)果，根據(jù)2.1節(jié)中的go.obo文件補充其父節(jié)點GO術(shù)語，補充之后的GO術(shù)語與訓(xùn)練樣本中的GO術(shù)語標(biāo)記空間一樣.

表5顯示了與CAFA蛋白質(zhì)功能預(yù)測方法的比較結(jié)果.其中↑代表評價指標(biāo)越大，性能越好，↓則相反；最優(yōu)的結(jié)果用粗體進行標(biāo)注.我們采用了不同的閾值(cutoff)來區(qū)分正負(fù)樣本，即預(yù)測值大于等于cutoff，為正樣本，否則為負(fù)樣本.實驗結(jié)果表明，在GPCRs的GO分子功能和生物學(xué)過程預(yù)測上，我們的方法均優(yōu)于這3種CAFA蛋白質(zhì)功能預(yù)測方法,如表5所示.例如，我們模型比最優(yōu)的CAFA蛋白質(zhì)功能預(yù)測方法Argot在GO分子功能上高0.265 3(AP值)和在GO生物學(xué)過程上高0.291 3(AP值).

Table 5 Comparison with CAFA Protein Function Prediction Methods表5 與CAFA蛋白質(zhì)功能預(yù)測方法的比較

Notes: “↑” indicates the larger the value, the better the performance; “↓” indicates the smaller the value, the better the performance; the best results on each evaluation criterion are highlighted in boldface.

3 總 結(jié)

本文提出了一種基于快速多示例多標(biāo)記學(xué)習(xí)MIMLfast的G蛋白偶聯(lián)受體的生物學(xué)功能預(yù)測方法，該方法采用了一種新的混合特征，它考慮了GPCRs結(jié)構(gòu)域的三聯(lián)氨基酸、氨基酸關(guān)聯(lián)、進化、二級結(jié)構(gòu)關(guān)聯(lián)、信號肽、無序殘基等多種信息.方法包括2個步驟：1)基于G蛋白偶聯(lián)受體的特征空間，學(xué)習(xí)出一個低維的子空間，該子空間為所有GO標(biāo)記所共享；2)為每個GO標(biāo)記根據(jù)子空間中相關(guān)標(biāo)記信息，學(xué)習(xí)出一個分類模型.實驗結(jié)果證明，我們的模型獲得了很好的性能，優(yōu)于目前最優(yōu)的多示例多標(biāo)記學(xué)習(xí)、多標(biāo)記學(xué)習(xí)的預(yù)測方法和CAFA蛋白質(zhì)功能預(yù)測方法.

WuJiansheng, born in 1979. PhD, associate professor. His main research interests include machine learning and bioinformatics.

FengQiaoyu, born in 1991. Master. Her main research interest is machine learning.

YuanJingzhou, born in 1997. Undergraduate in biomedical engineering. His main research interest is biomedical hardware.

HuHaifeng, born in 1973. PhD, associate professor. His main research interests include large-scale similarity search, wireless sensor networks, wireless networking and distributed systems.

ZhouJiate, born in 1997. Undergraduate in biomedical engineering. His main research interests include biomedical imaging and image processing.

GaoHao, born in 1997. Undergraduate in biomedical engineering. His main research interests include biomedical hardware and software developing.