Samtools中文说明书

samtools-1.2中文使用手册,英文原文发布于2015年12月15日。 本翻译为2016年4月20日第一稿,如有建议或疑问请联系 roden.luo [AT] gmail.com。

samtools-1.3中文使用手册正在筹备翻译中,欢迎参与。

NAME – 名字

samtools – 序列比对/图谱(SAM)格式的工具箱

SYNOPSIS – 概要

samtools view -bt ref_list.txt -o aln.bam aln.sam.gz

samtools sort -T /tmp/aln.sorted -o aln.sorted.bam aln.bam

samtools index aln.sorted.bam

samtools idxstats aln.sorted.bam

samtools view aln.sorted.bam chr2:20,100,000-20,200,000

samtools merge out.bam in1.bam in2.bam in3.bam

samtools faidx ref.fasta

samtools fixmate in.namesorted.sam out.bam

samtools mpileup -C50 -gf ref.fasta -r chr3:1,000-2,000 in1.bam in2.bam

samtools tview aln.sorted.bam ref.fasta

samtools flags PAIRED,UNMAP,MUNMAP

samtools bam2fq input.bam > output.fastq

DESCRIPTION – 描述

Samtools是一个用来处理BAM格式(SAM的二进制格式,译者注)的比对文件的工具箱。它能够输入和输出SAM(序列比对/图谱)格式的文件,对其进行排序、合并、建立索引等处理,还可以快速获取任意区域的reads。

Samtools被设计成以数据流的形式工作。它把标准输入作为输入(stdin)文件 `-‘ ,把标准输出作为输出(stdout)文件 `-‘ 。因此多个命令可以通过Unix的管道技术结合起来。Samtools总是把警告和错误消息输出到标准错误输出(stderr)。

Samtools也能够打开一个远程FTP或HTTP服务器上的BAM(SAM不行)文件,但是这个BAM文件名要以 `ftp://’ 或 `http://’ 开头。Samtools会在当前工作目录下检查相应的索引文件,如果不存在的话,它会把它下载下来。Samtools不会获取整个比对文件,除非你要求它这么做。

COMMANDS AND OPTIONS – 命令和参数

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

如果没有指定参数或者区域,这条命令会以SAM格式(不含头文件)打印输入文件(SAM,BAM或CRAM格式)里的所有比对到标准输出。

你可以在输入文件的文件名后面指定一个或多个以空格分隔的区域来限制输出,这样只会输出落在指定区域内的比对。要指定区域,需要一个对坐标排过序的并且建立了索引的输入文件(BAM或CRAM格式)。

-b, -C, -1, -u, -h, -H, 和 -c 参数会改变默认的没有头文件的SAM输出格式。 -o -U 参数可以设置输出文件的名字。

-t-T 参数用来提供额外的参考数据,当输入的SAM没有@SQ头文件的时候,必须需要提供这两个中的一个参数。无论何时,只要是输出CRAM格式,就必须提供 -T 参数。

-L, -r, -R, -q, -l, -m, -f, 和 -F 参数用来过滤比对,只有那些满足特定标准的比对会被输出。

-x, -B, 和 -s 参数用来修饰每个比对中的数据。

最后 -@ 参数可用来为压缩分配额外的线程,而 -? 参数可以得到一个很长的帮助信息。

REGIONS – 区域:
区域可以用如下方法制定:RNAME[:STARTPOS[-ENDPOS]] ,即:染色体名[:起始位点[-终止位点]],([ ] 中的内容表示为可选内容,可以指定,也可以不指定。译者注),并且,所有的坐标系统都是以1标记第一个碱基的(1-based)。

要非常注意的是:如果给定了多个区域,有些比对可能会被输出多次,因为它们和多个区域都匹配。

指定区域的例子:

`chr1′
把所有比对到参考序列中名为`chr1’(也就是头文件里标明的 @SQ SN:chr1)的比对全部输出。(chr一般为染色体的缩写,译者注)
`chr2:1000000′
chr2上从第1,000,000碱基开始的序列到这条染色体的结尾的区域。
`chr3:1000-2000′
chr3上从第1,000个碱基到第2,000个碱基之间的共1001bp的区域(包含两个端点)。

OPTIONS – 参数:

