RNA-Seq differential expression analysis: An extended review and a software tool RNA-Seq差异表达分析： 扩展评论和软件工具

RNA-Seq differential expression analysis: An extended review and a software tool   RNA-Seq差异表达分析： 扩展评论和软件工具

正确鉴定特定条件之间的差异表达基因（DEG）是理解表型变异的关键。高通量转录组测序（RNA-Seq）已成为这些研究的主要选择。

因此，用于RNA-Seq数据的差异表达分析的方法和软件的数量也迅速增加。 但是，对于最合适的管道还是没有达成共识用于从RNA-Seq数据鉴定差异表达基因的方案。这项工作对该主题进行了扩展审查，其中包括对六种绘图方法的评估读取，包括伪对齐和准映射以及来自RNA-Seq数据的九种差异表达分析方法。

使用qRT-PCR数据作为参考（金标准），基于真实RNA-Seq数据评估所采用的方法。作为结果的一部分，我们开发了一个软件，可以执行本工作中提供的所有分析，可在https://github.com/costasilvati/consexpression免费获取。结果表明，考虑到映射方法对最终DEG分析的影响很小采用的数据有一个带注释的参考基因组。关于采用的实验模型中，具有更一致结果的DEGs识别方法是limma + voom，NOIseq和DESeq2。此外，五种DEG识别方法的共识保证了DEG的列表具有很高的准确性，表明了这种不同的方法组合可以产生更合适的结果。共识选项也包括在内用于可用的软件。

介绍
高通量测序已成为测量表达水平的主要选择，即RNA-Seq [1]。RNA-Seq可以在不事先了解参考文献的情况下进行感兴趣的序列，允许广泛的应用，如：“从头”重建转录组（没有参考基因组），核苷酸变异的评估，评估甲基化模式[2]，仅举几例。RNA-seq技术比cDNA微阵列具有一些优势，例如高通过泳道和流通池的数据再现性水平，这减少了实验的技术重复数量。

此外，RNA-seq可以识别和量化同种型和未知转录本的表达[3]。关于越来越受欢迎高通量测序方法，下一代测序实验的成本

已大幅下降。然而，对于RNA-Seq的定性和定量分析的清晰理解尚未实现，特别是与较老的方法如cDNA微阵列相比时[4]。

一般而言，RNA-Seq技术对于涉及某些特定条件的差异表达分析非常有用[5]，其中通常采用五个步骤[6,7]。 首先，将RNA样品片段化成小的互补DNA序列（cDNA）和然后从高吞吐量平台排序。

其次，将小的生成序列映射到基因组或转录组。第三，估计每种基因或同种型的表达水平。第四，映射的数据被标准化，例如使用统计和机器

在学习方法中，鉴定了差异表达的基因（DEG）。 最后，最终从生物学背景中评估所产生数据的相关性[8]。随着RNA-Seq技术的日益普及，开发了许多软件和管道

这些数据的差异基因表达分析。来自RNA-Seq的差异基因表达分析方法可以分为两个主要子集：参数和非参数。参数方法捕获有关参数内数据的所有信息。

在这些情况下，可以通过观察采用的模型及其参数来预测未知数据的值。当参数化方法应用于差异基因表达假设，通常在归一化后，给定基因的每个表达值被映射到特定分布，例如泊松[9±11]或负二项[12±14]。

另一方面，非参数方法可以捕获关于数据分布的更多细节，即，不强加拟合的刚性模型。这是可能的，因为非参数模型考虑到不能从有限的参数集定义数据分布，因此关于数据的信息量可以随着其体积而增加。

关于RNA-Seq差异表达分析，一些工具如edgeR [13]和baySeq [11]采用负二项模型作为主要方法。

其他软件工具，如NOIseq [15]和SAMseq [16]，采用非参数方法。

