New upstream release, with among many changes a new python script.
[samtools.git] / samtools.1
1 .TH samtools 1 "11 July 2010" "samtools-0.1.8" "Bioinformatics tools"
2 .SH NAME
3 .PP
4 samtools - Utilities for the Sequence Alignment/Map (SAM) format
5 .SH SYNOPSIS
6 .PP
7 samtools view -bt ref_list.txt -o aln.bam aln.sam.gz
8 .PP
9 samtools sort aln.bam aln.sorted
10 .PP
11 samtools index aln.sorted.bam
12 .PP
13 samtools idxstats aln.sorted.bam
14 .PP
15 samtools view aln.sorted.bam chr2:20,100,000-20,200,000
16 .PP
17 samtools merge out.bam in1.bam in2.bam in3.bam
18 .PP
19 samtools faidx ref.fasta
20 .PP
21 samtools pileup -f ref.fasta aln.sorted.bam
22 .PP
23 samtools mpileup -f ref.fasta -r chr3:1,000-2,000 in1.bam in2.bam
24 .PP
25 samtools tview aln.sorted.bam ref.fasta
26
27 .SH DESCRIPTION
28 .PP
29 Samtools is a set of utilities that manipulate alignments in the BAM
30 format. It imports from and exports to the SAM (Sequence Alignment/Map)
31 format, does sorting, merging and indexing, and allows to retrieve reads
32 in any regions swiftly.
33
34 Samtools is designed to work on a stream. It regards an input file `-'
35 as the standard input (stdin) and an output file `-' as the standard
36 output (stdout). Several commands can thus be combined with Unix
37 pipes. Samtools always output warning and error messages to the standard
38 error output (stderr).
39
40 Samtools is also able to open a BAM (not SAM) file on a remote FTP or
41 HTTP server if the BAM file name starts with `ftp://' or `http://'.
42 Samtools checks the current working directory for the index file and
43 will download the index upon absence. Samtools does not retrieve the
44 entire alignment file unless it is asked to do so.
45
46 .SH COMMANDS AND OPTIONS
47
48 .TP 10
49 .B view
50 samtools view [-bhuHS] [-t in.refList] [-o output] [-f reqFlag] [-F
51 skipFlag] [-q minMapQ] [-l library] [-r readGroup] [-R rgFile] <in.bam>|<in.sam> [region1 [...]]
52
53 Extract/print all or sub alignments in SAM or BAM format. If no region
54 is specified, all the alignments will be printed; otherwise only
55 alignments overlapping the specified regions will be output. An
56 alignment may be given multiple times if it is overlapping several
57 regions. A region can be presented, for example, in the following
58 format: `chr2' (the whole chr2), `chr2:1000000' (region starting from
59 1,000,000bp) or `chr2:1,000,000-2,000,000' (region between 1,000,000 and
60 2,000,000bp including the end points). The coordinate is 1-based.
61
62 .B OPTIONS:
63 .RS
64 .TP 8
65 .B -b
66 Output in the BAM format.
67 .TP
68 .B -u
69 Output uncompressed BAM. This option saves time spent on
70 compression/decomprssion and is thus preferred when the output is piped
71 to another samtools command.
72 .TP
73 .B -h
74 Include the header in the output.
75 .TP
76 .B -H
77 Output the header only.
78 .TP
79 .B -S
80 Input is in SAM. If @SQ header lines are absent, the
81 .B `-t'
82 option is required.
83 .TP
84 .B -t FILE
85 This file is TAB-delimited. Each line must contain the reference name
86 and the length of the reference, one line for each distinct reference;
87 additional fields are ignored. This file also defines the order of the
88 reference sequences in sorting. If you run `samtools faidx <ref.fa>',
89 the resultant index file
90 .I <ref.fa>.fai
91 can be used as this
92 .I <in.ref_list>
93 file.
94 .TP
95 .B -o FILE
96 Output file [stdout]
97 .TP
98 .B -f INT
99 Only output alignments with all bits in INT present in the FLAG
100 field. INT can be in hex in the format of /^0x[0-9A-F]+/ [0]
101 .TP
102 .B -F INT
103 Skip alignments with bits present in INT [0]
104 .TP
105 .B -q INT
106 Skip alignments with MAPQ smaller than INT [0]
107 .TP
108 .B -l STR
109 Only output reads in library STR [null]
110 .TP
111 .B -r STR
112 Only output reads in read group STR [null]
113 .TP
114 .B -R FILE
115 Output reads in read groups listed in
116 .I FILE
117 [null]
118 .RE
119
120 .TP
121 .B tview
122 samtools tview <in.sorted.bam> [ref.fasta]
123
124 Text alignment viewer (based on the ncurses library). In the viewer,
125 press `?' for help and press `g' to check the alignment start from a
126 region in the format like `chr10:10,000,000' or `=10,000,000' when
127 viewing the same reference sequence.
128
129 .TP
130 .B pileup
131 samtools pileup [-f in.ref.fasta] [-t in.ref_list] [-l in.site_list]
132 [-iscgS2] [-T theta] [-N nHap] [-r pairDiffRate] <in.bam>|<in.sam>
133
134 Print the alignment in the pileup format. In the pileup format, each
135 line represents a genomic position, consisting of chromosome name,
136 coordinate, reference base, read bases, read qualities and alignment
137 mapping qualities. Information on match, mismatch, indel, strand,
138 mapping quality and start and end of a read are all encoded at the read
139 base column. At this column, a dot stands for a match to the reference
140 base on the forward strand, a comma for a match on the reverse strand,
141 `ACGTN' for a mismatch on the forward strand and `acgtn' for a mismatch
142 on the reverse strand. A pattern `\\+[0-9]+[ACGTNacgtn]+' indicates
143 there is an insertion between this reference position and the next
144 reference position. The length of the insertion is given by the integer
145 in the pattern, followed by the inserted sequence. Similarly, a pattern
146 `-[0-9]+[ACGTNacgtn]+' represents a deletion from the reference. The
147 deleted bases will be presented as `*' in the following lines. Also at
148 the read base column, a symbol `^' marks the start of a read segment
149 which is a contiguous subsequence on the read separated by `N/S/H' CIGAR
150 operations. The ASCII of the character following `^' minus 33 gives the
151 mapping quality. A symbol `$' marks the end of a read segment.
152
153 If option
154 .B -c
155 is applied, the consensus base, Phred-scaled consensus quality, SNP
156 quality (i.e. the Phred-scaled probability of the consensus being
157 identical to the reference) and root mean square (RMS) mapping quality
158 of the reads covering the site will be inserted between the `reference
159 base' and the `read bases' columns. An indel occupies an additional
160 line. Each indel line consists of chromosome name, coordinate, a star,
161 the genotype, consensus quality, SNP quality, RMS mapping quality, #
162 covering reads, the first alllele, the second allele, # reads supporting
163 the first allele, # reads supporting the second allele and # reads
164 containing indels different from the top two alleles.
165
166 The position of indels is offset by -1.
167
168 .B OPTIONS:
169 .RS
170 .TP 10
171 .B -s
172 Print the mapping quality as the last column. This option makes the
173 output easier to parse, although this format is not space efficient.
174 .TP
175 .B -S
176 The input file is in SAM.
177 .TP
178 .B -i
179 Only output pileup lines containing indels.
180 .TP
181 .B -f FILE
182 The reference sequence in the FASTA format. Index file
183 .I FILE.fai
184 will be created if
185 absent.
186 .TP
187 .B -M INT
188 Cap mapping quality at INT [60]
189 .TP
190 .B -m INT
191 Filter reads with flag containing bits in
192 .I
193 INT
194 [1796]
195 .TP
196 .B -d INT
197 Use the first
198 .I NUM
199 reads in the pileup for indel calling for speed up. Zero for unlimited. [0]
200 .TP
201 .B -t FILE
202 List of reference names ane sequence lengths, in the format described
203 for the
204 .B import
205 command. If this option is present, samtools assumes the input
206 .I <in.alignment>
207 is in SAM format; otherwise it assumes in BAM format.
208 .TP
209 .B -l FILE
210 List of sites at which pileup is output. This file is space
211 delimited. The first two columns are required to be chromosome and
212 1-based coordinate. Additional columns are ignored. It is
213 recommended to use option
214 .B -s
215 together with
216 .B -l
217 as in the default format we may not know the mapping quality.
218 .TP
219 .B -c
220 Call the consensus sequence using SOAPsnp consensus model. Options
221 .B -T,
222 .B -N,
223 .B -I
224 and
225 .B -r
226 are only effective when
227 .B -c
228 or
229 .B -g
230 is in use.
231 .TP
232 .B -g
233 Generate genotype likelihood in the binary GLFv3 format. This option
234 suppresses -c, -i and -s.
235 .TP
236 .B -T FLOAT
237 The theta parameter (error dependency coefficient) in the maq consensus
238 calling model [0.85]
239 .TP
240 .B -N INT
241 Number of haplotypes in the sample (>=2) [2]
242 .TP
243 .B -r FLOAT
244 Expected fraction of differences between a pair of haplotypes [0.001]
245 .TP
246 .B -I INT
247 Phred probability of an indel in sequencing/prep. [40]
248 .RE
249
250 .TP
251 .B mpileup
252 samtools mpileup [-r reg] [-f in.fa] in.bam [in2.bam [...]]
253
254 Generate pileup for multiple BAM files. Consensus calling is not
255 implemented.
256
257 .B OPTIONS:
258 .RS
259 .TP 8
260 .B -r STR
261 Only generate pileup in region
262 .I STR
263 [all sites]
264 .TP
265 .B -f FILE
266 The reference file [null]
267 .RE
268
269 .TP
270 .B reheader
271 samtools reheader <in.header.sam> <in.bam>
272
273 Replace the header in
274 .I in.bam
275 with the header in
276 .I in.header.sam.
277 This command is much faster than replacing the header with a
278 BAM->SAM->BAM conversion.
279
280 .TP
281 .B sort
282 samtools sort [-no] [-m maxMem] <in.bam> <out.prefix>
283
284 Sort alignments by leftmost coordinates. File
285 .I <out.prefix>.bam
286 will be created. This command may also create temporary files
287 .I <out.prefix>.%d.bam
288 when the whole alignment cannot be fitted into memory (controlled by
289 option -m).
290
291 .B OPTIONS:
292 .RS
293 .TP 8
294 .B -o
295 Output the final alignment to the standard output.
296 .TP
297 .B -n
298 Sort by read names rather than by chromosomal coordinates
299 .TP
300 .B -m INT
301 Approximately the maximum required memory. [500000000]
302 .RE
303
304 .TP
305 .B merge
306 samtools merge [-h inh.sam] [-nr] <out.bam> <in1.bam> <in2.bam> [...]
307
308 Merge multiple sorted alignments.
309 The header reference lists of all the input BAM files, and the @SQ headers of
310 .IR inh.sam ,
311 if any, must all refer to the same set of reference sequences.
312 The header reference list and (unless overridden by
313 .BR -h )
314 `@' headers of
315 .I in1.bam
316 will be copied to
317 .IR out.bam ,
318 and the headers of other files will be ignored.
319
320 .B OPTIONS:
321 .RS
322 .TP 8
323 .B -h FILE
324 Use the lines of
325 .I FILE
326 as `@' headers to be copied to
327 .IR out.bam ,
328 replacing any header lines that would otherwise be copied from
329 .IR in1.bam .
330 .RI ( FILE
331 is actually in SAM format, though any alignment records it may contain
332 are ignored.)
333 .TP
334 .B -r
335 Attach an RG tag to each alignment. The tag value is inferred from file names.
336 .TP
337 .B -n
338 The input alignments are sorted by read names rather than by chromosomal
339 coordinates
340 .RE
341
342 .TP
343 .B index
344 samtools index <aln.bam>
345
346 Index sorted alignment for fast random access. Index file
347 .I <aln.bam>.bai
348 will be created.
349
350 .TP
351 .B idxstats
352 samtools idxstats <aln.bam>
353
354 Retrieve and print stats in the index file. The output is TAB delimited
355 with each line consisting of reference sequence name, sequence length, #
356 mapped reads and # unmapped reads.
357
358 .TP
359 .B faidx
360 samtools faidx <ref.fasta> [region1 [...]]
361
362 Index reference sequence in the FASTA format or extract subsequence from
363 indexed reference sequence. If no region is specified,
364 .B faidx
365 will index the file and create
366 .I <ref.fasta>.fai
367 on the disk. If regions are speficified, the subsequences will be
368 retrieved and printed to stdout in the FASTA format. The input file can
369 be compressed in the
370 .B RAZF
371 format.
372
373 .TP
374 .B fixmate
375 samtools fixmate <in.nameSrt.bam> <out.bam>
376
377 Fill in mate coordinates, ISIZE and mate related flags from a
378 name-sorted alignment.
379
380 .TP
381 .B rmdup
382 samtools rmdup [-sS] <input.srt.bam> <out.bam>
383
384 Remove potential PCR duplicates: if multiple read pairs have identical
385 external coordinates, only retain the pair with highest mapping quality.
386 In the paired-end mode, this command
387 .B ONLY
388 works with FR orientation and requires ISIZE is correctly set. It does
389 not work for unpaired reads (e.g. two ends mapped to different
390 chromosomes or orphan reads).
391
392 .B OPTIONS:
393 .RS
394 .TP 8
395 .B -s
396 Remove duplicate for single-end reads. By default, the command works for
397 paired-end reads only.
398 .TP 8
399 .B -S
400 Treat paired-end reads and single-end reads.
401 .RE
402
403 .TP
404 .B calmd
405 samtools calmd [-eubS] <aln.bam> <ref.fasta>
406
407 Generate the MD tag. If the MD tag is already present, this command will
408 give a warning if the MD tag generated is different from the existing
409 tag. Output SAM by default.
410
411 .B OPTIONS:
412 .RS
413 .TP 8
414 .B -e
415 Convert a the read base to = if it is identical to the aligned reference
416 base. Indel caller does not support the = bases at the moment.
417 .TP
418 .B -u
419 Output uncompressed BAM
420 .TP
421 .B -b
422 Output compressed BAM
423 .TP
424 .B -S
425 The input is SAM with header lines
426 .RE
427
428 .SH SAM FORMAT
429
430 SAM is TAB-delimited. Apart from the header lines, which are started
431 with the `@' symbol, each alignment line consists of:
432
433 .TS
434 center box;
435 cb | cb | cb
436 n | l | l .
437 Col     Field   Description
438 _
439 1       QNAME   Query (pair) NAME
440 2       FLAG    bitwise FLAG
441 3       RNAME   Reference sequence NAME
442 4       POS     1-based leftmost POSition/coordinate of clipped sequence
443 5       MAPQ    MAPping Quality (Phred-scaled)
444 6       CIAGR   extended CIGAR string
445 7       MRNM    Mate Reference sequence NaMe (`=' if same as RNAME)
446 8       MPOS    1-based Mate POSistion
447 9       ISIZE   Inferred insert SIZE
448 10      SEQ     query SEQuence on the same strand as the reference
449 11      QUAL    query QUALity (ASCII-33 gives the Phred base quality)
450 12      OPT     variable OPTional fields in the format TAG:VTYPE:VALUE
451 .TE
452
453 .PP
454 Each bit in the FLAG field is defined as:
455
456 .TS
457 center box;
458 cb | cb | cb
459 l | c | l .
460 Flag    Chr     Description
461 _
462 0x0001  p       the read is paired in sequencing
463 0x0002  P       the read is mapped in a proper pair
464 0x0004  u       the query sequence itself is unmapped
465 0x0008  U       the mate is unmapped
466 0x0010  r       strand of the query (1 for reverse)
467 0x0020  R       strand of the mate
468 0x0040  1       the read is the first read in a pair
469 0x0080  2       the read is the second read in a pair
470 0x0100  s       the alignment is not primary
471 0x0200  f       the read fails platform/vendor quality checks
472 0x0400  d       the read is either a PCR or an optical duplicate
473 .TE
474
475 .SH LIMITATIONS
476 .PP
477 .IP o 2
478 Unaligned words used in bam_import.c, bam_endian.h, bam.c and bam_aux.c.
479 .IP o 2
480 In merging, the input files are required to have the same number of
481 reference sequences. The requirement can be relaxed. In addition,
482 merging does not reconstruct the header dictionaries
483 automatically. Endusers have to provide the correct header. Picard is
484 better at merging.
485 .IP o 2
486 Samtools paired-end rmdup does not work for unpaired reads (e.g. orphan
487 reads or ends mapped to different chromosomes). If this is a concern,
488 please use Picard's MarkDuplicate which correctly handles these cases,
489 although a little slower.
490
491 .SH AUTHOR
492 .PP
493 Heng Li from the Sanger Institute wrote the C version of samtools. Bob
494 Handsaker from the Broad Institute implemented the BGZF library and Jue
495 Ruan from Beijing Genomics Institute wrote the RAZF library. Various
496 people in the 1000 Genomes Project contributed to the SAM format
497 specification.
498
499 .SH SEE ALSO
500 .PP
501 Samtools website: <http://samtools.sourceforge.net>