基因共表达网络的构建及分析方法研究综述

汪涛等

摘 要：随着高通量生物实验技术的快速发展，特别是基因芯片和新一代测序技术的发展，全基因组范围内的基因表达数据呈爆炸式增长。利用网络生物学的方法对高通量基因表达数据进行分析和挖掘已经成为生物信息学重要的研究方向。对基因共表达网络的研究与分析从系统层面上加深了研究人员对生物系统的认识。本文综述了基因共表达网络的构建和分析的常用方法，主要包括基因相似性度量方法、阈值选择方法、拓扑分析方法、基因模块识别及其功能注释注释方法，并对一些常用的分析工具进行了分析总结。

关键词：基因共表达网络；基因模块；功能和拓扑分析

中图分类号：TP391 文献标识号：A 文章编号：2095-2163（2014）06-

Abstract： With the rapid development of high-throughput biological experiments technology， particularly the development of gene microarray and next generation sequencing technology， the genome wide gene expression data grow explosively. Network biology has become an important research direction in bioinformatics， which enhances the analysis of high throughput gene expression data. The study and analysis of gene co-expression network help the researchers to understand the biological system in system level better. This article reviewes the common methods on constructing and analyzing gene co-expression network， including the methods of measuring similarity of genes， selecting proper threshold， analyzing topological structure， detecting gene module and functional annotation on gene module. In the end， the paper summarizes some commonly used tools for analyzing gene co-expression network.

Keywords： Gene Co-expression Network； Gene Module； Functional and Topological Analysis

0 引 言

后基因组时代生命科学研究的一个主要目的是理清生物细胞内所有分子以及分子之间的联系，并且揭示分子之间相互作用以及如何决定细胞生命功能的内在机理[1]。随着系统生物学和复杂图理论的发展，分子生物网络的研究为探索复杂生命活动提供了有力工具。分子生物网络在系统层面反映了生物分子的相互作用关系，因而在相当过程中有助于研究者深入理解生物细胞中各种生物分子是如何相互作用、进而行使生物功能的完整处理实现过程。目前人们已经对各种类型的分子生物网络进行了广泛的研究，如基因共表达网络（gene co-expression network）、基因调控网络（gene regulatory network）、蛋白质相互作用网络（protein-protein interaction network）、代谢网络（metabolic network）等。这些探索已经从对单个生物分子研究层面上升到解析研究分子相互作用关系的系统研究层面，并且产生了丰硕的研究成果。

基因芯片技术以及新一代测序技术的应用，使得全基因组范围内基因表达数据得以快速累积。仅仅对单个基因功能水平的研究已经限制了人们探索生物细胞行使生命功能的视野和进程。利用系统生物学的方法构建基因共表达网络，从而由系统层面揭示基因之间的相互关系已经成为一个主要的研究方向。本文中，主要对基因共表达网络的构建方法以及常用的分析方法进行了关注及论述。

1基因共表达网络的构建

基因共表达网络大多是以基因间表达谱数据的相关性为基础而实现构建的。在基因共表达网络的表示中，经常使用图模型来描述基因之间的关系。图中的节点代表基因，边表示两个基因之间的共表达相互作用关系。基因共表达网络的构建主要分为两个步骤，第一是对所有基因进行相似性度量；第二是通过阈值的选择确定共表达网络的边。以下则对其展开具体分析。

1.1数据来源及表示

在分子生物学领域，基因表达谱是指基因表达活性的有效度量。从基因表达谱的数据来源上进行分类，常用于基因共表达网络构建的表达谱数据主要将划定为两类：一类是基因芯片（microarray）数据，另一类是RNA-seq数据。基因的表达谱数据可以用一个n*m的矩阵X=[xij]来表示，数学表述如公式（1）所示。其中，第i行数据xi （i = 1，…，n）对应一个基因的表达谱，矩阵中的列则能够反映在不同样本或者时间点下该基因的表达水平。

因相似性度量方法