一些方法基于转录物检测，其已被开发用于鉴定未知的转录物或同种型，并且还可用于鉴定DEG，例如EBSeq [17]和Cuffdiff2 [18]。如今，关于哪种方法最为重要，尚未达成共识适当的或哪种方法确保结果在稳健性，准确性方面的有效性和再现性。生物信息学研究中的这一主题仍在发展[5,19,20]。为了评估DEG的标准化和检测的统计方法以及文库制备对结果的影响，开发了一些研究工作[10]，通过考虑微生物来评估差异表达分析的方法，包括用于分析的映射方法[21]，并用模拟数据评估软件和管道[20,22]。特别是Rapaport等人。 [23]评估了实际数据集中差异表达分析软件的概要，考虑了分析的特征，如准确性，标准化，DEG的检测和没有检测到表达的条件。张等人。 [5]评估了重复次数，测序覆盖率和比较组的影响。郭等人。表明三种DEGs鉴定方法之间的排序可以产生更准确的鉴定[24]。李等人。 [8]评估标准化用于DEG检测的方法，表明两种标准化方法的结合导致更好的结果。Seyednasrollah等。 [25]介绍了八种软件工具的比较

实际数据中的DEG分析。Germain等人。 [26]介绍了有关RNA-Seq数据分析步骤的工作，比较了不同的转录物定位和定量方法，同时也为采用的方法比较提供了一个在线工具。

最近，Yu等人。 [27]提出了一个基于模拟的程序，采用负二项分布和广义线性模型（在基因水平）。该方法的主要目标是降低先前研究中报道的I型错误的高发率[17]，即

假阴性。Abedalrhman和Rueda [28]介绍了Zseq工具，指出了高通量测序数据分析中预处理步骤的重要性。 更具特异性，Zseq专注于改善成绩单的组装，用不同的预处理方法评估DEGs的结果。

另一方面，其他方法一直关注评估其他观点，例如RNA-Seq实验所需的生物复制数量，以及基于实验重复次数分析差异表达的最合适工具[29]。 。从Sahraeian等人的不同角度对RNA-seq数据进行全面而系统的分析。除了先前产生的RNA-Seq数据的表达分析外，al（2017）可以作出显着贡献[30]。与这些研究不同，我们评估了映射方法对差异基因表达分析结果的影响。

我们还通过不同的视角评估DEGs分析的方法，不仅指出了更好的方法。先前的研究及其结果表明，DEGs分析受许多因素的影响，如文库的制备和实验的结构。在这种情况下，我们用RNA-Seq数据分析了DEG鉴定中必需步骤的影响，并开发了一种软件，可以获得主要DEGs鉴定方法的结果。六个映射器之间的比较研究，包括一个伪对齐和一个通常用于差异表达研究的准映射工具，导致识别

这一步骤在DEG的分析和鉴定中的重要性。还采用了金标准qRT-PCR数据，以评估DEG鉴定工具的准确性，并指出其结果具有高可靠性的工具。这项工作的另一项贡献是评估DEGs识别方法的综合结果，

我们的工具允许执行五种不同的差异表达分析方法的共识，结果表明差异表达的基因具有更高的可靠性和准确性。在本研究中，我们用RNA-Seq数据对差异基因表达分析的主要方法进行了扩展回顾，评估了作图和量化方法的影响。在本研究中，我们采用了绘图软件Bowtie2 [31]，TopHat [32]，BWA [33]和STAR [34]。对于其他方法，例如伪对齐和准映射，我们采用了kallisto [35]和Salmon [36]。我们分析了差异表达分析软件，它代表了该领域的最新技术水平，例如baySeq [11]，DESeq [12]，DESeq2 [37]，EBSeq [17]，edgeR [13]，limma + voom [38]，NOIseq [15,39]和SAMseq [16]。将映射结果用作一些差异表达分析软件工具的输入，并将其结果与qRT-PCR进行比较[40]，从而验证与不同映射器相关的每个软件的准确性。

结果表明，NOIseq [15,39]，limma + voom [38]和DESeq2 [37]是考虑精度，准确度和灵敏度最平衡的软件。我们评估了不同方法中个体和综合方式的结果。

结果表明，一组软件可以比单个解决方案一起产生高精度和准确性。 最后，这项工作仍然很容易作为软件工具的贡献适用于差异基因表达分析的不同实验。该软件工具提供了一个集成执行，具有映射，映射计数（如果需要）和表达水平的量化，表明所采用的方法在识别DEG时的性质和准确性的特征。

材料和方法数据集

