Gene co-expression analysis for functional classification and gene-disease predictions
这篇文章是2018年在brief bioinform上发表的。最近研究基因共表达网络，看到这篇文章，有些内容上面还是不太理解，全文翻译了一下。和大家分享。
原文网址：Gene co-expression analysis for functional classification and gene-disease predictions - PubMed (nih.gov)
Abstract
基因共表达网络可用于将未知功能的基因与生物过程联系起来，对候选疾病基因进行优先排序，或辨别转录调节程序。随着最近转录组学和下一代测序技术的发展，从RNA测序数据中构建的共表达网络也能推断出非编码基因和剪接变体的功能和疾病关联。尽管基因共表达网络通常不提供关于因果关系的信息，但新兴的差异性共表达分析方法能够识别各种表型的调节基因。在此，我们介绍并指导研究人员进行（差异）共表达分析。我们概述了用于创建和分析由基因表达数据构建的共表达网络的方法和工具，并解释了这些方法和工具如何被用来识别在疾病中具有调节作用的基因。此外，我们讨论了其他数据类型与共表达网络的整合，并提供了共表达分析的未来前景。
Introduction
生物研究的一个关键目标是系统地确定活细胞内的所有分子以及它们如何相互作用。然而，许多基因的功能仍未被理解，这种情况随着最近许多新的非编码基因的识别而变得更加复杂[1]。随着高通量技术的发展，包括微阵列和RNA测序（RNA-seq），以及它们各自的数据分析方法，现在可以从系统的角度来确定一个基因的功能状态[2, 3] 。从全基因组基因表达推断基因功能和基因-疾病关联的一种方法是共表达网络分析（图1），这种方法构建了在一组样本中具有共同激活倾向的基因网络，随后对这个网络进行询问和分析。
图1
基因共表达网络可用于各种目的，包括候选疾病基因的优先排序（gene prioritization）、功能基因注释(functional gene annotation)（图1）和识别调节基因( identification of regulatory genes)。然而，共表达网络实际上只能识别相关性；它们表明哪些基因同时活跃，这往往表明它们在相同的生物过程中活跃，但通常不赋予因果关系的信息或区分调节基因和被调节基因。一种越来越多的超越传统共表达网络的方法是差异性共表达分析(differential co-expression analysis)[4-7]。这种方法可以识别在不同条件下具有不同共表达伙伴的基因，如疾病状态[4, 8-10]、组织类型[11]和发育阶段[12]，因为这些基因更可能是支撑表型差异的调节因子。通过整合数据类型，如蛋白质-蛋白质相互作用、甲基组数据、转录因子（TFs）与其靶点之间的相互作用，以及与共表达基因的序列主题分析，可以进一步研究此类基因的调节作用[13-15]。这有助于识别影响共表达模块的表达和组成的调控元素，如TFs、表达定量性状位点（eQTLs）和甲基化模式。
基因表达和调控可能具有高度的组织特异性(highly tissue-specific)，大多数与疾病相关的基因具有组织特异性表达异常[16, 17]。多个组织的表达数据越来越多，使得差异性共表达分析成为可能，它可以识别组织特异性的特征和共享的共表达特征[11]。这些组织特异性特征在组织特异性疾病中可能会被破坏，在汇总多个组织的分析中不会被发现。即使没有样本分类，也可以解决亚群特异性模块，这种方法在对不同癌症亚型进行分类以提供预后标记方面特别成功[18-20]。差异共表达分析也适用于分析亚群未知的数据集，如大规模单细胞RNA-seq数据[5, 12]。虽然差异共表达方法对噪声很敏感[21]，但随着RNA-seq数据数量和质量的增加，它们变得越来越有效。RNA-seq进一步允许共表达分析关注剪接变体和非编码RNA。
在这篇综述中，我们提供了关于什么是共表达网络的介绍和概述，然后是利用RNA-seq数据进行共表达分析的不同步骤指南。然后，我们描述了常用的和新出现的共表达分析的方法和工具，重点是差异性共表达分析，以确定疾病的调节基因。最后，我们讨论了共表达网络与其他类型数据的整合，例如推断调控过程，以及该领域的未来前景和剩余挑战。
提出问题：
- 什么是共表达网络？（co-expression network）
- 如何利用RNA-seq 数据进行共表达分析
- 常用的共表达分析工具和方法有哪些
- 共表达网络分析和其他类型的数据还有其他结合与应用吗？
Co-expression networks
共表达网络确定了哪些基因在一组样本中具有协调表达模式（coordinated expression pattern）的趋势。这种共表达网络可以表示为一个基因-基因相似性矩阵（gene–gene similarity matrix），可以用于下游分析（图1）。典型的共表达网络的构建和分析可以用以下三个步骤来描述。
第一步，根据每对基因之间的相关性度量（correlation measures）或相互信息（mutual information ）[22-24]，定义基因之间的个体关系。这些关系描述了所有样本中基因对的表达模式的相似性。不同的相关措施被用来构建网络，包括Pearson或Spearman的相关关系[25, 26]。另外，最小绝对误差回归[27]或贝叶斯方法[28]也可用于构建共表达网络。后两种方法还有一个好处，就是可以用来识别因果关系，在其他地方也有解释[29]。关于其他类型的相似性度量（ similarity measures）的讨论，我们参考了[30]。许多这些相似性指标也可以用来构建蛋白质-蛋白质相互作用网络（protein–protein interaction networks），在[31]中使用癌症数据进行了比较。
在第二步，共表达关联被用来构建一个网络，每个节点（node）代表一个基因，每条边（edge）代表共表达关系的存在（presence）和强度（strength）（图1）[32]。
在第三步，使用几种可用的聚类技术（available clustering techniques）之一确定模块（共同表达基因群）(groups of co-expressed genes)。共表达分析中的聚类是用来对多个样本中表达模式相似的基因进行分组，以产生共表达基因群（ groups of co-expressed genes ），而不仅仅是一对（pairs）。聚类方法的选择需要考虑，因为它可以大大影响分析的结果和意义。许多聚类方法可用，包括k-means聚类和层次聚类（hierarchal clustering），并在[33]中详细讨论。模块随后可以通过功能富集分析来解释，这是一种在基因列表中识别和排列代表性过强的功能类别的方法[34-36]。
在共表达分析中，考虑样本的异质性（heterogeneity）很重要。在由多个组织或条件构建的共表达网络中，组织特异性或条件特异性的共表达模块可能无法检测到，因为组织/条件特异性模块的相关信号被其他组织/条件中的相关性不足所稀释。然而，将共表达分析限制在特定的组织或条件下也减少了样本量，从而也降低了检测共享共表达模块的统计能力。因此，不区分组织或条件的方法应该用于识别共同的共表达模块，而比较不同条件或组织的差异性共表达将更有利于识别特定条件或组织的独特模块。
summary
典型的共表达网络的构建和分析：
- 第一步是寻找信息，根据不同基因的共表达的相似性寻求基因间的关系，比如说这个基因的表达模式是什么样的。
- 第二步将共表达的这部分基因构建网络，节点为基因，边为联系的强度和有无
- 第三步将构建的这些基因进行聚类，对表达模式相似的这部分基因进行分组。
Types of co-expression networks 共表达网络的类型
Signed and unsigned co-expression networks 有方向和无方向性的共表达网络
在一个基于相关性的共表达网络中，相关度的数值在-1（完全负相关）和1（完全正相关）之间。在无符号网络（unsigned network）中，使用绝对相关值，这意味着两个负相关的基因将被认为是共表达的。这导致负相关的基因组合在一起。因为这些基因很可能也与一组完全不同的基因正向共表达，这些基因也会被分组到同一个模块中，并破坏网络的结构。有符号的网络（signed network）通过在0和1之间缩放相关值来解决这个问题，因此，值<0.5表示负相关，值>0.5表示正相关。有符号的方法创建的网络中，有生物意义的模块（如代表特定生物过程的模块）被更好地分开[37]。因此，接近0的换算值（scaled value）表示负相关，当microRNAs （miRNA）被纳入网络时，这一特征可能特别有趣，因为已知这些基因主要通过下调其他基因来发挥其功能[38]。这对一些长基因间非编码RNAs（lincRNAs）也是如此[39]。
Weighted and un-weighted co-expression networks
在加权网络中，所有的基因都是相互连接的，这些连接具有0和1之间的连续权重值，表示基因之间共同调节的强度。在非加权网络中，基因对之间的相互作用是二进制的，即要么是0，要么是1，基因要么连接，要么不连接。一个非加权网络可以从一个加权网络中创建，例如，将所有相关度超过一定阈值的基因视为连接的，而将所有其他的基因视为不连接的。在本综述中，我们重点讨论加权网络，因为（到目前为止）它们比非加权网络产生了更多的稳健结果[40]。
Microarrays versus RNA-seq data微阵列与RNA-seq数据的比较
协同表达网络可以通过从微阵列或RNA-seq技术获得的基因表达数据来构建。RNA-seq的主要好处之一是它可以量化通常不用芯片测量的7万多个非编码RNA的表达[1]，包括最近被注释的lincRNA，其中许多被认为具有调节作用[41]，并在疾病中发挥作用[42，43]。因此，为了更好地了解驱动生物过程的调控机制，在分析中需要考虑非编码RNAs。
RNA-seq也有其他好处[35]。它能提高低丰度转录物的准确性[44]，对识别组织特异性表达具有更高的分辨率，并比芯片衍生的表达谱更好地区分密切相关的同源物的表达谱[45]。RNA-seq还可以区分不同剪接变体的表达[46, 47]，它们可以有不同的相互作用伙伴[48]和生物功能[49]。RNA-seq数据的共表达分析可以为这些剪接变体和lincRNAs分配假定的角色[2]，并确定它们可能发挥的作用的疾病[2]。在剪接变体水平上的联合表达分析的局限性是引入偏见，因为如果多个剪接变体共享同一个表达的外显子，就很难确定哪个剪接变体被表达。
作为RNA-seq对同工酶和外显子特异性表达水平测量的一个例子，外显子水平的表达被用来构建一个共表达网络[50, 51]。在基因共表达网络中，源于同一基因的不同转录本的表达通常是聚集的，这可能导致有偏见的共表达信号[50]。在共剪切网络中，在计算基因共表达相关性时，通过考虑基因内的外显子表达水平分布来解决这个问题。在生物学上，这意味着只有当两个基因的不同剪接变体显示出协调的表达时，它们的表达才会被认为是相关的。如果不是这种情况，即使基因的总体表达水平相关，也不认为它们是共表达的。这种方法发现了新的功能模块，而这些功能模块用传统的共表达网络是检测不到的[51]。此外，使用这种方法，包含多个外显子和转录物的基因在网络中获得了更多的相关位置[50]，这是一个令人欣慰的结果，因为剪接变体可能具有不同的功能，因此很可能与功能不同的伙伴共同表达，而共同剪接网络则说明了这一点。
一个不同的方法是根据映射到其不同外显子的读数的分布来确定源自同一基因的不同异构体的表达。SpliceNet使用了这种方法，它有效地将映射到与两个异构体共享的外显子的读数按比例划分为两个整体异构体各自的总表达[52]。这意味着，如果两个异构体，即异构体A和异构体B，只共享一个外显子X（一些读数映射到该外显子），但没有读数映射到异构体A的其他外显子，而一些读数映射到异构体B的外显子，那么所有映射到外显子X的读数就被分配到异构体B，导致异构体A被认为根本没有表达。虽然这个优雅的解决方案是用模拟验证的，但没有进行实验验证。
构建基于RNA-seq的共表达网络的最常见的方法是在RNA-seq数据分析中合并所有重叠的基因异构体，然后在基因水平上构建网络。然而，这种方法失去了关于同一基因所编码的不同转录本的信息。另外，可以构建基于转录物的共表达网络。这些网络的缺点是由于许多基因异构体和非编码RNA的存在，其规模急剧增加。由于共表达网络是方形矩阵，网络的大小随着所包含的基因数量的增加而呈二次方（n2）。由于人类基因组中有20万个被注释的转录物（根据Ensembl GRCh38.p5（人类）注释[53]），而只有2万个蛋白编码基因，因此产生的网络大小增加100倍，大大增加了分析所需的计算资源。解决这个问题的一个办法是，从数据的子集建立共表达网络块，并在分析的后期将这些块结合起来[54]。然而，我们建议用户谨慎使用块状聚类（block-wise clustering），因为它可能会影响后续模块检测分析的结果，而且目前还不清楚使用大量块状时这些分析的表现如何。
RNA-seq data for co-expression networks用于共表达网络的RNA-seq数据
RNA-seq分析需要多个步骤，包括从测序读数中获得表达量估计，数据归一化和质量控制。最近在[55]中回顾了从RNA-seq数据中获得可靠的表达量的不同工具和方法，这里就不做回顾了。
根据我们的经验，不同的归一化方法在共表达分析中引入了不同的偏差，通常是正相关的偏差。为了解决这些归一化的问题，新的方法正在不断产生。例如，最近发表的从RNA-Seq数据中提取模式和识别共表达基因（EPIG-seq）的方法，旨在计算RNA-seq样本间的基因相关性，不受样本间读深差异和RNA-seq衍生的表达矩阵中大量0值的影响[56]。在单细胞实验中，由于每个细胞的RNA数量较少，源于大量0值的偏差甚至更加明显。为分析单细胞RNA-seq数据，已经创建了专门的工具，并在[57]中进行了回顾。虽然有一些研究比较了RNA-seq数据的不同归一化方法[58]，但还需要更全面的比较研究，包括新的方法。
Minimum read depth and sample size required for co-expression analyses共表达分析所需的最小读深度和样本量
为了从RNA-seq数据中创建共表达网络，有人建议至少要有20个样本[21, 54]，而增加样本量会产生具有更高功能连接的网络[21, 59]。毫不奇怪，高质量的数据往往能产生更准确的共表达网络[21, 59]。因此，为数据质量控制设定截止阈值是非常重要的。RNA-seq样本的总读深度越高，表达测量的准确性就越高，特别是对于低表达的基因[21, 59]。对于RNA-seq数据，测序深度截止阈值通常是任意选择的。一些共表达研究使用了每个样品1000万读数的截止值[2, 21, 60]。有学者认为，如果从相同数量的样本构建的共表达网络，其质量与基于芯片的共表达网络相似[21]，但随着读数的减少，质量也会下降。映射读数的百分比是另一个经常被考虑的分界线，其中<70%或80%的读数映射(reads mapping)到基因组的样本被删除。Giorgi等人使用65个拟南芥样本的1200万条读数，但只应用了30%的映射临界值（cut-off threshold），结果基于RNA-seq的共表达网络与生物网络的相似度低于微阵列网络[61]。根据基因组注释的质量等因素，每个物种的截断阈值可能有所不同。随着更多、更高质量的数据的出现，更高的截断阈值可能是更好的。
为了确保网络的稳健性（robust），可以使用自举法（bootstrapping）[62]。这是通过使用数据中的随机样本集（random sets of samples）（一个样本可以是多个子集的一部分）重复构建网络，随后用来评估从整个数据集创建的网络的可重复性。将数据集随机化（例如，通过随机地将表达值重新分配给它们的基因/转录物标识符并重建网络）也可以帮助识别由于特定的偏差而随机发生的相关性，而不是生物相关的相互作用的结果[2]。
Clustering and network analysis
Identifying modules 识别模块
聚类是用来对在多个样本中具有类似表达模式的基因进行分组。由此产生的模块通常代表生物过程[63, 64]，并且可以是表型特定的[65]。
用于共表达分析的最广泛的聚类包是加权基因相关网络分析（WGCNA）[40]。这个易于使用的工具在由表达数据创建的相关网络上使用层次聚类法构建共表达模块[54]。层次聚类法迭代地将每个聚类划分为子聚类，形成一棵树，其分支代表共表达模块。然后，通过在一定高度上切割分支来定义模块（图1）。
WGCNA是第一个应用于RNA-seq数据的共表达工具；它有效地识别了表型和模块之间的生物相关关系[19, 66, 67]，表现与基于微阵列的分析相似。一项基于RNA-seq的对正常和衰竭小鼠心脏的共表达研究发现，许多lincRNAs存在于与衰竭小鼠心脏表型相关的集群中，表明这些非编码RNAs可能在这种疾病中发挥作用[67]。对瘦猪和肥胖猪的RNA-seq数据进行联合表达分析，发现了与肥胖有关的模块[66]，并发现肥胖、免疫系统和骨重塑之间存在联系，研究发现CCR1、MSR1和SPI1可能是这些过程中的调节器。WGCNA也被用来从单细胞RNA-seq数据中识别生物相关的关联。通过使用dGCNA[12]定义的保存检测，确定了人类和小鼠之间保守的胚胎移植前过程的调控机制和基因，这一功能后来被添加到该软件包中[68]。分别为人类和小鼠的不同发育阶段确定共表达模块。然后将每个阶段确定的模块在人和小鼠之间进行比较，发现小鼠卵母细胞形成的共表达模块与人的卵母细胞和单细胞阶段的共表达模块有强烈的重叠。这表明，人类和小鼠在早期发育中共享核心转录程序，但在后期阶段出现分歧[12]。
Identifying hub genes识别枢纽基因
通过聚类确定的共表达模块往往很大，因此，确定每个模块中哪个（些）基因最能解释其行为是很重要的。一个广泛使用的方法是识别共表达网络中高度连接的基因（枢纽基因）（hub-gene）。与其他节点相比，枢纽基因（hub-gene）经常与网络的功能更相关[69]。在生物网络中也是如此[32]，尽管数学推导表明，只有模内枢纽基因（相对于模间枢纽基因而言[64, 65]）才是如此。模内枢纽基因（ Intra-modular hubs）是网络中特定模块的中心，而模间枢纽基因（ inter-modular hubs）是整个网络的中心（图2）。为了识别枢纽基因，经常使用中心度量（ centrality measures），主要是 "间性中心度"（ ‘betweenness centrality’）。具有高中心度的基因作为通过网络的最短路径连接器是很重要的[70]。连通性（Connectivity）经常被用来衡量网络的稳健性，表明在剩余的基因断开连接之前需要从网络中移除多少个基因。识别共表达网络中的枢纽基因已经导致识别了几个在癌症[71, 72]、2型糖尿病[73]、慢性疲劳[74]、其他疾病[75, 76]和组织再生[77]中必不可少的基因。
图2
假设的网络解释了模间和模内枢纽和网络中心性。模块间枢纽具有很高的网络中心性，因为它是所有可能的节点对之间最短路径的要求。红线表示一对节点之间通过网络的最短路径的一个例子。模块内的枢纽（用橙色标记）是单个模块的中心，通常具有很高的生物相关性。
由于一个模块中通常有多个中心基因或差异连接的基因，因此并不总是清楚哪个是表型的最重要基因。也不能保证任何一个中心基因对表型有因果关系[78]。一项研究使用了1617个样本的心肌数据，发现在发育中的心肌和患病的心肌所共有的基因网络中上调的已知胎儿基因标志物并不是枢纽基因[79]。另一项沙门氏菌的共表达研究发现，枢纽基因对于生长、压力适应和毒力来说是可有可无的，这说明枢纽基因不一定是必不可少的[80]。
一个广泛使用的为模块附加生物学意义的方法是使用表1中描述的工具确定模块内基因之间的功能富集度（functional enrichment ）。假设共表达的基因在功能上是相关的（functionally related），富集的功能可以分配给同一共表达模块中注释不全( poorly annotated genes )的基因，这种方法通常被称为 "关联罪"( ‘guilt by association’ )（GBA）[121]。如果一个模块内有相当比例的基因与某种疾病相关，GBA方法也被广泛用于识别新的潜在疾病基因[26, 121-126]（图1）
当使用GBA方法时，重要的是要记住，模块中的每个基因不一定与它被富集的功能或疾病关联相关。因为共表达模块通常由大量的基因组成，任何功能过程或疾病相关基因组的过度代表很快就会变得具有统计学意义，这通常由欺骗性的低P值表示。对这些低P值的误解可能会导致错误的结论，即一个模块中的所有基因在特定过程或疾病中发挥重要作用。实际上，一个模块中与主要生物学功能有关的基因比例往往<20%[127]，模块与性状的相关性可能相对较低（相关性<0.5），即使有统计学意义[128]。
Regulatory network construction调控网络的构建
尽管有充分的证据表明共表达分析可以帮助识别在疾病和生物功能中起重要作用的基因，但从共表达网络中推断因果关系仍然很困难。ARACNE[23]和GENIE3[113]等工具试图从共表达网络中构建调节网络。ARACNE删除了基因之间的间接联系（即一个基因的伙伴之间的相关性比该基因本身的相关性更强），只留下那些预计具有调控作用的联系。GENIE3结合TF信息，通过确定最能解释其每个目标基因表达的TF表达模式来构建一个调控网络。GENIE3的一个局限性是，需要TF信息才能使其表现优于随机机会[113]。这些方法的性能已经与黄金标准进行了比较，黄金标准是由超过150项研究中实验验证的调控相互作用定义的。比较结果表明，只有在使用扰动实验数据构建网络的情况下，试图仅从共表达网络中推导出调控网络的方法才能可靠地区分真性和假性调控相互作用[129] 。这些工具与包括WGCNA在内的其他工具的比较表明，WGCNA和ARACNE在定义大肠杆菌的网络结构方面表现最好[130]，对于大肠杆菌来说，一个定义明确的调控网络被作为一个黄金标准[131]。
Differential co-expression analysis差异共表达分析
差异共表达分析可以确定生物学上重要的差异共表达模块，这些模块使用常规的共表达或差异表达分析是无法发现的。不同样本组之间差异共表达的基因更有可能是调节器（ regulators），因此很可能解释表型之间的差异[4, 8-10]。差异共表达分析已被用于识别健康和疾病样本之间[4, 8-10]或不同组织[11]、细胞类型[5]或物种[132, 133]之间差异的基础基因。下面，我们提供一个常用的和新出现的方法和工具的概述，分为两类：（1）识别预定义样本组（如条件、时间点或组织类型）之间的差异共表达的方法；（2）不需要关于样本组的先验知识，使用一种算法来识别样本中先验未知亚群的共表达集群的方法。
Differential co-expression analysis between sample groups样本组之间的差异性共表达分析
大多数差异性共表达分析依赖于差异性聚类；它们识别包含不同基因或在变化的条件或表型下表现不同的聚类。最经常使用的差异聚类分析程序是WGCNA[54]、DICER[4]和DiffCoEx[100]，它们都是首先识别在全部研究样本中共同表达的模块，也与其他程序进行过比较。然后，这些共表达模块可以与代表疾病状态或组织类型等的预定义样本亚群相关联。
WGCNA确定每个模块在每个样本亚群中的活性和重要性（图3A和3C）。对于每个模块，计算一个eigengene，它是最能描述该模块内所有基因在分析中包括的样本中的表达行为（以线性方式）的向量。然后，它通过识别与该模块eigengene表现相似的基因或作为模块内枢纽基因的基因（这些基因往往是重合的），优先确定这些模块中的哪些基因可能是与该模块相关的表型的基础。根据设计，DICER专门用于识别在样本组之间有不同关联的模块对，例如，在一个组中形成一个大的相互连接的模块，而在另一个组中形成几个小的模块（图3D）。DICER可能对时间序列实验特别有用，在这些实验中，共表达的变化是渐进的，例如细胞周期序列实验，其中的模块是特定于某一阶段的，并在各阶段之间的过渡期共表达。DiffCoEx关注的是与同一组基因有差异的共表达的模块。这种行为的最极端的情况是基因组以协调的方式从一组相关的基因 "跳 "到另一组（图3E）。在这种情况下，DiffCoEx会以类似的方式对'跳'的基因进行聚类。DINGO是一个较新的工具，其工作方式与DiffCoEx类似，它根据基因在特定样本子集（代表特定条件）中的表现与所有样本确定的基线共表达的不同程度进行分组[102]。这些是最有可能解释与两个不同网络相关的不同表型的基因。每种方法都通过设计检测特定的模块变化，但它们也可以检测它们没有专门设计的模块变化，并且在识别这些变化方面可能优于其他工具[130]。
图3
样本之间可能发生的基因共表达模式的变化。差异性共表达可以作为一个模块只存在于一个样本组中（A），作为模块结构的差异（B）或作为模块成员之间的相关强度的差异（C）而发生。此外，如果一个较大的相互联系的模块分裂成几个较小的模块（D），或者如果一组基因改变其相关伙伴['基因跳跃'（E）]，可以检测到不同的共表达。如果在差异共表达分析之前没有定义样本组，或者是未知的样本组，双聚类方法可以通过将样本同时分类到存在这些模块的组中来识别样本亚群的独特模块（F）。
一些研究使用差异性共表达网络分析来确定特定组织[11]或疾病状态下的独特网络[134]。公开的RNA-seq数据和GTEx和ENCODE等项目的迅速增加，产生了大规模的RNA-seq资料，使得不同组织内和不同组织间的共表达分析成为可能[11, 15]。GTEx项目收集并提供多个人体组织的表达数据，用于研究基因表达、调控及其与遗传变异的关系[135]。在一项比较GTEx数据集中35个组织的RNA-seq数据的研究中，根据每个组织的平均基因表达量，构建了一个组织层次结构。相关的组织，如来自不同脑区的组织，聚集在一起。这个层次结构被用来构建一个由特定组织的共表达网络衍生出来的单一的联合共表达网络--元网络。结果表明，在组织特异性网络中，具有该组织特有功能的TFs往往与组织特异性基因一起高表达。这些基因之间往往形成比其他基因更强的联系，但仍处于网络的外围（因此中心度较低），而组织特异性TFs则成为该模块的中心[11]。因此，组织特异性TFs可以通过识别组织特异性网络中共表达强度增加的模块（图3A和3C），并通过确定这些模块的中心枢纽而被发现。相反，那些不是TFs但具有组织特异性的基因应该可以通过识别这些模块中处于外围的基因来发现（图3B）。此外，一些TFs在不同组织中有不同的作用。这些TFs可望成为枢纽基因，在一个条件下是一个模块的中心，而在另一个条件下是另一个模块的中心。
差异连接的基因是那些在两个样本组之间有不同的共同表达伙伴的基因。这些基因似乎在两组之间观察到的表型差异中起着调节作用（图3D）[8-10]。例如，一项研究使用与DiffCoEx类似的方法，比较了肌肉增长的突变体牛的共表达与非突变体的共表达。通过识别差异表达最大的基因和与这些基因显示出最高差异联系的TFs[10]（图3D），确定了含有因果突变（肌肽）的TF。有趣的是，编码该TF的Mstn基因本身的表达几乎没有变化，这就提供了一个例子，说明差异性共表达分析可以发现单独的差异性表达分析所不能揭示的生物学上的重要发现。
Generalized Single Value Decomposition (GSVD)广义单值分解(GSVD)
广义单值分解（GSVD）是一种独特的差异性共表达分析，它依靠谱系分解（spectral decomposition）来确定共调基因的模块（modules of co-regulated genes）。这种方法的独特之处在于，它将样本和所有基因的表达总结为较少的变量，旨在用尽可能少的变量解释尽可能多的表达变化。在这里，我们重点关注将基因表达总结为主成分或 "小基因"（principal components or ‘genelets’），这个术语（ term）是在[104]中引入的，可以解释为与共表达模块相类似（analogy ），它代表了多个基因的部分表达。这些小基因的相对意义--描述小基因的信号在数据集中出现的程度（即小基因的表达程度）--可以在两个数据集中进行比较。如果显著性相似，该基因小体代表了两个数据集之间共享的共同表达模式，而显著性的差异则表明该共同表达模式是其中一个数据集所独有。高阶（HO）-GSVD是最近开发的，使用类似的方法在两个以上的数据矩阵之间进行比较[105]。
GSVD在2003年首次用于分析人类和芽殖酵母的微阵列表达数据，以确定这两个物种之间共同和独特的信息素和压力反应模式[104]。最近，HO-GSVD被证明能有效地识别对神经祖细胞的自我更新很重要的途径[136]。GSVD被证明可以识别多形性胶质母细胞瘤（一种脑瘤）的独特模式，这对预后很有用[137]。同样，在正常样本中活跃的小基因也被识别出来[138]。然后将这些小基因信号从癌症样本的总信号中去除，显示出癌症特定的特征[138]。这两项研究都表明，癌症特有的签名对癌症中重复的基因有很强的信号[137, 138]，这在癌症中是很常见的，这表明识别的轮廓反映了基因组中的致癌事件。
随着高质量表达数据成本的降低，差异性共表达方法越来越受欢迎，这并不令人惊讶。虽然这些方法尚未应用于RNA-seq数据，但最近从微阵列研究中得到的发现使其成为一个令人兴奋的前景。然而，由于这些方法对异常值很敏感，它们需要高质量的数据。
Differential co-expression without prior grouping没有事先分组的差异性共表达
检测数据亚群之间差异表达集群（subpopulations of data）的另一种方法是双聚类（ biclustering）。如果一个数据集包含几个生物学上不同但未知的样本组，双聚类可以识别出只在样本的一个子集中具有类似表达模式的基因，而不需要事先进行样本分类（图3F）。当这种信息不可用时，这一点特别有用，像那些使用Drop-seq系统[139]或inDrop[140]的大规模单细胞RNA-seq实验就可能是如此。
在临床研究（clinical study）中，通常可以预先定义健康和疾病样本组（ healthy and diseased samples）。然而，同一疾病可以通过不同的机制表现出来。这是癌症中常见的情况，不同的突变会导致共表达模式的不同改变，但表型相似（similar phenotype）[7]。双聚类（Biclustering）允许研究人员在难以预先定义生物相关样本组的情况下，对机制进行分解。为此，双聚类比其他共表达分析方法更有效[7]。
Cheng等人首先在共表达分析中使用了双簇法（biclustering）[141]，随后又开发和应用了更多的双簇法（由Pontes等人[106]回顾）。双聚类方法的选择取决于样本的数量（number of samples）和一些因素，如样本是否具有物种（species-）或组织特异性(tissue-specific )，以及所包括的样本是否构成疾病表型和/或不同的时间点。根据所使用的方法，双聚类方法在计算上可能具有挑战性[106]。应谨慎选择方法，因为不同的二聚类方法在同一数据集中会有不同的结果[142]。
最近，双聚类方法被应用于基于RNA-seq的表达数据。对来自蠕虫和果蝇几个发育阶段的表达数据进行分析，通过识别包含这两个物种之间不同发育阶段所特有的类似的正交基因组的双簇，导致识别出在发育过程中具有类似的、因而也是保守的功能的基因［132］。双聚类也被应用于单细胞RNA-seq数据[5]。由于双聚类同时对基因和样本进行分组，它能够同时识别细胞类型组和相应的基因模块，揭示了49种不同的细胞类型及其相应的细胞类型特异性基因模块，这些结果后来被实验验证所支持[5]。随着单细胞RNA-seq的出现，双聚类方法可能能够识别出存在于患病细胞而非健康细胞中的细胞类型特异性模块。
另一种双簇方法通过其在癌症样本特有的双簇中的存在，确定了乳腺癌中的miRNAs减控[7]。这些miRNAs被认为是诊断和治疗反应的标志物[7]。双聚类也被用来识别癌症患者亚群特有的紧密共表达的蛋白质编码基因组，这可以用来了解患者的预后，并进一步推进精准医疗方法[18, 20]。在另一个癌症数据集中，一个三维聚类方法（triclustering）被用来识别跨样本亚群和时间点共同表达的基因[6]。这种方法通过使用每个三簇样本之间的eigengene变化，确定三簇中在早期和晚期时间点的癌症样本之间差异表达的枢纽基因，有效地识别了乳腺癌细胞系中的几个已知乳腺癌基因[6]。
Comparison of differential co-expression analysis methods差异性共表达分析方法的比较
尽管对用于差异共表达分析的方法进行全面和无偏见的比较是可取的，但这些工具的性能可能取决于情况，在不同的物种、疾病状态甚至可能是数据集之间有所不同，因此难以确定每种情况下的最佳方法。最近，有人试图对10种差异性共表达算法进行比较，但结论是，由于缺乏验证这些方法结果的黄金标准基因集，仍然很难对这些算法进行评估[143]。本综述中描述的几个工具已经在介绍竞争性方法的出版物中进行了比较。基于差异表达模块的功能富集分析，DICER被认为比DiffCoEx和CoXpress[4]的表现更好。基于其在模拟数据中检测集群的能力，HO-GSVD的表现优于WGCNA和DiffCoEx[136]。虽然双聚类是一个强大的方法，但它的表现不一定比其他网络分析方法如WGCNA更好，正如在模拟数据上使用不同工具的比较所显示的那样[144]。然而，正如前面所讨论的，二聚类可以在不需要事先进行样本组分类的情况下进行。
尽管本综述中描述的许多工具和方法最初是为微阵列数据创建的，但它们也适用于RNA-seq数据。有一些RNA-seq特定的差异共表达分析方法，利用外显子或等位基因特异性表达信息或等位基因特异性表达效应，据报道比其他不考虑这些信息的工具表现更好[52, 145]。然而，目前还不清楚的是，如果其他方法提供相同的异构体特异性表达信息，这些差异性共表达方法是否也能表现得更好，这可以在差异性共表达分析之前确定。因此，仍然很难评估这些新工具是否比已经成熟的工具（如WGCNA）表现更好，WGCNA也可以用于同构体特异表达数据[50]。
由于本综述中描述的工具是根据不同的标准来创建模块的，所以比较中使用的措施是否代表了所有情况下的理想属性也是值得怀疑的。例如，DiffCoEx根据不同的共表达行为对基因进行分组，而WGCNA则是识别在多个样本和条件下共表达的模块。在同质化的数据集中，DiffCoEx可能会检测到较少和较小的模块，表明没有很多差异共表达的基因。当DiffCoEx与WGCNA在这种同质数据上的性能比较时，这将可能导致较低的富集分数，而这只是表明数据中发生的共表达伙伴变化不多。相比之下，专注于共享共表达模块的工具很可能会发现具有高富集分数的强相关模块，如果目标是识别调控模块，这可能并不相关。WGCNA已被广泛证明在许多不同的情况和不同的目的下表现良好[54]。然而，它需要关于样品条件的信息，以便将模块分配给条件。如果没有这种信息（如在大规模的单细胞RNA-seq实验中），或者研究人员希望在样本组中确定子组，那么双聚类是一种更合适的方法。
为了系统地评估不同工具和方法的性能，诸如DREAM4和DREAM5[146]等项目是非常有价值的。这些项目要求研究人员从模拟和体内的基准数据集中构建调控网络。由于这些挑战是预先定义的，它们允许研究人员以无偏见的方式测试其方法/工具。然而，这些挑战最后一次提出是在2010年，此后又开发了许多新的方法和工具。
Integrated network analysis
综合网络分析
实验验证往往集中在单一基因上。由于这些实验成本高、耗时长，对因果基因的高置信度预测是非常重要的。仅仅基于共表达的分析（还）不能提供这种程度的信心。因此，结合其他类型的数据信息可以帮助确定可能是表型基础的基因的优先次序。例如，可以使用描述哪些基因是TFs的信息来实现，正如GENIE3[113]对调控预测的做法。然而，仅仅关注TFs是不够的，通常需要整合多种数据类型，以提高所形成网络的准确性和实用性[13, 147]。
将不同层次的数据信息结合起来，可能会以多种方式产生新的生物学上可解释的关联。如果模块内的枢纽基因是TFs或TF的目标，这个TF就更有可能在被调查的表型中起因果作用[10]。如果在同一模块中存在多个全基因组关联研究（GWAS）的命中率，它们的累积存在可以大大促进疾病的发展[120, 163, 164]。共表达模块内基因的不同甲基化状态可以阐明疾病背后的甲基化模式[165]。如果多个基因受同一遗传变异体的调控（在反式eQTL效应下），就有可能通过识别驱动反式eQTL效应的顺式eQTL基因来确定负责网络的改变的基因（图4）。疾病相关基因组变异的反式调控下的基因有时与相应的疾病相关的过程或途径有功能上的联系，这一点得到了支持。这方面的好例子是IFN（干扰素）-α和补体途径，其中有几个基因处于系统性红斑狼疮相关变体的反调控之下，可能是通过顺式调控IKZF1[155]。将调控基因变异信息整合到共表达网络分析中，以顺式eQTLs作为因果锚，确定TYROBP是晚发性阿尔茨海默病患者最可能的因果因素，这一发现得到了观察的支持，该基因的突变已知会导致Nasu-Hakola病[128]。最后，拷贝数变异可以影响基因表达水平，包括这些信息可能有助于识别和/或解释疾病或性状中存在的共表达网络结构的改变[138]。
图4
整合多组学数据与共表达分析的策略。如果使用相关组织的特定表达数据构建网络，其信息量更大。基因组变异可以被映射到共表达网络中，方法是将暗示性的GWAS命中率与网络中的基因联系起来，或者首先确定对基因表达水平有影响的遗传变异（顺式和反式eQTLs），然后将其映射到共表达网络中。额外的数据层可以包括TFBS（基于结合图案或ChIP-seq/ChIP-chip实验）、miRNA靶点结合点（基于硅学预测或实验技术）和已建立的蛋白-蛋白相互作用。共表达网络可用于识别模块、中心基因和预测未知性状相关基因的功能。确定的模块可以通过富集分析来确定叠加的特征。此外，如果有相应性状的病例和对照的omics数据，可以通过额外的差异表达、共表达和甲基化分析来支持研究假设。 eQTL：表达定量性状位点；GWAS：全基因组关联研究；OMIM：人类在线孟德尔遗传；miRNA：微RNA；PPI：蛋白质-蛋白质相互作用；TF：转录因子；TFBS：TF结合位点。
总的来说，多种数据类型的整合可以提高预测结果的准确性[13, 147]。例如，通过整合肿瘤基因组序列和基因网络，发现了不同亚型癌症的独特模块[166]，这些模块可能对预后和识别基于药物的个性化治疗的假定目标很有用。本综述前面描述的一些工具可用于差异共表达分析，但也可应用于其他数据类型。在最初的DINGO出版物中，作者对mRNA表达、DNA拷贝数变化和甲基化数据进行了综合分析。通过叠加每种数据类型的差异网络并识别所有数据中存在的边缘，PI3K通路的一些基因被确定为多形性胶质母细胞瘤患者的重要角色[102]。这一途径是一个已经确立的治疗目标，支持了这一观点，即这是一个识别疾病研究相关目标的有效方法[167]。最近发表的一个工具CoRegNet允许在共表达分析中整合不同类型的数据，从不同的数据类型中确定基因的合作调节器[114]。另一个已建立的方法，cMonkey，通过计算不同数据类型的联合双簇成员概率，实现了类似的数据整合，即识别多个数据类型中的基因组[115]。
Future prospects未来展望
近年来，差异共表达分析已被越来越多地用于分析大型数据集。这可能是由于大规模基因表达谱分析的成本降低，特别是RNA-seq，样本量增加，以及来自扰动实验的组织特异性数据的更多可用性，这对于富有成效的差异共表达分析是必需的[103, 168]。同样，双聚类算法也受益于更大的样本量和更高的数据质量，如癌症亚型特有的共表达模块的鉴定[18, 20]。双聚类在单细胞RNA-seq数据上的作用已经通过对不同细胞类型的分类和对特定细胞类型中唯一共表达的基因簇的鉴定得到证明[5]。我们期待这些方法在未来得到更广泛的应用，因为它们受益于RNA-seq数据数量和质量的增加，这将使我们能够更准确地识别组织特异性和细胞类型特异性的疾病相关模块和调节器。
大规模的单细胞测序技术被越来越多地使用，首次使用这种技术的共表达研究发现了细胞类型特异性的共表达模块，这些模块在多细胞类型的共表达分析中是不会被发现的[5, 12]。因为后者代表了多种细胞类型的聚集信号，它们通常不能检测到不同实验组之间的细胞亚群的改变。在非细胞类型特异性数据的分析中，与衰老相关的细胞周期基因的表达下降，这一观察结果支持了这一点[169]。然而，来自单细胞实验的数据显示，这一观察是由高度表达细胞周期基因的G1/S细胞比例下降而不是整个细胞群的表达改变引起的[170]。
一个额外的前景是从RNA-seq数据中检测突变[171]。由于突变在不同的细胞中随着年龄的增长而积累，这些突变可以被用来识别细胞的起源。突变的积累已经被用来研究癌症的发展和转移的起源[172]。在大规模的单细胞RNA-seq实验中，突变可用于根据细胞的起源来分离细胞，或根据它们所携带的突变来分组[173]。携带相同突变的细胞可以调查共同表达模式，并且可以检测到具有特定突变的细胞所特有的模块。这可能允许将突变与表达模块直接联系起来，但有一个限制，即只有编码区的突变可以在RNA-seq数据中检测到。
尽管单细胞RNA-seq数据有许多令人兴奋的新可能性，但重要的挑战仍然存在。通常情况下，每个细胞被测序的读数数量较少，然后来自同一类型的多个细胞的信号被聚集起来，以获得细胞类型的特定基因表达谱。对于较罕见的细胞群，如干细胞，很难获得足够的数据，这也是目前对这些细胞类型分析的限制。此外，每个细胞的低读数导致稀疏的表达矩阵，目前用于典型RNA-seq分析的归一化方法并不适应。这些归一化方法通常还假设大多数基因在不同样本之间的表达没有变化，而在单细胞RNA-seq中，由于不同细胞之间的表达变化，情况并不一定如此。从单细胞获得高质量的RNA的困难进一步加剧了这种情况。这些问题和其他问题将在[174]中进一步讨论。
除了在单细胞RNA-seq中出现的归一化问题外，对大量RNA-seq数据进行归一化的最佳方法也仍然不清楚。广泛使用的片段/每千碱基百万读数（FPKM）归一化方法已经引起了争论[58]，尽管已经和正在创建替代方法，但每种方法都有其局限性。此外，根据我们的经验，在某些情况下，使用不同的映射工具会导致不同的结果。尽管已经对不同的工具和方法进行了一些比较[175]，但使用公共数据等进行大规模的比较，可以确定这种情况，并确定追求每个研究问题的最佳做法。
随着RNA-seq、基因组序列、ChIP-seq、甲基组和蛋白质组数据等不同数据类型的可用性增加，整合这些数据集以更准确地预测调控基因将成为可能。像GTEx[156]、表观基因组路线图[176]和ENCODE[15]这样的大型联盟的项目已经在产生来自多个组学水平的数据，以促进这些综合分析。为了识别调控关系，扰动数据是最好的，因为典型数据不能区分调控关系中的真假阳性[129, 168]。此外，调控关系可能具有高度的细胞类型、组织或发育阶段的特异性[129]。目前只有少数工具和方法可以用来研究多组学数据，而且现有的工具大多只整合了两层全息数据[177]。综合网络分析伴随着额外的数学挑战，而且最佳实践还远未建立。对这一主题的进一步研究是研究界非常感兴趣的，因为这将使人们更好地理解能够解释共表达模式和疾病机制的调控机制。对这些疾病机制和相应的共表达模式的更好理解将有助于确定适当的干预研究目标。
网友评论