-b
以BAM格式输出。
-C
以CRAM格式输出(要求有 -T 参数)。
-1
启用快速BAM压缩(隐含 -b)。
-u
输出未压缩的BAM。这个参数省下了压缩/解压缩的时间,因此更适合在运用管道技术把另一个samtools命令结合在一起时使用。
-h
在输出中包含头文件。
-H
只输出头文件。
-c
与输出比对不同,只对它们计数,并且输出总数。所有的过滤参数,比如 -f, -F, 和 -q, 都可以被考虑进来。
-?
输出长篇帮助文档并且立即退出。
-o FILE
输出到 FILE [stdout] ([ ]中表示默认值,这里表示如果使用了-o参数,但没有提供文件名,输出到标准输出。译者注)。
-U FILE
把那些没有被各种过滤参数所选中的比对,输出到 FILE里面。 当这个参数被使用是,所有的比对(如果指定了区域,则是所有和指定区域相交的比对)要么被写进输出文件,要么被写进这个文件,但绝不会同时被写进两个文件。
-t FILE
一个以制表符(tab)分隔的文件 FILE。 每一行必须含有参考序列的名字,放在第一列,和这个参考序列的长度,放在第二列,每一行代表一个不同的参考序列。第二列后面的任何字段都会被忽略。这个文件同时也定义了在排序时的参考序列的顺序。如果你运行`samtools faidx <ref.fa>’,那么它输出的索引文件 <ref.fa>.fai 可以被用作这里的 FILE
-T FILE
一个FASTA格式的参考序列文件 FILE,可以是被 bgzip 压缩过的,并且最好是被 samtools faidx建立过索引的。 如果索引文件不存在,我们会为你生成一个。
-L FILE
只输出那些和输入的BEDFILE重合的比对 [null](默认为空)。
-r STR
只输出 STR 这个read group里面的比对 [null]。
-R FILE
输出文件 FILE 里面列出来的read groups里的比对 [null]。
-q INT
当比对的MAPQ小于 INT 时,忽略这个比对 [0]。
-l STR
只输出 STR 这个文库(library)里面的比对 [null]。
-m INT
Only output alignments with number of CIGAR bases consuming query sequence ≥ INT [0]
-f INT
Only output alignments with all bits set in INT present in the FLAG field. INT can be specified in hex by beginning with `0x’ (i.e. /^0x[0-9A-F]+/) or in octal by beginning with `0′ (i.e. /^0[0-7]+/) [0].
-F INT
Do not output alignments with any bits set in INT present in the FLAG field. INT can be specified in hex by beginning with `0x’ (i.e. /^0x[0-9A-F]+/) or in octal by beginning with `0′ (i.e. /^0[0-7]+/) [0].
-x STR
Read tag to exclude from output (repeatable) [null]
-B
Collapse the backward CIGAR operation.
-s FLOAT
Integer part is used to seed the random number generator [0]. Part after the decimal point sets the fraction of templates/pairs to subsample [no subsampling].
-@ INT
进行BAM压缩操作时,除了主线程之外,额外分配的线程数 [0]。
-S
这个参数会被忽略,它只是为了和旧版本的samtools相兼容。之前,如果输入文件是SAM格式的话,必须添加这个参数,但是现在,软件会侦测输入的前几个字符,自动确定正确的格式。
tview
samtools tview [-p chr:pos] [-s STR] [-d display] <in.sorted.bam> [ref.fasta]

文字版的reads和参考序列的比对情况的查看器(基于 ncurses library )。在这个查看器中,按下 `?’ 来寻求帮助,按下 `g’ 可以从某个区域开始查看比对,这个区域用如下格式来表示:`chr10:10,000,000’,或者,当继续查看相同的参考序列时,可以使用 `=10,000,000’。

Options:

-d display
以Html格式(H)、Curses格式(C)或文本格式(T)输出。
-p chr:pos
直接到达这个位置。
-s STR
只显示 STR 这个样本或read group的比对。
mpileup
samtools mpileup [-EBugp] [-C capQcoef] [-r reg] [-f in.fa] [-l list] [-Q minBaseQ] [-q minMapQ] in.bam [in2.bam []]

对一个或多个BAM文件,生成VCF、BCF或者pileup(堆积文件)。比对记录会按照头文件中@RG行里面样本标识符(SM)来分组。如果样本标识符不存在,那么每一个输入文件都被视作一个独立的样本。

在(没有使用 -u -g 参数的)pileup格式中,每一行代表基因组上的一个位置,它由染色体名字,1-based坐标(该染色体的第一个碱基标为1),参考序列上的碱基,覆盖这个位点的reads数目,read上的碱基,碱基的质量和比对的质量构成。下列信息:匹配、不匹配、插入删除、正反向链、比对质量、是否是某read的起始或终止碱基,全部都编码在read上的碱基这里面了。在这一列中,一个点代表和参考序列正向链相同的一个匹配,一个逗号代表和反向链相同的一个匹配。一个 ‘>’ 或 ‘<‘ 代表a reference skip(抱歉,译者还不确定含义,不知该如何翻译),`ACGTN’代表正向链上的一个mismatch(不匹配),`acgtn’代表反向链上的一个mismatch。这样的一个样式,`\\+[0-9]+[ACGTNacgtn]+’,代表了在参考序列上这个位置和下个位置之间的一个insertion(插入)。Insertion的长度由样式中的数字给出,它后面是插入的序列。同样的,这个样式,`-[0-9]+[ACGTNacgtn]+’,代表一个deletion(删除)。在之后的行中,被删除的碱基会用 `*’ 来代表。同样在这一列中,一个 `^’ 符号表示是某read的起始点。 `^’的ASCII码减去33给出了比对质量。一个 `$’ 表示是某read的终止点。

Input Options – 输入参数:

-6, –illumina1.3+
假设质量是用Illumina 1.3+编码的话。
-A, –count-orphans
在识别变异时,不要跳过反常的read。
-b, –bam-list FILE
输入的BAM文件的地址的列表,一个文件一行。 [null]
-B, –no-BAQ
Disable probabilistic realignment for the computation of base alignment quality (BAQ). BAQ is the Phred-scaled probability of a read base being misaligned. Applying this option greatly helps to reduce false SNPs caused by misalignments.
-C, –adjust-MQ INT
Coefficient for downgrading mapping quality for reads containing excessive mismatches. Given a read with a phred-scaled probability q of being generated from the mapped position, the new mapping quality is about sqrt((INT-q)/INT)*INT. A zero value disables this functionality; if enabled, the recommended value for BWA is 50. [0]
-d, –max-depth INT
在一个位置上,对每个BAM输入文件,最多读取 INT 个reads。 [250]
-E, –redo-BAQ
在运行时重新计算BAQ,忽略存在的BQ标签。
-f, –fasta-ref FILE
faidx建立了索引的FASTA格式的参考序列文件,这个文件也可以是被bgzip压缩过的。 [null]
-G, –exclude-RG FILE
排除列在文件FILE里面的readgroups里面的reads。(每一个@RG-ID占一行)
-l, –positions FILE
BED或者是列了位置坐标的文件,这个文件含有一列区间或位点,这些地方是生产pileup或BCF的位置。如果是BED文件,坐标系统是0-based half-open(第一个碱基是0,区间是半开的)的。[null]
-q, -min-MQ INT
最小的比对质量,某个比对的质量超过该值时才会被使用。[0]
-Q, –min-BQ INT
最小的碱基质量,某个碱基的质量超过该值时才会被使用。[13]
-r, –region STR
只产生这个区域的pileup。需要建立了索引的BAM文件。如果和 -l 同时使用的话,会考虑两个要求的交集。STR[all sites]
-R, –ignore-RG
忽略RG标签,把一个BAM文件中所有的reads都视为一个样本。
–rf, –incl-flags STR|INT
需要的flag的标签:不在这里的mask bits(flag的值)的reads会被忽略 [null]
–ff, –excl-flags STR|INT
过滤flags:忽略有这些mask bits标记的reads [UNMAP,SECONDARY,QCFAIL,DUP]
-x, –ignore-overlaps
取消read-pair(双端测序的两个reads)的重叠的检查。

Output Options – 输出参数:

-o, –output FILE
把pileup或者VCF/BCF输出到FILE,而不是默认的标准输出里面。

(这个短参数同时被 –open-prob–output 使用。如果-o的参数含有任何非数字的字符,而不是由 + 或 – 符号开头,它会被视作是–output。通常文件的扩展名会考虑到这一点,但是如果要写入到一个完全由数字组成的文件名,使用 -o ./123–output 123 。)

-g, –BCF
计算基因型的概率,并输出到binary call format(BCF,二进制识别格式)。在v1.0版本的samtools中,它是BCF2,和BCF1格式是不兼容的,而以前的版本(0.1.x)输出的是BCF1。
-v, –VCF
计算基因型的概率,并输出到variant call format(VCF,变异识别格式)。输出是bgzip压缩后的VCF,除非设定了-u参数。

Output Options for mpileup format (without -g or -v) – mpileup格式的输出参数(没有 -g 或 -v):

-O, –output-BP
输出reads里面碱基的位置。
-s, –output-MQ
输出mapping的质量。

Output Options for VCF/BCF format (with -g or -v) – VCF/BCF格式的输出参数(有 -g 或 -v):

-D
输出每个样本的read depth(深度)[DEPRECATED – use -t DP instead](不再宜用,使用-t DP来替代它)
-S
输出每个样本的Phred标准的的链的bias(偏差)的P-value(不再宜用,使用-t SP来替代它)
-t, –output-tags LIST
输出逗号分隔的 FORMAT 和 INFO 标签的列表(大小写不敏感): DP (高质量碱基的数目,FORMAT), DV (高质量的非参考序列的碱基数,FORMAT), DPR (支持每个观察到的allele(等位基因)的高质量的碱基数,FORMAT), INFO/DPR (支持每个观察到的allele(等位基因)的高质量的碱基数,INFO), DP4 (高质量的 ref-forward, ref-reverse, alt-forward and alt-reverse (正向-参考序列,反向-参考序列,正向-突变序列,反向-突变序列)的碱基数目,FORMAT), SP (Phred标准的的链的bias(偏差)的P-value,FORMAT) [null]
-u, –uncompressed
输出未被压缩的VCF/BCF,在使用unix管道(pipe)时推荐这个参数。
-V
输出每个样本的非参考序列的reads数 (不再宜用,使用-t DV来替代它)

Options for SNP/INDEL Genotype Likelihood Computation (for -g or -v): 对SNP/INDEL基因型的概率计算的参数 (for -g or -v):

-e, –ext-prob INT
Phred标准的gap延伸的测序错误概率。降低INT值会产生更长的indels。 [20]
-F, –gap-frac FLOAT
最小的有gap的reads的比例。 [0.002]
-h, –tandem-qual INT
对同聚物错误建模的系数。考虑一个长度为 l 的同聚物而言,大小为 s 的indel的测序错误是如下建模的: INT * s / l. [100]
-I, –skip-indels
不进行INDEL的识别。
-L, –max-idepth INT
当每个样本的平均深度大于 INT 时,略过INDEL的识别。 [250]
-m, –min-ireads INT
识别候选indel时,最小的存在gap的reads数目 INT。 [1]
-o, –open-prob INT
由测序错误引入的gap的Phred-标准的概率。减少 INT 值会识别更多的indel。 [40]

(这个短参数同时被 –open-prob–output 使用。当 -o 的参数是可选的 + 或 – 号,跟随着数字0-9的话,它会被视作 –open-prob。)

-p, –per-sample-mF
Apply -m and -F thresholds per sample to increase sensitivity of calling. By default both options are applied to reads pooled from all samples.
-P, –platforms STR
逗号分隔的一列仪器平台(由@RG-PL决定的),候选indel会从这些reads中获取。我们建议从低indel错误率的测序平台,比如ILLUMINA,中去收集候选indel。 [all]
reheader
samtools reheader <in.header.sam> <in.bam>

in.header.sam 里面的头文件替换 in.bam 里的。这个命令比通过BAM→SAM→BAM转换来替换头文件要快很多。

cat
samtools cat [-h header.sam] [-o out.bam] <in1.bam> <in2.bam> [ … ]

串连 BAM 文件。每一个输入的BAM里面的 sequence dictionary 必须相同,虽然这个命令并不检查这一点。这个命令使用了和 reheader 类似的小技巧来实现了快速的BAM串连。

sort
samtools sort [-l level] [-m maxMem] [-o out.bam] [-O format] [-n] -T out.prefix [-@ threads] [in.bam]

对序列进行排序。默认通过最左边的坐标排序;当使用 -n 时,通过read的名字排序。 一个相应的 @HD-SO 排序顺序标签会被添加在头文件里面,或者如果必要的话会更新已经存在的标签。

排好序后默认输出到标准输出,如果使用了 -o 的话,就输出到指定的文件里面 (out.bam) 。当整个的序列数据不能放入内存(可以通过 -m 参数来控制)的时候,这个命令也会产生一些临时文件 out.prefix.%d.bam

Options:

-l INT
为最终的输出文件设置想要的压缩程度,从0(不压缩)或1(最快但是压缩得最少)到9(最大程度压缩但是写得最慢),和 gzip(1) 的压缩程度设置相似。

如果 -l 没有被使用,会使用默认的压缩程度。

-m INT
每个线程的最大内存,以 bytes 为单位设定或者添加一个 K, M, 或 G 后缀。 [768 MiB]
-n
根据read的名字(也就是 QNAME 那一列)来排序,而不是染色体坐标。
-o FILE
最后排好序后输出到 FILE 中,而不是标准输出。
-O FORMAT
最终输出 sam, bam, 或 cram

默认地,samtools 会基于 -o 文件的扩展名尝试选择一个格式;如果是到标准输出或者格式不能被推测出, -O 必须要被指定。

-T PREFIX
把临时文件写到 PREFIX.nnnn.bam. 里面。这个参数是必须的。
-@ INT
设置排序和压缩的线程数,该操作默认是单线程的。

为了和之前的脚本兼容,samtools sort 也支持以前的不太灵活的方式去指定临时的和最终的输出文件名:

samtools sort [-nof] [-m maxMem] in.bam out.prefix

排好序的BAM会被输出到 out.prefix.bam中,(或者由下面提到的 -o-f 决定),而且,任何临时文件都会被写进out.prefix.%d.bam

-o
把最后比对的结果输出到标准输出。
-f
使用 out.prefix 作为输出的全路径,并且不要添加 .bam 后缀。
-l, -m, -n, -@
以上面提到的相同的含义被接收。

最终这个会被删除的,您应该改去使用新的更灵活的调用方法。

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

合并多个排好序的比对文件,产生一个排好序的输出文件,这个文件含有所有的输入记录,并且保留了他们原来的顺序。

如果 -h 被指定了的话,输入文件里的@SQ头文件将会被合并到指定的头中,否则的话,他们会被合并成由输入的头文件组成的复合的头中。如果在合并@SQ行的过程中,这些行体现了输入文件排好序的坐标,而在序列当中又出现了一个矛盾(比如说,input1.bam 有一个@SQ是 a,b,c 而input2.bam是 b,a,c),那么输出的结果文件需要被重新按照坐标序列排序。

除非指定了 -c-p 标志,在合并 @RG 和 @PG 记录到输出的头文件的时候,任何和现有 IDs 重复的 IDs 将会连接一个后缀把它们和其他头文件中相似的记录区分开来,并且read记录也会被更新来反映这一点。

OPTIONS:

-1
使用 zlib 压缩等级1来压缩输出文件。
-b FILE
输入 BAM 文件的列表,每行一个文件。
-f
如果输出文件存在的话,强制覆盖它。
-h FILE
FILE 的`@’头复制到 out.bam 中,如果不指定该参数的话,会复制 out.bam 的头。( FILE 实际是 SAM 格式的,但是它里面可能存在的任何比对记录都会被忽略。)
-n
使输入的比对序列按照read的名字来排序,而不是染色体坐标。
-R STR
在由 STR 指定的特定区域中合并文件。 [null]
-r
给每个比对添加一个 RG 标签,标签的值是通过文件名来推测的。
-u
输出未被压缩的 BAM 。
-c
把 RG 标签和 colliding IDs 结合,而不是添加一个后缀去区分它们。
-p
把 PG 标签和 colliding IDs 结合,而不是添加一个后缀去区分它们。
index
samtools index [-bc] [-m INT] aln.bam|aln.cram

