胡鶴還， 孟 軍， 趙思遠(yuǎn)， 紀(jì)騰其

(大連理工大學(xué) 計(jì)算機(jī)科學(xué)與技術(shù)學(xué)院 遼寧 大連 116023)

0 引言

長非編碼RNA(iong non-coding RNA, lncRNA)是長度大于200個(gè)核苷酸(nt)的非編碼RNA，在植物生長、發(fā)育等進(jìn)程中發(fā)揮作用[1]。從是否編碼蛋白質(zhì)的角度，lncRNA起初被認(rèn)為不具備編碼能力而歸結(jié)為非編碼RNA(non-coding RNA, ncRNA)。然而，最近研究表明，部分lncRNA中含有不超過300 nt的短開放閱讀框(short open reading frames, sORFs)，具備編碼小肽的能力。這些lncRNA能夠在細(xì)胞質(zhì)中與核糖體結(jié)合。核糖體在lncRNA上不斷移動(dòng)，對(duì)長度不超過300 nt的sORFs進(jìn)行翻譯，從而形成一類長度小于100個(gè)氨基酸(amino acid, aa)的小肽。這類由lncRNA的sORFs編碼小肽(sORFs-encoded small peptides, SEPs)在植物生命活動(dòng)中發(fā)揮了調(diào)節(jié)作用[2]。sORFs以及SEPs的發(fā)現(xiàn)，使得mRNA和ncRNA的界限變得模糊(部分ncRNA同樣具有編碼能力)，同時(shí)有助于提升人們對(duì)基因組學(xué)的整體認(rèn)知。

目前識(shí)別SEPs的方法主要分為生物實(shí)驗(yàn)方法和計(jì)算方法兩類。生物實(shí)驗(yàn)方法一方面造價(jià)高、耗時(shí)長，另一方面不適用于大規(guī)模的鑒定。計(jì)算方法大多基于人類和動(dòng)物數(shù)據(jù)訓(xùn)練的機(jī)器學(xué)習(xí)模型，考慮到動(dòng)植物ncRNA之間由于生成過程中聚合酶的不同而存在差異[3]，因此，植物SEPs與動(dòng)物SEPs之間可能同樣存在一定的差異，傳統(tǒng)機(jī)器學(xué)習(xí)模型涉及過多的人工干預(yù)，未充分挖掘SEPs存在的深層特征。所以，采用深度學(xué)習(xí)方法挖掘植物SEPs，已成為發(fā)展趨勢(shì)。

本文采用生物信息學(xué)軟件sORF finder和ORF finder挖掘植物lncRNA中的sORFs，使用基因組學(xué)中生物序列的連續(xù)編碼方式，對(duì)sORFs序列進(jìn)行編碼后作為模型輸入，提出了一種結(jié)合多尺度卷積膠囊網(wǎng)絡(luò)(capsule network, CapsNet)的深度學(xué)習(xí)模型。該模型兼顧了卷積層充分提取初級(jí)局部特征與CapsNet提取高級(jí)特征并自動(dòng)進(jìn)行特征聚類的特性，使兩者達(dá)到互補(bǔ)，從而更好地實(shí)現(xiàn)對(duì)lncRNA中sORFs的分類預(yù)測。通過與單一、簡單融合的深度學(xué)習(xí)模型的比較，以及對(duì)多個(gè)物種數(shù)據(jù)集的測試，結(jié)果表明本文提出的模型具有良好的分類效果和泛化能力。

本文的主要貢獻(xiàn)如下：

1) 采用生物信息學(xué)軟件挖掘植物lncRNA中的sORFs，將植物lncRNA編碼小肽的預(yù)測問題轉(zhuǎn)化為判斷sORF是否具有編碼潛力的問題。

2) 基于邏輯推理的思想，將生物信息學(xué)軟件挖掘出的sORFs進(jìn)一步篩選，提升數(shù)據(jù)的可信度。

3) 提出了一種結(jié)合多尺度卷積膠囊網(wǎng)絡(luò)的深度學(xué)習(xí)模型來實(shí)現(xiàn)植物lncRNA編碼小肽的預(yù)測。