在基因共表达网络的构建和分析中，经常需要对两个基因进行表达相似性度量。基因间的相似性有多种度量方式，若从计算方法上进行分类，主要可以分为基于表达谱的相似性度量和基于拓扑结构的相似性度量。具体地，基于基因表达谱的相似性度量指标主要通过计算不同基因表达谱的线性或非线性相关系数而最终得到。其中，常见的线性相关性指标主要有皮尔森相关系数（Pearson Correlation Coefficient， Pcc）、斯皮尔曼相关系数（Spearman Correlation Coefficient），以及排除其他变量影响的偏相关系数（Partial Correlation Coefficient）等。而非线性相关性指标则可列举实例，诸如互信息（Mutual Information， MI）等。在基因模块识别的过程中，有些研究者又利用网络拓扑结构对基因相似性实现了进一步的度量。

1.3阈值选取方法

设定合理的阈值，将具有潜在生物功能的边保留下来，是基因共表达网络构建的重点及关键。以计算方法为依据，大致可以分为以下几类：基于人工设定的阈值选取方法、基于统计显著水平的阈值选取方法、基于对照试验的阈值选取方法、基于网络拓扑结构的阈值选取方法以及基于多种方法综合的阈值选取方法等。在此，即对各类方法做以详细的探讨与分析。

首先，基于人工设定的阈值选取方法是最简单、也是应用最广的方法。其中有两类模型使用最为广泛：

（1）选取固定的阈值t，将相似性大于阈值t的基因对保留下来，例如Tom C Freeman[4]等人的文章中就采取了这种策略。

（2）将所有基因对按照相似性系数进行排序，选择前百分之x（例如前1%）的基因对进行保留，Ala等人[5]选取前0.5%作为阈值，构建保守的共表达网络。但这种方法最大的缺点在于其具有的武断性，而并未考虑到生物网络本身的特性，也没有验证其生物统计显著水平。对应地，基于统计显著水平的阈值选取方法往往通过与随机网络进行比对，得到不同阈值下的显著性水平，由此而进行阈值选择。并且，显著性水平往往多采用p值（p-value）、q值（q-value）等来衡量。然而这种方法却无法保留那些相似性系数较低、但却具有生物意义的边。

其次，基于对照试验的阈值选取方法解决了计算方法带来的随机性。这种方法在基因芯片的实验阶段，通过在基因芯片中引入与目标物种亲缘关系较远的物种基因而形成对照组。理论上，对照组的基因由于未与任何RNA杂交而不曾产生任何信号，但是现实中随机杂交的情况却使得对照组将产生一定的噪音信号。而且，这种噪声信号可以作为参考，用于辅助目标网络的阈值选择。Voy等人[6]证明了利用这种方法保留下的边具有显著的生物统计意义。但是，该种方法却增加了实验成本，并且结果将十分依赖所选择的对照样本。

再有，基于网络拓扑结构的阈值选取方法考虑到了基因共表达网络所具有的特性，如无标度（scale-free）以及小世界（small world），并通过优化阈值的方法，使得网络的某些拓扑特性得以明显体现。Horvath等人[7]即提出了一种“软”阈值选取方法，利用网络达到无标度的拓扑结构而以此来确定阈值。Elo等人在文献[8]中则提出了一种基于聚集系数（Clustering Coefficient）的阈值选择策略。由于基因共表达网络具有小世界的特性，因此往往比随机网络具有更高的聚集系数。作者利用这种特点，将阈值选择的问题转化为具体的关于聚集系数C*的优化问题。Borate等人又在文章[9]中将基于最大团和基于谱图理论的阈值选择方法等进行了对比研究。当阈值从高到低下降时，网络中的最大团的数目是指数增长的。根据这一原理，基于最大团的方法在最大团的数目随阈值变化的趋势线中将选择一个关键拐点作为阈值。该选取拐点往往是最大团的数目突然增加2倍或者3倍时的阈值点。而与此类似，利用谱图理论的阈值选择方是基于网络的拉普拉斯矩阵的特征值和特征向量，来发现网络的基因模块，并会选择一个模块划分最优时的阈值。