为按坐标排好序的 BAM 或 CRAM 文件,用以快速随机访问。下列情况是需要这个索引的:用 region 参数限制了 samtools view , 和类似的针对特定感兴趣区域的命令。

对于一个 CRAM 文件, aln.cram ,名为 aln.cram.crai 的索引文件将会被创建。 对于一个 BAM 文件, aln.bam , 名为 aln.bam.baialn.bam.csi 的索引文件中的一个将会被创建,这由选择的索引格式决定。

Options:

-b
创建一个 BAI 索引,这是目前没用指定格式参数时的默认值。
-c
创建一个 CSI 索引。默认情况下,索引的最小间隔大小是 2^14 ,这和 BAI 使用的的固定值是一样的。
-m INT
创建一个 CSI 索引,使用 2^INT 作为最小间隔大小。
idxstats
samtools idxstats <aln.bam>

获取并且打印索引文件的统计信息。输出文件是制表符分隔的,每一行含有参考序列名,序列长度,比对上去的read数和没有比对上去的read数。

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

对FASTA格式的参考序列建立索引,或者从已建立索引的参考序列中获取子序列。 如果没有指定区间, faidx 会为这个文件在硬盘上创建索引文件 <ref.fasta>.fai 。如果指定了区间,相应的子序列会被获取并以FASTA格式打印到标准输出。输入文件可以是以 BGZF 格式的压缩文件。