这项工作采用了为微阵列质量控制（MAQC）项目[10,40]生成的实际数据集。使用Illumina的Genome Analyzer II获得数据集。该实验分析了两个生物样本：来自Ambion人类大脑的RNA和Stratagene的人类通用参考RNA，我们将在这项工作中分别称为Brain和UHR集[10]。我们只使用了使用PhiX Control的Brain和UHR样本。数据集可在NCBI Short-Read Archive（SRA）上获得，登记号为SRA010153。读数针对人类基因组/转录组，版本19（GRCh37.p13）进行了定位。作为MAQC项目的一部分，通过qRT-PCR分析了大约一千个基因[41]。qRT-PCR数据可在Gene Expression Omnibus获得，访问：GSE5350，平台GPL4097 [40]。

Ambion人脑和Stratagene通用人体样本也被用作该实验的生物学参考。我们将qRT-PCR数据视为评估DEGs鉴定方法的金标准。从RNA-Seq数据（ENSEMBL）的注释到qRT-PCR数据的转换由在线工具bioDBnet [42]进行，不包括重复的ID或同义词。该转换产生了997个独特qRT-PCR基因的列表。有关qRT-PCR基因列表的详细信息，请参阅S1表。序列比对和基因计数

采用的RNA-Seq数据集映射在人类基因组/转录组（hg19）中，同一版本的注释文件均来自GENCODE项目[43]。

转录组到基因组注释的转换由R包txImport [44]完成。对于映射和量化，使用了各种方法：拼接读取对准器，未拼写读取对齐器，伪对齐和准映射。对于拼接读取对准器方法，使用了TopHat软件（v.2.1.0）[18]，该软件应用了外显子优先方法。对于无符号读取对准器方法，使用了两个映射软件，BWA（v.0.7.12-r1039）[33]和Bowtie（v.2.2.6）[31]，它们应用了Burrows-Wheeler变换。

对于伪对齐方法，使用kallisto软件（v.0.43.1）[35]。对于准映射方法，使用Salmon（v0.8.2）软件[36]。对于映射执行，默认

采用各软件的参数。表1列出了采用的映射器。采用HTSeq软件（v.0.6.0）[12]生成计数矩阵，默认参数。

用于生成计数矩阵的采用的注释文件与映射中使用的相同。图1展示了目前的工作。通过每种采用的方法将RNA-Seq数据集（表示为图1中的ΒNCBI-SRA）映射到人类基因组（hg19），从而获得计数矩阵。矩阵用作所采用的差分表达方法的输入。为了评估映射软件对DEGs识别的影响，我们使用六个生成的计数矩阵分析了四个差分表达软件。 对于Salmon，STAR和kallisto，我们分析了两种差异表达软件。通过考虑以下DEG识别方法获得映射器性能：edgeR，DESeq，baySeq和NOIseq，到Tophat，Bowtie2和BWA，分析Salmon，STAR和kallisto分别进行了edgeR和NOISeq。 DESeq和baySeq只能使用计数数据运行（图1）。

差异表达

在这项工作中，我们比较了八种DEGs或转录本鉴定方法。在应用每个软件时，我们关注最常用的方法。因此，我们遵循手册中提供的指南，应用默认参数，包括每个软件的标准化方法。

 使用由HTSeq产生的相同计数矩阵进行所有差异表达分析。对于每个评估的映射器，生成计数矩阵或丰度矩阵，以这种方式对映射器工具评估DEG方法。表2总结了采用的DEG识别方法及其性质。

使用BWA，TopHat，Bowtie和STAR映射器绘制RNA-Seq数据。量化是从Salmon和kallisto工具获得的。每个映射器的计数表用作DEG识别方法（edgeR，DESeq，baySeq和NOISeq）的输入，从而为不同的映射器生成每种DEG识别方法的DEG列表。Salmon，STAR和kallisto的结果用作edgeR和NOISeq的输入。将结果与qRT-PCR（金标准）进行比较，以评估映射是否影响DEGs检测的性能。EBSeq，SAMSeq和limma + voom，DESeq2和侦探方法被添加到研究中，用于单独评估DEG识别工具，仅使用TopHat映射器的映射结果。如其用户指南所示，特殊侦探收到kallisto工具的量化输出。图1给出了概述

