馬 揚(yáng) 劉澤一 梁星星 程光權(quán) 陽方杰 成 清 劉 忠

(國防科技大學(xué)系統(tǒng)工程學(xué)院 長沙 410073)

基因組數(shù)據(jù)是目前最流行的生物信息數(shù)據(jù)類型之一，自人類基因組計(jì)劃開始，生物基因序列數(shù)據(jù)開始呈指數(shù)級增長.如何將有效的生物信息特征從海量高維的生物序列數(shù)據(jù)中提取出來，是目前生物信息領(lǐng)域最重要的問題之一.

近年來，基于基因序列的深度學(xué)習(xí)模型研究激增,主要原因之一是該類模型可以直接從原始數(shù)據(jù)進(jìn)行表示學(xué)習(xí)或特征學(xué)習(xí).深度學(xué)習(xí)可以自動(dòng)提取特征，省去了大量手工選擇特征的工作量[1].同時(shí)，深度學(xué)習(xí)模型也被證明可以有效識別局部基因序列在特定環(huán)境中的作用，并進(jìn)一步可以用來推斷基因的調(diào)控規(guī)則和表達(dá)模式.自編碼器[2]、多層感知機(jī)(multilayer perception, MLP)[3-4]、玻爾茲曼機(jī)(restricted Boltzmann machines, RBMs)[5]、深度信念網(wǎng)絡(luò)(deep belief networks, DBNs)[6]、卷積神經(jīng)網(wǎng)絡(luò)(convolutional neural network, CNN)、循環(huán)神經(jīng)網(wǎng)絡(luò)(recurrent neural network, RNN)[7]、長短記憶網(wǎng)絡(luò)(long short-term memory, LSTM)[8]等深度學(xué)習(xí)方法與基因組研究的結(jié)合促進(jìn)了對基因組學(xué)的深刻理解，并為精密醫(yī)學(xué)、藥物設(shè)計(jì)[9]等多個(gè)領(lǐng)域提供幫助.

Urda等人[10]使用深度學(xué)習(xí)方法在分析RNA序列基因表達(dá)數(shù)據(jù)方面優(yōu)于Lasso(least absolute shrinkage and selection operator).Chen等人[11]提出了一種3層前饋神經(jīng)網(wǎng)絡(luò)對基因表達(dá)進(jìn)行預(yù)測，其性能優(yōu)于線性回歸.Xie等人[12]的研究表明，他們基于MLP和SDAs的深度模型在SNP基因型預(yù)測基因表達(dá)量化方面優(yōu)于Lasso和隨機(jī)森林.曹明宇等人[13]使用了一種新標(biāo)注模式,將藥物實(shí)體及關(guān)系的聯(lián)合抽取轉(zhuǎn)化為端對端的序列標(biāo)注任務(wù).Karlic等人[14-17]在實(shí)驗(yàn)中提出了組蛋白修飾與基因調(diào)控之間的相關(guān)性，并且在之前的一些機(jī)器學(xué)習(xí)工作中已經(jīng)進(jìn)行了研究.Singh等人[18]提出了將MLP疊加在CNN之上的統(tǒng)一判別框架DeepChrome，在預(yù)測基因表達(dá)水平高低的二元分類任務(wù)中，平均AUC(area under curve)為0.8. DeepTarget和deepMirGene分別使用RNN和LSTM模型使用表達(dá)數(shù)據(jù)進(jìn)行mRNA目標(biāo)預(yù)測，所提方法的一個(gè)主要優(yōu)勢是不需要手工制作的特征集[19-20].與此同時(shí)，生成模型也被使用來捕捉基因序列中高階、隱藏的相關(guān)性，Way和Greene在癌癥基因組RNA序列數(shù)據(jù)上訓(xùn)練了一個(gè)變分自編碼模型來捕獲生物學(xué)相關(guān)特征[21-22].

然而，在實(shí)際應(yīng)用中，由于基因組序列之間存在著復(fù)雜的相互依賴和長期的相互作用，直接從基因組序列中學(xué)習(xí)特征非常耗時(shí)[23]，且成功的預(yù)測通常在很大程度上依賴于對生物學(xué)知識的正確利用[24].如何在更普遍的無監(jiān)督環(huán)境下構(gòu)建一個(gè)基因序列表示模型，使得模型學(xué)習(xí)到的基因表示可以更有效地表示基因特征是目前生物信息領(lǐng)域一個(gè)主要的研究問題.且目前基因序列的研究大都集中在基因表達(dá)、調(diào)控等領(lǐng)域上[25],對基因數(shù)據(jù)進(jìn)行低維表示也會(huì)結(jié)合基因交互表達(dá)網(wǎng)絡(luò)的信息[26]，直接對病毒基因序列進(jìn)行表示學(xué)習(xí)的研究并不多.

因此，本文針對病毒傳播網(wǎng)絡(luò)構(gòu)建了一種基于傳播網(wǎng)絡(luò)的基因序列表示模型(gene representation based on transmission network, GRTN)，將病毒的基因序列信息與傳播網(wǎng)絡(luò)的結(jié)構(gòu)信息相結(jié)合，進(jìn)而得到新的基因序列表示.模型首先把輸入序列編碼至特定長度的隱藏層以得到降維嵌入，根據(jù)網(wǎng)絡(luò)傳播的特點(diǎn)使用正向LSTM和反向LSTM分別對節(jié)點(diǎn)的圖上下文節(jié)點(diǎn)信息進(jìn)行聚合，然后使用注意力機(jī)制對其鄰居節(jié)點(diǎn)的序列信息進(jìn)行聚合，進(jìn)而得到目標(biāo)節(jié)點(diǎn)基因序列的新的表示.依據(jù)病毒傳播網(wǎng)絡(luò)中鄰居節(jié)點(diǎn)的基因序列相似性高于非鄰居節(jié)點(diǎn)這一特點(diǎn)，對基因序列表示模型進(jìn)行訓(xùn)練和優(yōu)化.