fixmate
samtools fixmate [-rpc] [-O format] in.nameSrt.bam out.bam

从一个按照名字排序后的比对文件中,添加mate的坐标, ISIZE ,和mate相关的flag。

OPTIONS:

-r
删除次要的和没有比对上的reads。
-p
不使用 FR 匹配检查。
-c
添加模板 cigar ct 标签。
-O FORMAT
Write the final output as sam, bam, or cram.

默认地, samtools 会尝试基于输出文件的扩展名选择一个格式;如果是输出到标准输出或者格式不能被推出,就必须要指定 -O 参数。

rmdup
samtools rmdup [-sS] <input.srt.bam> <out.bam>

去除可能的由PCR造成的重复序列:如果多个 read pair 有相同的坐标,只保留比对质量最高的那一对。在双端模式下,这个命令用于FR方向并且要求正确设置 ISIZE 。它不处理不配对的 reads (比如两端比对到不同的染色体上了或者孤reads(pair中只有一个read比对到基因组了))。

OPTIONS:

-s
去除重复的单端测序reads。默认地,这个命令只对双端测序reads工作。
-S
处理双端测序reads和单端测序reads。
calmd
samtools calmd [-EeubSr] [-C capQcoef] <aln.bam> <ref.fasta>

产生 MD 标签。如果 MD 标签已经存在,并且新产生的MD标签和已有的不一样的话,这个命令会发出一个警告。默认输出成 SAM 。

