|
1简介 序列比对意味着需要排列序列来找出它们相似的位点以及探寻这些位点具有多高的相似度，将读长比对并且比对到参考基因组或者转录组能够使得我们得知这个Read的来源。将读Read比对到基因组可得知其在基因组位置的信息，这样就可以接着用来寻找新的基因和转录本以及量化基因表达。如果你手头没有参考基因组或者仅仅想对已知的转录本进行定量，那么我们可以考虑将这些Read比对至转录组。 将序列比对到参考基因组是个非常具有挑战性的工作，这些序列数量一般可达百万并且长度还非常小，但是相比之下，基因组却是非常之大，并且重复序列和假基因的存在降低了匹配率。此外，序列比对还需要面对基因组变异和测序误差造成的压力。因为真核生物体染色体内还有内含子，所以RNA-Seq中序列比对到基因组上是不连续的，将剪接的序列穿过内含子，正确地定位到外显子-内含子边界是困难的，因为剪接位点的序列信号是确定的，内含子却横跨上千个碱基。 2 比对程序 现如今人们已经开发了数十种比对方式，提供各种各样的程序克服这些挑战。Fonseca等人对市面上存在的比对软件做了比较，分析了这些软件的优缺点并且在网上发布了他们研究的信息。 2.1 Bowtie 由于它的高速度以及低内存需求的特点，Bowtie是已经成为我们现在市面上最常用的比对方式之一，在这里我们来关注一下它的最新版本，Bowtie2，Bowtie2较于Bowtie1更适合较长的读长（从50bp到数千bp）的比对而且还可以运用于存在空位的比对。在早期的版本，Bowtie1对于短读长的比对是非常合适的，但是它不能允许空位。并且Bowtie2自身不能去寻找可变剪接的位点。 Bowtie2实现高速度以及低内存的目标是通过以FM-index的方式对参考基因组建立索引，这里的FM index是以BWT的方法为基础的，这对参考基因组进行了一次有规律的重排生成索引文件。为了加速比对，它通过多个种子序列比对来缩小比对范围。这是最开始的一步，bowtie2比对上的一些种子序列是不允许有空位或是N的碱基，但是用户个人是可以通过设定来控制种子序列的长度以及容错的数量。 Bowtie2含有两种比对模式，分别是局部比对和全局比对。全局比对需要所有的碱基都要进行比对，而局部比对则可以从一端或是双端跳过一些碱基以实现比对的最高得分。Bowtie2默认的方式四end-to-end, 这也是Tophat2所采用的比对模式。在这样的比对模式中，每次出现错配或是漏配都是会被罚分的，最理想的得分情况是0分。错配到高测序质量的碱基要比错配到低测序质量的碱基罚分要高，因为这有可能测序发生了错误。用户个人应该通过是否考虑碱基的质量值然后确定是否考虑罚分，而不是想要设置个人认为的参数来设置罚分多少。最后一点，我们可以挑选任意现成的参数值组合使用：非常快、快速、灵敏（默认）和非常敏感。 在默认情况下，bowtie2可以搜寻多个可能匹配的位点，并且查找出最适合的那个显示，如果存在多个质量一样好的位点，那么将会随机选择一个显示，并会显示其他潜在可能位点的数量。在这里如果比对的质量不行，那么很可能是原始Read的质量不够好。下面Bowtie2使用的实例将会展现如何将Read比对到基因组上。同样Bowtie2也适合于比对到转录组上，这个实例将会在后面进行详尽展示。 构建或下载参考索引 在Bowtie2的官网和iGenomes网站上，许多生物参考基因组的索引是公布的。当你下载index文件时，确保是适用于Bowtie2的,如果是这个index是适用于先前Bowtie的版本，那么这个index是无效的。同时你也可以通过如下的Bowtie2-build的命令来自己建立索引，此外，你还可能使用FASTA格式的基因组索引文件（后面介绍的GTF格式的文件），这也可以在前面介绍的网站中得到。如果下面你要利用HTSeq进行定量，那么选用Ensembl数据库的数据是个很好的选择，因为Ensembl数据库的GTF文件格式与HTSeq需要的格式十分契合，而从iGenomes下载的GTF格式的文件产生的结果适用于下面cuffdiff做差异表达，所以依据你是用的包来确定不同数据库的文件。 在我们的实例中，我们使用从Emsembl数据库下载的FASTA格式的基因组来建立索引： 打开网页（http://www./info/data/ftp/index.html），选择相应的生物以及其对应的DNA文件，在这里我们需要的是把所有染色体信息都包含的dna.toplevel.fa.gz，而不是那些会含有重复标记DNA的dna_rm和dna_sm。这里注意的是这个FASTA格式的文件含有序列拼接的补丁和单倍体序列，但没有染色体组的信息。如果你想利用染色体组的信息，那么你应该分别下载每一个染色体的FSATA格式的文件，在将他们组装成染色体组。 下载: wget ftp://ftp./pub/release-74/fasta/homo_sapiens/dna/Homo_sapiens.GRCh37.74.dna.toplevel.fa.gz 解压: gunzip Homosapiens.GRCh37.74.dna.toplevel.fa.gz 或是用Bowtie2-build建立索引： bowtie2-build –f Homo_sapiens.GRCh37.74.dna.toplevel.fa GRCh37.74 命令最后的这一部分(GRCh37.74)将会是你生成的6个.bt2索引文件的前缀 GRCh37.74.1.bt2 GRCh37.74.2.bt2 GRCh37.74.3.bt2 GRCh37.74.4.bt2 GRCh37.74.rev.1.bt2 GRCh37.74.rev.2.bt2 基因组的FASTA文件在Bowtie2比对中是不需要的，但是如果下面你需要使用Tophat2，那么就应该保留，并且对其进行重命名以便于下面比对到索引文件。 mv Homo_sapiens.GRCh37.74.dna.toplevel.fa GRCh37.74.fa 将读长比对到基因组 Bowtie2可以使用FASTQ以及FASTA类型的文件作为输入文件，下面展示的一个经典的单端测序的Read匹配的案例。如果你的输入文件有好几个，那么输入时用逗号隔开。 bowtie2 -q --phred64 -p 8 --no-unal -x GRCh37.74 -U reads1.fastq.gz -S reads1aligned.sam 在Bowtie2的使用过程中，参数（-x）可以设置bowtie2-build所生成的索引文件的前缀（GRCh37.74），参数（-U）后跟着的是输入文件，参数（-q）设置其输入文件格式为FASTQ（-q），此项为默认值。参数（-phred64）是输入的碱基质量等于ASCII码值-64，在最近的illumina pipiline中参数-phred33得以运用，代表输入的碱基质量等于ASCII码值-33.如果使用参数（-p 8）可设置八个线程，可以加快比对的速度。Bowtie2的使用过程也具备其他功能，可以查询手册来了解每个参数的使用方法。如下图我们可以看到，这里47.42％的Read比对了一次，35.26％的Read比对了一次以上，所以共有82.68％的Read比对上。 34232081 reads; of these: 34232081(100.00%)were unpaired; of these: 5928253(17.32%)aligned 0 times 16232369(47.42%)aligned exactly 1 time 12071459(35.26%)aligned > 1 times 82.68% overall alignment rate 双端测序的Read的比对方法与单端类似，不过就是有两个输入文件分别-1，-2即可。 bowtie2 –q –phred64 –p 8 --no-unal -x GRCh37.74 -1 reads1.fastq.gz -2 reads2.fastq.gz -S paired.sam 两个输入文件必须要以相同的顺序排序才能够进行匹配，如果同时需要比对多个文件，应该将它们用逗号分隔开来，多个参数可以在双端测序的读长比对中使用实现特定的功能，如查询最大Read的片段长度，或是这些映射方向以及距离正确的Read以及没有匹配上的Read是否允许存在。在下面的例子中我们可以发现，39.19％的片段方向以及距离都正确匹配了一次，26.83％匹配了一次以上，共有81.67％的匹配率。 34232081 reads; of these: 34232081(100.00%)were paired; of these: 11633330(33.98%)aligned concordantly 0 times 13415597(39.19%)aligned concordantly exactly 1 time 9183154(26.83%)aligned concordantly >1 times ---- 11633330 pairs aligned concordantly 0 times; of these: 1999775(17.19%)aligned discordantly 1 time ---- 9633555 pairs aligned 0 times concordantly or discordantly; of these: 19267110 mates make up the pairs; of these: 12546751(65.12%)aligned 0 times 4349286(22.57%)aligned exactly 1 time 2371073(12.31%)aligned >1 times 81.67% overall alignment rate Bowtie2的输出文件是SAM 格式的，这是序列比对中默认的输出格式，为了节省空间，SAM可以转变为二进制格式的BAM文件，这种方法会在下面内容进行介绍。 2.2TopHat-预告（下次讲解了） 相对速度快并且占用内存小的TopHat在RNA-Seq中常常用于跨内含子比对，这里我们来关注一下TopHat2，TopHat2的安装是依赖Bowtie2的（当然Bowtie1也是可行的），TopHat2适合于75bp以上Read的比对。TopHat2是一个多步骤比对的程序，如果基因组的注释文件存在的话，那么它会先将Read比对到转录组上。这就大大提高了比对的准确性，还可以避免序列匹配到假基因上，加快了整个比对的进程。TopHat2容错率很低，并不会因为没有比对上就而截短Read从而去比对上，所以低质量的碱基比对的效果可能不是那么好。此外，TopHat2可以察觉基因组的易位，将Read比对到潜在的融合转录子上。 TopHat2比对有三个主要的过程，转录组的比对（步骤一），基因组的比对（步骤二），剪接比对（步骤三到六，如图显示的一样），双端测序的Read首先每端数据要分别单独进行比对，然后考虑比对的片段长度以及方向将单独比对结果合并在一起形成双端比对结果。 |
|
|