这项工作中提出的管道。采用的方法（使用TopHat映射器）确定的DEG用于评估结果部分中的性能统计数据。

·baySeq [11]：使用贝叶斯经验方法估计每组模型的后验概率，定义每个元组的差分表达模式。

·DESeq [12]：基于负二项分布，方差和均值受局部回归约束。

·EBSeq [17]：开发的主要目的是鉴定差异表达的同种型，它在鉴定DEG方面也很稳健。它类似于baySeq [11]，它采用贝叶斯经验方法。

edgeR [13]：Poisson超分散模型用于解释技术和生物学问题
变异。应用贝叶斯经验方法来调节过度分散对转录物的程度。

·limma + voom [38]：基于线性模型，最初开发用于分析微阵列数据，目前扩展用于RNA-Seq分析。 limma用户指南建议使用与使用voom转换相关的edgeR包的TMM标准化，它基本上将标准化计数转换为对数基数2并估计均值 - 方差关系以确定每个观察的权重
最初由线性模型[45]制作。

·NOIseq [15,39]：自适应数据和非参数，根据经验模拟计数数据中的噪声，并允许数据分析而无需复制。

·SAMseq：[16]：非参数方法，对不同的测序计数进行重新采样
深处。它可以应用于具有定量结果，两类或多类的数据。

·DESeq2 [37]：DESeq2首先建立一个具有观察计数的模型。其次，它使用与原始DESeq相同的方法拟合，或者分为两步：找到使似然最大的参数值，称为最大似然估计。
然后，它获取所有基因值并将这些值移向平均值。

 DESeq2使用贝叶斯定理来指导每个基因的运动量：如果基因的信息很低，其值接近平均值，如果基因的信息很高，则其值移动很少。因此，移动的值可用于评估不同的基因组以及应用阈值;

·侦探 sleuth[35]：
 sleuth工作流程首先过滤低丰度转录本，然后应用两个标准化，然后对每个转录本的模型进行参数估计。这使得生物方差的正则化有助于跨样品的转录物丰度变化，并最终导致每个转录物的总体方差估计。

·

结果与讨论

阅读参考基因组中的映射
为了评估映射方法，人类基因组在数据集小节中描述。

为了评估基因组作图工具对DEG分析的影响，
所有采用默认参数的绘图软件均采用。

每个映射的计数矩阵由HTSeq包[12]，通过htseq-count函数，使用基因组注释文件和默认参数生成。

每个映射工具的计数矩阵用作DEG检测方法的输入。

关于qRT-PCR数据，采用默认方法（Benjamini＆Hochberg）的GEO2R工具无法识别DEG。

它仅被认为是DEGs的转录物，log2FC> =±2且P±值<= 0.05。

完整的DEG列表可在S1表中找到。

我们比较了RNA-Seq（baySeq，edgeR，DESeq和NOIseq）中鉴定的DEG与qRT-PCR的DEGs。

 在图2和表3中可以观察到DEG集中在映射器之间的交叉点上，表明即使随着映射方法的改变，这些方法也保持了识别行为。

通过这种方式，我们观察到映射器对最终结果的影响很小。

 在图2A和2C中，可以观察到正确识别的DEG的数量（一致）
与qRT-PCR相比，与识别方法相比，与采用的映射器相关。

baySeq和NOIseq方法获得了少量未识别的DEG，并且该量没有随着不同的映射器而改变。为了评估不使用映射的方法，以及其他量化读数的策略，我们使用Salmon和STAR比较edgeR和NOISeq结果来量化

映射基因。表3显示了来自不同RNA-Seq作图方法的DEGs鉴定数。它只被认为是NOISeq和edgeR，因为baySeq和DESeq无法接收不同整数值的输入。结果再一次表明，差异表达分析更多地受到DEGs鉴定方法的影响，而不是所采用的读数定位或定量方法。S2表提供了有关使用不同映射器的每种DEGs识别方法的性能的更多细节。 DEG识别方法如上一节所述，映射器对表达式分析的影响很小。

