|
0
|
1 .TH samtools 1 "15 March 2013" "samtools-0.1.19" "Bioinformatics tools"
|
|
|
2 .SH NAME
|
|
|
3 .PP
|
|
|
4 samtools - Utilities for the Sequence Alignment/Map (SAM) format
|
|
|
5
|
|
|
6 bcftools - Utilities for the Binary Call Format (BCF) and VCF
|
|
|
7 .SH SYNOPSIS
|
|
|
8 .PP
|
|
|
9 samtools view -bt ref_list.txt -o aln.bam aln.sam.gz
|
|
|
10 .PP
|
|
|
11 samtools sort aln.bam aln.sorted
|
|
|
12 .PP
|
|
|
13 samtools index aln.sorted.bam
|
|
|
14 .PP
|
|
|
15 samtools idxstats aln.sorted.bam
|
|
|
16 .PP
|
|
|
17 samtools view aln.sorted.bam chr2:20,100,000-20,200,000
|
|
|
18 .PP
|
|
|
19 samtools merge out.bam in1.bam in2.bam in3.bam
|
|
|
20 .PP
|
|
|
21 samtools faidx ref.fasta
|
|
|
22 .PP
|
|
|
23 samtools pileup -vcf ref.fasta aln.sorted.bam
|
|
|
24 .PP
|
|
|
25 samtools mpileup -C50 -gf ref.fasta -r chr3:1,000-2,000 in1.bam in2.bam
|
|
|
26 .PP
|
|
|
27 samtools tview aln.sorted.bam ref.fasta
|
|
|
28 .PP
|
|
|
29 bcftools index in.bcf
|
|
|
30 .PP
|
|
|
31 bcftools view in.bcf chr2:100-200 > out.vcf
|
|
|
32 .PP
|
|
|
33 bcftools view -Nvm0.99 in.bcf > out.vcf 2> out.afs
|
|
|
34
|
|
|
35 .SH DESCRIPTION
|
|
|
36 .PP
|
|
|
37 Samtools is a set of utilities that manipulate alignments in the BAM
|
|
|
38 format. It imports from and exports to the SAM (Sequence Alignment/Map)
|
|
|
39 format, does sorting, merging and indexing, and allows to retrieve reads
|
|
|
40 in any regions swiftly.
|
|
|
41
|
|
|
42 Samtools is designed to work on a stream. It regards an input file `-'
|
|
|
43 as the standard input (stdin) and an output file `-' as the standard
|
|
|
44 output (stdout). Several commands can thus be combined with Unix
|
|
|
45 pipes. Samtools always output warning and error messages to the standard
|
|
|
46 error output (stderr).
|
|
|
47
|
|
|
48 Samtools is also able to open a BAM (not SAM) file on a remote FTP or
|
|
|
49 HTTP server if the BAM file name starts with `ftp://' or `http://'.
|
|
|
50 Samtools checks the current working directory for the index file and
|
|
|
51 will download the index upon absence. Samtools does not retrieve the
|
|
|
52 entire alignment file unless it is asked to do so.
|
|
|
53
|
|
|
54 .SH SAMTOOLS COMMANDS AND OPTIONS
|
|
|
55
|
|
|
56 .TP 10
|
|
|
57 .B view
|
|
|
58 samtools view [-bchuHS] [-t in.refList] [-o output] [-f reqFlag] [-F
|
|
|
59 skipFlag] [-q minMapQ] [-l library] [-r readGroup] [-R rgFile] <in.bam>|<in.sam> [region1 [...]]
|
|
|
60
|
|
|
61 Extract/print all or sub alignments in SAM or BAM format. If no region
|
|
|
62 is specified, all the alignments will be printed; otherwise only
|
|
|
63 alignments overlapping the specified regions will be output. An
|
|
|
64 alignment may be given multiple times if it is overlapping several
|
|
|
65 regions. A region can be presented, for example, in the following
|
|
|
66 format: `chr2' (the whole chr2), `chr2:1000000' (region starting from
|
|
|
67 1,000,000bp) or `chr2:1,000,000-2,000,000' (region between 1,000,000 and
|
|
|
68 2,000,000bp including the end points). The coordinate is 1-based.
|
|
|
69
|
|
|
70 .B OPTIONS:
|
|
|
71 .RS
|
|
|
72 .TP 10
|
|
|
73 .B -b
|
|
|
74 Output in the BAM format.
|
|
|
75 .TP
|
|
|
76 .BI -f \ INT
|
|
|
77 Only output alignments with all bits in INT present in the FLAG
|
|
|
78 field. INT can be in hex in the format of /^0x[0-9A-F]+/ [0]
|
|
|
79 .TP
|
|
|
80 .BI -F \ INT
|
|
|
81 Skip alignments with bits present in INT [0]
|
|
|
82 .TP
|
|
|
83 .B -h
|
|
|
84 Include the header in the output.
|
|
|
85 .TP
|
|
|
86 .B -H
|
|
|
87 Output the header only.
|
|
|
88 .TP
|
|
|
89 .BI -l \ STR
|
|
|
90 Only output reads in library STR [null]
|
|
|
91 .TP
|
|
|
92 .BI -o \ FILE
|
|
|
93 Output file [stdout]
|
|
|
94 .TP
|
|
|
95 .BI -q \ INT
|
|
|
96 Skip alignments with MAPQ smaller than INT [0]
|
|
|
97 .TP
|
|
|
98 .BI -r \ STR
|
|
|
99 Only output reads in read group STR [null]
|
|
|
100 .TP
|
|
|
101 .BI -R \ FILE
|
|
|
102 Output reads in read groups listed in
|
|
|
103 .I FILE
|
|
|
104 [null]
|
|
|
105 .TP
|
|
|
106 .BI -s \ FLOAT
|
|
|
107 Fraction of templates/pairs to subsample; the integer part is treated as the
|
|
|
108 seed for the random number generator [-1]
|
|
|
109 .TP
|
|
|
110 .B -S
|
|
|
111 Input is in SAM. If @SQ header lines are absent, the
|
|
|
112 .B `-t'
|
|
|
113 option is required.
|
|
|
114 .TP
|
|
|
115 .B -c
|
|
|
116 Instead of printing the alignments, only count them and print the
|
|
|
117 total number. All filter options, such as
|
|
|
118 .B `-f',
|
|
|
119 .B `-F'
|
|
|
120 and
|
|
|
121 .B `-q'
|
|
|
122 , are taken into account.
|
|
|
123 .TP
|
|
|
124 .BI -t \ FILE
|
|
|
125 This file is TAB-delimited. Each line must contain the reference name
|
|
|
126 and the length of the reference, one line for each distinct reference;
|
|
|
127 additional fields are ignored. This file also defines the order of the
|
|
|
128 reference sequences in sorting. If you run `samtools faidx <ref.fa>',
|
|
|
129 the resultant index file
|
|
|
130 .I <ref.fa>.fai
|
|
|
131 can be used as this
|
|
|
132 .I <in.ref_list>
|
|
|
133 file.
|
|
|
134 .TP
|
|
|
135 .B -u
|
|
|
136 Output uncompressed BAM. This option saves time spent on
|
|
|
137 compression/decomprssion and is thus preferred when the output is piped
|
|
|
138 to another samtools command.
|
|
|
139 .RE
|
|
|
140
|
|
|
141 .TP
|
|
|
142 .B tview
|
|
|
143 samtools tview
|
|
|
144 .RB [ \-p
|
|
|
145 .IR chr:pos ]
|
|
|
146 .RB [ \-s
|
|
|
147 .IR STR ]
|
|
|
148 .RB [ \-d
|
|
|
149 .IR display ]
|
|
|
150 .RI <in.sorted.bam>
|
|
|
151 .RI [ref.fasta]
|
|
|
152
|
|
|
153 Text alignment viewer (based on the ncurses library). In the viewer,
|
|
|
154 press `?' for help and press `g' to check the alignment start from a
|
|
|
155 region in the format like `chr10:10,000,000' or `=10,000,000' when
|
|
|
156 viewing the same reference sequence.
|
|
|
157
|
|
|
158 .B Options:
|
|
|
159 .RS
|
|
|
160 .TP 14
|
|
|
161 .BI -d \ display
|
|
|
162 Output as (H)tml or (C)urses or (T)ext
|
|
|
163 .TP
|
|
|
164 .BI -p \ chr:pos
|
|
|
165 Go directly to this position
|
|
|
166 .TP
|
|
|
167 .BI -s \ STR
|
|
|
168 Display only reads from this sample or read group
|
|
|
169 .RE
|
|
|
170
|
|
|
171 .TP
|
|
|
172 .B mpileup
|
|
|
173 samtools mpileup
|
|
|
174 .RB [ \-EBugp ]
|
|
|
175 .RB [ \-C
|
|
|
176 .IR capQcoef ]
|
|
|
177 .RB [ \-r
|
|
|
178 .IR reg ]
|
|
|
179 .RB [ \-f
|
|
|
180 .IR in.fa ]
|
|
|
181 .RB [ \-l
|
|
|
182 .IR list ]
|
|
|
183 .RB [ \-M
|
|
|
184 .IR capMapQ ]
|
|
|
185 .RB [ \-Q
|
|
|
186 .IR minBaseQ ]
|
|
|
187 .RB [ \-q
|
|
|
188 .IR minMapQ ]
|
|
|
189 .I in.bam
|
|
|
190 .RI [ in2.bam
|
|
|
191 .RI [ ... ]]
|
|
|
192
|
|
|
193 Generate BCF or pileup for one or multiple BAM files. Alignment records
|
|
|
194 are grouped by sample identifiers in @RG header lines. If sample
|
|
|
195 identifiers are absent, each input file is regarded as one sample.
|
|
|
196
|
|
|
197 In the pileup format (without
|
|
|
198 .BR -u or -g ),
|
|
|
199 each
|
|
|
200 line represents a genomic position, consisting of chromosome name,
|
|
|
201 coordinate, reference base, read bases, read qualities and alignment
|
|
|
202 mapping qualities. Information on match, mismatch, indel, strand,
|
|
|
203 mapping quality and start and end of a read are all encoded at the read
|
|
|
204 base column. At this column, a dot stands for a match to the reference
|
|
|
205 base on the forward strand, a comma for a match on the reverse strand,
|
|
|
206 a '>' or '<' for a reference skip, `ACGTN' for a mismatch on the forward
|
|
|
207 strand and `acgtn' for a mismatch on the reverse strand. A pattern
|
|
|
208 `\\+[0-9]+[ACGTNacgtn]+' indicates there is an insertion between this
|
|
|
209 reference position and the next reference position. The length of the
|
|
|
210 insertion is given by the integer in the pattern, followed by the
|
|
|
211 inserted sequence. Similarly, a pattern `-[0-9]+[ACGTNacgtn]+'
|
|
|
212 represents a deletion from the reference. The deleted bases will be
|
|
|
213 presented as `*' in the following lines. Also at the read base column, a
|
|
|
214 symbol `^' marks the start of a read. The ASCII of the character
|
|
|
215 following `^' minus 33 gives the mapping quality. A symbol `$' marks the
|
|
|
216 end of a read segment.
|
|
|
217
|
|
|
218 .B Input Options:
|
|
|
219 .RS
|
|
|
220 .TP 10
|
|
|
221 .B -6
|
|
|
222 Assume the quality is in the Illumina 1.3+ encoding.
|
|
|
223 .B -A
|
|
|
224 Do not skip anomalous read pairs in variant calling.
|
|
|
225 .TP
|
|
|
226 .B -B
|
|
|
227 Disable probabilistic realignment for the computation of base alignment
|
|
|
228 quality (BAQ). BAQ is the Phred-scaled probability of a read base being
|
|
|
229 misaligned. Applying this option greatly helps to reduce false SNPs
|
|
|
230 caused by misalignments.
|
|
|
231 .TP
|
|
|
232 .BI -b \ FILE
|
|
|
233 List of input BAM files, one file per line [null]
|
|
|
234 .TP
|
|
|
235 .BI -C \ INT
|
|
|
236 Coefficient for downgrading mapping quality for reads containing
|
|
|
237 excessive mismatches. Given a read with a phred-scaled probability q of
|
|
|
238 being generated from the mapped position, the new mapping quality is
|
|
|
239 about sqrt((INT-q)/INT)*INT. A zero value disables this
|
|
|
240 functionality; if enabled, the recommended value for BWA is 50. [0]
|
|
|
241 .TP
|
|
|
242 .BI -d \ INT
|
|
|
243 At a position, read maximally
|
|
|
244 .I INT
|
|
|
245 reads per input BAM. [250]
|
|
|
246 .TP
|
|
|
247 .B -E
|
|
|
248 Extended BAQ computation. This option helps sensitivity especially for MNPs, but may hurt
|
|
|
249 specificity a little bit.
|
|
|
250 .TP
|
|
|
251 .BI -f \ FILE
|
|
|
252 The
|
|
|
253 .BR faidx -indexed
|
|
|
254 reference file in the FASTA format. The file can be optionally compressed by
|
|
|
255 .BR razip .
|
|
|
256 [null]
|
|
|
257 .TP
|
|
|
258 .BI -l \ FILE
|
|
|
259 BED or position list file containing a list of regions or sites where pileup or BCF should be generated [null]
|
|
|
260 .TP
|
|
|
261 .BI -q \ INT
|
|
|
262 Minimum mapping quality for an alignment to be used [0]
|
|
|
263 .TP
|
|
|
264 .BI -Q \ INT
|
|
|
265 Minimum base quality for a base to be considered [13]
|
|
|
266 .TP
|
|
|
267 .BI -r \ STR
|
|
|
268 Only generate pileup in region
|
|
|
269 .I STR
|
|
|
270 [all sites]
|
|
|
271 .TP
|
|
|
272 .B Output Options:
|
|
|
273
|
|
|
274 .TP
|
|
|
275 .B -D
|
|
|
276 Output per-sample read depth
|
|
|
277 .TP
|
|
|
278 .B -g
|
|
|
279 Compute genotype likelihoods and output them in the binary call format (BCF).
|
|
|
280 .TP
|
|
|
281 .B -S
|
|
|
282 Output per-sample Phred-scaled strand bias P-value
|
|
|
283 .TP
|
|
|
284 .B -u
|
|
|
285 Similar to
|
|
|
286 .B -g
|
|
|
287 except that the output is uncompressed BCF, which is preferred for piping.
|
|
|
288
|
|
|
289 .TP
|
|
|
290 .B Options for Genotype Likelihood Computation (for -g or -u):
|
|
|
291
|
|
|
292 .TP
|
|
|
293 .BI -e \ INT
|
|
|
294 Phred-scaled gap extension sequencing error probability. Reducing
|
|
|
295 .I INT
|
|
|
296 leads to longer indels. [20]
|
|
|
297 .TP
|
|
|
298 .BI -h \ INT
|
|
|
299 Coefficient for modeling homopolymer errors. Given an
|
|
|
300 .IR l -long
|
|
|
301 homopolymer
|
|
|
302 run, the sequencing error of an indel of size
|
|
|
303 .I s
|
|
|
304 is modeled as
|
|
|
305 .IR INT * s / l .
|
|
|
306 [100]
|
|
|
307 .TP
|
|
|
308 .B -I
|
|
|
309 Do not perform INDEL calling
|
|
|
310 .TP
|
|
|
311 .BI -L \ INT
|
|
|
312 Skip INDEL calling if the average per-sample depth is above
|
|
|
313 .IR INT .
|
|
|
314 [250]
|
|
|
315 .TP
|
|
|
316 .BI -o \ INT
|
|
|
317 Phred-scaled gap open sequencing error probability. Reducing
|
|
|
318 .I INT
|
|
|
319 leads to more indel calls. [40]
|
|
|
320 .TP
|
|
|
321 .BI -p
|
|
|
322 Apply -m and -F thresholds per sample to increase sensitivity of calling.
|
|
|
323 By default both options are applied to reads pooled from all samples.
|
|
|
324 .TP
|
|
|
325 .BI -P \ STR
|
|
|
326 Comma dilimited list of platforms (determined by
|
|
|
327 .BR @RG-PL )
|
|
|
328 from which indel candidates are obtained. It is recommended to collect
|
|
|
329 indel candidates from sequencing technologies that have low indel error
|
|
|
330 rate such as ILLUMINA. [all]
|
|
|
331 .RE
|
|
|
332
|
|
|
333 .TP
|
|
|
334 .B reheader
|
|
|
335 samtools reheader <in.header.sam> <in.bam>
|
|
|
336
|
|
|
337 Replace the header in
|
|
|
338 .I in.bam
|
|
|
339 with the header in
|
|
|
340 .I in.header.sam.
|
|
|
341 This command is much faster than replacing the header with a
|
|
|
342 BAM->SAM->BAM conversion.
|
|
|
343
|
|
|
344 .TP
|
|
|
345 .B cat
|
|
|
346 samtools cat [-h header.sam] [-o out.bam] <in1.bam> <in2.bam> [ ... ]
|
|
|
347
|
|
|
348 Concatenate BAMs. The sequence dictionary of each input BAM must be identical,
|
|
|
349 although this command does not check this. This command uses a similar trick
|
|
|
350 to
|
|
|
351 .B reheader
|
|
|
352 which enables fast BAM concatenation.
|
|
|
353
|
|
|
354 .TP
|
|
|
355 .B sort
|
|
|
356 samtools sort [-nof] [-m maxMem] <in.bam> <out.prefix>
|
|
|
357
|
|
|
358 Sort alignments by leftmost coordinates. File
|
|
|
359 .I <out.prefix>.bam
|
|
|
360 will be created. This command may also create temporary files
|
|
|
361 .I <out.prefix>.%d.bam
|
|
|
362 when the whole alignment cannot be fitted into memory (controlled by
|
|
|
363 option -m).
|
|
|
364
|
|
|
365 .B OPTIONS:
|
|
|
366 .RS
|
|
|
367 .TP 8
|
|
|
368 .B -o
|
|
|
369 Output the final alignment to the standard output.
|
|
|
370 .TP
|
|
|
371 .B -n
|
|
|
372 Sort by read names rather than by chromosomal coordinates
|
|
|
373 .TP
|
|
|
374 .B -f
|
|
|
375 Use
|
|
|
376 .I <out.prefix>
|
|
|
377 as the full output path and do not append
|
|
|
378 .I .bam
|
|
|
379 suffix.
|
|
|
380 .TP
|
|
|
381 .BI -m \ INT
|
|
|
382 Approximately the maximum required memory. [500000000]
|
|
|
383 .RE
|
|
|
384
|
|
|
385 .TP
|
|
|
386 .B merge
|
|
|
387 samtools merge [-nur1f] [-h inh.sam] [-R reg] <out.bam> <in1.bam> <in2.bam> [...]
|
|
|
388
|
|
|
389 Merge multiple sorted alignments.
|
|
|
390 The header reference lists of all the input BAM files, and the @SQ headers of
|
|
|
391 .IR inh.sam ,
|
|
|
392 if any, must all refer to the same set of reference sequences.
|
|
|
393 The header reference list and (unless overridden by
|
|
|
394 .BR -h )
|
|
|
395 `@' headers of
|
|
|
396 .I in1.bam
|
|
|
397 will be copied to
|
|
|
398 .IR out.bam ,
|
|
|
399 and the headers of other files will be ignored.
|
|
|
400
|
|
|
401 .B OPTIONS:
|
|
|
402 .RS
|
|
|
403 .TP 8
|
|
|
404 .B -1
|
|
|
405 Use zlib compression level 1 to comrpess the output
|
|
|
406 .TP
|
|
|
407 .B -f
|
|
|
408 Force to overwrite the output file if present.
|
|
|
409 .TP 8
|
|
|
410 .BI -h \ FILE
|
|
|
411 Use the lines of
|
|
|
412 .I FILE
|
|
|
413 as `@' headers to be copied to
|
|
|
414 .IR out.bam ,
|
|
|
415 replacing any header lines that would otherwise be copied from
|
|
|
416 .IR in1.bam .
|
|
|
417 .RI ( FILE
|
|
|
418 is actually in SAM format, though any alignment records it may contain
|
|
|
419 are ignored.)
|
|
|
420 .TP
|
|
|
421 .B -n
|
|
|
422 The input alignments are sorted by read names rather than by chromosomal
|
|
|
423 coordinates
|
|
|
424 .TP
|
|
|
425 .BI -R \ STR
|
|
|
426 Merge files in the specified region indicated by
|
|
|
427 .I STR
|
|
|
428 [null]
|
|
|
429 .TP
|
|
|
430 .B -r
|
|
|
431 Attach an RG tag to each alignment. The tag value is inferred from file names.
|
|
|
432 .TP
|
|
|
433 .B -u
|
|
|
434 Uncompressed BAM output
|
|
|
435 .RE
|
|
|
436
|
|
|
437 .TP
|
|
|
438 .B index
|
|
|
439 samtools index <aln.bam>
|
|
|
440
|
|
|
441 Index sorted alignment for fast random access. Index file
|
|
|
442 .I <aln.bam>.bai
|
|
|
443 will be created.
|
|
|
444
|
|
|
445 .TP
|
|
|
446 .B idxstats
|
|
|
447 samtools idxstats <aln.bam>
|
|
|
448
|
|
|
449 Retrieve and print stats in the index file. The output is TAB delimited
|
|
|
450 with each line consisting of reference sequence name, sequence length, #
|
|
|
451 mapped reads and # unmapped reads.
|
|
|
452
|
|
|
453 .TP
|
|
|
454 .B faidx
|
|
|
455 samtools faidx <ref.fasta> [region1 [...]]
|
|
|
456
|
|
|
457 Index reference sequence in the FASTA format or extract subsequence from
|
|
|
458 indexed reference sequence. If no region is specified,
|
|
|
459 .B faidx
|
|
|
460 will index the file and create
|
|
|
461 .I <ref.fasta>.fai
|
|
|
462 on the disk. If regions are speficified, the subsequences will be
|
|
|
463 retrieved and printed to stdout in the FASTA format. The input file can
|
|
|
464 be compressed in the
|
|
|
465 .B RAZF
|
|
|
466 format.
|
|
|
467
|
|
|
468 .TP
|
|
|
469 .B fixmate
|
|
|
470 samtools fixmate <in.nameSrt.bam> <out.bam>
|
|
|
471
|
|
|
472 Fill in mate coordinates, ISIZE and mate related flags from a
|
|
|
473 name-sorted alignment.
|
|
|
474
|
|
|
475 .TP
|
|
|
476 .B rmdup
|
|
|
477 samtools rmdup [-sS] <input.srt.bam> <out.bam>
|
|
|
478
|
|
|
479 Remove potential PCR duplicates: if multiple read pairs have identical
|
|
|
480 external coordinates, only retain the pair with highest mapping quality.
|
|
|
481 In the paired-end mode, this command
|
|
|
482 .B ONLY
|
|
|
483 works with FR orientation and requires ISIZE is correctly set. It does
|
|
|
484 not work for unpaired reads (e.g. two ends mapped to different
|
|
|
485 chromosomes or orphan reads).
|
|
|
486
|
|
|
487 .B OPTIONS:
|
|
|
488 .RS
|
|
|
489 .TP 8
|
|
|
490 .B -s
|
|
|
491 Remove duplicate for single-end reads. By default, the command works for
|
|
|
492 paired-end reads only.
|
|
|
493 .TP 8
|
|
|
494 .B -S
|
|
|
495 Treat paired-end reads and single-end reads.
|
|
|
496 .RE
|
|
|
497
|
|
|
498 .TP
|
|
|
499 .B calmd
|
|
|
500 samtools calmd [-EeubSr] [-C capQcoef] <aln.bam> <ref.fasta>
|
|
|
501
|
|
|
502 Generate the MD tag. If the MD tag is already present, this command will
|
|
|
503 give a warning if the MD tag generated is different from the existing
|
|
|
504 tag. Output SAM by default.
|
|
|
505
|
|
|
506 .B OPTIONS:
|
|
|
507 .RS
|
|
|
508 .TP 8
|
|
|
509 .B -A
|
|
|
510 When used jointly with
|
|
|
511 .B -r
|
|
|
512 this option overwrites the original base quality.
|
|
|
513 .TP 8
|
|
|
514 .B -e
|
|
|
515 Convert a the read base to = if it is identical to the aligned reference
|
|
|
516 base. Indel caller does not support the = bases at the moment.
|
|
|
517 .TP
|
|
|
518 .B -u
|
|
|
519 Output uncompressed BAM
|
|
|
520 .TP
|
|
|
521 .B -b
|
|
|
522 Output compressed BAM
|
|
|
523 .TP
|
|
|
524 .B -S
|
|
|
525 The input is SAM with header lines
|
|
|
526 .TP
|
|
|
527 .BI -C \ INT
|
|
|
528 Coefficient to cap mapping quality of poorly mapped reads. See the
|
|
|
529 .B pileup
|
|
|
530 command for details. [0]
|
|
|
531 .TP
|
|
|
532 .B -r
|
|
|
533 Compute the BQ tag (without -A) or cap base quality by BAQ (with -A).
|
|
|
534 .TP
|
|
|
535 .B -E
|
|
|
536 Extended BAQ calculation. This option trades specificity for sensitivity, though the
|
|
|
537 effect is minor.
|
|
|
538 .RE
|
|
|
539
|
|
|
540 .TP
|
|
|
541 .B targetcut
|
|
|
542 samtools targetcut [-Q minBaseQ] [-i inPenalty] [-0 em0] [-1 em1] [-2 em2] [-f ref] <in.bam>
|
|
|
543
|
|
|
544 This command identifies target regions by examining the continuity of read depth, computes
|
|
|
545 haploid consensus sequences of targets and outputs a SAM with each sequence corresponding
|
|
|
546 to a target. When option
|
|
|
547 .B -f
|
|
|
548 is in use, BAQ will be applied. This command is
|
|
|
549 .B only
|
|
|
550 designed for cutting fosmid clones from fosmid pool sequencing [Ref. Kitzman et al. (2010)].
|
|
|
551 .RE
|
|
|
552
|
|
|
553 .TP
|
|
|
554 .B phase
|
|
|
555 samtools phase [-AF] [-k len] [-b prefix] [-q minLOD] [-Q minBaseQ] <in.bam>
|
|
|
556
|
|
|
557 Call and phase heterozygous SNPs.
|
|
|
558 .B OPTIONS:
|
|
|
559 .RS
|
|
|
560 .TP 8
|
|
|
561 .B -A
|
|
|
562 Drop reads with ambiguous phase.
|
|
|
563 .TP 8
|
|
|
564 .BI -b \ STR
|
|
|
565 Prefix of BAM output. When this option is in use, phase-0 reads will be saved in file
|
|
|
566 .BR STR .0.bam
|
|
|
567 and phase-1 reads in
|
|
|
568 .BR STR .1.bam.
|
|
|
569 Phase unknown reads will be randomly allocated to one of the two files. Chimeric reads
|
|
|
570 with switch errors will be saved in
|
|
|
571 .BR STR .chimeric.bam.
|
|
|
572 [null]
|
|
|
573 .TP
|
|
|
574 .B -F
|
|
|
575 Do not attempt to fix chimeric reads.
|
|
|
576 .TP
|
|
|
577 .BI -k \ INT
|
|
|
578 Maximum length for local phasing. [13]
|
|
|
579 .TP
|
|
|
580 .BI -q \ INT
|
|
|
581 Minimum Phred-scaled LOD to call a heterozygote. [40]
|
|
|
582 .TP
|
|
|
583 .BI -Q \ INT
|
|
|
584 Minimum base quality to be used in het calling. [13]
|
|
|
585 .RE
|
|
|
586
|
|
|
587 .SH BCFTOOLS COMMANDS AND OPTIONS
|
|
|
588
|
|
|
589 .TP 10
|
|
|
590 .B view
|
|
|
591 .B bcftools view
|
|
|
592 .RB [ \-AbFGNQSucgv ]
|
|
|
593 .RB [ \-D
|
|
|
594 .IR seqDict ]
|
|
|
595 .RB [ \-l
|
|
|
596 .IR listLoci ]
|
|
|
597 .RB [ \-s
|
|
|
598 .IR listSample ]
|
|
|
599 .RB [ \-i
|
|
|
600 .IR gapSNPratio ]
|
|
|
601 .RB [ \-t
|
|
|
602 .IR mutRate ]
|
|
|
603 .RB [ \-p
|
|
|
604 .IR varThres ]
|
|
|
605 .RB [ \-m
|
|
|
606 .IR varThres ]
|
|
|
607 .RB [ \-P
|
|
|
608 .IR prior ]
|
|
|
609 .RB [ \-1
|
|
|
610 .IR nGroup1 ]
|
|
|
611 .RB [ \-d
|
|
|
612 .IR minFrac ]
|
|
|
613 .RB [ \-U
|
|
|
614 .IR nPerm ]
|
|
|
615 .RB [ \-X
|
|
|
616 .IR permThres ]
|
|
|
617 .RB [ \-T
|
|
|
618 .IR trioType ]
|
|
|
619 .I in.bcf
|
|
|
620 .RI [ region ]
|
|
|
621
|
|
|
622 Convert between BCF and VCF, call variant candidates and estimate allele
|
|
|
623 frequencies.
|
|
|
624
|
|
|
625 .RS
|
|
|
626 .TP
|
|
|
627 .B Input/Output Options:
|
|
|
628 .TP 10
|
|
|
629 .B -A
|
|
|
630 Retain all possible alternate alleles at variant sites. By default, the view
|
|
|
631 command discards unlikely alleles.
|
|
|
632 .TP 10
|
|
|
633 .B -b
|
|
|
634 Output in the BCF format. The default is VCF.
|
|
|
635 .TP
|
|
|
636 .BI -D \ FILE
|
|
|
637 Sequence dictionary (list of chromosome names) for VCF->BCF conversion [null]
|
|
|
638 .TP
|
|
|
639 .B -F
|
|
|
640 Indicate PL is generated by r921 or before (ordering is different).
|
|
|
641 .TP
|
|
|
642 .B -G
|
|
|
643 Suppress all individual genotype information.
|
|
|
644 .TP
|
|
|
645 .BI -l \ FILE
|
|
|
646 List of sites at which information are outputted [all sites]
|
|
|
647 .TP
|
|
|
648 .B -N
|
|
|
649 Skip sites where the REF field is not A/C/G/T
|
|
|
650 .TP
|
|
|
651 .B -Q
|
|
|
652 Output the QCALL likelihood format
|
|
|
653 .TP
|
|
|
654 .BI -s \ FILE
|
|
|
655 List of samples to use. The first column in the input gives the sample names
|
|
|
656 and the second gives the ploidy, which can only be 1 or 2. When the 2nd column
|
|
|
657 is absent, the sample ploidy is assumed to be 2. In the output, the ordering of
|
|
|
658 samples will be identical to the one in
|
|
|
659 .IR FILE .
|
|
|
660 [null]
|
|
|
661 .TP
|
|
|
662 .B -S
|
|
|
663 The input is VCF instead of BCF.
|
|
|
664 .TP
|
|
|
665 .B -u
|
|
|
666 Uncompressed BCF output (force -b).
|
|
|
667 .TP
|
|
|
668 .B Consensus/Variant Calling Options:
|
|
|
669 .TP 10
|
|
|
670 .B -c
|
|
|
671 Call variants using Bayesian inference. This option automatically invokes option
|
|
|
672 .BR -e .
|
|
|
673 .TP
|
|
|
674 .BI -d \ FLOAT
|
|
|
675 When
|
|
|
676 .B -v
|
|
|
677 is in use, skip loci where the fraction of samples covered by reads is below FLOAT. [0]
|
|
|
678 .TP
|
|
|
679 .B -e
|
|
|
680 Perform max-likelihood inference only, including estimating the site allele frequency,
|
|
|
681 testing Hardy-Weinberg equlibrium and testing associations with LRT.
|
|
|
682 .TP
|
|
|
683 .B -g
|
|
|
684 Call per-sample genotypes at variant sites (force -c)
|
|
|
685 .TP
|
|
|
686 .BI -i \ FLOAT
|
|
|
687 Ratio of INDEL-to-SNP mutation rate [0.15]
|
|
|
688 .TP
|
|
|
689 .BI -m \ FLOAT
|
|
|
690 New model for improved multiallelic and rare-variant calling. Another
|
|
|
691 ALT allele is accepted if P(chi^2) of LRT exceeds the FLOAT threshold. The
|
|
|
692 parameter seems robust and the actual value usually does not affect the results
|
|
|
693 much; a good value to use is 0.99. This is the recommended calling method. [0]
|
|
|
694 .TP
|
|
|
695 .BI -p \ FLOAT
|
|
|
696 A site is considered to be a variant if P(ref|D)<FLOAT [0.5]
|
|
|
697 .TP
|
|
|
698 .BI -P \ STR
|
|
|
699 Prior or initial allele frequency spectrum. If STR can be
|
|
|
700 .IR full ,
|
|
|
701 .IR cond2 ,
|
|
|
702 .I flat
|
|
|
703 or the file consisting of error output from a previous variant calling
|
|
|
704 run.
|
|
|
705 .TP
|
|
|
706 .BI -t \ FLOAT
|
|
|
707 Scaled muttion rate for variant calling [0.001]
|
|
|
708 .TP
|
|
|
709 .BI -T \ STR
|
|
|
710 Enable pair/trio calling. For trio calling, option
|
|
|
711 .B -s
|
|
|
712 is usually needed to be applied to configure the trio members and their ordering.
|
|
|
713 In the file supplied to the option
|
|
|
714 .BR -s ,
|
|
|
715 the first sample must be the child, the second the father and the third the mother.
|
|
|
716 The valid values of
|
|
|
717 .I STR
|
|
|
718 are `pair', `trioauto', `trioxd' and `trioxs', where `pair' calls differences between two input samples, and `trioxd' (`trioxs') specifies that the input
|
|
|
719 is from the X chromosome non-PAR regions and the child is a female (male). [null]
|
|
|
720 .TP
|
|
|
721 .B -v
|
|
|
722 Output variant sites only (force -c)
|
|
|
723 .TP
|
|
|
724 .B Contrast Calling and Association Test Options:
|
|
|
725 .TP
|
|
|
726 .BI -1 \ INT
|
|
|
727 Number of group-1 samples. This option is used for dividing the samples into
|
|
|
728 two groups for contrast SNP calling or association test.
|
|
|
729 When this option is in use, the following VCF INFO will be outputted:
|
|
|
730 PC2, PCHI2 and QCHI2. [0]
|
|
|
731 .TP
|
|
|
732 .BI -U \ INT
|
|
|
733 Number of permutations for association test (effective only with
|
|
|
734 .BR -1 )
|
|
|
735 [0]
|
|
|
736 .TP
|
|
|
737 .BI -X \ FLOAT
|
|
|
738 Only perform permutations for P(chi^2)<FLOAT (effective only with
|
|
|
739 .BR -U )
|
|
|
740 [0.01]
|
|
|
741 .RE
|
|
|
742
|
|
|
743 .TP
|
|
|
744 .B index
|
|
|
745 .B bcftools index
|
|
|
746 .I in.bcf
|
|
|
747
|
|
|
748 Index sorted BCF for random access.
|
|
|
749 .RE
|
|
|
750
|
|
|
751 .TP
|
|
|
752 .B cat
|
|
|
753 .B bcftools cat
|
|
|
754 .I in1.bcf
|
|
|
755 .RI [ "in2.bcf " [ ... "]]]"
|
|
|
756
|
|
|
757 Concatenate BCF files. The input files are required to be sorted and
|
|
|
758 have identical samples appearing in the same order.
|
|
|
759 .RE
|
|
|
760 .SH SAM FORMAT
|
|
|
761
|
|
|
762 Sequence Alignment/Map (SAM) format is TAB-delimited. Apart from the header lines, which are started
|
|
|
763 with the `@' symbol, each alignment line consists of:
|
|
|
764
|
|
|
765 .TS
|
|
|
766 center box;
|
|
|
767 cb | cb | cb
|
|
|
768 n | l | l .
|
|
|
769 Col Field Description
|
|
|
770 _
|
|
|
771 1 QNAME Query template/pair NAME
|
|
|
772 2 FLAG bitwise FLAG
|
|
|
773 3 RNAME Reference sequence NAME
|
|
|
774 4 POS 1-based leftmost POSition/coordinate of clipped sequence
|
|
|
775 5 MAPQ MAPping Quality (Phred-scaled)
|
|
|
776 6 CIAGR extended CIGAR string
|
|
|
777 7 MRNM Mate Reference sequence NaMe (`=' if same as RNAME)
|
|
|
778 8 MPOS 1-based Mate POSistion
|
|
|
779 9 TLEN inferred Template LENgth (insert size)
|
|
|
780 10 SEQ query SEQuence on the same strand as the reference
|
|
|
781 11 QUAL query QUALity (ASCII-33 gives the Phred base quality)
|
|
|
782 12+ OPT variable OPTional fields in the format TAG:VTYPE:VALUE
|
|
|
783 .TE
|
|
|
784
|
|
|
785 .PP
|
|
|
786 Each bit in the FLAG field is defined as:
|
|
|
787
|
|
|
788 .TS
|
|
|
789 center box;
|
|
|
790 cb | cb | cb
|
|
|
791 l | c | l .
|
|
|
792 Flag Chr Description
|
|
|
793 _
|
|
|
794 0x0001 p the read is paired in sequencing
|
|
|
795 0x0002 P the read is mapped in a proper pair
|
|
|
796 0x0004 u the query sequence itself is unmapped
|
|
|
797 0x0008 U the mate is unmapped
|
|
|
798 0x0010 r strand of the query (1 for reverse)
|
|
|
799 0x0020 R strand of the mate
|
|
|
800 0x0040 1 the read is the first read in a pair
|
|
|
801 0x0080 2 the read is the second read in a pair
|
|
|
802 0x0100 s the alignment is not primary
|
|
|
803 0x0200 f the read fails platform/vendor quality checks
|
|
|
804 0x0400 d the read is either a PCR or an optical duplicate
|
|
|
805 .TE
|
|
|
806
|
|
|
807 where the second column gives the string representation of the FLAG field.
|
|
|
808
|
|
|
809 .SH VCF FORMAT
|
|
|
810
|
|
|
811 The Variant Call Format (VCF) is a TAB-delimited format with each data line consists of the following fields:
|
|
|
812 .TS
|
|
|
813 center box;
|
|
|
814 cb | cb | cb
|
|
|
815 n | l | l .
|
|
|
816 Col Field Description
|
|
|
817 _
|
|
|
818 1 CHROM CHROMosome name
|
|
|
819 2 POS the left-most POSition of the variant
|
|
|
820 3 ID unique variant IDentifier
|
|
|
821 4 REF the REFerence allele
|
|
|
822 5 ALT the ALTernate allele(s), separated by comma
|
|
|
823 6 QUAL variant/reference QUALity
|
|
|
824 7 FILTER FILTers applied
|
|
|
825 8 INFO INFOrmation related to the variant, separated by semi-colon
|
|
|
826 9 FORMAT FORMAT of the genotype fields, separated by colon (optional)
|
|
|
827 10+ SAMPLE SAMPLE genotypes and per-sample information (optional)
|
|
|
828 .TE
|
|
|
829
|
|
|
830 .PP
|
|
|
831 The following table gives the
|
|
|
832 .B INFO
|
|
|
833 tags used by samtools and bcftools.
|
|
|
834
|
|
|
835 .TS
|
|
|
836 center box;
|
|
|
837 cb | cb | cb
|
|
|
838 l | l | l .
|
|
|
839 Tag Format Description
|
|
|
840 _
|
|
|
841 AF1 double Max-likelihood estimate of the site allele frequency (AF) of the first ALT allele
|
|
|
842 DP int Raw read depth (without quality filtering)
|
|
|
843 DP4 int[4] # high-quality reference forward bases, ref reverse, alternate for and alt rev bases
|
|
|
844 FQ int Consensus quality. Positive: sample genotypes different; negative: otherwise
|
|
|
845 MQ int Root-Mean-Square mapping quality of covering reads
|
|
|
846 PC2 int[2] Phred probability of AF in group1 samples being larger (,smaller) than in group2
|
|
|
847 PCHI2 double Posterior weighted chi^2 P-value between group1 and group2 samples
|
|
|
848 PV4 double[4] P-value for strand bias, baseQ bias, mapQ bias and tail distance bias
|
|
|
849 QCHI2 int Phred-scaled PCHI2
|
|
|
850 RP int # permutations yielding a smaller PCHI2
|
|
|
851 CLR int Phred log ratio of genotype likelihoods with and without the trio/pair constraint
|
|
|
852 UGT string Most probable genotype configuration without the trio constraint
|
|
|
853 CGT string Most probable configuration with the trio constraint
|
|
|
854 VDB float Tests variant positions within reads. Intended for filtering RNA-seq artifacts around splice sites
|
|
|
855 RPB float Mann-Whitney rank-sum test for tail distance bias
|
|
|
856 HWE float Hardy-Weinberg equilibrium test, Wigginton et al., PMID: 15789306
|
|
|
857 .TE
|
|
|
858
|
|
|
859 .SH EXAMPLES
|
|
|
860 .IP o 2
|
|
|
861 Import SAM to BAM when
|
|
|
862 .B @SQ
|
|
|
863 lines are present in the header:
|
|
|
864
|
|
|
865 samtools view -bS aln.sam > aln.bam
|
|
|
866
|
|
|
867 If
|
|
|
868 .B @SQ
|
|
|
869 lines are absent:
|
|
|
870
|
|
|
871 samtools faidx ref.fa
|
|
|
872 samtools view -bt ref.fa.fai aln.sam > aln.bam
|
|
|
873
|
|
|
874 where
|
|
|
875 .I ref.fa.fai
|
|
|
876 is generated automatically by the
|
|
|
877 .B faidx
|
|
|
878 command.
|
|
|
879
|
|
|
880 .IP o 2
|
|
|
881 Attach the
|
|
|
882 .B RG
|
|
|
883 tag while merging sorted alignments:
|
|
|
884
|
|
|
885 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
|
|
|
886 samtools merge -rh rg.txt merged.bam ga.bam 454.bam
|
|
|
887
|
|
|
888 The value in a
|
|
|
889 .B RG
|
|
|
890 tag is determined by the file name the read is coming from. In this
|
|
|
891 example, in the
|
|
|
892 .IR merged.bam ,
|
|
|
893 reads from
|
|
|
894 .I ga.bam
|
|
|
895 will be attached
|
|
|
896 .IR RG:Z:ga ,
|
|
|
897 while reads from
|
|
|
898 .I 454.bam
|
|
|
899 will be attached
|
|
|
900 .IR RG:Z:454 .
|
|
|
901
|
|
|
902 .IP o 2
|
|
|
903 Call SNPs and short INDELs for one diploid individual:
|
|
|
904
|
|
|
905 samtools mpileup -ugf ref.fa aln.bam | bcftools view -bvcg - > var.raw.bcf
|
|
|
906 bcftools view var.raw.bcf | vcfutils.pl varFilter -D 100 > var.flt.vcf
|
|
|
907
|
|
|
908 The
|
|
|
909 .B -D
|
|
|
910 option of varFilter controls the maximum read depth, which should be
|
|
|
911 adjusted to about twice the average read depth. One may consider to add
|
|
|
912 .B -C50
|
|
|
913 to
|
|
|
914 .B mpileup
|
|
|
915 if mapping quality is overestimated for reads containing excessive
|
|
|
916 mismatches. Applying this option usually helps
|
|
|
917 .B BWA-short
|
|
|
918 but may not other mappers.
|
|
|
919
|
|
|
920 .IP o 2
|
|
|
921 Generate the consensus sequence for one diploid individual:
|
|
|
922
|
|
|
923 samtools mpileup -uf ref.fa aln.bam | bcftools view -cg - | vcfutils.pl vcf2fq > cns.fq
|
|
|
924
|
|
|
925 .IP o 2
|
|
|
926 Call somatic mutations from a pair of samples:
|
|
|
927
|
|
|
928 samtools mpileup -DSuf ref.fa aln.bam | bcftools view -bvcgT pair - > var.bcf
|
|
|
929
|
|
|
930 In the output INFO field,
|
|
|
931 .I CLR
|
|
|
932 gives the Phred-log ratio between the likelihood by treating the
|
|
|
933 two samples independently, and the likelihood by requiring the genotype to be identical.
|
|
|
934 This
|
|
|
935 .I CLR
|
|
|
936 is effectively a score measuring the confidence of somatic calls. The higher the better.
|
|
|
937
|
|
|
938 .IP o 2
|
|
|
939 Call de novo and somatic mutations from a family trio:
|
|
|
940
|
|
|
941 samtools mpileup -DSuf ref.fa aln.bam | bcftools view -bvcgT pair -s samples.txt - > var.bcf
|
|
|
942
|
|
|
943 File
|
|
|
944 .I samples.txt
|
|
|
945 should consist of three lines specifying the member and order of samples (in the order of child-father-mother).
|
|
|
946 Similarly,
|
|
|
947 .I CLR
|
|
|
948 gives the Phred-log likelihood ratio with and without the trio constraint.
|
|
|
949 .I UGT
|
|
|
950 shows the most likely genotype configuration without the trio constraint, and
|
|
|
951 .I CGT
|
|
|
952 gives the most likely genotype configuration satisfying the trio constraint.
|
|
|
953
|
|
|
954 .IP o 2
|
|
|
955 Phase one individual:
|
|
|
956
|
|
|
957 samtools calmd -AEur aln.bam ref.fa | samtools phase -b prefix - > phase.out
|
|
|
958
|
|
|
959 The
|
|
|
960 .B calmd
|
|
|
961 command is used to reduce false heterozygotes around INDELs.
|
|
|
962
|
|
|
963 .IP o 2
|
|
|
964 Call SNPs and short indels for multiple diploid individuals:
|
|
|
965
|
|
|
966 samtools mpileup -P ILLUMINA -ugf ref.fa *.bam | bcftools view -bcvg - > var.raw.bcf
|
|
|
967 bcftools view var.raw.bcf | vcfutils.pl varFilter -D 2000 > var.flt.vcf
|
|
|
968
|
|
|
969 Individuals are identified from the
|
|
|
970 .B SM
|
|
|
971 tags in the
|
|
|
972 .B @RG
|
|
|
973 header lines. Individuals can be pooled in one alignment file; one
|
|
|
974 individual can also be separated into multiple files. The
|
|
|
975 .B -P
|
|
|
976 option specifies that indel candidates should be collected only from
|
|
|
977 read groups with the
|
|
|
978 .B @RG-PL
|
|
|
979 tag set to
|
|
|
980 .IR ILLUMINA .
|
|
|
981 Collecting indel candidates from reads sequenced by an indel-prone
|
|
|
982 technology may affect the performance of indel calling.
|
|
|
983
|
|
|
984 Note that there is a new calling model which can be invoked by
|
|
|
985
|
|
|
986 bcftools view -m0.99 ...
|
|
|
987
|
|
|
988 which fixes some severe limitations of the default method.
|
|
|
989
|
|
|
990 For filtering, best results seem to be achieved by first applying the
|
|
|
991 .IR SnpGap
|
|
|
992 filter and then applying some machine learning approach
|
|
|
993
|
|
|
994 vcf-annotate -f SnpGap=n
|
|
|
995 vcf filter ...
|
|
|
996
|
|
|
997 Both can be found in the
|
|
|
998 .B vcftools
|
|
|
999 and
|
|
|
1000 .B htslib
|
|
|
1001 package (links below).
|
|
|
1002
|
|
|
1003 .IP o 2
|
|
|
1004 Derive the allele frequency spectrum (AFS) on a list of sites from multiple individuals:
|
|
|
1005
|
|
|
1006 samtools mpileup -Igf ref.fa *.bam > all.bcf
|
|
|
1007 bcftools view -bl sites.list all.bcf > sites.bcf
|
|
|
1008 bcftools view -cGP cond2 sites.bcf > /dev/null 2> sites.1.afs
|
|
|
1009 bcftools view -cGP sites.1.afs sites.bcf > /dev/null 2> sites.2.afs
|
|
|
1010 bcftools view -cGP sites.2.afs sites.bcf > /dev/null 2> sites.3.afs
|
|
|
1011 ......
|
|
|
1012
|
|
|
1013 where
|
|
|
1014 .I sites.list
|
|
|
1015 contains the list of sites with each line consisting of the reference
|
|
|
1016 sequence name and position. The following
|
|
|
1017 .B bcftools
|
|
|
1018 commands estimate AFS by EM.
|
|
|
1019
|
|
|
1020 .IP o 2
|
|
|
1021 Dump BAQ applied alignment for other SNP callers:
|
|
|
1022
|
|
|
1023 samtools calmd -bAr aln.bam > aln.baq.bam
|
|
|
1024
|
|
|
1025 It adds and corrects the
|
|
|
1026 .B NM
|
|
|
1027 and
|
|
|
1028 .B MD
|
|
|
1029 tags at the same time. The
|
|
|
1030 .B calmd
|
|
|
1031 command also comes with the
|
|
|
1032 .B -C
|
|
|
1033 option, the same as the one in
|
|
|
1034 .B pileup
|
|
|
1035 and
|
|
|
1036 .BR mpileup .
|
|
|
1037 Apply if it helps.
|
|
|
1038
|
|
|
1039 .SH LIMITATIONS
|
|
|
1040 .PP
|
|
|
1041 .IP o 2
|
|
|
1042 Unaligned words used in bam_import.c, bam_endian.h, bam.c and bam_aux.c.
|
|
|
1043 .IP o 2
|
|
|
1044 Samtools paired-end rmdup does not work for unpaired reads (e.g. orphan
|
|
|
1045 reads or ends mapped to different chromosomes). If this is a concern,
|
|
|
1046 please use Picard's MarkDuplicate which correctly handles these cases,
|
|
|
1047 although a little slower.
|
|
|
1048
|
|
|
1049 .SH AUTHOR
|
|
|
1050 .PP
|
|
|
1051 Heng Li from the Sanger Institute wrote the C version of samtools. Bob
|
|
|
1052 Handsaker from the Broad Institute implemented the BGZF library and Jue
|
|
|
1053 Ruan from Beijing Genomics Institute wrote the RAZF library. John
|
|
|
1054 Marshall and Petr Danecek contribute to the source code and various
|
|
|
1055 people from the 1000 Genomes Project have contributed to the SAM format
|
|
|
1056 specification.
|
|
|
1057
|
|
|
1058 .SH SEE ALSO
|
|
|
1059 .PP
|
|
|
1060 Samtools website: <http://samtools.sourceforge.net>
|
|
|
1061 .br
|
|
|
1062 Samtools latest source: <https://github.com/samtools/samtools>
|
|
|
1063 .br
|
|
|
1064 VCFtools website with stable link to VCF specification: <http://vcftools.sourceforge.net>
|
|
|
1065 .br
|
|
|
1066 HTSlib website: <https://github.com/samtools/htslib>
|