1 相關(guān)工作

目前，對(duì)SEPs的研究大多采用生物實(shí)驗(yàn)和計(jì)算方法。生物實(shí)驗(yàn)方法主要有核糖體檢測、質(zhì)譜分析、肽組學(xué)分析等[4]，研究人員采用上述生物學(xué)實(shí)驗(yàn)方法，鑒定出了一系列SEPs。隨著研究的深入，已有研究者采用上述方法對(duì)SEPs進(jìn)行了大規(guī)模分析。Fesenko等[5]針對(duì)小立碗蘚以及其他10種植物(擬南芥、玉米和卷柏等)，通過質(zhì)譜分析技術(shù)得到具有高可信度的SEPs。然而，生物實(shí)驗(yàn)耗費(fèi)大量人力物力，不適合進(jìn)行大規(guī)模的SEPs鑒定。

與生物實(shí)驗(yàn)方法相比，計(jì)算方法節(jié)省了大量的時(shí)間和成本。識(shí)別sORFs和SEPs的生物信息學(xué)工具主要使用機(jī)器學(xué)習(xí)方法。例如，Hanada等[6]基于編碼序列(CDS)與非編碼序列(NCDS)之間的六聚體組成偏差，采用perl語言編寫了適用于擬南芥等11種生物的工具包sORF finder。Zhu等[7]開發(fā)了一種專門用于鑒定SEPs的基于機(jī)器學(xué)習(xí)方法的工具M(jìn)iPepid，該工具通過提取人類氨基酸序列對(duì)應(yīng)RNA序列，以及人類非編碼RNA中sORFs序列的4-mer特征，并使用邏輯回歸模型提出了工具M(jìn)iPEPid。Tong等[8]分別提取了人類和動(dòng)物的小編碼RNA和小非編碼RNA的序列與理化特征，并依托支持向量機(jī)(SVM)構(gòu)建了工具CPPred。Zhang等[9]使用了與CPPred相同的數(shù)據(jù)集，提出了一種基于CNN的RNA編碼潛力預(yù)測工具DeepCPP，表明深度學(xué)習(xí)方法挖掘SEPs已逐步進(jìn)入公眾視野。然而，已有的生物信息學(xué)工具大多挖掘人類和動(dòng)物SEPs，由于動(dòng)植物SEPs之間可能存在一定的差異，因此目前的工具不能挖掘具有高可信度的植物SEPs。

2017年Sabour等[10]在NIPS會(huì)議中提出了膠囊網(wǎng)絡(luò)(CapsNet)結(jié)構(gòu)。CapsNet采用了一種新的“向量進(jìn)，向量出”的傳遞方案，使用向量代替標(biāo)量表示對(duì)應(yīng)特征，從而使得輸出向量的每個(gè)維度代表“特征的特征”。CapsNet結(jié)構(gòu)有兩個(gè)重要的創(chuàng)新：Squash壓縮激活函數(shù)和動(dòng)態(tài)路由。Squash壓縮激活函數(shù)將向量的模長壓縮到0～1之間，從而使得向量模長能夠代表對(duì)應(yīng)特征的概率；動(dòng)態(tài)路由通過強(qiáng)化相似特征、弱化離群特征，從而做到特征聚類。

為了增加特征多樣性，減少關(guān)鍵信息的損失，本文結(jié)合CNN與CapsNet的優(yōu)勢(shì)，提出了一種結(jié)合多尺度卷積膠囊網(wǎng)絡(luò)的深度學(xué)習(xí)模型，用三種不同尺度的卷積核分別提取特征，并將提取結(jié)果分別轉(zhuǎn)化成三個(gè)不同尺度的膠囊矩陣。多尺度卷積膠囊網(wǎng)絡(luò)不僅能夠更加充分地提取特征，還能夠考慮sORFs序列中數(shù)據(jù)之間的相關(guān)性，從而更好地實(shí)現(xiàn)lncRNA編碼小肽的預(yù)測，并為相應(yīng)的生物學(xué)實(shí)驗(yàn)打下了基礎(chǔ)。