通过这种方式，仅考虑TopHat映射结果开发了所有后续分析。在本工作的这一步，我们分析了以下软件工具的结果：limma + voom [38]，EBSeq [17]，SAMseq [15]，DESeq2 [37]和sleuth [35]。有关评估工具的更多详细信息，请参阅差分表达式小节。

 我们比较了9种工具差异表达的基因和qRT-PCR指示的DEGs。软件按照每本手册中的定义执行，并且通过每个工具的手册指示的限制，工具列出的基因被认为是差异表达的。基于每种方法结果与qRT-PCR之间的匹配来评估所采用的DEG鉴定方法的性能。表4列出了每种采用方法的性能。可以注意到EBSeq，SAMseq和DESeq方法虽然使用不同的DEG识别方法具有相似的行为，但呈现低TPR（真阳性率）和低ACC（准确度）。DESeq的性能可以通过以下事实来证明：小样本（每个条件两个样本）可以获得更好的结果，如[22]中所示。SAMseq的结果很大程度上受样本大小和重复次数的影响。SAMseq能够对最相关的DEG进行排序，不过它是结果产生许多误报[22,25]。NOIseq，DESeq2和limma + voom方法表现良好，具有较高的TPR和ACC率。limma + voom工具已经在之前的工作中被指出是DEG排名和两个以上样本分析的更好结果之一[22]。 NOIseq和DESeq2工具显示了一致的结果，表明这些方法适用于具有大量样品和注释基因组的实验。

整合DEG识别方法

对DEG识别方法的个别评估清楚地表明，每种方法都产生了截然不同的结果。

此外，一些方法有更好的结果和更多的方法
样本，而其他人表现出其结果的变化受其他特征的影响，如测序深度和异常高计数的异常值。

为了验证每种DEG识别方法的各个结果之间的兼容性并确定可能的性能改进，我们通过将所采用的方法集成到这项工作中来评估结果。

我们通过整合九种方法评估了结果中的性能，因此对于通过x方法鉴定为差异表达的每种基因，其中x是已经鉴定每种DEG的方法的数量。

将DEG鉴定方法的每种组合的结果与来自qRT-PCR的金标准进行比较。

为了评估DEG集成方法的性能，我们验证了性能更好的方法组合。图3显示了从1到9°的积分
识别方法。 可以注意到差异没有一致性

表达了九种方法整合的成绩单。从九种评估方法来看，

对于165个转录物，发生8个同时适应症的频率，qRT-PCR也表示差异表达。但是，当观察qRT-PCR指示的DEG数时，有可能观察到8种方法的整合未能鉴定出qRT-PCR所示的大量基因。为了确定具有更有效的DEGs指示的方法组合，并且在误差最小的情况下，我们评估了每个方法的DEGs识别性能。方法的子集：九，八，七，六，五，四，三，二和一。

每个子集的性能结果如表5所示。正如所料，每个子集的性能表明，一起考虑更多方法往往会提高准确性并降低错误率。正如在基因网络推断的背景下报道的那样，集体知识或数据整合可以产生比个体结果更好的结果[48,49]。基于这一原理，我们发现五种方法的整合可以获得比任何其他测试子集更高的TPR和SPC值。为了识别每个基数的DEG识别方法的最佳组合（1,2，...，9），我们采用ROC（接收器工作特性）曲线[50]，一种标准模式识别工具。图4给出了DEG识别方法共识的更好组合。可以注意到，五种方法的组合在所有测试组合中提供了最有效的解决方案。六种方法的共识导致FPR略有改善，但TPR也出现下降。5种DEG识别方法的共识提出了最佳的综合结果，具有更高的SPC和TPR值，从而得到高精度的结果。图5显示了通过增加DEGs识别方法的整合而相关的TPR和SPC值的演变。

包含方法带来了相当大的特异性增益（SPC），但是从六种方法的整合来看，TPR值经历了相当大的下降。该结果表明了https://github.com/costasilvati/consexpression上可用的软件的默认值，用户可以在执行默认共识的方法（五种方法）之间进行选择， 从而获得SPC和TPR之间的最佳平衡。应用这种方法的另一种可能性是改变定义所需方法数量共识的方法数量，同时考虑到变更的成本，其他绩效衡量标准，以及在这项工作中只使用其中一种采用的方法。表6概述了根据qRT-PCR正确鉴定DEG的方法组。关于来自qRT-PCR的413°，其中19种未通过任何方法鉴定。当将一个指示视为差异表达时（可接受任何工具），可以观察到它不可能达到指示的413个基因定量RT-PCR。

