比对软件STAR的简单使用

第一次听说START这款比对软件是因为其是ENCODE计划的御用软件,ENCODE计划(ENCyclopedia Of DNA Elements)又称人类基因组DNA元件百科全书计划,是2003年在人类基因组计划完成之后紧接着的又一个大型国际科研项目。

第二次听说则的由于Gaining comprehensive biological insight into the transcriptome by performing a broad-spectrum RNA-seq analysis这篇发表于2017年的文章,主要是针对转录组各个分析流程的不同分析工具的比较,里面针对mRNA的比对方法总结了基于参考基因组的三款比对软件:TopHat,STAR和HASAT2。其中讲到STAT相比较其他两款软件有较高的唯一比对率;STAR会将没有paired mapping上的reads都剔除,避免single reads比对到基因组上;并且STAR对lower-quality(包括more soft-clipped和错配碱基)比对有较高的容忍度

第三次听说也是由于恰好需要使用GATK对RNA-Seq Call Variants,因而在GATK刚好查到一篇教程Calling variants in RNAseq GATK 将reads比对至Reference上是采用STAR的STAR 2-pass模式,所以为了学习该教程,必须先学习如何使用STAR了

STAR的下载及安装

先进入START官网https://github.com/alexdobin/STAR

下载STAR,无须编译即可使用

wget https://github.com/alexdobin/STAR/archive/2.5.3a.tar.gz
tar -xzf 2.5.3a.tar.gz
cd STAR-2.5.3a

STAR的使用

  1. 作为一款比对软件,建index肯定是必不可少的一步

    STAR --runThreadN 6 --runMode genomeGenerate \
    --genomeDir ~/reference/index/STAR/mm10/ \
    --genomeFastaFiles ~/reference/genome/mm10/GRCm38.p5.genome.fa \ 
    --sjdbGTFfile ~/annotation/mm10/gencode.vM13.annotation.gtf \
    --sjdbOverhang 100
    

    这个命令参数也很好理解:
    --runThreadN :设置线程数
    --genomeDir :index输出的路径
    --genomeFastaFiles :参考基因组序列
    --sjdbGTFfile :参考基因组注释文件
    --sjdbOverhang :这个是reads长度的最大值减1,默认是100

  2. 然后进行比对

    STAR --runThreadN 20 --genomeDir ~/reference/index/STAR/mm10/ \
    --readFilesIn SRR3589959_1.fastq SRR3589959_2.fastq \
    --outSAMtype BAM SortedByCoordinate \
    --outFileNamePrefix ./SRR3589959
    

    --readFilesIn :paired reads文件
    --outSAMtype :表示输出默认排序的bam文件,类似于samtools sort(还有–outSAMtype BAM Unsorted和–outSAMtype BAM Unsorted SortedByCoordinate)
    --outFileNamePrefix :输出文件路径即前缀

  3. 结果文件:

    SRR3589959Aligned.sortedByCoord.out.bam
    SRR3589959Log.final.out
    SRR3589959Log.out
    SRR3589959Log.progress.out
    SRR3589959SJ.out.tab

    可以通过samtools view SRR3589959Aligned.sortedByCoord.out.bam |less -S来查看对应文件的每列信息

    前面12列一般也是规范的sam格式,最后一列attributes信息的话,STAR默认是输出NH HI AS nM attributes,这里需要注意的是HI,其表示多重比对的reads的起始位置,默认是以1开始算的,但是如果下游分析需要用到Cufflinks or StringTie的话,需要用–outSAMattrIHstart设置为0比对软件STAR的使用—高通量测序数据处理学习记录(一)

    SRR3589959SJ.out.tab则是Splice junctions的一些信息,其中需要注意的是:对于junction的位置信息,STAR则是按照intron的起始和终止位置来定,而其他的一些软件则是按照exon的位置来决定的;至于每列代表的含义可以看mannul,很好理解

STAR 2-pass mode

为了发现更加灵敏的new junction,STAR建议使用2-pass mode,其能增加检测到的new junction数目,使得更多的splices reads能mapping到new junction。因此STAR先用一般参数做一遍mapping,收集检测到的junction信息,然后利用这已经annotated junction来做第二次mapping

  1. STAR对于2-pass mode有新旧两种方式,比如original 2-pass 方法:

    首先做一遍常规的比对,结果中会生成一个SJ.out.tab文件,如上面所提到的SRR3589959SJ.out.tab。然后用--sjdbFileChrStartEnd参数将所有样品的SJ.out.tab文件作为输入的annotated junction进行第二次建index

    STAR --runThreadN 20 --runMode genomeGenerate 
    --genomeDir  ~/reference/index/STAR/mm10/index_2-pass/ \
    --genomeFastaFiles ~/reference/genome/mm10/GRCm38.p5.genome.fa \ 
    --sjdbGTFfile ~/annotation/mm10/gencode.vM13.annotation.gtf \
    --sjdbFileChrStartEnd SRR3589959SJ.out.tab SRR3589960SJ.out.tab SRR3589961SJ.out.tab SRR3589962SJ.out.tab \
    --sjdbOverhang 100
    

    然后用第二次建立的index再一次对每个样品进行STAR比对,以SRR3589959为例

    STAR --runThreadN 20 --genomeDir ~/reference/index/STAR/mm10/index_2-pass/ \
    --readFilesIn SRR3589959_1.fastq SRR3589959_2.fastq \
    --outSAMtype BAM SortedByCoordinate \
    --outFileNamePrefix ./SRR3589959_2-pass
    
  2. 上述方法original方法适用于多样本和单个样本的处理,但是如果是per-sample(单个样本?)的2-pass mapping,可以直接用--twopassMode Basic参数将第两步mapping中的make index省去,直接再mapping

    STAR --runThreadN 20 --genomeDir ~/reference/index/STAR/mm10/ \
    --twopassMode Basic \
    --readFilesIn SRR3589959_1.fastq SRR3589959_2.fastq \
    --outSAMtype BAM SortedByCoordinate \
    --outFileNamePrefix ./SRR3589959
    

    这个比常规的结果还多2个临时产生的文件夹(SRR3589959_STARgenome,SRR3589959_STARpass1)

    至于bam文件则是跟上述的original 2-pass

STAR还有其他一些不太常用的参数,可以参看manual,Download后即可查看

本文出自于http://www.bioinfo-scrounger.com转载请注明出处