2 數(shù)據(jù)預(yù)處理

2.1 數(shù)據(jù)集構(gòu)建

從GreeNC數(shù)據(jù)庫(http:∥greenc.sequentiabiotech.com/wiki/Main_Page)下載擬南芥(Arabidopsisthaliana)、大豆(Glycinemax)以及苔蘚(Physcomitrellapatens)的lncRNA數(shù)據(jù)。分別使用生物信息學(xué)軟件sORF finder(http:∥hanadb01.bio.kyutech.ac.jp/sORFfinder/)、ORF finder[11](https:∥www.ncbi.nlm.nih.gov/orffinder/)獲取lncRNA的sORFs。ORF finder是一個(gè)圖形的序列分析工具，能夠分析并找到序列中的sORFs。sORF finder基于編碼序列(CDS)與非編碼序列(NCDS)之間的六聚體組成偏差，進(jìn)而識(shí)別序列中具有編碼能力的sORFs。取兩種工具識(shí)別結(jié)果的交集與差集，通過CD-HIT[12]本地化工具，以0.8作為閾值，進(jìn)行去冗余處理，得到候選sORFs序列，其中兩種工具結(jié)果的交集作為正集，差集作為負(fù)集。獲取候選sORFs的流程如圖1所示。

圖1 候選sORFs的獲取流程

2.2 序列編碼

對(duì)于候選的sORFs序列，采用p-nts編碼方式[13]，將每p個(gè)連續(xù)的堿基作為一個(gè)子序列，子序列間不重疊。由于連續(xù)三個(gè)相鄰的堿基構(gòu)成一個(gè)密碼子，所以本文的p取3，即對(duì)正、負(fù)樣本中所有序列進(jìn)行分詞處理后，統(tǒng)計(jì)得出一個(gè)大小為43=64的生物單詞表。按照單詞在生物序列中出現(xiàn)的頻率，從大到小進(jìn)行編碼。由于sORFs的長度不超過300，所以可將sORFs序列嵌入到一個(gè)100維向量中。例如，當(dāng)輸入序列S=(GAGGCCGTT……ACTCTATGT)時(shí)，根據(jù)上述編碼方式，每三個(gè)連續(xù)堿基視為一個(gè)單詞，再按詞頻大小進(jìn)行編碼，即可將S編碼視為一個(gè)固定長度的向量SC=(8,55,11,…,37,58,14)。向量SC為模型的最終輸入格式。

3 特征分析與邏輯推理

3.1 特征分析

為了探究CDS與候選正集，NCDS中sORFs與候選負(fù)集在序列組成和理化特性方面存在的差異性，針對(duì)sORFs序列和氨基酸序列分別采用不同的特征編碼方式進(jìn)行分析。對(duì)于sORFs序列，提取其k-mer特征[14]；對(duì)于氨基酸序列[15]，分別提取其188D、雙氨基酸組成(di-peptide composition,DPC)特征。然后基于奇異值分解(singular value decomposition, SVD)、主成分分析(principal component analysis, PCA)、t分布隨機(jī)相鄰嵌入(t-distributed stochastic neighbor embedding, t-SNE)、核主成分分析(kernel principal component analysis, KPCA)四種特征降維方法，將特征降為二維。上述3類特征編碼方式具體如下。

1)k-mer特征描述序列組成信息，每個(gè)k-mer為序列中相鄰的連續(xù)k個(gè)核苷酸。本文設(shè)定k為4，共256維，其中每一維度代表對(duì)應(yīng)k-mer的頻率fk-mer：

(1)

2) 188D特征綜合考慮了氨基酸組成(amino acid composition,AAC)特征，氨基酸類別組成過濾和分布特征(composition transition distribution,CTD)，共188維。其中：前20維描述20種氨基酸的組成頻率fAAC；后168維描述8大類氨基酸的理化特性fCTD，

(2)

3) DPC特征描述了氨基酸序列中雙核苷酸的組成信息，共400維，其中每一維度代表對(duì)應(yīng)雙核苷酸的頻率fDPC，