該模型學(xué)習(xí)到的基因序列表示可以作為基因特征或先驗(yàn)信息用于與基因序列相關(guān)的任務(wù)，為了驗(yàn)證基因序列表示模型的有效性，本文設(shè)計(jì)一套在病毒傳播網(wǎng)絡(luò)中發(fā)現(xiàn)未采樣感染者的任務(wù)流程.并且通過在模擬的病毒傳播網(wǎng)絡(luò)、真實(shí)的澳大利亞新冠病毒傳播網(wǎng)絡(luò)和真實(shí)的艾滋病毒傳播網(wǎng)絡(luò)上進(jìn)行了發(fā)現(xiàn)未采樣感染者的實(shí)驗(yàn)，驗(yàn)證了本文提出的基因序列表示模型的有效性.并通過與各類基因表示方法的對比實(shí)驗(yàn)表明模型的高效性.

1 問題描述

病毒傳播網(wǎng)絡(luò)是一種復(fù)雜網(wǎng)絡(luò)，該網(wǎng)絡(luò)可以表示為G=(V,E,A).其中，V表示網(wǎng)絡(luò)中節(jié)點(diǎn)集合，每個(gè)節(jié)點(diǎn)表示一個(gè)病毒感染者；E表示連邊集合，每條邊(vi,vj)∈E表示節(jié)點(diǎn)vi將疾病傳染給節(jié)點(diǎn)vj，由于疾病傳播是有方向的，故該網(wǎng)絡(luò)為有向圖；A表示節(jié)點(diǎn)屬性集合，在G中，ai∈A表示節(jié)點(diǎn)i的病毒基因序列.

給定一個(gè)網(wǎng)絡(luò)G，基于網(wǎng)絡(luò)的基因序列表示的核心任務(wù)是將網(wǎng)絡(luò)的特征結(jié)合到基因序列中，從而為網(wǎng)絡(luò)中每一個(gè)節(jié)點(diǎn)的基因序列學(xué)習(xí)一個(gè)統(tǒng)一的低維度(β)的表示向量，從而將網(wǎng)絡(luò)中每個(gè)節(jié)點(diǎn)的屬性維度由|A|映射到β(β?|A|)，使其能夠捕捉網(wǎng)絡(luò)中的結(jié)構(gòu)和聚合鄰居節(jié)點(diǎn)的特征[27-28].

2 模型描述

本節(jié)介紹基于傳播網(wǎng)絡(luò)的基因序列表示模型，如圖1所示，該模型包括病毒網(wǎng)絡(luò)中節(jié)點(diǎn)病毒序列信息編碼、基于網(wǎng)絡(luò)傳播的結(jié)構(gòu)提取、基于圖上下文信息的特征表示和基于注意力的節(jié)點(diǎn)信息聚合.

Fig. 1 Gene sequence representation model based on transmission network圖1 基于傳播網(wǎng)絡(luò)的基因序列表示模型

作為病毒傳播網(wǎng)絡(luò)的輸入，節(jié)點(diǎn)的基因序列首先經(jīng)過獨(dú)熱編碼和維度規(guī)約轉(zhuǎn)化為適合模型輸入的數(shù)據(jù)；然后基于網(wǎng)絡(luò)傳播結(jié)構(gòu)對每個(gè)節(jié)點(diǎn)抽取子圖并獲得子圖中的前向節(jié)點(diǎn)信息、后向節(jié)點(diǎn)信息和自身表示信息；基于注意力機(jī)制，模型將前向信息、后向信息和節(jié)點(diǎn)本身的表示信息進(jìn)行聚合，得到節(jié)點(diǎn)的聚合表示；為了能對網(wǎng)絡(luò)節(jié)點(diǎn)的病毒序列信息進(jìn)行更好的表達(dá)，通過設(shè)置基于圖上下文信息的損失函數(shù)，對整個(gè)模型經(jīng)過優(yōu)化與訓(xùn)練，最終得到節(jié)點(diǎn)的表示向量.

2.1 病毒基因序列編碼

不同病毒的基因序列長度迥異，但目前已知的大多數(shù)病毒序列中，存在著大量的無效基因序列，而且在對基因序列進(jìn)行測序時(shí)，往往也會(huì)存在部分序列信息缺失.同時(shí)，病毒在傳播中會(huì)不斷發(fā)生變異，不同網(wǎng)絡(luò)節(jié)點(diǎn)的病毒基因序列都存在一定差異.但同種病毒傳播網(wǎng)絡(luò)中，所有節(jié)點(diǎn)的基因序列大部分是相同的.針對病毒傳播網(wǎng)絡(luò)中基因序列的這種特征，本文采用獨(dú)熱編碼，該編碼對于不同節(jié)點(diǎn)的相同基因片段可以保證其一致性，同時(shí)可以有效地表達(dá)病毒基因序列的不同.

病毒在特定的傳播網(wǎng)絡(luò)中，存在基因序列相同的片段，本文通過對這些基因片段的刪減，對輸入數(shù)據(jù)維度進(jìn)行了維度歸約：若該網(wǎng)絡(luò)中所有節(jié)點(diǎn)的某一位基因值相同，則將該位基因從輸入數(shù)據(jù)中舍棄.通過本節(jié)操作，網(wǎng)絡(luò)中任一節(jié)點(diǎn)v的病毒基因序列g(shù)v在規(guī)約后被編碼為相同長度的0/1編碼序列xv，作為網(wǎng)絡(luò)中節(jié)點(diǎn)的屬性向量.

2.2 網(wǎng)絡(luò)傳播結(jié)構(gòu)抽取

網(wǎng)絡(luò)傳播結(jié)構(gòu)抽取旨在抽取病毒傳播網(wǎng)絡(luò)中與每個(gè)節(jié)點(diǎn)關(guān)聯(lián)度高的鄰居節(jié)點(diǎn)集合，從而在保證計(jì)算效率的同時(shí)，使模型從網(wǎng)絡(luò)中獲取節(jié)點(diǎn)的拓?fù)湫畔?在該抽取的子圖中，目標(biāo)節(jié)點(diǎn)v以外的節(jié)點(diǎn)可以劃分為2類：前向信息節(jié)點(diǎn)集合Γf(v)和后向信息節(jié)點(diǎn)集合Γb(v):

Γf(v)=F(v)∪F(i),i∈F(v),

(1)