OPTIONS:

-A
当和 -r 一起使用时,这个选项会重写原始的碱基质量。
-e
如果某个碱基和它比对上的参考序列是一致的,就把它转换成 = 。现在的插入删除的识别软件并不支持 = 的碱基。
-u
输出没有压缩的 BAM 文件。
-b
输出压缩的 BAM 文件。
-S
输入文件是带头文件的 SAM 。
-C INT
用来对低比对质量的reads加帽的系数。更多细节请查看 pileup 命令。 [0]
-r
计算 BQ 标签(没有 -A) 或者对碱基质量加帽(有 -A)。
-E
扩展版的 BAQ 计算,这个参数削弱了特异性而增加了敏感性,虽然只有小的影响。
targetcut
samtools targetcut [-Q minBaseQ] [-i inPenalty] [-0 em0] [-1 em1] [-2 em2] [-f ref] <in.bam>

这个命令鉴别目标区域,它通过检查 read 深度的连续性,来计算目标区域的单倍体的一致序列,并为每个相应目标区域输出一个序列到 SAM 文件中。当指定了参数 -f 时,会使用 BAQ 。这个命令的设计,只是为了 cutting fosmid clones from fosmid pool sequencing [Ref. Kitzman et al. (2010)].

phase
samtools phase [-AF] [-k len] [-b prefix] [-q minLOD] [-Q minBaseQ] <in.bam>