(3)

3.2 邏輯推理

本文采用的邏輯推理框架如圖2所示。

圖2 邏輯推理結(jié)構(gòu)

反繹學(xué)習(xí)[16]由機(jī)器學(xué)習(xí)和邏輯推理兩部分構(gòu)成。給定一組數(shù)據(jù)，機(jī)器學(xué)習(xí)部分中的初始分類器給出偽標(biāo)簽，導(dǎo)致偽事實(shí)，然后基于最小化與知識(shí)庫的不一致性，通過邏輯推理修正偽標(biāo)簽，替換原分類器，重復(fù)迭代直到偽標(biāo)簽與知識(shí)庫完全一致為止。本文將候選數(shù)據(jù)集的正負(fù)標(biāo)簽視為偽標(biāo)簽，將CDS、NCDS中sORFs的特征視為知識(shí)庫，對(duì)候選正負(fù)集進(jìn)行篩選。具體做法如下。

1) 將候選正集與CDS的特征、候選負(fù)集與NCDS中sORFs的特征進(jìn)行比對(duì)。

2) 獲取CDS、NCDS中sORFs的特征的橫、縱坐標(biāo)范圍，然后將其視為知識(shí)庫，篩選掉不在CDS特征坐標(biāo)范圍內(nèi)的候選正集、不在NCDS特征坐標(biāo)范圍內(nèi)的候選負(fù)集。

4 多尺度卷積膠囊網(wǎng)絡(luò)模型

4.1 嵌入階段

嵌入層的作用是將輸入序列映射成卷積層易于處理的矩陣向量的形式，方便卷積層充分提取特征。主要工作是將輸入序列的每個(gè)數(shù)字映射成一個(gè)1×n維的向量，這樣輸入序列被映射成m×n維的矩陣形式，其中：m代表序列長度；n代表嵌入維度。本文中序列長度m為100，嵌入維度n為64，即嵌入階段每條序列可映射為100×64的矩陣，作為多尺度卷積膠囊階段的輸入。嵌入層使用Keras庫的Embedding()方法，參數(shù)output_dim設(shè)置為64，input_length設(shè)置為100。

4.2 多尺度卷積膠囊階段

卷積神經(jīng)網(wǎng)絡(luò)(convolutional neural networks, CNN)是深度學(xué)習(xí)中一種被廣泛應(yīng)用的網(wǎng)絡(luò)[17]，主要由卷積層、池化層、全連接層構(gòu)成。卷積層采用卷積核沿著矩陣向量以指定滑動(dòng)步長逐個(gè)進(jìn)行卷積計(jì)算，具體計(jì)算公式為

(4)

由于CNN中的池化層會(huì)導(dǎo)致關(guān)鍵信息的丟失，同時(shí)也會(huì)忽略整體與部分之間的關(guān)聯(lián)，因此采用膠囊網(wǎng)絡(luò)代替池化層。研究表明，模型性能與特征多樣性呈正相關(guān)，即特征種類越豐富，模型性能越好。由于固定尺度的卷積核只能提取某一局部特征，難以捕捉到不同局部特征，從而忽略某些潛在信息。為克服這一缺陷，引入多尺度卷積核和多個(gè)膠囊層代替單尺度卷積核和單一膠囊，從而避免特征提取不充分導(dǎo)致有效信息丟失。實(shí)驗(yàn)中采用三個(gè)不同尺度的卷積核對(duì)嵌入層的輸出分別進(jìn)行卷積操作，每個(gè)卷積操作可捕捉序列不同位置的局部特征，然后對(duì)每個(gè)卷積操作的特征映射矩陣，采用膠囊網(wǎng)絡(luò)將其轉(zhuǎn)換為三個(gè)膠囊矩陣，最后將膠囊矩陣進(jìn)行串聯(lián),作為多尺度卷積膠囊層的輸出。