Γb(v)=S(v)∪S(i),i∈S(v),

(2)

其中,F(v)表示節(jié)點(diǎn)v的父節(jié)點(diǎn)，S(v)表示節(jié)點(diǎn)v的子節(jié)點(diǎn).前向節(jié)點(diǎn)是位于病毒傳播鏈中目標(biāo)節(jié)點(diǎn)的祖先節(jié)點(diǎn)，即目標(biāo)節(jié)點(diǎn)的父節(jié)點(diǎn)及其祖父節(jié)點(diǎn)；后向節(jié)點(diǎn)是位于病毒傳播鏈中目標(biāo)節(jié)點(diǎn)的子孫節(jié)點(diǎn)，即目標(biāo)節(jié)點(diǎn)的子節(jié)點(diǎn)和孫子節(jié)點(diǎn).由于病毒基因序列傳播變異較快，間隔多代后基因序列之間存在較大差異.此外，研究發(fā)現(xiàn)[29-31]，網(wǎng)絡(luò)中距離超過3的節(jié)點(diǎn)對之間的互相作用非常小，在本文模型中，選取與目標(biāo)節(jié)點(diǎn)距離小于3的節(jié)點(diǎn)，構(gòu)成目標(biāo)節(jié)點(diǎn)的鄰居集合.如圖1所示，節(jié)點(diǎn)子圖抽取得到的網(wǎng)絡(luò)中，中心節(jié)點(diǎn)到所有節(jié)點(diǎn)的距離均小于3，由于每個(gè)節(jié)點(diǎn)只有一個(gè)感染源，所以前向信息節(jié)點(diǎn)集合的規(guī)模較小(0～2個(gè))，后向信息節(jié)點(diǎn)集合的規(guī)模相對更多.通過本節(jié)操作，模型中每個(gè)節(jié)點(diǎn)v的結(jié)構(gòu)信息只和抽取子圖中節(jié)點(diǎn)Γf(v)和Γb(v)相關(guān).

2.3 基于圖上下文的信息表示

基于圖上下文的信息表示旨在將每個(gè)子圖中的信息聚合在一起，子圖中的節(jié)點(diǎn)稱為目標(biāo)節(jié)點(diǎn)的上下文鄰居：節(jié)點(diǎn)的上文鄰居指沿著網(wǎng)絡(luò)中連邊經(jīng)過一定長度隨機(jī)游走能到達(dá)目標(biāo)節(jié)點(diǎn)的節(jié)點(diǎn)集合，節(jié)點(diǎn)的下文鄰居指從目標(biāo)節(jié)點(diǎn)出發(fā)，經(jīng)過指定長度的隨機(jī)游走能到達(dá)的節(jié)點(diǎn)集合.由于節(jié)點(diǎn)之間存在傳播關(guān)系，本文采用LSTM將傳播序列信息聚合在一起，對每個(gè)節(jié)點(diǎn)v，在計(jì)算前向表示向量時(shí)，直接采用LSTM將前向信息聚合在一起，并用LSTM最后的隱藏層輸出為前向表示向量；在計(jì)算后向表示時(shí)，在每個(gè)分支上，采用和前向表示類似的LSTM模型，用LSTM最后的隱藏層作為該分支的表示向量，并用所有分支向量的均值作為后向表示向量.對于自身向量，通過多層全連接層將輸入數(shù)據(jù)轉(zhuǎn)換為需要的維度.通過充分利用子圖中所有節(jié)點(diǎn)的信息，可以獲得前向表示Ef、后向表示Eb和自身表示Em三個(gè)表示向量:

Ef(v)=LSTMf(BN(FCf(xi))),i∈Γf(v),

(3)

(4)

Em(v)=BN(FC(xv)),

(5)

其中，BN表示批標(biāo)準(zhǔn)化操作，F(xiàn)C表示全連接層，3個(gè)表示向量的維度相同，記為do.通過該操作，模型獲得了每個(gè)目標(biāo)節(jié)點(diǎn)在網(wǎng)絡(luò)子圖中抽取的結(jié)構(gòu)信息的表示.

2.4 基于注意力的節(jié)點(diǎn)信息聚合

注意力機(jī)制已經(jīng)成為許多序列學(xué)習(xí)任務(wù)中的一項(xiàng)標(biāo)準(zhǔn)配置[32]，可以為不同節(jié)點(diǎn)分配不同權(quán)重，訓(xùn)練時(shí)依賴于成對的相鄰節(jié)點(diǎn)，而不依賴具體的網(wǎng)絡(luò)結(jié)構(gòu).本文采用注意力機(jī)制學(xué)習(xí)節(jié)點(diǎn)前向、自身和后向向量Ek，k∈{f,m,b}的聚合權(quán)重，進(jìn)而獲得每個(gè)節(jié)點(diǎn)的信息聚合E.

E=AGGREGATE(Ef,Em,Eb),

(6)

其中

AGGREGATE(Ef,Em,Eb)=
αfEf+αmEm+αbEb,

(7)

αk,k∈{f,m,b}為注意力權(quán)重系數(shù)，本文通過一個(gè)共享的注意力向量a∈2d×1來反映向量Ek對節(jié)點(diǎn)自身向量Em的重要性:

emk=aT(Em⊕Ek),

(8)

其中⊕為拼接操作，a為可學(xué)習(xí)的注意力參數(shù).相應(yīng)的注意力系數(shù)αk可以表示為

(9)

其中σ為激活函數(shù)LeakyRelu.

為了減少數(shù)據(jù)的方差，提升表示結(jié)果的穩(wěn)定性，將該模式推廣到多頭注意力形式：注意力機(jī)制被獨(dú)立地重復(fù)nh次(每個(gè)頭的參數(shù)初始化不同)，并將得到的nh個(gè)表示向量拼接作為最終向量:

E=E1⊕E2⊕…⊕Enh,

(10)

2.5 目標(biāo)優(yōu)化與模型訓(xùn)練