识别并且 phase (分相位)杂合SNPs。

OPTIONS:

-A
不处理那些 phase (相位)不清晰的 reads 。 Drop reads with ambiguous phase.
-b STR
Prefix of BAM output. When this option is in use, phase-0 reads will be saved in file STR.0.bam and phase-1 reads in STR.1.bam. Phase unknown reads will be randomly allocated to one of the two files. Chimeric reads with switch errors will be saved in STR.chimeric.bam. [null]
-F
不要尝试修复嵌合型的 reads 。
-k INT
局部 phasing 的最大长度。 [13]
-q INT
用以识别杂合子的最小的 Phred 标准的 LOD 值。 [40]
-Q INT
用于识别杂合的最小的碱基质量。 [13]
flags
samtools flags INT|STR[,…]

文本式和数字式的标志法之间的转换。

FLAGS:

0x1 PAIRED .. 双端(多段)测序技术
0x2 PROPER_PAIR .. 每短序列都被 aligner 合适的比对上了
0x4 UNMAP .. 片段未比对成功
0x8 MUNMAP .. 模板中的下一个片段没有被比对上
0x10 REVERSE .. SEQ序列是反向互补的
0x20 MREVERSE .. SEQ序列的下一段是反向互补的
0x40 READ1 .. 模板中的第一个片段
0x80 READ2 .. 模板中的最后一个片段
0x100 SECONDARY .. 次要比对
0x200 QCFAIL .. 没有通过质量控制
0x400 DUP .. PCR或者光学的重复
0x800 SUPPLEMENTARY .. 增补的比对
bam2fq
samtools bam2fq [-nO] [-s <outSE.fq>] <in.bam>

