陶斯涵 , 丁彥蕊,2

1(江蘇省媒體設(shè)計與軟件技術(shù)重點實驗室(江南大學(xué)),江蘇 無錫 214122)

2(工業(yè)生物技術(shù)教育部重點實驗室(江南大學(xué)),江蘇 無錫 214122)

生物網(wǎng)絡(luò),例如蛋白質(zhì)-蛋白質(zhì)相互作用網(wǎng)絡(luò)、代謝網(wǎng)絡(luò)、殘基相互作用網(wǎng)絡(luò)、基因表達網(wǎng)絡(luò)等[1],已經(jīng)成為大數(shù)據(jù)時代生命科學(xué)相關(guān)研究的重要數(shù)據(jù)資源.這使得生物網(wǎng)絡(luò)比對在近年來成為研究代謝、結(jié)構(gòu)、功能和進化的有效的方法.通過生物網(wǎng)絡(luò)比對的方法可以發(fā)現(xiàn)兩個或兩個以上相互作用網(wǎng)絡(luò)在拓撲和功能上的相似區(qū)域,用于研究生物分子的結(jié)構(gòu)和功能,分析生物的進化和演變.

通常,功能相似的蛋白質(zhì)分子具有相似的空間結(jié)構(gòu),而結(jié)構(gòu)上局部的差異可能會導(dǎo)致其性質(zhì)的不同,如蛋白質(zhì)的熱穩(wěn)定性、親水性、疏水性、耐酸性、耐堿性等[2].殘基相互作用網(wǎng)絡(luò)對于從系統(tǒng)角度研究蛋白質(zhì)空間結(jié)構(gòu)和蛋白質(zhì)性質(zhì)、功能的關(guān)系有著至關(guān)重要的作用.而殘基相互作用網(wǎng)絡(luò)比對,對于研究蛋白質(zhì)的分子基礎(chǔ)和空間結(jié)構(gòu)非常重要,是探究蛋白質(zhì)空間結(jié)構(gòu)與蛋白質(zhì)性質(zhì)、功能異同關(guān)系的有效方法.

目前,絕大多數(shù)的網(wǎng)絡(luò)比對方法是針對蛋白質(zhì)-蛋白質(zhì)相互作用網(wǎng)絡(luò)(簡稱PPI 網(wǎng)絡(luò))提出的.PPI 網(wǎng)絡(luò)比對可以分為全局網(wǎng)絡(luò)比對和局部網(wǎng)絡(luò)比對.局部網(wǎng)絡(luò)比對旨在找到小的、高度保守的子網(wǎng)絡(luò),而不考慮比較網(wǎng)絡(luò)的整體相似性[1],例如PathBLAST[3]、MaWISH[4]和AlignMCL[5]等.PathBLAST 是最早的局部網(wǎng)絡(luò)比對算法.它通過蛋白質(zhì)-蛋白質(zhì)相互作用網(wǎng)絡(luò)中保守路徑和比對蛋白質(zhì)之間同源性的概率來搜索高分比對.MaWISH 模擬生物的復(fù)制和刪除,并利用網(wǎng)絡(luò)的加權(quán)邊找到最大權(quán)重誘導(dǎo)子圖.AlignMCL 將單個比對圖中的許多蛋白質(zhì)相互作用網(wǎng)絡(luò)合并,并在隨后將其挖掘用于識別不同物種中的保守子網(wǎng).全局網(wǎng)絡(luò)比對旨在最大化網(wǎng)絡(luò)之間的整體匹配.它產(chǎn)生一對一的節(jié)點映射,較小網(wǎng)絡(luò)中的每個節(jié)點都被映射到較大網(wǎng)絡(luò)中的一個唯一的節(jié)點[1],例如IsoRank[6]、GRAAL[7]、MI-GRAAL[8]、HubAlign[9]、MAGNA[10]和ModuleAlign[11]等.IsoRank 是第一種全局網(wǎng)絡(luò)比對算法,將PPI 網(wǎng)絡(luò)中蛋白質(zhì)的易匹配序列和鄰域拓撲結(jié)構(gòu)比作特征值問題,并求得網(wǎng)絡(luò)的最佳匹配.GRAAL 基于所有鄰居節(jié)點的拓撲相似性,利用貪心算法找到最優(yōu)比對.MI-GRAAL 整合了網(wǎng)絡(luò)節(jié)點之間的任何數(shù)量和類型的相似性標準,并決定了相似性標準間的權(quán)重,利用最大權(quán)重雙邊圖找出最優(yōu)比對.HubAlign 使用啟發(fā)式算法,從網(wǎng)絡(luò)拓撲信息的方面對蛋白質(zhì)的拓撲結(jié)構(gòu)和功能的重要性進行評估,并優(yōu)先比對拓撲結(jié)構(gòu)重要性高的蛋白質(zhì).MAGNA 是一種基于PPI 網(wǎng)絡(luò)拓撲信息的全局網(wǎng)絡(luò)比對算法.它利用遺傳算法框架,將兩個“父代”比對通過交叉函數(shù)“交叉”產(chǎn)生優(yōu)秀的“子代”比對,保留并進入下一代,直到達到停止條件,輸出最優(yōu)的比對.ModuleAlign 利用局部信息來定義模塊的同源性分數(shù),基于參與相同模塊的功能相關(guān)蛋白質(zhì)的分層聚類,并采用迭代方案找到兩個網(wǎng)絡(luò)之間的比對.目前,蛋白質(zhì)-蛋白質(zhì)相互作用網(wǎng)絡(luò)比對主要應(yīng)用于功能預(yù)測和系統(tǒng)發(fā)生分析等相關(guān)研究.