為了客觀衡量病毒傳播網(wǎng)絡(luò)的基因序列表示的效果，本文采用圖上下文損失作為損失函數(shù).圖上下文損失可以看作交叉熵的一種變形：對每個(gè)節(jié)點(diǎn)，網(wǎng)絡(luò)中的所有其他節(jié)點(diǎn)，根據(jù)距離該節(jié)點(diǎn)的距離被劃分為正鄰居節(jié)點(diǎn)和負(fù)鄰居節(jié)點(diǎn)，因此其核心是使目標(biāo)節(jié)點(diǎn)的向量表示與采樣正鄰居的相似性更高，并與采樣負(fù)鄰居的差異性更大，基于此的優(yōu)化目標(biāo)被定義為

(11)

其中Cv表示從節(jié)點(diǎn)v開始的隨機(jī)游走能到達(dá)的節(jié)點(diǎn)集合，p(vc|v;Θ)被定義為softmax形式的概率:

(12)

式(12)通過表示向量內(nèi)積大小衡量節(jié)點(diǎn)對的相似性，通過負(fù)采樣[33]技術(shù)，可以將式(12)中softmax函數(shù)進(jìn)行對數(shù)近似為

(13)

其中,S表示負(fù)采樣數(shù)目，P(vc)表示負(fù)采樣的概率分布，在本模型中，S=1從而將目標(biāo)函數(shù)簡約為交叉熵的形式[34]

lnθ(EvEvc)+lnθ(-EvEvs),

(14)

從而將損失函數(shù)轉(zhuǎn)化定義為

(15)

算法1.基因序列表示算法.

輸入：病毒傳播網(wǎng)絡(luò)G=(V,E,A)、節(jié)點(diǎn)基因序列g(shù);

輸出：網(wǎng)絡(luò)中每個(gè)節(jié)點(diǎn)的基因表示E.

① 通過病毒基因序列編碼將g轉(zhuǎn)化為x；

② for eachv∈V：

③ 分別按照式(1)～(2)抽取節(jié)點(diǎn)v的前向信息節(jié)點(diǎn)集合Γf(v)和后向信息節(jié)點(diǎn)集合Γb(v)；

④ 分別按照式(3)～(5)計(jì)算前向表示Ef、后向表示Eb和自身表示Em三個(gè)表示向量；

⑤ 按照式(8)計(jì)算每個(gè)節(jié)點(diǎn)的注意力系數(shù)αk,k∈{f,m,b}，并獲得每個(gè)節(jié)點(diǎn)聚合表示E；

⑦ end for

⑧ 根據(jù)式(14)進(jìn)行模型訓(xùn)練和優(yōu)化；

⑨ 循環(huán)步驟②～⑧直至到達(dá)結(jié)束條件.

3 實(shí)驗(yàn)設(shè)計(jì)

本節(jié)設(shè)計(jì)一種病毒傳播網(wǎng)絡(luò)中發(fā)現(xiàn)未采樣感染者的實(shí)驗(yàn)流程，通過在仿真病毒傳播網(wǎng)絡(luò)、真實(shí)的新冠病毒傳播網(wǎng)絡(luò)和艾滋病毒傳播網(wǎng)絡(luò)上發(fā)現(xiàn)未采樣感染者來驗(yàn)證模型的有效性，并對實(shí)驗(yàn)結(jié)果進(jìn)行分析.

3.1 未采樣感染者發(fā)現(xiàn)流程

實(shí)際情況下，受限于醫(yī)療水平、檢測條件、疾病癥狀的不同，檢測到的病毒傳播網(wǎng)絡(luò)往往會(huì)存在缺失節(jié)點(diǎn)(未采樣感染者)，其網(wǎng)絡(luò)規(guī)模小于實(shí)際的傳播網(wǎng)絡(luò).如圖2所示，虛線節(jié)點(diǎn)為未采樣感染者，由于這些病患節(jié)點(diǎn)在實(shí)際網(wǎng)絡(luò)構(gòu)建時(shí)未被發(fā)現(xiàn)，因此網(wǎng)絡(luò)的實(shí)際連邊關(guān)系存在異常，原本不應(yīng)有傳播關(guān)系的父節(jié)點(diǎn)和子節(jié)點(diǎn)被認(rèn)為存在直接傳播關(guān)系(存在連邊)，如右圖虛線所示.

Fig. 2 The un-sampled infections detection圖2 未采樣感染者示意圖

基于本文基因序列表示模型，我們設(shè)計(jì)2個(gè)未采樣感染者發(fā)現(xiàn)流程：

① 通過基于病毒傳播網(wǎng)絡(luò)的基因表示模型，獲取網(wǎng)絡(luò)中每個(gè)節(jié)點(diǎn)的基因表示E；

② 將發(fā)現(xiàn)傳播網(wǎng)絡(luò)中未采樣節(jié)點(diǎn)任務(wù)轉(zhuǎn)化為連邊得分排序任務(wù)，即通過計(jì)算兩節(jié)點(diǎn)的表示的相似度，評價(jià)這條連邊的穩(wěn)固性.對于節(jié)點(diǎn)對i，j，其連邊相似度為

simij=EiEj,

(16)

在所有存在連邊的排序中，得分最低的k條連邊所對應(yīng)的節(jié)點(diǎn)對被認(rèn)為是存在未采樣感染者的節(jié)點(diǎn)組合，進(jìn)而在實(shí)際的醫(yī)學(xué)調(diào)查診斷相應(yīng)環(huán)節(jié)進(jìn)行排查.

3.2 數(shù)據(jù)描述

本文采用仿真和真實(shí)的病毒傳播網(wǎng)絡(luò)數(shù)據(jù)集進(jìn)行實(shí)驗(yàn).

仿真數(shù)據(jù)集是根據(jù)病毒傳播和變異的規(guī)律，生成傳播網(wǎng)絡(luò)和網(wǎng)絡(luò)中每個(gè)節(jié)點(diǎn)的基因序列，其中，每次仿真基因序列的變異是基于其父節(jié)點(diǎn)的基因序列，該網(wǎng)絡(luò)的生成規(guī)則為：

1) 給定病毒基因序列長度L，每位基因取值數(shù)量為n，假定每位基因的變異是互相獨(dú)立的，病毒基因在每一次傳播時(shí)序列發(fā)生變異的概率為p；