本文的卷積核尺寸分別為3×64、6×64、9×64，即卷積層的輸出維度分別為98×64、95×64、92×64；然后將每個(gè)特征映射矩陣分別轉(zhuǎn)化為8、12、16個(gè)膠囊，每個(gè)膠囊的維度為16，即膠囊矩陣的維度分別為8×16、12×16、16×16；最后采用Concatenate()函數(shù)將所有膠囊矩陣串聯(lián)，形成108×16的膠囊矩陣作為次級(jí)膠囊層的輸入。

4.3 分類階段

膠囊之間通過動(dòng)態(tài)路由算法進(jìn)行連接，低一級(jí)膠囊轉(zhuǎn)換成更高一級(jí)的膠囊，從而自動(dòng)進(jìn)行特征聚類，更好地表達(dá)高級(jí)特征。其具體實(shí)現(xiàn)如圖3所示。

圖3 動(dòng)態(tài)路由算法流程圖

(5)

(6)

(7)

利用Squash非線性激活函數(shù)處理輸入向量sj，得到高級(jí)膠囊層的輸出向量vj，

(8)

(9)

動(dòng)態(tài)路由機(jī)制是一個(gè)迭代算法，當(dāng)?shù)螖?shù)等于指定值時(shí)迭代終止，即bij停止更新。此時(shí)得到的vj即為高級(jí)膠囊層的最終輸出向量。

本文將多尺度卷積膠囊階段生成的108×16的膠囊矩陣視為108個(gè)16維膠囊，基于動(dòng)態(tài)路由算法將其轉(zhuǎn)換成10個(gè)16維膠囊，即10×16的膠囊矩陣。利用Flatten層處理該膠囊矩陣得到160維向量，再添加參數(shù)為0.4的Dropout層防止過擬合，最后使用參數(shù)為1的Dense層得出一個(gè)具體數(shù)字，并使用sigmoid()函數(shù)將其映射在[0,1]之間，即得出預(yù)測標(biāo)簽。

4.4 模型實(shí)現(xiàn)

輸入序列(Input)經(jīng)過p-nts編碼(Coded)形成100維向量。首先利用嵌入層(Embedding)將輸入向量映射成一個(gè)100×64的矩陣向量，方便進(jìn)行卷積操作；然后通過多尺度卷積膠囊層(Multi-scale Convolution and CapsNet)進(jìn)行操作，輸出108×16的膠囊矩陣，再基于動(dòng)態(tài)路由算法(Dynamic Routing)形成10×16的膠囊矩陣，最后經(jīng)分類得出預(yù)測結(jié)果，模型的整體結(jié)構(gòu)如圖4所示。

圖4 模型整體結(jié)構(gòu)

5 實(shí)驗(yàn)與結(jié)果

5.1 驗(yàn)證方法與評(píng)價(jià)指標(biāo)

采用苔蘚數(shù)據(jù)作為訓(xùn)練集，擬南芥和大豆作為獨(dú)立測試集。通過與現(xiàn)有深度學(xué)習(xí)方法比較，驗(yàn)證提出的模型在sORFs編碼小肽預(yù)測方面的性能與泛化能力。

實(shí)驗(yàn)采用5折交叉驗(yàn)證法來驗(yàn)證模型的性能。使用準(zhǔn)確率(Accuracy，ACC)、精確率(Precision，P)、召回率(Recall，R)以及F1值(F1_score，F(xiàn)1)對(duì)模型性能進(jìn)行評(píng)估。評(píng)價(jià)指標(biāo)的定義如下。

1) 準(zhǔn)確率(ACC)表示被正確預(yù)測的樣本所占比例，

(10)

2) 召回率(R)表示實(shí)際的正樣本在被正確預(yù)測的樣本中所占比例，

(11)

3) 精確率(P)表示被預(yù)測為正的樣本中，實(shí)際的正集樣本所占比，

(12)

4)F1值(F1)為召回率與精確率的加權(quán)平均值，

(13)

其中：TP表示正確分類的正樣本；FN表示被分為負(fù)樣本的正樣本；TN表示正確分類的負(fù)樣本；FP表示被分為正樣本的負(fù)樣本。

5.2 數(shù)據(jù)集

