0%

Samtools常用命令的总结

  1. flags

    1 0x1 这序列是PE双端测序
    2 0x2 这序列和参考序列完全匹配,没有错配和缺失
    4 0x4 这序列没有mapping到参考序列上
    8 0x8 这序列的mate序列没有mapping到参考序列上
    16 0x10 这序列比对到参考序列的负链上
    32 0x20 这序列的mate序列比对到参考序列的负链上
    64 0x40 这序列是read1
    128 0x80 这序列是read2
    256 0x100 这序列不是主要的比对,因为序列可能比对到参考序列的多个位置上
    512 0x200 这序列没有通过QC
    1024 0x400 这序列是PCR重复序列
    2048 0x800 这序列是补充比对

  2. view

     samtools view [options] in.sam|in.bam|in.cram [region...]

    -f 提取 ## -f 4 提取出没有mapping上的reads
    -F 过滤 ## -F 4 过滤掉没有mapping上的reads,也就是说提取出mapping上的reads
    -u 输出格式为未压缩的bam格式
    -q 过滤掉MAPQ值低某个阈值 ## -q 1 过滤掉MAPQ值低于1的情况
    -h 设定输出的SAM文件带有header
    -b 输出格式设定为BAM
    -S 输入格式为SAM

    提取比对到参考序列的结果:

     samtools view -bF 4 tmp.bam > tmp_F.bam

    提取双端序列都比对到参考序列(4+8)的结果:

     samtools view -bF 12 tmp.bam > tmp_F.bam

    提取比对到chr1的结果

     samtools view -b tmp.bam chr1 > tmp_chr1.bam

    注:With no options or regions specified, prints all alignments in the specified input alignment file (in SAM, BAM, or CRAM format) to standard output in SAM format (with no header),也就是说,没有设定输出格式的话,默认是输出SAM格式,并且是没有header的SAM

  3. index

     samtools index [-bc] [-m INT] aln.bam|aln.cram [out.index]

    -b 创建一个bai索引,默认设定这个参数(如果在命令中没这个参数)

    建索引(必须是已经使用默认排序后的):

     samtools index tmp.bam
  4. sort

     samtools sort [-l level] [-m maxMem] [-o out.bam] [-O format] [-n] [-t tag] [-T tmpprefix] [-@ threads] [in.sam|in.bam|in.cram]

    -m 设置内存使用大小,默认是500,000,000(现在支持K,M,G等缩写)
    -n Sort by read names (i.e., the QNAME field) rather than by chromosomal coordinates(似乎一般也是使用默认排序,即Sort alignments by leftmost coordinates,因为index需要默认排序...) -@ 设置线程数
    -O 输出的格式(sam,bam,cram),默认是bam

    使用默认排序:

     sort -@ 5 tmp.bam >tmp.sorted.bam
  5. merge

     samtools merge [-nur1f] [-h inh.sam] [-R reg] [-b <list>] <out.bam> <in1.bam> [<in2.bam> <in3.bam> ... <inN.bam>]

    -b 一个bam文件一行的一个bam list文件
    -n The input alignments are sorted by read names rather than by chromosomal coordinates
    -h Use the lines of FILE as `@' headers to be copied to out.bam, replacing any header lines that would otherwise be copied from in1.bam. (FILE is actually in SAM format, though any alignment records it may contain are ignored.)
    -c 当多个输入文件包含相同的@RG头ID时,只保留第一个到合并后输出的文件。当合并多个相同样本的不同文件时,非常有用
    -p 与-c参数类似,对于要合并的每一个文件中的@PG ID只保留第一个文件中的@PG

    merge前必须是已经sort的文件,如果只是单纯的merge:

     samtools merge tmp1.bam tmp2.bam
  6. mpileup

     samtools mpileup [-EBugp] [-C capQcoef] [-r reg] [-f in.fa] [-l list] [-Q minBaseQ] [-q minMapQ] in.bam [in2.bam [...]]

    从官方说明:Generate VCF, BCF or pileup for one or multiple BAM files可看出,可以用来和bcftools搭配Call SNPs

    最常用的几个参数:
    >-f The faidx-indexed reference file in the FASTA format(有索引(faidx)文件的参考序列)
    >-l BED or position list file containing a list of regions or sites where pileup or BCF should be generated(bed格式的文件,如果需要只处理特定位点的bam文件的话)
    >-r Only generate pileup in region 搭配-l使用,比如可以指定染色体
    >-g Compute genotype likelihoods and output them in the binary call format (BCF).(输出bcf格式文件)
    >-u Generate uncompressed VCF/BCF output(如果后面接管道符的话,必须使用这个指定不进行压缩)

    搭配bcftools使用:

     samtools mpileup -ugf <ref.fa> <sample1.bam>| bcftools call -vmO z -o <study.vcf.gz>
  7. tview

     samtools tview [-p chr:pos] [-s STR] [-d display] <in.sorted.bam> [ref.fasta]

    显示reads比对到基因组的情况,类似于基因组浏览器

  8. faidx

     samtools faidx <ref.fasta> [region1 [...]]

    给参考序列建索引,或者从已建索引的参考序列中提取一定位置范围内的序列

  9. depth

     samtools depth [options] [in1.sam|in1.bam|in1.cram [in2.sam|in2.bam|in2.cram] [...]]

    计算bam/sam文件每个位点的测序深度

  10. flagstat

    samtools flagstat in.sam|in.bam|in.cram

    统计bam文件中reads的比对情况,如多少reads比对上等信息

samtools官网手册还介绍了其他好多的功能,具体可参见:
http://www.htslib.org/doc/samtools.html

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