2) 給定病毒傳播能力值R，每一代傳播人數(shù)在以R為中心服從泊松分布(為計(jì)算方便，模擬生成網(wǎng)絡(luò)中傳染人數(shù)在0～6之間對稱分布).假定網(wǎng)絡(luò)中所有節(jié)點(diǎn)起始于“零號病人”，即有一個(gè)同一父節(jié)點(diǎn)，則根據(jù)R值，零號病人生成子節(jié)點(diǎn)，每個(gè)子節(jié)點(diǎn)的基因序列的變異是獨(dú)立的；

3) 重復(fù)步驟2)，直至傳播網(wǎng)絡(luò)節(jié)點(diǎn)數(shù)達(dá)到給定網(wǎng)絡(luò)規(guī)模N.

基于上述規(guī)則，本文生成病毒傳播仿真網(wǎng)絡(luò)，相關(guān)參數(shù)如表1所示:

Table 1 Parameters of Simulating Network表1 仿真網(wǎng)絡(luò)參數(shù)設(shè)置

真實(shí)數(shù)據(jù)集包括Aussie數(shù)據(jù)集和HIV數(shù)據(jù)集.Aussie數(shù)據(jù)集是部分澳大利亞新型冠狀病毒患者地基因序列，病人信息已經(jīng)匿名化處理，該數(shù)據(jù)集基因序列數(shù)據(jù)和標(biāo)準(zhǔn)參考序列比對出現(xiàn)多處缺失，序列質(zhì)量不高.我們使用TransPhylo對該基因序列集進(jìn)行分析，TransPhylo是一個(gè)可以通過基因數(shù)據(jù)構(gòu)建病毒傳播網(wǎng)絡(luò)的工具[35].分析得到病毒傳播網(wǎng)絡(luò)如圖3所示，該網(wǎng)絡(luò)中共有1 032個(gè)節(jié)點(diǎn)，且約有10%的未采樣節(jié)點(diǎn).網(wǎng)絡(luò)的每個(gè)節(jié)點(diǎn)代表一名新冠病毒患者，節(jié)點(diǎn)的有向連邊表示新冠病毒傳播關(guān)系，每個(gè)節(jié)點(diǎn)的屬性為其感染的新冠病毒的基因序列，基因長度為29 915.基因的每個(gè)位置的取值范圍為{A,T,C,G,-}，其中-表示該位基因缺失.如圖3所示，其中節(jié)點(diǎn)顏色越深，表示該節(jié)點(diǎn)距離零號病人越近.

Fig. 3 Aussie: 2019-nCoV transmission network圖3 Aussie新冠病毒傳播網(wǎng)絡(luò)

HIV數(shù)據(jù)集是加拿大阿爾伯塔省南部部分HIV患者的病毒序列信息，并且包含經(jīng)TransPhylo分析且經(jīng)醫(yī)院驗(yàn)證的病毒傳播網(wǎng)絡(luò)信息.該網(wǎng)絡(luò)節(jié)點(diǎn)數(shù)量是140，病毒基因序列長度是1 522，該網(wǎng)絡(luò)中未采樣節(jié)點(diǎn)數(shù)約為41，即該傳播網(wǎng)絡(luò)包含較多未采樣節(jié)點(diǎn)，網(wǎng)絡(luò)如圖4所示.病毒傳播信息采集中，往往無法獲取0號病人的信息，該HIV數(shù)據(jù)集中同樣缺失0號病人基因序列信息.本文參考Lauren等人[36]研究推斷與0號病人傳播關(guān)系最接近的病毒序列，隨機(jī)改變約30%的序列信息后作為0號病人的病毒序列.

Fig. 4 HIV: HIV transmission network圖4 HIV病毒傳播網(wǎng)絡(luò)

3.3 驗(yàn)證集數(shù)據(jù)生成

為了驗(yàn)證基因表示模型的有效性，我們會(huì)在仿真和真實(shí)傳播網(wǎng)絡(luò)中刪除一些節(jié)點(diǎn)，這些刪除的節(jié)點(diǎn)即未采樣感染者，然后生成一個(gè)存在缺失節(jié)點(diǎn)的傳播網(wǎng)絡(luò)，存在缺失節(jié)點(diǎn)的連邊將作為驗(yàn)證集，之后使用模型在該網(wǎng)絡(luò)上進(jìn)行缺失節(jié)點(diǎn)預(yù)測.通過模型預(yù)測的結(jié)果和驗(yàn)證集的比較，進(jìn)而計(jì)算模型在缺失節(jié)點(diǎn)還原結(jié)果的準(zhǔn)確度.其中缺失網(wǎng)絡(luò)的生成規(guī)則為：

1) 生成指定數(shù)目缺失節(jié)點(diǎn)，檢查缺失節(jié)點(diǎn)合理性：若生成節(jié)點(diǎn)處于網(wǎng)絡(luò)兩端(沒有父節(jié)點(diǎn)或沒有子節(jié)點(diǎn))，移除該節(jié)點(diǎn)；

2) 將缺失節(jié)點(diǎn)的父節(jié)點(diǎn)與子節(jié)點(diǎn)連接，生成新的連邊；

3) 重復(fù)步驟1)～2)，直至缺失節(jié)點(diǎn)數(shù)目達(dá)到給定值.

3.4 參數(shù)設(shè)置

對于基因編碼序列屬性，模型形式采用獨(dú)熱編碼；模型訓(xùn)練代數(shù)為1 000，初始學(xué)習(xí)速率lr=0.01，每進(jìn)行100代訓(xùn)練，將學(xué)習(xí)率減半.表示向量維度為do=128，隱藏層維度為dh=1024，L2正則化參數(shù)β=0.001，注意力頭數(shù)nh=4，隨機(jī)數(shù)種子為1，隨機(jī)游走長度Lw=10，隨機(jī)游走重啟概率pr=20%.

3.5 對比算法

為了驗(yàn)證GRGC模型對基因序列表示的高效性，本文采用5種算法作為基準(zhǔn)對比算法，包括不對基因序列降維處理直接用于計(jì)算的Direct方法與可以對高維基因序列進(jìn)行降維及特征提取的PCA和Auto Encoder(AE)方法.此外，由于本文的方法是基于病毒傳播網(wǎng)絡(luò)構(gòu)建的，因此也選用了GCN和SDNE兩種網(wǎng)絡(luò)表示方法作為對比.