從GreeNC數(shù)據(jù)庫中下載相關(guān)數(shù)據(jù)并進(jìn)行數(shù)據(jù)預(yù)處理，然后通過邏輯推理過程篩選數(shù)據(jù)集。為保證正、負(fù)樣本均衡，從負(fù)集樣本庫中隨機(jī)抽取與正集相同數(shù)目的樣本作為負(fù)集。各物種具體數(shù)據(jù)如表1所示。

表1 各物種數(shù)據(jù)集

5.3 基于不同方法的分類結(jié)果

為驗(yàn)證提出方法的有效性和優(yōu)勢(shì)，實(shí)驗(yàn)采用苔蘚數(shù)據(jù)集，按照p-nts編碼方式對(duì)sORFs序列進(jìn)行編碼，并使用卷積神經(jīng)網(wǎng)絡(luò)(CNN)、長短期記憶網(wǎng)絡(luò)(Bi-LSTM)、膠囊網(wǎng)絡(luò)(CapsNet)、CNN+CapsNet、CNN+Bi-LSTM方法進(jìn)行對(duì)比實(shí)驗(yàn)，5折交叉驗(yàn)證的實(shí)驗(yàn)結(jié)果如表2。表中性能最好的用粗體表示。

表2 基于不同方法的分類結(jié)果

從表2可看出，較之單一深度學(xué)習(xí)方法，提出方法在準(zhǔn)確率、精確率、召回率和F1值4個(gè)指標(biāo)性能上都具有明顯的優(yōu)勢(shì)。在準(zhǔn)確率上比CNN、Bi-LSTM和CapsNet分別提升了5.2%、3.0%、4.2%，說明提出方法在預(yù)測sORFs是否具有編碼能力方面具有良好的分類能力；同時(shí)，與單尺度模型CNN+CapsNet相比，提出模型能夠提取不同局部特征，充分挖掘潛在信息；與CNN+Bi-LSTM模型相比，提出的模型既能提取豐富的特征，又能自動(dòng)進(jìn)行特征聚類，輸出更具有表達(dá)力的特征向量。在準(zhǔn)確率上比CNN+CapsNet、CNN+Bi-LSTM分別提升了2.1%、2.0%。

5.4 基于不同物種的分類結(jié)果

為證明提出方法的泛化能力，選用擬南芥、大豆數(shù)據(jù)集作為獨(dú)立測試集，使用CNN+CapsNet、CNN+Bi-LSTM方法進(jìn)行對(duì)比實(shí)驗(yàn)。獨(dú)立測試的實(shí)驗(yàn)結(jié)果如表3。表中性能最好的用粗體表示。從表3結(jié)果可看出，在擬南芥、大豆兩個(gè)物種的獨(dú)立測試集上，提出方法的性能與CNN+CapsNet和CNN+Bi-LSTM方法相比都有一定的提升。提出方法在預(yù)測擬南芥、大豆的sORFs是否具有編碼能力的性能較好，表明模型具有良好泛化能力。

表3 基于不同物種的分類結(jié)果

6 結(jié)束語

本文提出一種多尺度卷積膠囊網(wǎng)絡(luò)的深度學(xué)習(xí)模型，使用多尺度卷積核捕捉不同局部特征，然后使用多層膠囊網(wǎng)絡(luò)捕捉深層次特征并自動(dòng)進(jìn)行特征聚類，從而實(shí)現(xiàn)lncRNA編碼小肽的分類預(yù)測。實(shí)驗(yàn)結(jié)果表明，提出模型與傳統(tǒng)深度學(xué)習(xí)和單尺度模型對(duì)比，取得了最好的分類效果。此外，模型在多物種獨(dú)立測試集上也取得了良好的分類效果，驗(yàn)證了提出模型具有健壯的泛化能力。未來將對(duì)植物lncRNA編碼小肽開展更加深入細(xì)致的研究，如直接使用模型對(duì)輸入的植物lncRNA進(jìn)行分析，識(shí)別出其中的編碼小肽區(qū)域，進(jìn)而進(jìn)行生物學(xué)驗(yàn)證。