眾所周知,影響蛋白質(zhì)相似性的因素是多方面的,其中最主要的是蛋白質(zhì)序列的相似性和空間結(jié)構(gòu)的相似性.蛋白質(zhì)序列相似的分子往往具有相似的結(jié)構(gòu)或相似的功能.但序列相似的分子在高級結(jié)構(gòu)和功能上并不具有必然的相似性[12].因此,將蛋白質(zhì)的三維結(jié)構(gòu)編碼為殘基相互作用網(wǎng)絡(luò),對網(wǎng)絡(luò)進行比對,是從系統(tǒng)角度分析蛋白質(zhì)功能與結(jié)構(gòu)及序列的關(guān)系的主要途徑.

然而,目前一般采用PPI 網(wǎng)絡(luò)比對算法對殘基相互作用網(wǎng)絡(luò)進行比對,還沒有針對殘基相互作用網(wǎng)絡(luò)的比對算法,這就忽視了氨基酸殘基本身的信息.因此,在對上述PPI 網(wǎng)絡(luò)比對的各種算法分析的基礎(chǔ)上,本文基于MAGNA 網(wǎng)絡(luò)比對方法,將蛋白質(zhì)的序列信息即殘基匹配度引入到優(yōu)化函數(shù)中,提出了針對殘基相互作用網(wǎng)絡(luò)比對的SI-MAGNA 算法.同時,與其他現(xiàn)有的基于網(wǎng)絡(luò)拓撲信息的經(jīng)典比對方法(GRAAL[7]、MI-GRAAL[8]、MAGNA[10]、CytoGEDEVO[13])相比較,邊正確性表明,SI-MAGNA 算法在殘基相互作用網(wǎng)絡(luò)比對方面優(yōu)于其他的方法.最后,對不同耐熱性的同源蛋白質(zhì)的殘基相互作用網(wǎng)絡(luò)使用SI-MAGNA 方法進行網(wǎng)絡(luò)比對和分析,探索蛋白質(zhì)結(jié)構(gòu)對其熱穩(wěn)定性的影響.

本文第1 節(jié)對SI-MAGNA 方法的算法框架和原理進行總結(jié).第2 節(jié)構(gòu)建一系列蛋白質(zhì)的殘基相互作用網(wǎng)絡(luò)對,并使用SI-MAGNA 算法進行比對.其中,第2.1 節(jié)中總結(jié)拓撲信息-序列信息權(quán)重α對比對結(jié)果的影響程度,并與MAGNA 方法進行比較;第2.2 節(jié)中使用SI-MAGNA 方法與其他現(xiàn)有的基于網(wǎng)絡(luò)拓撲信息的經(jīng)典比對方法(GRAAL[7]、MI-GRAAL[8]、CytoGEDEVO[13])相比較,以驗證SI-MAGNA 方法的優(yōu)越性;第2.3 節(jié)中介紹SI-MAGNA 方法的應(yīng)用.第3 節(jié)總結(jié)全文.

1 原理與方法

1.1 殘基相互作用網(wǎng)絡(luò)的構(gòu)建

構(gòu)成蛋白質(zhì)序列的基本單位是氨基酸殘基,簡稱殘基.將殘基相互作用網(wǎng)絡(luò)定義為一個無加權(quán)圖,節(jié)點表示殘基,邊表示不同殘基間的相互作用.為了判定相互作用的存在性,將殘基的Cα的坐標作為殘基的位置,計算不同殘基間的距離,當殘基間的距離小于6.5? 時,認為相互作用存在[14-16].假定兩個殘基相互作用網(wǎng)絡(luò)分別為G1=(V1,E1)和G2=(V2,E2),V1和V2表示節(jié)點的集合,E1和E2表示邊的集合.令m=|V1|,n=|V2|,|V1|≤|V2|.殘基相互作用網(wǎng)絡(luò)G1到G2的比對定義為f:v1→v2.節(jié)點集V1中每一個節(jié)點到V2中的節(jié)點都存在一對一的映射.

1.2 SI-MAGNA方法思想

MAGNA 方法是一種基于PPI 網(wǎng)絡(luò)的兩兩全局比對算法.它只利用網(wǎng)絡(luò)的拓撲信息產(chǎn)生比對.它基于遺傳算法框架,使用交叉函數(shù)使兩個“父代”比對產(chǎn)生一個優(yōu)秀的“子代”比對,并迭代計算直到達到停止條件,獲得比對結(jié)果[10].