1) Direct方法.直接使用基因序列進(jìn)行計(jì)算，每個(gè)節(jié)點(diǎn)的基因序列編碼后，不對其進(jìn)行壓縮和特征提取，直接將其用于發(fā)現(xiàn)缺失節(jié)點(diǎn)任務(wù).

2) AE方法.自編碼器[37]由編碼器和解碼器2部分組成，其中編碼器能將輸入屬性壓縮成潛在的空間表征，解碼器能從潛在空間表征重構(gòu)數(shù)據(jù)的輸入.自編碼器常被用于降維或特征學(xué)習(xí)，本文首先使用AE模型對基因序列進(jìn)行訓(xùn)練，然后使用訓(xùn)練好的AE模型的中間層即降維后的序列表示，作為發(fā)現(xiàn)缺失節(jié)點(diǎn)任務(wù)的輸入.

3) 主成分分析方法(principal component analysis， PCA)[38].是一種使用最廣泛的數(shù)據(jù)降維算法.PCA的主要思想是通過正交變化將可能存在相關(guān)性的N維特征數(shù)據(jù)映射到線性不相關(guān)的k維上，在原有N維特征的基礎(chǔ)上重新構(gòu)造k維特征，進(jìn)而實(shí)現(xiàn)對原始數(shù)據(jù)的高效和精確的表示，轉(zhuǎn)換后的k維特征則被稱為主成分.本文使用基因序列的主成分作為發(fā)現(xiàn)缺失節(jié)點(diǎn)任務(wù)的輸入.

4) 圖卷積神經(jīng)網(wǎng)絡(luò)方法(graph convolutional neural network, GCN)[39].通過低通濾波器產(chǎn)生節(jié)點(diǎn)嵌入，將每個(gè)節(jié)點(diǎn)的信息傳遞到圖中的鄰居節(jié)點(diǎn)上去，通過卷積運(yùn)算使得圖能夠逐層傳遞信息.本文構(gòu)建了一個(gè)2層GCN模型(隱藏層維度為256，輸出層維度為128)，并用上下文損失函數(shù)訓(xùn)練一個(gè)無監(jiān)督GCN模型，通過輸出層的基因表示，計(jì)算節(jié)點(diǎn)對相似度并發(fā)現(xiàn)缺失節(jié)點(diǎn).

5) 結(jié)構(gòu)網(wǎng)絡(luò)嵌入方法(structural deep network embedding， SDNE)[40].使用深層神經(jīng)網(wǎng)絡(luò)對節(jié)點(diǎn)表示間的非線性進(jìn)行建模，整個(gè)模型可以被分為2個(gè)部分：一個(gè)是由Laplace矩陣監(jiān)督的建模第1級相似度的模塊，另一個(gè)是由無監(jiān)督的深層自編碼器對第2級相似度關(guān)系進(jìn)行建模.本文采用3層神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)(2層隱藏層維度分別為500和300，輸出層維度為128).

3.6 評價(jià)指標(biāo)

對所有對比算法和本文方法，本文通過AUC,RS和k-RS三個(gè)指標(biāo)評價(jià)其效果.

鏈路預(yù)測中常用AUC指標(biāo)[41]衡量算法的效果：每次隨機(jī)地從驗(yàn)證集中選取1條存在的邊和1條不存在的邊并比較這2條邊的指標(biāo)值，如果存在的邊的分值大于不存在的邊的分值，則記為1分；如果2個(gè)邊的分?jǐn)?shù)值相等，則記為0.5分；若存在邊的分值小于不存在邊的分值，則記為0分.獨(dú)立重復(fù)地比較N次，如果有n′次驗(yàn)證集中的邊的分?jǐn)?shù)值大于不存在的邊，有n″次驗(yàn)證集中邊分?jǐn)?shù)值等于不存在邊的分?jǐn)?shù)，則AUC指標(biāo)計(jì)算結(jié)果為

(17)

AUC的指標(biāo)值越接近1，說明獲得的表示向量在連邊預(yù)測上的能力越強(qiáng).

RS即排序分(ranking score)，考慮了驗(yàn)證集中的邊在模型預(yù)測結(jié)果中排序的位置.令H=U-ET為未知邊的集合(相當(dāng)于驗(yàn)證集中的邊和不存在的邊的集合)，re表示e∈EP在排序中的排名.那么這條驗(yàn)證邊的排序分為

(18)

該數(shù)據(jù)集上的總的排序分為

(19)

根據(jù)目標(biāo)，本文希望存在未采樣感染者的連邊對的排名更靠后，因此，RS指標(biāo)越大，表明該表示向量在預(yù)測任務(wù)中表現(xiàn)越好.

k-RS統(tǒng)計(jì)RS得分前k條邊中驗(yàn)證集邊所占的比例.如果模型預(yù)測結(jié)果中有m條邊在驗(yàn)證集中，k-RS定義為

(20)

k-RS指標(biāo)側(cè)重于評價(jià)驗(yàn)證集中特征最明顯的那部分邊，在實(shí)際應(yīng)用時(shí)，我們只會(huì)關(guān)心排名特征最明顯的部分(排名最靠前的部分)，從而在準(zhǔn)確率和覆蓋率中進(jìn)行權(quán)衡，k-RS指標(biāo)不再關(guān)注本身特征就不明顯的邊，從而更好地反映和評價(jià)算法在真實(shí)使用時(shí)的效用.

4 實(shí)驗(yàn)分析

本文分別在仿真病毒傳播網(wǎng)絡(luò)和Aussie新冠病毒傳播網(wǎng)絡(luò)、HIV病毒傳播網(wǎng)絡(luò)進(jìn)行實(shí)驗(yàn).為了使實(shí)驗(yàn)更符合真實(shí)場景，即現(xiàn)實(shí)中獲取的病毒傳播網(wǎng)絡(luò)往往存在未采樣感染者，因此我們首先在仿真或者真實(shí)數(shù)據(jù)集的基礎(chǔ)上生成存在一定節(jié)點(diǎn)缺失的網(wǎng)絡(luò)，然后在該存在缺失節(jié)點(diǎn)的網(wǎng)絡(luò)上訓(xùn)練基因表示模型.