把 bam 转换成 FASTQ 格式。

OPTIONS:

-n
默认情况下,根据 BAM_READ1 或 BAM_READ2 的标志,‘/1’ 或者 ‘/2’ 会被相应的添加在read名字的末尾。 使用 -n 会保留原有的read的名字。
-O
如果可行的话,同标准质量字符串相比,优先选择使用 OQ 标签里的质量值.
-s FILE
把单个的 reads 以 FASTQ 格式写入 FILE 中,而不是输出他们。
help, –help
列出samtools可用的命令,并显示他们的简要用法。 如果同时指定了某个命令,如 samtools help view ,那么这个命令的详细用法将会被展示出来。
–version
显示 samtools 的版本号和版权信息,同时也列出 samtools 所使用的重要的库。
–version-only
以机器可读的格式显示 samtools 的完整版本号。

REFERENCE SEQUENCES – 参考序列

对于读和写,CRAM格式都需要使用参考序列。

当读取一个 CRAM 文件的时候,会去查询 @SQ 头文件,以寻找参考序列的MD5sum (M5: 标签)和本地参考序列文件(UR: 标签)。注意,基于http://ftp:// 的网址是不可以在 UR: 字段使用的,但是本地的 fasta 文件名(不管有没有添加 file:// )是可以的。

当创建一个 CRAM 文件的时候,@SQ 头文件也会被读取然后识别参考序列,但是 M5: 和 UR: 标签可能不存在。 在这种情况下,可以使用 samtools view 的 -T-t 参数去分别指定 fasta 或 fasta.fai 文件名(提供的 .fasta.fai 文件同样需要有 .fasta 文件做支撑)。

获取参考序列的搜索顺序是:

使用由命令行参数(比如 -T)指定的任何文件。

通过 REF_CACHE 环境变量查找 MD5。

查询 REF_PATH 环境变量中每个元素的 MD5 。

查找列在 UR 头标签里的本地文件。

ENVIRONMENT VARIABLES – 环境变量

REF_PATH
由冒号分隔(在windows上是分号)的一列位置,用于去搜索由MD5sum指定的序列。它可以是 一列目录或者网址。注意,如果是网址的话,http:// 和 ftp:// 以及可能出现的端口号中的冒 号将不会被视为网址的一部分,而不会被当做 PATH 字段的分隔符。对于网址而言, %s 会被读入的 MD5 替换。

如果没有指定REF_PATH,那么将会使用默认值:http://www.ebi.ac.uk/ena/cram/md5/%s; 如果 REF_CACHE 也没有指定的话,那么将会被设置成 $XDG_CACHE_HOME/hts-ref/%2s/%2s/%s 。 如果 $XDG_CACHE_HOME 没有设置的话, 会使用$HOME/.cache (如果没有找到 home 文件夹的话, 会使用本地系统的临时文件夹)。

REF_CACHE
这可被定义为一个单一的文件夹,用作本地的序列缓存区。一旦下载了一个参考序列,它将会被储存在 由REF_CACHE指定的位置。当读取一个参考序列时,会先在这个文件夹下搜索,然后再去 REF_PATH 搜索。 为了避免同一个文件夹下存入了大量的文件,可能会使用 %nums and %s 标记去创建一个由 MD5sum 中 num 字符构成的路径。例如:/local/ref_cache/%2s/%2s/%s 这串字符 将会创建一个嵌套两层的文件夹,最内层的文件夹里面将会存入 md5sum 最后的28个字符。

REF_CACHE 会被先搜索,然后才会尝试通过 REF_PATH 搜索列表加载。 如果 REF_PATH 没有被设定,那么 REF_PATH 和 REF_CACHE 会被自动设置(参见上文)。 但是如果 REF_PATH 被设定了而 REF_CACHE 没有被设定,那么就不会使用本地缓存。