1.2.1 交叉函數(shù)

交叉函數(shù)是MAGNA 方法的核心[10].它主要由以下3 個步驟組成(如圖1 所示).交叉函數(shù)通過節(jié)點標號對應(yīng)節(jié)點標號的形式來表示任何比對f(即一個對應(yīng)的排列σ).建立圖Γn,它具有節(jié)點集Sn和邊集En,其中,Sn為所有排列σ的集合,當且僅當排列σ和τ鄰接時,En表示排列σ和τ之間的邊的集合.排列σ和τ交叉(表示為σ?τ)表示圖Γn中從σ到τ最短路徑的近中點.因此,將兩排列交叉產(chǎn)生的子比對定義為兩個父比對之間的“中間”比對,子比對預(yù)計將繼承其每個父比對的大約1/2.SI-MAGNA 算法采用相同的交叉函數(shù).

Fig.1 Main operation of the cross function圖1 交叉函數(shù)的主要操作

1.2.2 基于遺傳算法的框架

遺傳算法是一種啟發(fā)式搜索算法.它借鑒生物界適者生存的進化規(guī)律.MAGNA 方法基于遺傳算法,迭代搜索以獲得最優(yōu)的比對[10].SI-MAGNA 算法采用相同的搜索策略.

通過隨機比對以獲得給定規(guī)模p的初始種群P0,種群中的成員即比對.為了防止種群規(guī)模的無限增長,只有“適者”能夠生存到下一代.對于每一代種群P,通過適應(yīng)度函數(shù)F,將種群中的成員以其比對質(zhì)量從高到低排序,設(shè)定精英率e,將種群中占比為e的“高質(zhì)量”成員保留作為子代成員直接加入下一代種群中.

接下來,使用Roulette 選擇算法選擇種群P中的成員進行“交叉”產(chǎn)生新的子代,以補足下一代中剩余的部分,所選擇成員的概率與成員的適應(yīng)度成正比.種群中成員被選中的概率可表示為

隨著算法的進行,種群中成員的適應(yīng)度(比對質(zhì)量)逐漸增加,直到達到停止條件,輸出最優(yōu)比對.

1.3 適應(yīng)度函數(shù)

適應(yīng)度函數(shù)F作為網(wǎng)絡(luò)比對的優(yōu)化條件,在算法中起到至關(guān)重要的作用.考慮殘基相互作用網(wǎng)絡(luò)的特點,這里將適應(yīng)度函數(shù)F定義為

其中,TopoScore(f)表示拓撲信息(詳見第1.3.1 節(jié)),SeqScore(f)表示序列信息(詳見第1.3.2 節(jié)),權(quán)重α用來調(diào)整拓撲信息和序列信息對比對的影響作用.拓撲信息-序列信息權(quán)重α的取值范圍是[0,1],當α等于1 時,表示只引入拓撲信息,而不考慮序列信息;當α等于0 時,表示只引入序列信息,而不考慮拓撲信息.

1.3.1 拓撲信息TopoScore(f)

拓撲信息TopoScore(f)使用對稱子結(jié)構(gòu)得分(S3)[10].由網(wǎng)絡(luò)G1(V1,E1)和G2(V2,E2)的比對f:v1→v2,設(shè)G2[Y]為點集為Y的G2的子網(wǎng),f(V1)={f(v)∈V2:v∈V1},f(E1)={(f(u),f(v))∈E2:(u,v)∈E1}.并將保守邊定義為由通過f比對的兩個網(wǎng)絡(luò)的兩條邊組成(如圖2 所示),即當G1中的節(jié)點u,v通過f分別比對上G2中的節(jié)點u′,v′時,那么邊(u,v)和邊(u′,v′)構(gòu)成一條保守邊.對稱子結(jié)構(gòu)得分(S3)表示保守邊的數(shù)量占網(wǎng)絡(luò)G1和G2[f(V1)]疊加的復(fù)合圖邊的數(shù)量的比例,它既懲罰了比對從密集區(qū)域映射到稀疏區(qū)域,又懲罰了從稀疏區(qū)域映射到密集區(qū)域[10],表示為

Fig.2 Illustration of conserved edges圖2 保守邊示意圖

1.3.2 序列信息SeqScore(f)

序列信息SeqScore(f)使用基于BLOSUM 矩陣[12]的殘基匹配度矩陣進行評價.本文基于BLOSUM 矩陣是因為該矩陣基于蛋白質(zhì)進化的星狀模型(即忽略物種近端和遠端的關(guān)系)和區(qū)塊中的保守位置與置換關(guān)系進行計分,這對于發(fā)現(xiàn)同源蛋白質(zhì)中的保守區(qū)域有非常重要的作用.由于蛋白質(zhì)的功能是由序列和特定的空間結(jié)構(gòu)等因素共同決定,在殘基相互作用網(wǎng)絡(luò)比對中,只利用網(wǎng)絡(luò)的拓撲信息無法獲得在網(wǎng)絡(luò)拓撲和一級結(jié)構(gòu)方面更準確更合理的比對結(jié)果.因此,本文引入蛋白質(zhì)的序列信息,并定義殘基匹配度矩陣,加入到網(wǎng)絡(luò)比對的適應(yīng)度函數(shù)中.