最后，基于多种方法综合的阈值选取方法通常是融合了上述方法的优点，从而用多种指标来选取阈值。例如，Langston等人[10]利用了本体距离、统计显著性以及多种图的拓扑特性来进行阈值选择。

2基因共表达网络的拓扑分析方法

生物网络的快速发展表明分子网络遵从着某些普适性的规则[1]。这些规则通常会在网络的拓扑特性中得到相应体现。因此，对有关基因共表达网络的拓扑展开分析即是从系统层面了解基因共表达关系的重要工具。根据度量拓扑结构中的基因个数的不同，可以将拓扑分析方法大致分为两类，也就是：全局网络拓扑分析和网络中心性分析。

在研究进展中，为了从全局角度刻画网络的拓扑结构，一些学者即利用图论的的方法定义了许多衡量网络全局性拓扑结构的度量指标，常见的主要包括平均度（average degree）和度分布（degree distribution）、聚集系数（clustering coefficient）、平均路长（average path length）、直径（diameter）等。这些全局性的度量指标能够很好地反应出基因共表达网络的三种全局拓扑特性，分别是：无标度分布（scale-free distribution）、小世界特性（small world property）、功能模块网络（functional modular network）。另外，基因共表达网络中不同节点在网络中的重要性是不同的，通常用中心性（Centrality）来度量网络节点的地位。中心性往往体现在网络的拓扑结构上，度很高的节点或者是起着关键连接作用的节点都可能在某些生物途径中发挥着重要的作用，因此人们往往根据所研究问题的不同来定义节点的中心性。具体地，较为常用的中心性度量指标主要有度中心性（degree centrality）、亲近中心性（closeness centrality）、居间中心性（betweenness centrality）、特征向量中心性（eigenvector centrality）等。

CentiBiN[11]就是一款专门用于生物网络中心性的计算和可视化软件，其中集成了17种无向图的中心性度量方法（针对有向图则汇总了15种中心性度量方法），并且可以计算网络直径、平均路长等全局性网络拓扑指标。另外，诸如Cytoscape[12]、Pajek[13]、Visone[14]、VisANT[15]等复杂网络分析和可视化软件工具，虽然不是专门为了计算网络拓扑特性而产生的，但对于常用的网络拓扑特性分析却都能提供良好的支持。

3基因共表达网络的模块分析

3.1 基因模块识别

基因功能模块识别是基因共表达网络分析中的最重要方法之一。研究证明，在基因共表达网络中致密的连通子图往往具有特定的生物学功能。基因共表达网络模块识别算法大多来源于对蛋白质相互作用网络的研究中。Junzhong Ji等人[16]对蛋白质网络中的模块识别方法即做了详细的综述。在本文中，仅对目前常用于基因模块识别的四种典型算法进行了深入分析，如基于层次的模块识别方法（WGCNA[7]）、基于密度的模块识别方法（MCODE[17]）、基于流模拟的模块识别方法（MCL[18]）以及基于划分的模块识别方法（Qcut[19]）。下面依然对其展开综合性论述。

首先，WGCNA[7]（Weighted Gene Co-expression Network Analysis）是以基于相关系数构建基因共表达网络的代表性方法，其中集成了多种网络分析方法，主要具有网络构建、功能模块探测、基因选择、拓扑特性计算、数据模拟、可视化以及与其他软件交互等功能。WGCNA使用了基于拓扑结构交叠的度量指标，利用无先验知识的聚类策略来探测基因的功能模块，实际默认使用的则是层次聚类方法——hclust[20]。层次聚类得到的系统树图的分支对应着可能的功能模块，并且可以使用多种剪枝的方法如固定高度的剪枝方法或者两种动态剪枝方法[21]来确定功能模块。作为多种生物网络构建的代表性方法，WGCNA现已成功地应用到多种研究的开展和实现当中。