在仿真數(shù)據(jù)實(shí)驗(yàn)中，我們首先將模型在生成的病毒傳播網(wǎng)絡(luò)上進(jìn)行訓(xùn)練，然后使用節(jié)點(diǎn)缺失比例為10%和20%的傳播網(wǎng)絡(luò)對模型進(jìn)行測試.節(jié)點(diǎn)缺失網(wǎng)絡(luò)按照3.3節(jié)缺失網(wǎng)絡(luò)生成規(guī)則生成，即選取未采樣節(jié)點(diǎn)并將其在數(shù)據(jù)集中刪除.然后通過分析模型成功找出缺失節(jié)點(diǎn)的正確率來評價(jià)算法的表現(xiàn).在Aussie數(shù)據(jù)上，測試網(wǎng)絡(luò)的缺失節(jié)點(diǎn)比例分別為10%和20%.在HIV傳播網(wǎng)絡(luò)實(shí)驗(yàn)中，由于病毒傳播網(wǎng)絡(luò)規(guī)模較小，10%的節(jié)點(diǎn)缺失總數(shù)是14，較小的驗(yàn)證集會(huì)導(dǎo)致測試結(jié)果波動(dòng)較大，因此在HIV實(shí)驗(yàn)中我們使用全網(wǎng)絡(luò)進(jìn)行模型訓(xùn)練，在30%節(jié)點(diǎn)缺失的網(wǎng)絡(luò)上進(jìn)行測試.

4.1 結(jié)果分析

本節(jié)實(shí)驗(yàn)提到的k-RS指標(biāo)中，百分?jǐn)?shù)值表示排名靠前相應(yīng)比例的連邊數(shù)(25%表示RS得分排名前25%的連邊)，Top 20表示最可能包含缺失節(jié)點(diǎn)的20條連邊.表中實(shí)驗(yàn)結(jié)果均是10次生成缺失節(jié)點(diǎn)實(shí)驗(yàn)的平均值，分值越高代表算法預(yù)測存在缺失節(jié)點(diǎn)的正確率越高.

從表2和表3可以看出，在仿真數(shù)據(jù)中，采用Direct方法的基因表示在進(jìn)行缺失節(jié)點(diǎn)發(fā)現(xiàn)任務(wù)中的表現(xiàn)最好，這是因?yàn)榉抡鏀?shù)據(jù)的基因序列生成及變異是嚴(yán)格按照規(guī)則生成，因此仿真基因序列中無效序列較少，對基因序列進(jìn)行降維特征提取反而丟失了重要的特征.但由實(shí)驗(yàn)結(jié)果可見，GRTN方法在仿真數(shù)據(jù)集上的表現(xiàn)依然優(yōu)于Direct以外的對比方法.

Table 2 Results on Simulating Network with 10% Nodes Missing

Table 3 Results on Simulating Network with 20% Nodes Missing

為探究本文方法在實(shí)際應(yīng)用中的價(jià)值，我們也在真實(shí)新冠病毒網(wǎng)絡(luò)中進(jìn)行了實(shí)驗(yàn).分析表4和表5可得，在真實(shí)的新冠病毒病毒基因序列數(shù)據(jù)集中，采用Direct進(jìn)行基因表示的實(shí)驗(yàn)結(jié)果要遠(yuǎn)低于本文提出的GRTN方法.而GRTN方法通過聚合節(jié)點(diǎn)本身的基因序列信息和傳播網(wǎng)絡(luò)的特征信息，在多數(shù)驗(yàn)證指標(biāo)上取得了最優(yōu)的成績，在Top-5和RS等指標(biāo)上結(jié)果略差于GCN方法，但其綜合表現(xiàn)在所有方法中最好.

Table 4 Results on Aussie with 10% Nodes Missing表4 缺失10%節(jié)點(diǎn)Aussie中各指標(biāo)表現(xiàn) %

Table 5 Results on Aussie with 20% Nodes Missing表5 缺失20%節(jié)點(diǎn)Aussie中各指標(biāo)表現(xiàn) %

為探究本文模型的適用性，本文選取HIV病毒傳播網(wǎng)絡(luò)進(jìn)行試驗(yàn).因?yàn)镠IV病毒基因網(wǎng)絡(luò)節(jié)點(diǎn)數(shù)目相對較少，只關(guān)注排名靠前的節(jié)點(diǎn)的指標(biāo)時(shí)(比如k-RS指標(biāo)的前3%)，實(shí)驗(yàn)結(jié)果波動(dòng)極大，由表6中結(jié)果可得，雖然HIV的病毒序列特征與COVID-19病毒的極為不同，即HIV的序列長度要遠(yuǎn)遠(yuǎn)小于COVID-19，且該網(wǎng)絡(luò)本身即包含大量缺失節(jié)點(diǎn)，但從實(shí)驗(yàn)結(jié)果可以看出，基于病毒傳播網(wǎng)絡(luò)的基因表示模型(GRTN)在處理缺失節(jié)點(diǎn)發(fā)現(xiàn)的任務(wù)上，表現(xiàn)要遠(yuǎn)遠(yuǎn)好于其他對比方法.

Table 6 Results on HIV with 30% Nodes Missing表6 缺失30%節(jié)點(diǎn)HIV數(shù)據(jù)集中各指標(biāo)表現(xiàn) %

4.2 模型分析

4.1節(jié)實(shí)驗(yàn)驗(yàn)證了GRTN模型的有效性，本節(jié)從表示聚合分模塊能力和注意力使用等多個(gè)方面分析模型的表現(xiàn)效果.

4.2.1 表示聚合分析

在GRTN模型中，節(jié)點(diǎn)的基因表示由Ef,Eb和Em經(jīng)過注意力聚合得到，在本節(jié)，用GRTN-f表示去掉目標(biāo)節(jié)點(diǎn)上文節(jié)點(diǎn)信息聚合的模型，GRTN-b表示去掉下文節(jié)點(diǎn)信息聚合的模型，用GRTN-cut表示對輸入基因序列不進(jìn)行維度約減的操作，用GRTN-2n表示不考慮目標(biāo)節(jié)點(diǎn)距離為2的鄰居信息的模型.綜合上述多個(gè)角度，對方法進(jìn)行消融實(shí)驗(yàn).

