NGS项目一：RNA-Seq数据的Workflow

nature protocol
Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks

从原始的数据开始，进行reads回帖，到拼接转录本，计算表达量，分析差异表达，最后可视化分析结果。

    TopHat是一个把reads回帖到基因组上的工具。首先用Bowtie把reads回帖到基因组上，然后通过拼接，我们就可以在基因组上看到一些reads堆叠起来的区域，称为consensus,这些consensus可能是一个真的外显子，也有可能是几个外显子拼在一起的，或者一些别的情况。我们知道，经典的剪切位点一般都有GT和AG这样的序列标志，在consensus的边界和内部，TopHat会去找这样的剪切位点，并且得到他们可能的组合。然后对于那些没有被Bowtie贴到基因组上的reads,TopHat会对他们建立索引，去和这些可能的剪切位点比对，这样就把跨越剪切位点的reads准确地贴到基因组上。

    一些比较重要的命令行选项。

   

    关于插入片段长度的选项：在RNA-Seq中，会把mRNA打断成小的片段，然后对片段长度进行iding筛选后拿去测序，如果选择的片段长度是300bp，两端各测序75bp的reads,中间的插入片段长度就应该设为150bp.

    下面是设置插入片段长度的标准差，如果选择的片段长度比较集中，这个值可以设置的小一些，反之应该设置得大一些。

 

   

-G选项是提供哦呢一个已有的注释文件。如果你分析的基因组被注释得比较好了，最好能够提供这个文件，
这时TopHat就会先把reads往转录组上贴，没有贴到转录组上的再往基因组上贴，最后把结果合并起来。我们知道大多数的转录组都是比基因组小得多的，而且junctionreads可以直接贴到转录本上，所以这样回帖的效力和准确度都可以

得到提高。 

      标准的Illumina平台是不分链的，我们无法知道配对的reads哪个方向和转录本一致，哪个和转录本反向互补。对于分链的数据，也有两种情况，在firststrand这种分链方法中，第二个read和转录本方向一致，第一个read和转录本反向互补，在另一种fr-secondstrand分链方法中，就刚好反过来了。所以在分析的时候一定要弄清楚自己的数据有没有分链，是怎么分链的。

    下面是一个模拟的RNA-Seq数据集，双端测序，有两种处理，每种处理有3个重复，这里C代表处理，R代表重复，下面用C1R1进行演示

    首先，要有参考序列fasta文件，也就通常说的基因组序列。  TopHat是利用Bowtie2回帖reads,我们首先需要建立Bowtie2的索引文件：bowtie2-build genome.fa(基因组文件) genome  （注意程序和文件所在目录）

    我们还需要reads的fastq文件，双端测序的数据，两个fastq文件分别以下划线1和2这样的形式结尾。在实际分析中，需要对拿到的数据进行质量评估和过滤等依稀类预处理工作，这些工作都是非常重要的。需要准备注释文件，当然它不是必须的。它可以是GTF或者GFF3格式的文件，对于注释得比较好的基因组，在UCSC可以下载。

    准备好后就可以运行tophat了，-p是线程数，-G是注释文件，-o是输出文件夹，选项之后就是参考序列的索引，最后是两个reads的fastq文件。

    看里面生成的文件，align-summary文件爱你，这个文件是reads回帖的一些统计信息。90%以上的回帖比例就非常好了，当然百分之70以上一个可以接受的范围。bam文件详细记录了reads回帖到基因组上的情况，由于这是一个二进制的文件，我们需要用samtools查看它。

    Cufflinks是一套拼接转录本，计算表达量，计算差异表达的工具。尽可能拼接处最优可能的转录本的结构，并且估计它的表达量。

  -G是提供一个注释文件，并且告诉Cufflinks不要去拼接新的转录本，只能用注释文件里提供的转录本。

  -g也是提供一个注释文件，但是Cufflinks会在这些已知转录本的指导下，拼接新的转录本。

-u是告诉Cufflinks用更准确的方法去处理贴到多个位点上的reads，如果没有-u,Cufflinks只会把这些reads简单地平均分配。比如一个read贴到了10个位置，那么每个位置分得十分之一。加-u后会先进行平均分配，然后按照这10个位置各自的表达量，计算read被分配到每个位置的概率。实际上Cufflinks会用EM算法进行迭代，计算在观察到当前数据的情况下，最优可能的reads分配。

library type和TopHat里面差不多。