为了方便生成更多的 REF_CACHE 文件夹下的内容,我们在Samtools软件中提供了一个 脚本,misc/seq_ref_populate.pl。它会为一个 fasta 文件或一个文件夹下所有的 fasta 文件生成 MD5sum 命名的文件。

EXAMPLES – 示例

@SQ 行存在于头文件中时,把 SAM 转化为 BAM 格式:

samtools view -bS aln.sam > aln.bam

如果 @SQ 行不存在:

samtools faidx ref.fa
samtools view -bt ref.fa.fai aln.sam > aln.bam

ref.fa.fai 是由 faidx 命令自动生成的。

使用本地的参考序列,把 BAM 文件转换成 CRAM 文件。

samtools view -C -T ref.fa aln.bam > aln.cram

在合并排序后的比对序列的时候,添加 RG 标签:

perl -e 'print "@RG\\tID:ga\\tSM:hs\\tLB:ga\\tPL:Illumina\\n@RG\\tID:454\\tSM:hs\\tLB:454\\tPL:454\\n"' > rg.txt
samtools merge -rh rg.txt merged.bam ga.bam 454.bam

RG 标签的值是由包含这个read的文件名决定的。在这个例子中,在 merged.bam 里面, 来源于 ga.bam 的reads会被标记成 RG:Z:ga ,而来源于 454.bam 的reads则会被 标记成 RG:Z:454

识别 SNP 和短 INDEL (插入删除):

samtools mpileup -uf ref.fa aln.bam | bcftools call -mv > var.raw.vcf
bcftools filter -s LowQual -e '%QUAL<20 || DP>100' var.raw.vcf  > var.flt.vcf

bcftools filter 命令会标记低质量的位点,和read深度超过一个限度的位点, 这个限度应该被调整为大约两倍平均测序深度(更大的read深度通常意味着富含虚假信息的问题区域)。 对于含有过量mismatch的reads而言,它们mapping的质量如果被高估了的话,分析员可以考虑 给 mpileup 添加 -C50 。应用这个参数通常对 BWA-short 有所 帮助,对别的比对软件可能没有帮助。

个体是由 @RG header行的 SM 标签来确定的。不同个体可以放到一个比对文件中; 一个个体也可以被分到多个文件中。参数 -P 指明候选indel应该只从 @RG-PL 标签为 ILLUMINA 的read groups中收集。从容易出现indel的测序技术测出来的reads中收集候选 indel可能会影响indel识别的性能。

为一个二倍体个体产生保守序列

samtools mpileup -uf ref.fa aln.bam | bcftools call -c | vcfutils.pl vcf2fq > cns.fq

为一个个体解单倍体型(Phase one individual):

samtools calmd -AEur aln.bam ref.fa | samtools phase -b prefix - > phase.out

命令 calmd 用来减少INDELs附近的假杂合。

Dump BAQ applied alignment for other SNP callers:

samtools calmd -bAr aln.bam > aln.baq.bam

该命令同时添加并矫正 NMMD 标签。calmd 命令也有 -C 参数, 同 pileupmpileup 里面的是一样的。如果有帮助,可以使用它。

LIMITATIONS – 限制

在 bam_import.c, bam_endian.h, bam.c 和 bam_aux.c 中,单词没有对齐。

Samtools 对双端测序的 rmdup 命令对非配对的reads不适用。(例如:孤reads或者两条map到不同染 色体上的双端测序的reads).如果要考虑这个的话,请使用 Picard 的 MarkDuplicates ,虽然有一点慢, 但它能够正确地处理这些情况。

AUTHOR – 作者

桑格研究所的 Heng Li 写了 C 版本的 samtools 。Broad研究院的 Bob Handsaker 实现了 BGZF 库。 John Marshall 和 Petr Danecek 为源码做出了贡献,并且,千人基因组计划中的许多人为 SAM 格式标准 做出了贡献。

SEE ALSO – 另请参阅

bcftools(1), sam(5)

Samtools 网站: <http://samtools.sourceforge.net>
Samtools 最新源码: <https://github.com/samtools/samtools>
SAM/BAM,CRAM,VCF/BCF 的文件格式标准: <http://samtools.github.io/hts-specs>
HTSlib 网站: <https://github.com/samtools/htslib>
Bcftools 网站: <http://samtools.github.io/bcftools>