從表7結(jié)果可知，本文模型通過結(jié)合前向、后向和自身的表示信息表現(xiàn)效果最佳，缺失前向表示的信息對模型的效果影響更大，后向節(jié)點(diǎn)表示信息對模型表現(xiàn)也有較大影響，單純使用節(jié)點(diǎn)自身表示的退化模型，效果較差.當(dāng)只聚合目標(biāo)節(jié)點(diǎn)的子節(jié)點(diǎn)和父節(jié)點(diǎn)信息時(shí)，模型效果較差，這說明只使用一階鄰居的信息進(jìn)行基因表示是不充分的.當(dāng)不對輸入數(shù)據(jù)維度進(jìn)行約減時(shí)，模型效果有明顯的下降，這說明模型對基因序列降維后進(jìn)行了有效特征提取，相比不降維的全序列數(shù)據(jù)，可以去掉基因序列數(shù)據(jù)中無用的干擾信息.而對比同樣是對數(shù)據(jù)進(jìn)行降維提取特征的PCA和AE，說明基于傳播網(wǎng)絡(luò)上下文的基因表示模型效果更優(yōu).

Table 7 Ablationstudy on Aussie with 20% Nodes Missing表7 缺失20%節(jié)點(diǎn)Aussie中各模塊消融實(shí)驗(yàn) %

4.2.2 注意力權(quán)重分析

在Aussie中，本文提取節(jié)點(diǎn)的Ef、Eb和Em向量的注意力權(quán)重，統(tǒng)計(jì)如表8所示.對于處于網(wǎng)絡(luò)不同位置、屬性信息不同的節(jié)點(diǎn)，其注意力系數(shù)有明顯差異.從整體影響力來說，前向注意力的系數(shù)權(quán)重比較高，這和4.2.1節(jié)中前向節(jié)點(diǎn)影響大的結(jié)論一致.圖5展示了網(wǎng)絡(luò)中節(jié)點(diǎn)7和節(jié)點(diǎn)621的局部結(jié)構(gòu)，以節(jié)點(diǎn)7為例，該節(jié)點(diǎn)處于網(wǎng)絡(luò)中的末端位置，沒有后向節(jié)點(diǎn)，其屬性結(jié)果很大程度上和其父節(jié)點(diǎn)相關(guān)，故αf高，雖然模型給出了αb的系數(shù)，但是該節(jié)點(diǎn)的Eb為0；節(jié)點(diǎn)621同時(shí)存在子節(jié)點(diǎn)和父節(jié)點(diǎn)，注意力系數(shù)和三者皆有關(guān)系.

Table 8 Values of Forward, Backward and Self Attention表8 前向、后向和自身注意力系數(shù)

Fig. 5 Nodes attention analysis圖5 節(jié)點(diǎn)注意力分析

4.2.3 多頭注意力分析

為了探討不同注意力機(jī)制使用對模型效果的影響，列出采用不同注意力機(jī)制的模型在Aussie中的表現(xiàn)結(jié)果，如圖6所示.本節(jié)使用的不同注意力機(jī)制有：

1) average-4.采用4頭注意力機(jī)制，每個(gè)注意力的表示維度為128，將4個(gè)注意力的表示向量輸出在每個(gè)維度上取平均值.

2) concat-2.采用2頭注意力機(jī)制，每個(gè)注意力的表示維度為128，將2個(gè)注意力的表示向量輸出拼接在一起.

3) concat-4.采用4頭注意力機(jī)制，每個(gè)注意力的表示維度為128，將4個(gè)注意力的表示向量輸出拼接在一起.

4) single.單頭注意力表示，輸出的表示向量維度為128.

Fig. 6 Model performances with different attention mechanisms on Aussie with 20% nodes missing圖6 采用不同注意力機(jī)制的算法在Aussie上預(yù)測準(zhǔn)確率表現(xiàn)(20%節(jié)點(diǎn)缺失)

從圖6分析可知，采用單頭注意力的模型表示結(jié)果的隨機(jī)性較大，采用拼接操作和取均值操作均可以提高模型的效果，其中采用拼接操作的效果更優(yōu)，并且隨著拼接的注意力模塊的增加，模型表現(xiàn)有一定提高.

5 結(jié) 論

本文基于病毒傳播網(wǎng)絡(luò)構(gòu)建了一種新的基因序列表示方法，在對節(jié)點(diǎn)病毒序列進(jìn)行編碼后，使用注意力機(jī)制對節(jié)點(diǎn)及其鄰居節(jié)點(diǎn)的序列信息進(jìn)行聚合，進(jìn)而得到目標(biāo)節(jié)點(diǎn)病毒基因序列的表示.為了驗(yàn)證基因序列表示模型的有效性，本文設(shè)計(jì)一套在病毒傳播網(wǎng)絡(luò)中發(fā)現(xiàn)未采樣感染者的任務(wù)流程，可以在實(shí)際的醫(yī)學(xué)調(diào)查診斷的相應(yīng)環(huán)節(jié)進(jìn)行輔助排查.本文用網(wǎng)絡(luò)缺失點(diǎn)發(fā)現(xiàn)任務(wù)對得到的基因序列表示進(jìn)行評估，模型在模擬病毒傳播網(wǎng)絡(luò)、真實(shí)的澳大利亞新冠病毒傳播網(wǎng)絡(luò)和艾滋病毒網(wǎng)絡(luò)上的實(shí)驗(yàn)結(jié)果驗(yàn)證了模型的有效性.

作者貢獻(xiàn)聲明：馬揚(yáng)和劉澤一提出論文理論方法、思路，進(jìn)行實(shí)驗(yàn)及論文撰寫，為本文同等貢獻(xiàn)第一作者；程光權(quán)和劉忠指導(dǎo)論文思路與撰寫，為本文共同通信作者；梁星星、陽方杰、成清對論文提出修改建議.