其次，MCODE[17]（Molecular Complex Detection）算法是一种基于图论（或网络密度）的网络模块发现算法。算法共分为三个步骤：网络节点加权、模块预测以及模块优化处理。分别来说，网络节点加权是指根据节点所在的、最高的k-core的密度来为网络中的所有节点赋予一个权值。在第二步模块的预测中，首先选取一个具有最高权值的节点作为种子节点，而后依次向外扩增，再将那些权值在阈值VWP之上的周边节点依次纳入模块中，直到没有节点再可包含进入该模块为止。上述操作后，就将选择剩余节点中权值最高的作为种子节点进行同样的操作。而在第三步的模块优化处理中，算法会将那些不包含2-core的模块删除。并且用户可以指定是否进行“fluff”和“haircut”操作来处理模块边缘节点。AllegroMCODE[22]是一款基于MCODE算法的Cytoscape[12]插件，可通过GPU进行加速，多将用于高效地挖掘基因模块。

3.2 模块的功能分析方法及工具

为了发现基因模块的功能，往往需要对基因模块进行模块富集分析（Modular enrichment analysis， MEA）[23]。模块富集分析多数时候是通过计算基因本体注释数据库[24]（Gene Ontology， GO ）或通路注释数据库（如KEGG[25]）中的术语（term）在每个模块上的富集程度，再根据假设检验中的p值（p-value）来衡量富集到模块上的功能的显著性。计算p值的经典统计检验方法主要包括卡方检验（Chi-square test）、Fisher精确检验（Fisher's exact test）、累计超几何检验（cumulative hypergeometric test）等。模块富集分析不仅继承了单一富集分析（singular enrichment analysis， SEA）的一些特点（如在预选感兴趣的基因集后，再计算各基因与注释数据库中的术语的富集得分），而且模块富集分析在计算p值时还考虑到了模块基因-基因之间的关系以及对应的术语-术语关系。这种术语-术语间关系的着重关注将可提高功能富集的敏感性和准确度。常见的能够基于GO进行模块功能富集工具可概略给出为：GO：：TermFinder[26]、Ontologizer[27]、 topGO[28]、 GENECODIS[29]、 ADGO[30]、GoToolBox[31]、DAVID[32]等。

4 结束语

基因共表达网络的分析方法已成为高通量生物数据分析强有力的工具，并已经广泛应用于生物学研究中。但是目前的构建和分析方法还存在着一些不足，例如基因共表达网络的构建方法繁多，时下还缺乏可靠的评价系统；模块功能富集方法对连接度较低的节点的识别能力还有待进一步提高等。因此总地来说，对于基因共表达网络的构建和分析也依然是目前一项颇具挑战性的研究工作。

参考文献：

[1] BARAbASI A L，OLTVAI Z N. Network biology： understanding the cell's functional organization[J]. Nature Reviews Genetics， 2004，5（2）： 101-113.

[2] RAVASZ E， et al. Hierarchical organization of modularity in metabolic networks[J]. science， 2002， 297（5586）：1551-1555.

[3] YIP A M， HORVATH S. Gene network interconnectedness and the generalized topological overlap measure[J]. BMC bioinformatics， 2007，8（1）： 22.

[4] FREEMAN T C， et al. Construction， visualisation， and clustering of transcription networks from microarray expression data[J]. PLoS computational biology， 2007， 3（10）： e206.

[5] ALA U， et al. Prediction of human disease genes by human-mouse conserved coexpression analysis[J]. PLoS computational biology， 2008， 4（3）： e1000043.

[6] VOY B H， et al. Extracting gene networks for low-dose radiation using graph theoretical algorithms[J]. PLoS computational biology， 2006， 2（7）： e89.

[7] LANGFELDER P， HORVATH S. WGCNA： an R package for weighted correlation network analysis[J]. BMC bioinformatics， 2008， 9（1）： 559.

[8] ELO L L， et al. Systematic construction of gene coexpression networks with applications to human T helper cell differentiation process[J]. Bioinformatics， 2007， 23（16）： 2096-2103.

[9] BORATE B R. Comparative Analysis of Thresholding Algorithms for Microarray-derived Gene Correlation Matrices， 2008.