这里的bam文件就是刚才TopHat运行的结果。

Cuffmerge

当我们使用Cufflinks处理多个数据之后，我们需要将其转录本数据整合为一个全面的转录本集合，Cuffmerge是一个将Cufflinks生成的gtf文件融合为一个更加全面的转录本注释结果的工具。如下图所示，图中的6个转录本被整合为一个转录本集合。同时我们可以利用基因组注释文件，获得更加准确可靠的结果。合并后的转录本集合为计算每个基因和转录本的表达量提供了一个统一的基础。

-g参数 指向参考GTF文件。

    -p参数决定线程数

    -s参考指向基因组DNA序列。如果是一个文件夹，，每个contig是一个fasta文件。如果是一个fasta文件，则所有的contigs都需要在里面。

    最后一项是一个列表，内容包括经过Cufflinks拼接的转录本的文件路径。

流程如下：

首先，我们需要用cat命令创建一个所有拼接的转录本的文件路径列表。

  cat> assemblies.txt

  复制上图的列表

  然后运行cuffmerge,运行后的结果存储在merge_asm这个文件夹里面。其文件夹内包含一个Logs文件夹以及一个.gtf文件，也就是我们经过整合的转录本文件。

cuffdiff

    当我们利用Cufflinks获得了拼接的转录本时，我们就可以计算不同样品中转录本的表达量。计算的简单原理在于测序深度和外显子长度一定时，Read的数量与对应的转录本数量成正比。通过对Reads进行计数计算转录本的表达量。同时cuffdiff可以计算不同条件下转录本表达水平的显著性差异。

    -u命令指cuffdiff对回帖的基因组中多个位置的read进行一个初步的估计，然后加权分配到各个基因组位置。而不是简单的平均分配，其功能与Cufflinks中的u命令相同。

    -L为每个样品标上名称。

    接下来Cuffmerge产生的gtf文件，Cuffdiff需要它提供的注释进行初始转录产物和可变剪切等定量分析。最后是TopHat产生的bam文件，如果一个样品中有多个实验重复，那么我们需要提供哦呢bam文件列表，文件名之间以逗号隔开。   

    运行之后，cuffdiff输出的文件在diff_out目录之下。其中包括一些按类别统计的表达水平结果，如果有相同的转录起始位点，或具有相同的编码区的转录本的表达水平，我么你可以利用他们进行下一步的分析。

CummeRbund

当我们对数据进行分析之后，我们希望将结果可视化，CummeRbund作为一个分析结果数据的软件包，在作图方面具有出色的能力。我们也可以在R环境下输入以下命令进行下载和安装。

常见的绘图方式有：密度分布图，闪点图，火山图，柱形图

这些命令的输入都是cuffdiff的输出文件，经过处理后的结果，我们会在之后进行详细的介绍。这是为不同条件的样本标记名称。
以上参考：
http://blog.sina.com.cn/s/blog_670445240101jl08.html
 以下自行总结：

TopHat把readsaligns to genome and发现剪接位点

Cufflinks将这些readsmap到基因组上组装转录本,是分析多种情况下的reads基因表达差异工具。它定量每一个样本的表达梁，这个过程使用非常严格的统计学算法来过滤背景值与认为片段。Cufflinks可以通过参考基因组的注释而不是组装需类来定量转录组的值。对于多个样本，推荐使用Cufflinks来定量。在做RNA-seq差异表达分析之前，很有必要先导入数据在对数据进行组装。因此一种简单的办法就是将所有样本对比数据导入，然后再用Cufflinks组装。但一般不支持该策略。首先，随着序列测序深度提高，组装所需的计算时间增加。第二，当把所有样本数据导入时，很多基因剪切类型变得复杂，Cufflinks在组装时出现的亚型变多，导致组装变得复杂。推荐的方法师将每一个样本独自进行组装，再把结果融合在一起，即先使用cufflinks,再使用cuffmerge。

cuffmerge能够将Cufflinks组装的文件与参考注释文件进行比较，并帮助从已知基因中分类出新的基因，但由于组装过程中存在难度与错误，因此当发现新的基因时最后要辅助传统实验手段进行验证。Cuffdiff是计算两个或多个样本的表达量与用统计检验的方法评估每一个观察到的变化的可信度。统计模型是基于每个转录本序列的数目与它的丰度成比例，但由于测序与文库准备的技术，生物可变性会存在差异。