殘基匹配度矩陣是通過統(tǒng)計蛋白質(zhì)序列的替換率而得到的氨基酸矩陣,由蛋白質(zhì)序列塊比對推導(dǎo)得出.其大致步驟如下.

首先,消除相似度小于指定閾值的序列,計算數(shù)據(jù)中每個氨基酸組合發(fā)生的可能性和該組合預(yù)期發(fā)生的可能性比率Log-odds,表示為

其中,P(O)表示觀察的可能性,P(E)表示預(yù)期的可能性.

然后,基于此計算殘基匹配度矩陣,表示為

其中,ij是氨基酸i與j在同源序列中相互替換的概率;qi和qj是氨基酸出現(xiàn)在任意蛋白質(zhì)序列中的概率;λ是一個尺度參數(shù),使每個得分更易取整.

使用殘基匹配度矩陣作為序列信息計算網(wǎng)絡(luò)G1中的任意節(jié)點vi與網(wǎng)絡(luò)G2中任意節(jié)點vj相匹配時的替換率,并作為匹配度并用矩陣存儲.同時,為了保證適應(yīng)度函數(shù)的合理性,其數(shù)量級需要與拓撲信息一致.這里采用min-max 方法對殘基匹配度矩陣進行標準化,具體表示為

其中,x表示某兩個殘基的匹配度,x*表示標準化后某兩個殘基的匹配度,max 為矩陣中的最大值,min 為矩陣中的最小值.此時,序列信息與拓撲信息處于同一數(shù)量級,適合進行綜合比對評價.

綜上所述,算法步驟如下,并可由算法流程圖表示(如圖3 所示).

Fig.3 SI-MAGNA algorithmic framework圖3 SI-MAGNA 算法流程圖

(1)輸入網(wǎng)絡(luò)G1,G2及相關(guān)參數(shù):代數(shù)N、精英率e、種群規(guī)模p.

(2)隨機產(chǎn)生種群規(guī)模p的比對初始種群P0.

(3)以初始種群P0作為父代種群.

(4)設(shè)置代數(shù)計數(shù)器n=1.

(5)通過適應(yīng)度函數(shù)F,計算父代種群P中成員的適應(yīng)度,并進行排序.

(6)判斷并保留適應(yīng)度高的父代種群成員,保留比例為精英率e,即

(7)通過Roulette 選擇算法和交叉函數(shù)產(chǎn)生適應(yīng)度較高的比對成員,補足剩余的部分,即

(8)將步驟(5)和步驟(6)產(chǎn)生的比對成員組成子代成員P′.

(9)當n達到代數(shù)N時,終止循環(huán).

(10)輸出網(wǎng)絡(luò)比對結(jié)果.

1.4 網(wǎng)絡(luò)比對算法的性能評估

目前,邊正確性(edge correctness,簡稱EC)廣泛應(yīng)用于網(wǎng)絡(luò)比對算法的性能評估,通常通過計算EC 來評價網(wǎng)絡(luò)比對結(jié)果的優(yōu)劣[7,8,10,13].EC 表示比對中保守邊的數(shù)量占網(wǎng)絡(luò)G1邊的數(shù)量的比例.當一個輸入網(wǎng)絡(luò)G1與另一網(wǎng)絡(luò)G2是同構(gòu)的時,它取得最高值100%.EC 可表示為

2 結(jié)果與討論