[10] LANGSTON M A， et al. Innovative computational methods for transcriptomic data analysis： A case study in the use of FPT for practical algorithm design and implementation[J]. The Computer Journal， 2008， 51（1）： 26-38.

[11] KOSCHUTZKI D. CentiBiN Version 1.4. 2， in， 2006，CentiBiN Version， 2006，1（2）： 2004-2006.

[12] SHANNON P， et al. Cytoscape： a software environment for integrated models of biomolecular interaction networks[J]. Genome research， 2003， 13（11）： 2498-2504.

[13] BATAGELJ V， MRVAR A. Pajek-program for large network analysis[J]. Connections， 1998， 21（2）： 47-57.

[14] BAUR M， et al. Visone Software for visual social network analysis. in Graph Drawing. Springer，2002.

[15] HU Z， et al. VisANT： data-integrating visual framework for biological networks and modules[J]. Nucleic acids research， 2005， 33（suppl 2）： W352-W357.

[16] JI J， et al. Survey： Functional module detection from protein-protein interaction networks. Knowledge and Data Engineering[J]. IEEE Transactions on， 2014， 26（2）： 261-277.

[17] BADER G D， HOGUE C W. An automated method for finding molecular complexes in large protein interaction networks[J]. BMC Bioinformatics， 2003， 4： 2.

[18] HWANG W， et al. A novel functional module detection algorithm for protein-protein interaction networks[J]. Algorithms for Molecular Biology， 2006， 1（24）.

[19] RUAN J， ZHANG W. Identifying network communities with a high resolution[J]. Physical Review E， 2008， 77（1）： 016104.

[20] Kaufman L， Rousseeuw P J. Finding groups in data： an introduction to cluster analysis. John Wiley & Sons，2009，34.

[21] LANGFELDER P， ZHANG B， HORVATH S. Defining clusters from a hierarchical cluster tree： the Dynamic Tree Cut package for R[J]. Bioinformatics， 2008， 24（5）： 719-720.

[22] YOON J， JUNG W H. A GPU-accelerated bioinformatics application for large-scale protein interaction networks. APBC poster presentation， 2011.

[23] HUANG D W， SHERMAN B T， LEMPICKI R A. Bioinformatics enrichment tools： paths toward the comprehensive functional analysis of large gene lists[J]. Nucleic acids research， 2009， 37（1）： 1-13.

[24] GO T. Gene Ontology： tool for the unification of biology[J]. AMERICA N， editor. Nature Genetic， 2000， 25： 25-29.

[25] KANEHISA M， GOTO S. KEGG： kyoto encyclopedia of genes and genomes. Nucleic acids research， 2000， 28（1）： 27-30.

[26] BOYLE E I， et al. GO：： TermFinder—open source software for accessing Gene Ontology information and finding significantly enriched Gene Ontology terms associated with a list of genes[J]. Bioinformatics， 2004， 20（18）： 3710-3715.

[27] BAUER S， et al. Ontologizer 2.0—a multifunctional tool for GO term enrichment analysis and data exploration[J]. Bioinformatics， 2008， 24（14）： 1650-1651.

[28] Alexa A， Rahnenfuhrer J. topGO： enrichment analysis for gene ontology. R package version 2.8， 2010.

[29] CARMONA-SAEZ P， et al.， GENECODIS： a web-based tool for finding significant concurrent annotations in gene lists[J]. Genome biology， 2007， 8（1）： R3.

[30] NAM D， et al. ADGO： analysis of differentially expressed gene sets using composite GO annotation[J]. Bioinformatics， 2006， 22（18）： 2249-2253.

[31] MARTIN D， et al. GOToolBox： functional analysis of gene datasets based on Gene Ontology[J]. Genome biology， 2004， 5（12）： R101.

[32] ALVORD G， et al. The DAVID Gene Functional Classification Tool： a novel biological module-centric algorithm to functionally analyze large gene lists[J]. Genome Biol， 2007，8： R183.