为了确定哪组方法具有最佳共识，评估很重要每种方法在聚合结果中的表现如何，尤其是在五种适应症的组中。表7显示了汇总结果中每种方法的频率。比较表7和表6中的结果，我们可以观察到343种方法通过五种方法的共识指出，最正确指出的方法（几乎所有适应症）都是baySeq [11]，DESeq2 [37]，limma + voom [38]和NOISeq [15,39]。关于五种方法的共识，baySeq方法表明所有DEGs呈现

在五个共识结果中。DESeq2表示97.6％，limma + voom方法表明其中96.5％，NOISeq表示95.9％。对于baySeq的分析，

有必要定义模型集合，每个模型是样本细分为组，假设同一组中的样本共享底层分布的相同参数。在DESeq2方法中，为观察到的计数创建模型，使用贝叶斯定理拟合该模型以指导每个基因的运动。在NOISeq方法中，成绩单如果两个条件之间的log2与两个相应条件之间的差值的比值可能高于噪声，则差分表达。通过比较相同条件下的所有重复对来获得噪声分布。在limma + voom方法中，读取计数被转换为每百万计数的log2（logCPM），并且平均方差比用精确权重建模。总之，baySeq方法倾向于更高的FP值，如表4所示，表明五种方法的100％DEGs共识。参数共享样本这种方法的一组，减轻了同一组基因的变异，从而导致该方法的正确性更大的可能性。另一方面，NOISeq，DESeq2和limma + voom方法以均衡的方式执行与正确DEGs识别的关系，从而使结果具有高可靠性，这只能证明'未通过qRT-PCR鉴定，未鉴定出3.8％DEG。关于edgeR结果，我们可以证实其TPR倾向于较低的可靠性，呈现由qRT-PCR指示的正确鉴定DEG的81.3％。

结论
这项工作提出了关于鉴定差异表达基因（DEG）或转录本的方法的扩展审查。

 我们评估了六种映射方法的影响，包括一种伪对齐和一种准映射，
DEGs鉴定和这些方法整合的九种主要方法，以便从结果中产生共识。

通过比较来自相同测试转录物的参考qRT-PCR的各自结果来进行所采用方法的评估。

我们已经确定，绘图工具对最终结果的影响很小，表明DEGs鉴定方法是RNA-Seq数据中差异表达分析的主要选择。

对于评估的实验条件，我们没有在评估的方法中确定在所有性能测量中获得最佳结果的工具。 NOIseq，DESeq2和limma + vomm方法分别以95％，95％和93％的特异性和80％，84％和81％的真阳性率呈现最佳个体结果。

关于DEG识别方法的整合，我们发现五种方法的结合提高了识别的灵敏度，并提供了更可靠的结果。

使用的五种方法综合产生了91％的特异性和83％的特异性
真阳性率，因此表明五种方法的共识比单独的解决方案更好地平衡。

最后，本研究还有助于在https://github.com/costasilvati/consexpression上提供免费软件，
它实现了所提出的分析，可以很容易地用于复制这项工作，以及分析其他RNA-Seq数据源。

支持信息S1表。 qRT-PCR分析。通过qRT-PCR指示差异表达的基因。
（PDF）S2表。

映射分析。不同映射器的每种DEG识别方法的性能。（PDF）

致谢
我们感谢审稿人在以前版本的稿件中提供了深刻见解。

作者贡献
概念化：Douglas Domingues，Fabricio Martins Lopes。

正式分析：Douglas Domingues，Fabricio Martins Lopes。

资金收购：Fabricio Martins Lopes。

调查：Fabricio Martins Lopes。
方法论：Juliana Costa-Silva，Douglas Domingues，Fabricio Martins Lopes。

项目管理：Fabricio Martins Lopes。

软件：Juliana Costa-Silva，Fabricio Martins Lopes。

监督：Douglas Domingues，Fabricio Martins Lopes。

验证：Juliana Costa-Silva，Fabricio Martins Lopes。

写作原始草稿：Juliana Costa-Silva，Fabricio Martins Lopes。

写作评论与编辑：Douglas Domingues，Fabricio Martins Lopes。