在本文中,為了驗證SIMAGNA 網(wǎng)絡(luò)比對方法,從RCSB PDB 數(shù)據(jù)庫(http://www.rcsb.org/pdb/home/home.do)下載了不同蛋白質(zhì)的三維結(jié)構(gòu)信息和序列信息,構(gòu)建殘基相互作用網(wǎng)絡(luò)進行比對.在這里,主要對以下9 組網(wǎng)絡(luò)對數(shù)據(jù)進行比對分析,見表1,其中,每組網(wǎng)絡(luò)對的序列相似性由BLAST 序列比對算法(https://blast.cbi.nlm.nih.gov/Blast.cgi)兩兩比對得出.

Table 1 Information of residue interaction network for different proteins表1 不同蛋白質(zhì)的殘基相互作用網(wǎng)絡(luò)的基本信息

為了評估引入了蛋白質(zhì)序列信息(殘基匹配度)的SI-MAGNA 方法以及其他網(wǎng)絡(luò)比對方法(GRAAL[7]、MI-GRAAL[8]、MAGNA[10]和CytoGEDEVO[13]),選擇受到廣泛認可的EC 作為評估網(wǎng)絡(luò)比對質(zhì)量的標準[7,8,10,13].

2.1 基于SI-MAGNA方法的殘基相互作用網(wǎng)絡(luò)比對

本文選擇了來自不同物種的具有相似功能或結(jié)構(gòu)的蛋白質(zhì)進行殘基相互作用網(wǎng)絡(luò)的比對,以驗證SIMAGNA 網(wǎng)絡(luò)比對方法.實驗設(shè)置初始種群P0為15 000,精英率e為0.5.

拓撲信息-序列信息權(quán)重α是影響比對結(jié)果的重要因子,它的取值范圍是[0,1].在實驗中,設(shè)定代數(shù)N為2000,α取值步長為0.1.上述9 組殘基相互作用網(wǎng)絡(luò)比對組實驗結(jié)果如圖4 所示,圖中橫軸為拓撲信息-序列信息權(quán)重α取值,縱軸為EC 值,每條折線分別表示各組比對.

Fig.4 Influence of TopoScore-SeqScore weight α on the SI-MAGNA alignment results圖4 拓撲-序列信息權(quán)重α對SI-MAGNA 方法下比對結(jié)果的影響

從實驗中發(fā)現(xiàn),當α=0 時,表示只引入網(wǎng)絡(luò)的序列信息到適應(yīng)度函數(shù)中進行比對,其EC 接近0,表明網(wǎng)絡(luò)比對效果非常差;同樣地,當α=1 時,表示只引入網(wǎng)絡(luò)的拓撲信息到適應(yīng)度函數(shù)中進行比對,其邊EC 也很低.這說明在殘基相互作用網(wǎng)絡(luò)的比對中,僅依靠網(wǎng)絡(luò)的拓撲信息或序列信息均不能獲得很好的比對結(jié)果.當α={0.1,0.2,…,0.9}時,表示分別以一定比例引入拓撲信息和序列信息到適應(yīng)度函數(shù)中,發(fā)現(xiàn)比對結(jié)果的EC 有較大的提升;同時,比對得到的共同保守子圖規(guī)模也更大.這說明盡管網(wǎng)絡(luò)拓撲信息是網(wǎng)絡(luò)比對優(yōu)化的基礎(chǔ),但在殘基相互作用網(wǎng)絡(luò)比對中,僅用拓撲信息進行比對優(yōu)化不夠全面,序列信息對于殘基相互作用網(wǎng)絡(luò)而言也十分重要,將它引入比對優(yōu)化的適應(yīng)度函數(shù)中,能夠明顯地提升比對質(zhì)量,獲得更優(yōu)的結(jié)果.

為了探究改進后的方法在殘基相互作用網(wǎng)絡(luò)比對方面的可靠性,在原方法探討所得的基礎(chǔ)上對代數(shù)N進行了進一步實驗,同時以原方法作為參照.在實驗中,拓撲信息-序列信息權(quán)重α取各組在N=2000 時產(chǎn)生最優(yōu)比對時的值,代數(shù)N取值為2000、4000、6000、8000、10000、120000.上述9 組殘基相互作用網(wǎng)絡(luò)對比對實驗結(jié)果如圖5 所示,圖中橫軸為代數(shù)N,縱軸為EC,實線表示SI-MAGNA 算法,虛線表示MAGNA 算法,各組網(wǎng)絡(luò)對分別為:(a)1V8I&1MP2;(b)9AME&1WVO;(c)1TUX&1E0W;(d)1XNB&1XND;(e)1XXN&1M4W;(f)1YNA&1PVX;(g)3QMM&1I6W;(h)3D2B&1I6W;(i)3QMM&3D2B.從實驗中發(fā)現(xiàn),對于改進后的方法SI-MAGNA,代數(shù)N的取值變化對于比對結(jié)果的影響很小,比對結(jié)果的EC 保持在較高范圍內(nèi)浮動.而對于MAGNA 方法,代數(shù)N的取值變化對于比對結(jié)果的影響較大,比對結(jié)果的EC 在大多數(shù)情況下是在逐漸增加后進入平穩(wěn),也會在很少情況下出現(xiàn)比對結(jié)果的EC 保持在相對高的范圍內(nèi)浮動.這說明在殘基相互作用網(wǎng)絡(luò)比對方面,MAGNA 方法在大多數(shù)情況下比SI-MAGNA 需要更多的代數(shù)才能得到較好的比對結(jié)果,且最終比對結(jié)果的EC 往往不如改進后的方法SI-MAGNA 高.此外,由于某些殘基相互作用網(wǎng)絡(luò)對的網(wǎng)絡(luò)規(guī)模較小,會在較少情況下出現(xiàn)比對結(jié)果比較好(EC 值較高)且在小范圍浮動的情況.圖中結(jié)果說明,SI-MAGNA 方法表現(xiàn)更加穩(wěn)定,在較少的代數(shù)N時,已經(jīng)能夠獲得較好且相對穩(wěn)定的比對結(jié)果.并且相對于原方法,改進后的SI-MAGNA 方法不但取得了更優(yōu)的比對結(jié)果,而且因為更少的代數(shù)而提高了比對的效率.

Fig.5 Influence of generation N on the SI-MAGNA alignment results圖5 代數(shù)N 對SI-MAGNA 方法下比對結(jié)果的影響

2.2 與其他網(wǎng)絡(luò)比對方法對比

為了驗證引入蛋白質(zhì)序列信息的SI-MAGNA 算法在殘基相互作用網(wǎng)絡(luò)比對方面的優(yōu)越性,選擇基于網(wǎng)絡(luò)拓撲信息的其他網(wǎng)絡(luò)比對方法與之比較.由于目前絕大多數(shù)的生物網(wǎng)絡(luò)比對方法是針對蛋白質(zhì)-蛋白質(zhì)相互作用網(wǎng)絡(luò)的,這些方法中的許多方法融入了蛋白質(zhì)節(jié)點的生物信息(如BLAST-E 值、基因本體論正確性[17]等),因此,這類方法不適用于殘基相互作用網(wǎng)絡(luò)的比對.為了保證比對實驗的合理性,目前一般采用僅基于網(wǎng)絡(luò)拓撲信息的PPI 網(wǎng)絡(luò)比對算法對殘基相互作用網(wǎng)絡(luò)進行比對,選擇的方法是GRAAL、MI-GRAAL 和CytoGEDEVO方法.4 種比對方法的基本設(shè)置如下.

(1)SI-MAGNA 方法.在實驗中,設(shè)置初始種群P0為15 000,精英率e為0.5,代數(shù)N為2 000,拓撲信息-序列信息權(quán)重α在[0,1]中取值,步長為0.1,并在結(jié)果中擇優(yōu).

(2)GRAAL 方法.GRAAL 方法是基于graphlet 度特征相似性的網(wǎng)絡(luò)比對算法.它基于鄰居節(jié)點的拓撲相似性尋找最優(yōu)比對,適用于任何網(wǎng)絡(luò)的比對.為了保證算法比較的合理性,在實驗中,節(jié)點-邊權(quán)重α在[0,1]中取值,步長為0.1,其他使用默認設(shè)置,并選擇最優(yōu)比對結(jié)果.

(3)MI-GRAAL 方法.MI-GRAAL 方法整合了網(wǎng)絡(luò)節(jié)點之間多種類型的相似性標準,并決定它們的權(quán)重,以此找到最優(yōu)比對.為了保證算法比較的合理性,在實驗中,剔除默認方法中基于蛋白質(zhì)生物信息的相似性標準(BLAST-E 值),選擇基于網(wǎng)絡(luò)拓撲性的相似性標準的組合(graphlet 度特征、度、聚類系數(shù)、節(jié)點離心率和介數(shù)得分),使方法完全基于網(wǎng)絡(luò)的拓撲相似性進行比對,并選擇各相似性標準組合中的最優(yōu)比對結(jié)果.

(4)CytoGEDEVO 方法.GEDEVO 方法[18]基于進化算法,使用圖編輯距離作為優(yōu)化模型來找到最佳比對.CytoGEDEVO 方法是基于GEDEVO 方法的CytoScape 軟件上的擴展.為了保證算法比較的合理性,在實驗中,設(shè)置迭代次數(shù)為2 000,其他使用默認設(shè)置.

上述9 組殘基相互作用網(wǎng)絡(luò)對比對實驗結(jié)果如表2 和圖6 所示(圖6 中,橫軸表示各組殘基相互作用網(wǎng)絡(luò)對,縱軸為EC,每組柱狀圖從左到右分別表示SI-MAGNA、GRAAL、MI-GRAAL、CytoGEDEVO 的比對結(jié)果的EC).從實驗中發(fā)現(xiàn),SI-MAGNA 方法在各組網(wǎng)絡(luò)對中均獲得更優(yōu)的比對結(jié)果,具有更高的EC.這證實了SI-MAGNA 方法在殘基相互作用網(wǎng)絡(luò)比對方面的表現(xiàn)更為出色,優(yōu)于另外3 種拓撲網(wǎng)絡(luò)比對方法.

Table 2 Residue interaction network alignment results with different algorithms表2 不同方法下殘基相互作用網(wǎng)絡(luò)的比對結(jié)果

Fig.6 Residue interaction network alignment results (EC)with different algorithms圖6 不同方法下的殘基相互作用網(wǎng)絡(luò)比對結(jié)果(EC)

2.3 SI-MAGNA方法的應(yīng)用

為了探究蛋白質(zhì)空間結(jié)構(gòu)與熱穩(wěn)定性的聯(lián)系[19-22],本文構(gòu)建了一系列的殘基相互作用網(wǎng)絡(luò),利用SI-MAGNA 方法進行比對和分析.通過殘基相互作用網(wǎng)絡(luò)比對,發(fā)現(xiàn)殘基相互作用網(wǎng)絡(luò)中節(jié)點和邊的匹配,尋找相似的子結(jié)構(gòu)(共同保守子圖).殘基相互作用網(wǎng)絡(luò)對的基本信息和比對結(jié)果見表3.

Table 3 Information and SI-MAGNA alignment results of different residue interaction networks表3 不同網(wǎng)絡(luò)的基本信息和比對結(jié)果

使用CytoScape(http://www.cytoscape.org/)軟件對SI-MAGNA 方法產(chǎn)生的比對結(jié)果進行進一步處理分析,結(jié)果如圖7 所示,圖中各組網(wǎng)絡(luò)對分別為:(a)1V8I(上左)、1MP2(上右)和兩者比對的共同保守子圖(下);(b)9AME&1WVO(上右)和兩者比對的共同保守子圖(下);(c)1TUX(上左)、1E0W(上右)和兩者比對的共同保守子圖(下);(d)1XNB(上左)1XND(上右)和兩者比對的共同保守子圖(下);(e)1XXN(上左)1M4W(上右)和兩者比對的共同保守子圖(下);(f)1YNA(上左)1PVX(上右)和兩者比對的共同保守子圖(下);(g)3QMM(上左)1I6W(上右)和兩者比對的共同保守子圖(下);(h)3D2B(上左)1I6W(上右)和兩者比對的共同保守子圖(下);(i)3QMM(上左)3D2B(上右)和兩者比對的共同保守子圖(下).在圖7 每一組比對中,上方左右兩個殘基相互作用網(wǎng)絡(luò)分別代表比對的殘基相互作用網(wǎng)絡(luò)對G1和G2;下方一個殘基相互作用網(wǎng)絡(luò)為比對發(fā)現(xiàn)的子結(jié)構(gòu)(共同保守子圖),其中,節(jié)點表示為xi=yj,xi∈Vi,yj∈V2.這說明在比對中,節(jié)點xi和yj相匹配,邊表示比對中的保守邊.通過比較可以發(fā)現(xiàn)兩個源網(wǎng)絡(luò)與比對產(chǎn)生的子結(jié)構(gòu)(共同保守子圖)之間的異同,相似性高的部分通常對應(yīng)于實現(xiàn)相似蛋白質(zhì)功能的重要區(qū)域,而產(chǎn)生差異的部分則有很大的可能是蛋白質(zhì)性質(zhì)(如熱穩(wěn)定性)產(chǎn)生差異的原因.

以表3 中b 組比對9AME 和1WVO 為例,其與其余8 組網(wǎng)絡(luò)對相比較.盡管9AME 和1WVO 的序列相似性不高(39%),但通過殘基相互作用網(wǎng)絡(luò)比對發(fā)現(xiàn),兩者空間結(jié)構(gòu)相似性很高,殘基相互作用網(wǎng)絡(luò)比對的EC 達到85.2%,且兩者在各自生物體中實現(xiàn)相似功能.這也為“蛋白質(zhì)的功能由序列信息和特定空間結(jié)構(gòu)共同決定”[12]提供依據(jù).來自高緯度溫帶海域魚類的III 型抗凍蛋白質(zhì)異構(gòu)體9AME 和來自人類的唾液酸合成酶AFL 結(jié)構(gòu)1WVO 為兩種生物的同源蛋白質(zhì),在各自生物體中實現(xiàn)非常相似的功能(抗凍蛋白).9AME 的最適反應(yīng)溫度為273K,1WVO 的最適反應(yīng)溫度為293K～310K,兩者在穩(wěn)定性和活性方面具有不同的溫度依賴性[23].為了了解兩者網(wǎng)絡(luò)結(jié)構(gòu)的異同,對網(wǎng)絡(luò)比對發(fā)現(xiàn)的子結(jié)構(gòu)(共同保守子圖)和二級結(jié)構(gòu)相對應(yīng)進行標記,如圖8 所示,其中,黃色標記的節(jié)點表示兩源網(wǎng)絡(luò)匹配上的α-螺旋,橙色標記的節(jié)點表示兩源網(wǎng)絡(luò)匹配上的β-折疊,紫色標記的節(jié)點表示兩源網(wǎng)絡(luò)匹配上的310-螺旋,棕色標記的節(jié)點表示兩源網(wǎng)絡(luò)匹配上的β-橋,藍色標記的節(jié)點表示兩源網(wǎng)絡(luò)匹配上的彎曲,青色標記的節(jié)點表示兩源網(wǎng)絡(luò)匹配上的氫鍵轉(zhuǎn)折,綠色標記的節(jié)點表示9AME 中獨有的1 個310-螺旋.通過殘基相互作用網(wǎng)絡(luò)比對發(fā)現(xiàn),兩者的結(jié)構(gòu)均主要由1 個α-螺旋(9AME 中為殘基37-40,1WVO 中為殘基41-44)、2 個310-螺旋(9AME 中為殘基19-21,57-59,1WVO 中為殘基23-25,61-63)和2 個β-折疊(9AME 中為殘基4-7,22-25,1WVO 中為殘基8-11,26-29)組成,此外,兩者均包含一些β-橋和彎曲結(jié)構(gòu).同時,9AME 中還獨有1個310-螺旋(殘基34-36)結(jié)構(gòu),而在1WVO 中,無二級結(jié)構(gòu)相同的結(jié)構(gòu)與之匹配.根據(jù)殘基相互作用網(wǎng)絡(luò)比對結(jié)果,使用PyMOL(https://pymol.org/2/)軟件進一步對9AME 和1WVO 的序列和二級結(jié)構(gòu)圖進行標記,這能夠更直觀地對照和分析兩者在序列信息和空間結(jié)構(gòu)上的異同,如圖9 所示,其中,黃色標記的序列和二級結(jié)構(gòu)表示兩源網(wǎng)絡(luò)匹配上的α-螺旋,橙色標記的序列和二級結(jié)構(gòu)表示兩源網(wǎng)絡(luò)匹配上的β-折疊,紫色標記的序列和二級結(jié)構(gòu)表示兩源網(wǎng)絡(luò)匹配上的310-螺旋,棕色標記的序列和二級結(jié)構(gòu)表示兩源網(wǎng)絡(luò)匹配上的β-橋,藍色標記的序列和二級結(jié)構(gòu)表示兩源網(wǎng)絡(luò)匹配上的彎曲,青色標記的序列和二級結(jié)構(gòu)表示兩源網(wǎng)絡(luò)匹配上的氫鍵轉(zhuǎn)折,綠色標記的序列和二級結(jié)構(gòu)表示9AME 中獨有的1 個310-螺旋.

Fig.7 Common conserved subgraph by SI-MAGNA alignment algorithm圖7 SI-MAGNA 網(wǎng)絡(luò)比對方法產(chǎn)生的共同保守子圖

Fig.8 Common conserved subgraph of the 9AME&1WVO by SI-MAGNA alignment algorithm圖8 SI-MAGNA 方法發(fā)現(xiàn)9AME&1WVO 的共同保守子圖

Fig.9 Illustration of sequence and secondary structure of the 9AME&1WVO圖9 9AME&1WVO 的序列和二級結(jié)構(gòu)對照圖

從圖9 中可以發(fā)現(xiàn),9AME 和1WVO 構(gòu)成相同二級結(jié)構(gòu)的序列相似性很低,但兩者主要的二級結(jié)構(gòu)卻非常相似.因此可以推斷,盡管形成這些特定結(jié)構(gòu)的殘基有所不同,但兩者共同包含的重要結(jié)構(gòu)(1 個α-螺旋、2 個310-螺旋和2 個β-折疊)以及這些結(jié)構(gòu)間的重要相互作用與兩者在各自生物體中實現(xiàn)類似的功能具有極大的關(guān)系[24-25],可以從這部分的重要殘基及其與周圍結(jié)構(gòu)間的重要相互作用入手,進一步探討實現(xiàn)抗凍蛋白功能的結(jié)構(gòu).9AME 中獨有的1 個310-螺旋結(jié)構(gòu)以及它與周圍結(jié)構(gòu)間的重要相互作用,可能是導(dǎo)致兩者在穩(wěn)定性和活性方面具有不同溫度依賴性的原因[26],因此可以從這部分的重要殘基切入,從結(jié)構(gòu)的角度進一步研究熱穩(wěn)定性產(chǎn)生差異的根源.

3 結(jié)果與討論

針對殘基相互作用網(wǎng)絡(luò)比對,對蛋白質(zhì)-蛋白質(zhì)相互作用網(wǎng)絡(luò)比對算法MAGNA 進行了改進,在其優(yōu)化函數(shù)中引入了蛋白質(zhì)的序列信息(即殘基匹配度),并探討拓撲信息和序列信息對網(wǎng)路比對的影響程度,使改進的SI-MAGNA 算法更加適用于殘基相互作用網(wǎng)絡(luò)的比對.與此同時,通過實驗證實,SI-MAGNA 方法在殘基相互作用網(wǎng)絡(luò)比對方面的表現(xiàn)是優(yōu)秀且穩(wěn)定的,并且優(yōu)于現(xiàn)有的基于網(wǎng)絡(luò)拓撲的經(jīng)典網(wǎng)絡(luò)比對方法.此外,構(gòu)建了同源蛋白質(zhì)的殘基相互作用網(wǎng)絡(luò),使用SI-MAGNA 方法進行網(wǎng)絡(luò)比對和分析,探索蛋白質(zhì)結(jié)構(gòu)對蛋白質(zhì)性質(zhì)、功能的影響.

殘基相互作用網(wǎng)絡(luò)比對方法將會成為研究蛋白質(zhì)的空間結(jié)構(gòu)、性質(zhì)和功能的重要工具.通過殘基相互作用網(wǎng)絡(luò)比對方法進行比對分析,解釋蛋白質(zhì)結(jié)構(gòu)的形成機理、探究蛋白質(zhì)結(jié)構(gòu)-功能關(guān)系,并可以將其運用于分子設(shè)計、分子篩選、藥物設(shè)計等諸多領(lǐng)域.