基因調控網絡的重構可以挖掘出基因之間潛在的調控作用關系，幫助我們更深刻地理解復雜的調控機制，在系統生物學中，這已經成為人們研究的一個熱點問題。現在已經產生了大量推斷網絡構建的模擬理論和計算方法，其中一些信息論方法在計算中容易遺漏一些調控關系，這些調控關系可能代表基因間真實的調控關系。為了克服這個缺點和進一步提高網絡重構的精度，本文提出了一種新的基因調控網絡構建方法，它通過重采樣和加權平均策略，利用條件互信息來進行基因調控網絡重構。算法首先從基因表達數據中重采樣得到一系列子數據，然后通過條件互信息方法構建一系列子數據的網絡，最后利用這些子網絡加權平均得到最終的基因調控網絡。本文的算法在DREAM3模擬數據和SOS真實網絡數據上進行驗證，和其它一些流行的方法相比，本文提出的網絡重構算法取得了較高的精確度和準確度。
引用本文： 劉飛. 重采樣條件互信息構建基因調控網絡算法研究. 生物醫學工程學雜志, 2016, 33(5): 985-990. doi: 10.7507/1001-5515.20160158 復制
0 引言
基因調控網絡(gene regulatory networks, GRNs)描述了細胞內基因之間復雜的相互作用關系，可以識別生命現象并揭示生命活動的本質，它已成為生物信息學領域一個重要的研究話題。現在國內外已經產生了大量的GRNs重構模型和算法[1-4]，各種方法都有其應用的局限性和優缺點[5-7]。一種很流行的方法是基于動態貝葉斯網絡來推測基因之間相互作用的關系[8-10]，這種方法計算性能良好，但所處理的網絡規模因其計算時間復雜度而受到限制。此外，常微分方程是基于迭代的方法來構建GRNs[11-12]，它是一種相對簡單的模型，可以定量地預測GRNs的行為。回歸模型算法也可以用來推測調控基因和靶基因之間的相互作用關系，如LASSO(the least absolute shrinkage and selectionator operator)算法[13-14]、彈性網絡回歸算法[15-16]和貝葉斯平均模型算法[17]等。信息論也被用來構建基因調控網絡[18-20]，它不僅可以計算基因間的非線性關系，而且對基因樣本個數要求較低，可以克服基因表達數據高維低樣本的缺點。但是信息論方法在計算中會遺漏一些調控關系，而這些調控關系可能代表基因間真實的相互作用關系。為了克服這個缺點以及進一步提高網絡重構的精度，本文提出了一種新的基因調控網絡構建的方法。算法首先從基因表達數據重采樣得到一系列子數據，然后通過條件互信息(conditional mutual information, CMI)方法構建一系列子數據的網絡，最后利用這些子網絡加權平均得到最終的基因調控網絡。
本文算法的流程如圖 1所示，它是一種基于多個子網絡的集成算法，它的優點就是能恢復早期計算過程中誤刪除的“真實基因調控關系”，以此來提高基因調控網絡構建的精度。首先對原始基因表達數據進行刀切法(jackknife)重采樣[21]，得到一系列基因表達數據的子數據集；然后在各個子數據集上依次使用條件互信息推斷出對應的基因調控子網絡；最終使用算術平均策略融合各個子網絡中邊的支持度，根據與閾值的關系得出最終基因調控網絡拓撲結構。
圖1
算法流程圖
Figure1.
Flow chart of the method
1 數據集和理論方法
1.1 數據集
為了檢測本文網絡構建算法的精確性，我們對不同的實驗數據集進行網絡構建，然后與實驗數據集的真實網絡進行比對，以此來評價本文算法的優劣。基因表達數據采用微陣列基因芯片技術，它運用計算機芯片把核酸固定在平面載體上形成檢測器件，將待測樣本標記后同芯片進行雜交，芯片上的探針分子和樣本中的標記分子特異性結合后，用激光共聚焦熒光掃描儀截取信息，然后處理分析信息得到信號值。信號值代表了結合在探針上的待測樣本中特定大分子的信息。對基因芯片數據加以整理和分析，并發展對應的數學模型和算法來闡明生物現象，是近幾年系統生物學研究的熱點。基因表達譜如表 1所示，行代表不同的基因名稱，列代表每個基因在不同條件下的表達值。通過基因芯片數據分析基因之間的調控關系，構建復雜系統的基因調控網絡，從整體和系統的角度來研究生命的規律。
表1
基因表達值
Table1.
Gene expression profiles
本文采用不同類型的數據對所提出的算法進行檢測評估，這些數據包括計算機模擬數據和真實的生物分子網絡數據。模擬數據來源于廣泛采用的DREAM(Dialogue for Reverse Engineering Assessments and Methods)競賽數據[6]，該競賽數據中包含基因表達數據和標準網絡，它是經實驗驗證了的酵母(Yeast)和大腸桿菌(Escherichia coli)調控網絡。這里選取了DREAM3的兩個子網絡對應的基因表達數據進行實驗，這兩個網絡分別含有10、50個頂點(基因)，以及10、77條真實作用邊(基因間調控關系)。對于真實基因表達數據，我們采用大腸桿菌SOS DNA修復網絡的真實干擾實驗數據[22-23]，該標準網絡含有9個節點和24條相互作用邊。
1.2 信息論相關定義
基于信息論的互信息已經廣泛應用于基因調控網絡的重構，而且取得了較好的效果。它是衡量兩個基因變量X和Y之間的相關性，現給定一個基因變量X，它的信息熵H(x)描述如下：
| $H\left( X \right) = - \sum\limits_{x \in X} {p\left( x \right)\log p\left( x \right)} $ |
其中p(x)表示基因變量X為x時的概率值，基因變量X和Y的聯合熵H(X, Y)可以被定義為：
| $H\left( {X,Y} \right) = - \sum\limits_{x \in X,y \in Y} {p\left( {x,y} \right)\log p\left( {x,y} \right)} $ |
其中p(x, y)表示基因變量X和Y分別為x和y時的聯合概率值，基因變量X和Y的互信息值I(X, Y)定義如下：
| $I\left( {X,Y} \right) = - \sum\limits_{x \in X,y \in Y} {p\left( {x,y} \right)\log \frac{{p\left( {x,y} \right)}}{{p\left( x \right)p\left( y \right)}}} $ |
由公式(2)和(3)互信息值I(X, Y)被重寫為：
| $I\left( {X,Y} \right) = H\left( X \right) + H\left( Y \right) - H\left( {X,Y} \right)$ |
這里H(X, Y)表示基因變量X和Y的聯合熵，互信息值的大小表示基因相互作用關系的強弱，當互信息值非常小或者接近于零時表示兩個基因變量相互獨立。條件互信息公式定義如下：
| $I\left( {X,Y|Z} \right) = H\left( {X,Z} \right) + H\left( {Y,Z} \right) - H\left( Z \right) - H\left( {X,Y,Z} \right)$ |
其中H(X, Y)表示基因變量X，Y和Z的聯合熵，這里我們的熵用高斯核概率密度函數來估計[24]，它的計算公式下：
| $\begin{gathered} P\left( {{X_i}} \right) = \frac{1}{N}\sum\limits_{j = 1}^N {\frac{1}{{{{\left( {2\pi } \right)}^{n/2}}{{\left| C \right|}^{n/{\text{2}}}}}}} \hfill \\ \exp \left( { - \frac{1}{2}{{\left( {{X_j} - {X_i}} \right)}^T}{C^{ - 1}}\left( {{X_j} - {X_i}} \right)} \right) \hfill \\ \end{gathered} $ |
其中C表示基因變量X的協方差矩陣，｜C｜表示矩陣C的行列式，N是每個基因的樣本個數，n是基因變量的個數。由公式(1)和(6)得出基因變量X的信息熵為：
| $H\left( X \right) = \log \left[ {{{\left( {2\pi e} \right)}^{n/2}}{{\left| C \right|}^{1/2}}} \right] = \frac{1}{2}\log {\left( {2\pi e} \right)^n}\left| C \right|$ |
由公式(7)，互信息公式(4)被重寫為：
| $I\left( {X,Y} \right) = \frac{1}{2}\log \frac{{\left| {C\left( X \right)} \right| \cdot \left| {C\left( Y \right)} \right|}}{{\left| {C\left( {X,Y} \right)} \right|}}$ |
同理，條件互信息公式(5)被重寫為：
| $I\left( {X,Y,Z} \right) = \frac{1}{2}\log \frac{{\left| {C\left( {X,Z} \right)} \right| \cdot \left| {C\left( {Y,Z} \right)} \right|}}{{\left| {C\left( Z \right)} \right| \cdot \left| {C\left( {X,Y,Z} \right)} \right|}}$ |
1.3 重采樣構建網絡
矩陣G=[gim]N×M, (i=1, 2, …, N; m=1, 2, …, M)表示基因表達數據，其中i代表基因，m代表樣本或者采樣條件，gim表示第i個基因在第m個采樣條件下的基因表達值。Xi=[gi1, gi2, …, giM]則表示矩陣G中的第i行，由第i個基因在M個不同采樣條件下所得到的基因表達值組成。Cm=[g1m, g2m, …, gNm]T則表示矩陣G中的第m列，由在N個基因在第m個采樣條件下所得到的基因表達值組成。重采樣過程就是依次刪除G矩陣中的一列Cm，則剩余的M-1列則組成了一個新的矩陣G′=[gil]N×(M-1)。經過重采樣過程之后，我們可以獲得原始基因表達值的M個數據子集。
分別在M個基因表達子數據集上用互信息理論得到對應的M個子基因調控網絡，使用Lijm表示第m個子網絡中節點i和j之間邊的條件互信息值，使用算術平均融合策略計算一致網絡中節點i和j之間邊的條件互信息值Lij：
| ${L_{ij}} = \frac{1}{M}\sum\limits_{m = 1}^M {L_{ij}^m} $ |
根據一致網絡中邊的Lij值與所設定的閾值α的大小關系，可以確定在最終網絡中節點i和j之間是否存在連接關系。M個子網絡的相應邊經過算術平均融合策略處理之后，如果計算得到的互信息值Lij大于所設定的閾值，則認為在最終的一致網絡中節點i和j之間存在著連接，反之則沒有連接。經過上述過程就可以重構出最終的基因調控網絡。
2 實驗仿真與結果討論
2.1 評價標準
網絡預測的結果用如下指標來衡量：真實網絡的邊數(true edges, TE)，網絡預測出的邊數(predicted edges, PE)，預測出正確的邊數(predicted correct edges, PCE)，預測出錯誤的邊數(predicted wrong edges, PWE)。同時，按如下公式計算：正確率=PCE/PE，錯誤率=PWE/PE，相對正確率=PCE/TE，以此來衡量構建算法的優劣性。
為了驗證本文算法的有效性，我們和已有的具有較好性能的算法進行比較。在本文中我們選取了3種常用的GRNs構建算法：線性規劃算法(linear programming, LP)[25]、回歸模型方法(LASSO)[26]和基于信息論的方法(algorithm for the reconstruction of accurate cellular networks, ARACNE)[27]。LP算法是一種線性規劃模型，結構較為簡單，計算速度快。LASSO是基于回歸的一種GRNs構建算法，在回歸系數的絕對值之和小于一個常數的約束條件下，使殘差平方和最小化。這兩種模型可以優化算法的獎懲項來控制網絡的稀疏度，并且取得較好的效果，但是缺點是精度過分依賴模型的參數。ARACNE算法是一種基于互信息的基因調控網絡構建算法，能夠處理高維低樣本數據，具有較高的性能，但不能區分直接和間接調控關系而導致較高的假陽性。
2.2 模擬數據的網絡重構
首先，我們采用DREAM3競賽數據中的酵母基因表達數據來構建基因調控網絡。我們對10個基因的酵母數據網絡進行分析。如圖 2(a)所示是10個基因節點的真實網絡，包含10條邊。如圖 2(b)所示是本文算法推測出的基因調控網絡，推測出10條正確的邊，推測出1條錯誤的邊(G2-G4)，預測邊的相對正確率為100%，網絡預測的正確率高達90.91%。如表 2所示詳細列出了本文算法和其它算法在性能上的比較結果，從表中可以看出本文算法在預測出正確的邊數、預測的正確率、相對正確率、錯誤的邊數以及預測錯誤率方面都明顯好于LP算法、LASSO算法和ARANCE算法，我們算法預測的正確率分別比LP算法、LASSO算法和ARANCE算法高出49.24%、36.36%和15.91%，這些實驗結果充分說明本文算法的性能明顯高于其它算法。
圖2
構建DREAM10基因調控網絡圖
(a)真實網絡；(b)預測網絡
Figure2. Constructing GRNs on DREAM10 dataset(a) real network; (b) predicted network
表2
不同方法構建DREAM10節點基因調控網絡的性能比較
Table2.
Comparison of DREAM10 inferring GRNs with different methods
其次，我們對50個基因的酵母數據網絡進行仿真分析。如圖 3(a)所示是50個基因節點的真實網絡，包含了77條基因調控關系(即網絡中的邊)。如圖 3(b)所示是本文算法推測出的調控網絡，實線代表真實的調控關系，虛線代表預測錯誤的調控關系，推測出39條正確的邊，網絡的預測相對正確率為50.65%。如表 3所示，本文算法預測出正確的邊數、預測的正確率、相對正確率、錯誤的邊數以及預測錯誤率都明顯好于LP算法、LASSO算法和ARANCE算法，我們算法預測的正確率分別比LP算法、LASSO算法和ARANCE算法高出22.45%、18.09%和9.53%，這些實驗結果說明本文算法能更好地預測網絡。
圖3
構建DREAM50基因調控網絡圖
(a)真實網絡；(b)預測網絡
Figure3. Constructing GRNs on DREAM50 dataset(a) real network; (b) predicted network
表3
不同方法構建DREAM50節點基因調控網絡的性能比較
Table3.
Comparison of DREAM50 inferring GRNs with different methods
2.3 真實數據的網絡重構
同樣我們還在大腸桿菌SOS DNA修復網絡實驗數據集上進行了算法性能驗證，這是一個被實驗驗證了的真實基因調控網絡，它包含9個基因節點和24條真實的相互作用關系。如圖 4(a)所示是真實的SOS DNA修復網絡；如圖 4(b)所示是我們的算法推測出的調控網絡，推測出18條正確的邊和5條錯誤的邊，正確率為78.26%。如表 4所示，本文算法在預測正確的邊數、預測正確率以及相對正確率等方面比LP算法分別高出6、38.26%和25%；此外我們預測出錯誤的邊只有5條，預測錯誤率為21.74%。為了說明本文算法的時間效率，我們在同樣的電腦上(i5-2400 CPU, 4GB RAM)計算出幾種算法的運行時間。我們的算法在運行時間上比LP、LASSO和ARACNE算法分別快6.53、6.49和3.89 s，這進一步說明了本文算法的時效性。
圖4
構建大腸桿菌SOS DNA修復網絡圖
(a)真實網絡；(b)預測網絡
Figure4. Constructing GRNs on SOS dataset(a) real network; (b) predicted network
表4
不同方法構建SOS DNA修復網絡的性能比較
Table4.
Comparison of SOS DNA GRNs with different methods
3 總結
本文提出了一種新型的基因調控網絡構建算法，它以信息論為基礎，通過重采樣和加權平均策略，利用條件互信息來進行基因調控網絡重構。算法首先從基因表達數據中重采樣得到一系列子數據集，然后通過條件互信息方法構建出一系列子數據集的子網絡，最后利用這些子網絡加權平均得到最終的基因調控網絡。本算法可以克服信息論構建GRNs遺漏真陽性作用關系的缺點，從而提高網絡重構的精度。我們的算法在DREAM3模擬數據和SOS DNA修復網絡數據上進行驗證，和其它一些流行的方法相比，我們網絡重構算法取得了較高的精確度和準確度。
0 引言
基因調控網絡(gene regulatory networks, GRNs)描述了細胞內基因之間復雜的相互作用關系，可以識別生命現象并揭示生命活動的本質，它已成為生物信息學領域一個重要的研究話題。現在國內外已經產生了大量的GRNs重構模型和算法[1-4]，各種方法都有其應用的局限性和優缺點[5-7]。一種很流行的方法是基于動態貝葉斯網絡來推測基因之間相互作用的關系[8-10]，這種方法計算性能良好，但所處理的網絡規模因其計算時間復雜度而受到限制。此外，常微分方程是基于迭代的方法來構建GRNs[11-12]，它是一種相對簡單的模型，可以定量地預測GRNs的行為。回歸模型算法也可以用來推測調控基因和靶基因之間的相互作用關系，如LASSO(the least absolute shrinkage and selectionator operator)算法[13-14]、彈性網絡回歸算法[15-16]和貝葉斯平均模型算法[17]等。信息論也被用來構建基因調控網絡[18-20]，它不僅可以計算基因間的非線性關系，而且對基因樣本個數要求較低，可以克服基因表達數據高維低樣本的缺點。但是信息論方法在計算中會遺漏一些調控關系，而這些調控關系可能代表基因間真實的相互作用關系。為了克服這個缺點以及進一步提高網絡重構的精度，本文提出了一種新的基因調控網絡構建的方法。算法首先從基因表達數據重采樣得到一系列子數據，然后通過條件互信息(conditional mutual information, CMI)方法構建一系列子數據的網絡，最后利用這些子網絡加權平均得到最終的基因調控網絡。
本文算法的流程如圖 1所示，它是一種基于多個子網絡的集成算法，它的優點就是能恢復早期計算過程中誤刪除的“真實基因調控關系”，以此來提高基因調控網絡構建的精度。首先對原始基因表達數據進行刀切法(jackknife)重采樣[21]，得到一系列基因表達數據的子數據集；然后在各個子數據集上依次使用條件互信息推斷出對應的基因調控子網絡；最終使用算術平均策略融合各個子網絡中邊的支持度，根據與閾值的關系得出最終基因調控網絡拓撲結構。
圖1
算法流程圖
Figure1.
Flow chart of the method
1 數據集和理論方法
1.1 數據集
為了檢測本文網絡構建算法的精確性，我們對不同的實驗數據集進行網絡構建，然后與實驗數據集的真實網絡進行比對，以此來評價本文算法的優劣。基因表達數據采用微陣列基因芯片技術，它運用計算機芯片把核酸固定在平面載體上形成檢測器件，將待測樣本標記后同芯片進行雜交，芯片上的探針分子和樣本中的標記分子特異性結合后，用激光共聚焦熒光掃描儀截取信息，然后處理分析信息得到信號值。信號值代表了結合在探針上的待測樣本中特定大分子的信息。對基因芯片數據加以整理和分析，并發展對應的數學模型和算法來闡明生物現象，是近幾年系統生物學研究的熱點。基因表達譜如表 1所示，行代表不同的基因名稱，列代表每個基因在不同條件下的表達值。通過基因芯片數據分析基因之間的調控關系，構建復雜系統的基因調控網絡，從整體和系統的角度來研究生命的規律。
表1
基因表達值
Table1.
Gene expression profiles
本文采用不同類型的數據對所提出的算法進行檢測評估，這些數據包括計算機模擬數據和真實的生物分子網絡數據。模擬數據來源于廣泛采用的DREAM(Dialogue for Reverse Engineering Assessments and Methods)競賽數據[6]，該競賽數據中包含基因表達數據和標準網絡，它是經實驗驗證了的酵母(Yeast)和大腸桿菌(Escherichia coli)調控網絡。這里選取了DREAM3的兩個子網絡對應的基因表達數據進行實驗，這兩個網絡分別含有10、50個頂點(基因)，以及10、77條真實作用邊(基因間調控關系)。對于真實基因表達數據，我們采用大腸桿菌SOS DNA修復網絡的真實干擾實驗數據[22-23]，該標準網絡含有9個節點和24條相互作用邊。
1.2 信息論相關定義
基于信息論的互信息已經廣泛應用于基因調控網絡的重構，而且取得了較好的效果。它是衡量兩個基因變量X和Y之間的相關性，現給定一個基因變量X，它的信息熵H(x)描述如下：
| $H\left( X \right) = - \sum\limits_{x \in X} {p\left( x \right)\log p\left( x \right)} $ |
其中p(x)表示基因變量X為x時的概率值，基因變量X和Y的聯合熵H(X, Y)可以被定義為：
| $H\left( {X,Y} \right) = - \sum\limits_{x \in X,y \in Y} {p\left( {x,y} \right)\log p\left( {x,y} \right)} $ |
其中p(x, y)表示基因變量X和Y分別為x和y時的聯合概率值，基因變量X和Y的互信息值I(X, Y)定義如下：
| $I\left( {X,Y} \right) = - \sum\limits_{x \in X,y \in Y} {p\left( {x,y} \right)\log \frac{{p\left( {x,y} \right)}}{{p\left( x \right)p\left( y \right)}}} $ |
由公式(2)和(3)互信息值I(X, Y)被重寫為：
| $I\left( {X,Y} \right) = H\left( X \right) + H\left( Y \right) - H\left( {X,Y} \right)$ |
這里H(X, Y)表示基因變量X和Y的聯合熵，互信息值的大小表示基因相互作用關系的強弱，當互信息值非常小或者接近于零時表示兩個基因變量相互獨立。條件互信息公式定義如下：
| $I\left( {X,Y|Z} \right) = H\left( {X,Z} \right) + H\left( {Y,Z} \right) - H\left( Z \right) - H\left( {X,Y,Z} \right)$ |
其中H(X, Y)表示基因變量X，Y和Z的聯合熵，這里我們的熵用高斯核概率密度函數來估計[24]，它的計算公式下：
| $\begin{gathered} P\left( {{X_i}} \right) = \frac{1}{N}\sum\limits_{j = 1}^N {\frac{1}{{{{\left( {2\pi } \right)}^{n/2}}{{\left| C \right|}^{n/{\text{2}}}}}}} \hfill \\ \exp \left( { - \frac{1}{2}{{\left( {{X_j} - {X_i}} \right)}^T}{C^{ - 1}}\left( {{X_j} - {X_i}} \right)} \right) \hfill \\ \end{gathered} $ |
其中C表示基因變量X的協方差矩陣，｜C｜表示矩陣C的行列式，N是每個基因的樣本個數，n是基因變量的個數。由公式(1)和(6)得出基因變量X的信息熵為：
| $H\left( X \right) = \log \left[ {{{\left( {2\pi e} \right)}^{n/2}}{{\left| C \right|}^{1/2}}} \right] = \frac{1}{2}\log {\left( {2\pi e} \right)^n}\left| C \right|$ |
由公式(7)，互信息公式(4)被重寫為：
| $I\left( {X,Y} \right) = \frac{1}{2}\log \frac{{\left| {C\left( X \right)} \right| \cdot \left| {C\left( Y \right)} \right|}}{{\left| {C\left( {X,Y} \right)} \right|}}$ |
同理，條件互信息公式(5)被重寫為：
| $I\left( {X,Y,Z} \right) = \frac{1}{2}\log \frac{{\left| {C\left( {X,Z} \right)} \right| \cdot \left| {C\left( {Y,Z} \right)} \right|}}{{\left| {C\left( Z \right)} \right| \cdot \left| {C\left( {X,Y,Z} \right)} \right|}}$ |
1.3 重采樣構建網絡
矩陣G=[gim]N×M, (i=1, 2, …, N; m=1, 2, …, M)表示基因表達數據，其中i代表基因，m代表樣本或者采樣條件，gim表示第i個基因在第m個采樣條件下的基因表達值。Xi=[gi1, gi2, …, giM]則表示矩陣G中的第i行，由第i個基因在M個不同采樣條件下所得到的基因表達值組成。Cm=[g1m, g2m, …, gNm]T則表示矩陣G中的第m列，由在N個基因在第m個采樣條件下所得到的基因表達值組成。重采樣過程就是依次刪除G矩陣中的一列Cm，則剩余的M-1列則組成了一個新的矩陣G′=[gil]N×(M-1)。經過重采樣過程之后，我們可以獲得原始基因表達值的M個數據子集。
分別在M個基因表達子數據集上用互信息理論得到對應的M個子基因調控網絡，使用Lijm表示第m個子網絡中節點i和j之間邊的條件互信息值，使用算術平均融合策略計算一致網絡中節點i和j之間邊的條件互信息值Lij：
| ${L_{ij}} = \frac{1}{M}\sum\limits_{m = 1}^M {L_{ij}^m} $ |
根據一致網絡中邊的Lij值與所設定的閾值α的大小關系，可以確定在最終網絡中節點i和j之間是否存在連接關系。M個子網絡的相應邊經過算術平均融合策略處理之后，如果計算得到的互信息值Lij大于所設定的閾值，則認為在最終的一致網絡中節點i和j之間存在著連接，反之則沒有連接。經過上述過程就可以重構出最終的基因調控網絡。
2 實驗仿真與結果討論
2.1 評價標準
網絡預測的結果用如下指標來衡量：真實網絡的邊數(true edges, TE)，網絡預測出的邊數(predicted edges, PE)，預測出正確的邊數(predicted correct edges, PCE)，預測出錯誤的邊數(predicted wrong edges, PWE)。同時，按如下公式計算：正確率=PCE/PE，錯誤率=PWE/PE，相對正確率=PCE/TE，以此來衡量構建算法的優劣性。
為了驗證本文算法的有效性，我們和已有的具有較好性能的算法進行比較。在本文中我們選取了3種常用的GRNs構建算法：線性規劃算法(linear programming, LP)[25]、回歸模型方法(LASSO)[26]和基于信息論的方法(algorithm for the reconstruction of accurate cellular networks, ARACNE)[27]。LP算法是一種線性規劃模型，結構較為簡單，計算速度快。LASSO是基于回歸的一種GRNs構建算法，在回歸系數的絕對值之和小于一個常數的約束條件下，使殘差平方和最小化。這兩種模型可以優化算法的獎懲項來控制網絡的稀疏度，并且取得較好的效果，但是缺點是精度過分依賴模型的參數。ARACNE算法是一種基于互信息的基因調控網絡構建算法，能夠處理高維低樣本數據，具有較高的性能，但不能區分直接和間接調控關系而導致較高的假陽性。
2.2 模擬數據的網絡重構
首先，我們采用DREAM3競賽數據中的酵母基因表達數據來構建基因調控網絡。我們對10個基因的酵母數據網絡進行分析。如圖 2(a)所示是10個基因節點的真實網絡，包含10條邊。如圖 2(b)所示是本文算法推測出的基因調控網絡，推測出10條正確的邊，推測出1條錯誤的邊(G2-G4)，預測邊的相對正確率為100%，網絡預測的正確率高達90.91%。如表 2所示詳細列出了本文算法和其它算法在性能上的比較結果，從表中可以看出本文算法在預測出正確的邊數、預測的正確率、相對正確率、錯誤的邊數以及預測錯誤率方面都明顯好于LP算法、LASSO算法和ARANCE算法，我們算法預測的正確率分別比LP算法、LASSO算法和ARANCE算法高出49.24%、36.36%和15.91%，這些實驗結果充分說明本文算法的性能明顯高于其它算法。
圖2
構建DREAM10基因調控網絡圖
(a)真實網絡；(b)預測網絡
Figure2. Constructing GRNs on DREAM10 dataset(a) real network; (b) predicted network
表2
不同方法構建DREAM10節點基因調控網絡的性能比較
Table2.
Comparison of DREAM10 inferring GRNs with different methods
其次，我們對50個基因的酵母數據網絡進行仿真分析。如圖 3(a)所示是50個基因節點的真實網絡，包含了77條基因調控關系(即網絡中的邊)。如圖 3(b)所示是本文算法推測出的調控網絡，實線代表真實的調控關系，虛線代表預測錯誤的調控關系，推測出39條正確的邊，網絡的預測相對正確率為50.65%。如表 3所示，本文算法預測出正確的邊數、預測的正確率、相對正確率、錯誤的邊數以及預測錯誤率都明顯好于LP算法、LASSO算法和ARANCE算法，我們算法預測的正確率分別比LP算法、LASSO算法和ARANCE算法高出22.45%、18.09%和9.53%，這些實驗結果說明本文算法能更好地預測網絡。
圖3
構建DREAM50基因調控網絡圖
(a)真實網絡；(b)預測網絡
Figure3. Constructing GRNs on DREAM50 dataset(a) real network; (b) predicted network
表3
不同方法構建DREAM50節點基因調控網絡的性能比較
Table3.
Comparison of DREAM50 inferring GRNs with different methods
2.3 真實數據的網絡重構
同樣我們還在大腸桿菌SOS DNA修復網絡實驗數據集上進行了算法性能驗證，這是一個被實驗驗證了的真實基因調控網絡，它包含9個基因節點和24條真實的相互作用關系。如圖 4(a)所示是真實的SOS DNA修復網絡；如圖 4(b)所示是我們的算法推測出的調控網絡，推測出18條正確的邊和5條錯誤的邊，正確率為78.26%。如表 4所示，本文算法在預測正確的邊數、預測正確率以及相對正確率等方面比LP算法分別高出6、38.26%和25%；此外我們預測出錯誤的邊只有5條，預測錯誤率為21.74%。為了說明本文算法的時間效率，我們在同樣的電腦上(i5-2400 CPU, 4GB RAM)計算出幾種算法的運行時間。我們的算法在運行時間上比LP、LASSO和ARACNE算法分別快6.53、6.49和3.89 s，這進一步說明了本文算法的時效性。
圖4
構建大腸桿菌SOS DNA修復網絡圖
(a)真實網絡；(b)預測網絡
Figure4. Constructing GRNs on SOS dataset(a) real network; (b) predicted network
表4
不同方法構建SOS DNA修復網絡的性能比較
Table4.
Comparison of SOS DNA GRNs with different methods
3 總結
本文提出了一種新型的基因調控網絡構建算法，它以信息論為基礎，通過重采樣和加權平均策略，利用條件互信息來進行基因調控網絡重構。算法首先從基因表達數據中重采樣得到一系列子數據集，然后通過條件互信息方法構建出一系列子數據集的子網絡，最后利用這些子網絡加權平均得到最終的基因調控網絡。本算法可以克服信息論構建GRNs遺漏真陽性作用關系的缺點，從而提高網絡重構的精度。我們的算法在DREAM3模擬數據和SOS DNA修復網絡數據上進行驗證，和其它一些流行的方法相比，我們網絡重構算法取得了較高的精確度和準確度。
