Package Bio :: Package Sequencing :: Package Applications :: Module _samtools
[hide private]
[frames] | no frames]

Source Code for Module Bio.Sequencing.Applications._samtools

  1  """Command line wrapper for samtools""" 
  2  # Last Checked with samtools [0.1.20 and 1.2] 
  3  # TODO samtools 1.x has additional options over 0.x which 
  4  # are missing from this wrapper 
  5   
  6  from __future__ import print_function 
  7   
  8  from Bio.Application import _Option, _Argument, _Switch 
  9  from Bio.Application import AbstractCommandline, _ArgumentList 
 10  from Bio.Application import _StaticArgument 
 11   
 12   
13 -class SamtoolsViewCommandline(AbstractCommandline):
14 15 """Command line wrapper for samtools view. 16 17 Extract/print all or sub alignments in SAM or BAM format, equivalent to:: 18 19 $ samtools view [-bchuHS] [-t in.refList] [-o output] [-f reqFlag] 20 [-F skipFlag] [-q minMapQ] [-l library] [-r readGroup] 21 [-R rgFile] <in.bam>|<in.sam> [region1 [...]] 22 23 See http://samtools.sourceforge.net/samtools.shtml for more details 24 25 Example 26 ------- 27 28 >>> from Bio.Sequencing.Applications import SamtoolsViewCommandline 29 >>> input_file = "/path/to/sam_or_bam_file" 30 >>> samtools_view_cmd = SamtoolsViewCommandline(input_file=input_file) 31 >>> print(samtools_view_cmd) 32 samtools view /path/to/sam_or_bam_file 33 34 """
35 - def __init__(self, cmd="samtools", **kwargs):
36 self.program_name = cmd 37 self.parameters = [ 38 _StaticArgument("view"), 39 _Switch(["-b", "b"], "Output in the BAM format"), 40 _Switch(["-c", "c"], 41 """Instead of printing the alignments, only count them and 42 print the total number. 43 44 All filter options, such as '-f', '-F' and '-q', 45 are taken into account"""), 46 _Switch(["-h", "h"], "Include the header in the output"), 47 _Switch(["-u", "u"], 48 """Output uncompressed BAM. 49 50 This option saves time spent on compression/decompression 51 and is thus preferred when the output is piped to 52 another samtools command"""), 53 _Switch(["-H", "H"], "Output the header only"), 54 _Switch(["-S", "S"], 55 """Input is in SAM. 56 If @SQ header lines are absent, 57 the '-t' option is required."""), 58 _Option(["-t", "t"], 59 """This file is TAB-delimited. 60 Each line must contain the reference name and the 61 length of the reference, one line for each 62 distinct reference; additional fields are ignored. 63 64 This file also defines the order of the reference 65 sequences in sorting. 66 If you run 'samtools faidx <ref.fa>', 67 the resultant index file <ref.fa>.fai can be used 68 as this <in.ref_list> file.""", 69 filename=True, equate=False, 70 checker_function=lambda x: isinstance(x, str)), 71 _Option(["-o", "o"], "Output file", 72 filename=True, equate=False, 73 checker_function=lambda x: isinstance(x, str)), 74 _Option(["-f", "f"], 75 """Only output alignments with all bits in 76 INT present in the FLAG field""", 77 equate=False, 78 checker_function=lambda x: isinstance(x, int)), 79 _Option(["-F", "F"], 80 "Skip alignments with bits present in INT", 81 equate=False, 82 checker_function=lambda x: isinstance(x, int)), 83 _Option(["-q", "q"], 84 "Skip alignments with MAPQ smaller than INT", 85 equate=False, 86 checker_function=lambda x: isinstance(x, int)), 87 _Option(["-r", "r"], 88 "Only output reads in read group STR", 89 equate=False, 90 checker_function=lambda x: isinstance(x, str)), 91 _Option(["-R", "R"], 92 "Output reads in read groups listed in FILE", 93 filename=True, equate=False, 94 checker_function=lambda x: isinstance(x, str)), 95 _Option(["-l", "l"], 96 "Only output reads in library STR", 97 equate=False, 98 checker_function=lambda x: isinstance(x, str)), 99 _Switch(["-1", "fast_bam"], 100 "Use zlib compression level 1 to compress the output"), 101 _Argument(["input", "input_file"], 102 "Input File Name", filename=True, is_required=True), 103 _Argument(["region"], "Region", is_required=False), 104 ] 105 AbstractCommandline.__init__(self, cmd, **kwargs)
106 107
108 -class SamtoolsMpileupCommandline(AbstractCommandline):
109 """Command line wrapper for samtools mpileup. 110 111 Generate BCF or pileup for one or multiple BAM files, equivalent to:: 112 113 $ samtools mpileup [-EBug] [-C capQcoef] [-r reg] [-f in.fa] 114 [-l list] [-M capMapQ] [-Q minBaseQ] 115 [-q minMapQ] in.bam [in2.bam [...]] 116 117 See http://samtools.sourceforge.net/samtools.shtml for more details 118 119 Example: 120 121 >>> from Bio.Sequencing.Applications import SamtoolsMpileupCommandline 122 >>> input = ["/path/to/sam_or_bam_file"] 123 >>> samtools_mpileup_cmd = SamtoolsMpileupCommandline(input_file=input) 124 >>> print(samtools_mpileup_cmd) 125 samtools mpileup /path/to/sam_or_bam_file 126 127 """
128 - def __init__(self, cmd="samtools", **kwargs):
129 self.program_name = cmd 130 self.parameters = [ 131 _StaticArgument("mpileup"), 132 _Switch(["-E", "E"], 133 """Extended BAQ computation. 134 This option helps sensitivity especially 135 for MNPs, but may hurt specificity a little bit"""), 136 _Switch(["-B", "B"], 137 """Disable probabilistic realignment for the 138 computation of base alignment quality (BAQ). 139 140 BAQ is the Phred-scaled probability of a read base being 141 misaligned. 142 Applying this option greatly helps to reduce false SNPs 143 caused by misalignments"""), 144 _Switch(["-g", "g"], 145 """Compute genotype likelihoods and output them in the 146 binary call format (BCF)"""), 147 _Switch(["-u", "u"], 148 """Similar to -g except that the output is 149 uncompressed BCF, which is preferred for piping"""), 150 _Option(["-C", "C"], 151 """Coefficient for downgrading mapping quality for 152 reads containing excessive mismatches. 153 154 Given a read with a phred-scaled probability q of 155 being generated from the mapped position, 156 the new mapping quality is about sqrt((INT-q)/INT)*INT. 157 A zero value disables this functionality; 158 if enabled, the recommended value for BWA is 50""", 159 equate=False, 160 checker_function=lambda x: isinstance(x, int)), 161 _Option(["-r", "r"], 162 "Only generate pileup in region STR", 163 equate=False, 164 checker_function=lambda x: isinstance(x, str)), 165 _Option(["-f", "f"], 166 """The faidx-indexed reference file in the FASTA format. 167 168 The file can be optionally compressed by razip""", 169 filename=True, equate=False, 170 checker_function=lambda x: isinstance(x, str)), 171 _Option(["-l", "l"], 172 """BED or position list file containing a list of regions 173 or sites where pileup or BCF should be generated""", 174 filename=True, equate=False, 175 checker_function=lambda x: isinstance(x, str)), 176 _Option(["-M", "M"], 177 "Cap Mapping Quality at M", 178 equate=False, 179 checker_function=lambda x: isinstance(x, int)), 180 _Option(["-q", "q"], 181 "Minimum mapping quality for an alignment to be used", 182 equate=False, 183 checker_function=lambda x: isinstance(x, int)), 184 _Option(["-Q", "Q"], 185 "Minimum base quality for a base to be considered", 186 equate=False, 187 checker_function=lambda x: isinstance(x, int)), 188 _Switch(["-6", "illumina_13"], 189 "Assume the quality is in the Illumina 1.3+ encoding"), 190 _Switch(["-A", "A"], 191 "Do not skip anomalous read pairs in variant calling."), 192 _Option(["-b", "b"], 193 "List of input BAM files, one file per line", 194 filename=True, equate=False, 195 checker_function=lambda x: isinstance(x, str)), 196 _Option(["-d", "d"], 197 "At a position, read maximally INT reads per input BAM", 198 equate=False, 199 checker_function=lambda x: isinstance(x, int)), 200 _Switch(["-D", "D"], "Output per-sample read depth"), 201 _Switch(["-S", "S"], """Output per-sample Phred-scaled 202 strand bias P-value"""), 203 _Option(["-e", "e"], 204 """Phred-scaled gap extension sequencing error probability. 205 206 Reducing INT leads to longer indels""", 207 equate=False, 208 checker_function=lambda x: isinstance(x, int)), 209 _Option(["-h", "h"], 210 """Coefficient for modeling homopolymer errors. 211 212 Given an l-long homopolymer run, the sequencing error 213 of an indel of size s is modeled as INT*s/l""", 214 equate=False, 215 checker_function=lambda x: isinstance(x, int)), 216 _Switch(["-I", "I"], "Do not perform INDEL calling"), 217 _Option(["-L", "L"], 218 """Skip INDEL calling if the average per-sample 219 depth is above INT""", 220 equate=False, 221 checker_function=lambda x: isinstance(x, int)), 222 _Option(["-o", "o"], 223 """Phred-scaled gap open sequencing error probability. 224 225 Reducing INT leads to more indel calls.""", 226 equate=False, 227 checker_function=lambda x: isinstance(x, int)), 228 _Option(["-p", "p"], 229 """Comma delimited list of platforms (determined by @RG-PL) 230 from which indel candidates are obtained. 231 232 It is recommended to collect indel candidates from 233 sequencing technologies that have low indel error rate 234 such as ILLUMINA""", 235 equate=False, 236 checker_function=lambda x: isinstance(x, str)), 237 _ArgumentList(["input_file"], 238 "Input File for generating mpileup", 239 filename=True, is_required=True), 240 241 ] 242 AbstractCommandline.__init__(self, cmd, **kwargs)
243 244
245 -class SamtoolsReheaderCommandline(AbstractCommandline):
246 """Command line wrapper for samtools reheader. 247 248 Replace the header in in.bam with the header 249 in in.header.sam, equivalent to:: 250 251 $ samtools reheader <in.header.sam> <in.bam> 252 253 See http://samtools.sourceforge.net/samtools.shtml for more details 254 255 Example: 256 257 >>> from Bio.Sequencing.Applications import SamtoolsReheaderCommandline 258 >>> input_header = "/path/to/header_sam_file" 259 >>> input_bam = "/path/to/input_bam_file" 260 >>> samtools_reheader_cmd = SamtoolsReheaderCommandline(\ 261 input_header=input_header,\ 262 input_bam=input_bam) 263 >>> print(samtools_reheader_cmd) 264 samtools reheader /path/to/header_sam_file /path/to/input_bam_file 265 266 """
267 - def __init__(self, cmd="samtools", **kwargs):
268 self.program_name = cmd 269 self.parameters = [ 270 _StaticArgument("reheader"), 271 _Argument(["input_header", "header_sam", "sam_file"], 272 "Sam file with header", 273 filename=True, is_required=True), 274 _Argument(["input_bam", "input_file", "bam_file"], 275 "BAM file for writing header to", 276 filename=True, is_required=True) 277 ] 278 AbstractCommandline.__init__(self, cmd, **kwargs)
279 280
281 -class SamtoolsCatCommandline(AbstractCommandline):
282 """Command line wrapper for samtools cat. 283 284 Concatenate BAMs, equivalent to:: 285 286 $ samtools cat [-h header.sam] [-o out.bam] <in1.bam> <in2.bam> [ ... ] 287 288 See http://samtools.sourceforge.net/samtools.shtml for more details 289 290 Example: 291 292 >>> from Bio.Sequencing.Applications import SamtoolsCatCommandline 293 >>> input_bam1 = "/path/to/input_bam1" 294 >>> input_bam2 = "/path/to/input_bam2" 295 >>> input_bams = [input_bam1, input_bam2] 296 >>> samtools_cat_cmd = SamtoolsCatCommandline(input_bam=input_bams) 297 >>> print(samtools_cat_cmd) 298 samtools cat /path/to/input_bam1 /path/to/input_bam2 299 300 """
301 - def __init__(self, cmd="samtools", **kwargs):
302 self.program_name = cmd 303 self.parameters = [ 304 _StaticArgument("cat"), 305 _Option(["-h", "h"], "Header SAM file", 306 filename=True, equate=False, 307 checker_function=lambda x: isinstance(x, str)), 308 _Option(["-o", "o"], "Output SAM file", 309 filename=True, equate=False, 310 checker_function=lambda x: isinstance(x, str)), 311 _ArgumentList(["input", "input_bam", "bams"], "Input BAM files", 312 filename=True, is_required=True) 313 314 ] 315 AbstractCommandline.__init__(self, cmd, **kwargs)
316 317
318 -class SamtoolsVersion0xSortCommandline(AbstractCommandline):
319 """Command line wrapper for samtools version 0.1.x sort. 320 321 Concatenate BAMs, equivalent to:: 322 323 $ samtools sort [-no] [-m maxMem] <in.bam> <out.prefix> 324 325 See http://samtools.sourceforge.net/samtools.shtml for more details 326 327 Example 328 ------- 329 330 >>> from Bio.Sequencing.Applications import SamtoolsVersion0xSortCommandline 331 >>> input_bam = "/path/to/input_bam" 332 >>> out_prefix = "/path/to/out_prefix" 333 >>> samtools_sort_cmd = SamtoolsVersion0xSortCommandline(input=input_bam, out_prefix=out_prefix) 334 >>> print(samtools_sort_cmd) 335 samtools sort /path/to/input_bam /path/to/out_prefix 336 337 """
338 - def __init__(self, cmd="samtools", **kwargs):
339 self.program_name = cmd 340 341 # options for version samtools 0.0.19 342 self.parameters = [ 343 _StaticArgument("sort"), 344 _Switch(["-o", "o"], """Output the final alignment 345 to the standard output"""), 346 _Switch(["-n", "n"], """Sort by read names rather 347 than by chromosomal coordinates"""), 348 _Option(["-m", "m"], "Approximately the maximum required memory", 349 equate=False, 350 checker_function=lambda x: isinstance(x, int)), 351 _Argument(["input"], "Input BAM file", 352 filename=True, is_required=True), 353 _Argument(["out_prefix"], "Output prefix", 354 filename=True, is_required=True), 355 ] 356 AbstractCommandline.__init__(self, cmd, **kwargs)
357 358
359 -class SamtoolsVersion1xSortCommandline(AbstractCommandline):
360 """Command line wrapper for samtools version 1.3.x sort. 361 362 Concatenate BAMs, equivalent to:: 363 364 $ samtools sort [-n] [-T FREFIX] [-o file] [-I INT] [-m maxMem] <in.bam> 365 366 See http://samtools.sourceforge.net/samtools.shtml for more details 367 368 Example 369 ------- 370 371 >>> from Bio.Sequencing.Applications import SamtoolsVersion1xSortCommandline 372 >>> input_bam = "/path/to/input_bam" 373 >>> FREFIX = "/path/to/out_prefix" 374 >>> file_name = "/path/to/out_file" 375 >>> samtools_sort_cmd = SamtoolsVersion1xSortCommandline(input=input_bam, T=FREFIX, o=file_name) 376 >>> print(samtools_sort_cmd) 377 samtools sort -o /path/to/out_file -T /path/to/out_prefix /path/to/input_bam 378 379 """
380 - def __init__(self, cmd="samtools", **kwargs):
381 self.program_name = cmd 382 383 # options for version samtools 1.3.1 384 self.parameters = [ 385 _StaticArgument("sort"), 386 _Switch(["-n", "n"], """Sort by read names rather 387 than by chromosomal coordinates"""), 388 _Option(["-o", "o"], """(file) Write the final sorted output to FILE, 389 rather than to standard output""", 390 equate=False, checker_function=lambda x: isinstance(x, str)), 391 _Option(["-O", "O"], """(FORMAT) Write the final output as sam, bam, or cram""", 392 equate=False, checker_function=lambda x: isinstance(x, str)), 393 _Option(["-T", "T"], """(PREFIX) Write temporary files to PREFIX.nnnn.bam, or if the specified PREFIX 394 is an existing directory, to PREFIX/samtools.mmm.mmm.tmp.nnnn.bam, 395 where mmm is unique to this invocation of the sort command""", 396 equate=False, checker_function=lambda x: isinstance(x, str)), 397 _Option(["-I", "I"], """(INT) Set the desired compression level for the final output file, 398 ranging from 0 (uncompressed) or 1 (fastest but minimal compression) 399 to 9 (best compression but slowest to write), similarly to gzip(1)'s compression level setting.""", 400 equate=False, checker_function=lambda x: isinstance(x, str)), 401 _Option(["-m", "m"], "Approximately the maximum required memory", 402 equate=False, 403 checker_function=lambda x: isinstance(x, int)), 404 _Argument(["input"], "Input SAM/BAM/CRAM file", 405 filename=True, is_required=True), 406 ] 407 AbstractCommandline.__init__(self, cmd, **kwargs)
408 409
410 -class SamtoolsMergeCommandline(AbstractCommandline):
411 """Command line wrapper for samtools merge. 412 413 Merge multiple sorted alignments, equivalent to:: 414 415 $ samtools merge [-nur1f] [-h inh.sam] [-R reg] 416 <out.bam> <in1.bam> <in2.bam> [...] 417 418 See http://samtools.sourceforge.net/samtools.shtml for more details 419 420 Example 421 ------- 422 423 >>> from Bio.Sequencing.Applications import SamtoolsMergeCommandline 424 >>> out_bam = "/path/to/out_bam" 425 >>> in_bam = ["/path/to/input_bam1", "/path/to/input_bam2"] 426 >>> samtools_merge_cmd = SamtoolsMergeCommandline(\ 427 out_bam=out_bam,\ 428 input_bam=in_bam) 429 >>> print(samtools_merge_cmd) 430 samtools merge /path/to/out_bam /path/to/input_bam1 /path/to/input_bam2 431 432 """
433 - def __init__(self, cmd="samtools", **kwargs):
434 self.program_name = cmd 435 self.parameters = [ 436 _StaticArgument("merge"), 437 _Switch(["-n", "n"], 438 """The input alignments are sorted by read names 439 rather than by chromosomal coordinates"""), 440 _Switch(["-r", "r"], """Attach an RG tag to each alignment. 441 The tag value is inferred from file names"""), 442 _Switch(["-u", "u"], "Uncompressed BAM output"), 443 _Switch(["-1", "fast_bam"], """Use zlib compression level 1 444 to compress the output"""), 445 _Switch(["-f", "f"], """Force to overwrite the 446 output file if present"""), 447 _Option(["-h", "h"], """Use the lines of FILE as '@' 448 headers to be copied to out.bam""", 449 filename=True, equate=False, 450 checker_function=lambda x: isinstance(x, str)), 451 _Option(["-R", "R"], 452 "Merge files in the specified region indicated by STR", 453 equate=False, 454 checker_function=lambda x: isinstance(x, str)), 455 _Argument(["output_bam", "out_bam", "out", "output"], 456 "Output BAM file", 457 filename=True, is_required=True), 458 _ArgumentList(["input_bam", "in_bam", "input", "bam"], 459 "Input BAM", 460 filename=True, is_required=True), 461 ] 462 AbstractCommandline.__init__(self, cmd, **kwargs)
463 464
465 -class SamtoolsIndexCommandline(AbstractCommandline):
466 """Command line wrapper for samtools index. 467 468 Index sorted alignment for fast random access, equivalent to:: 469 470 $ samtools index <aln.bam> 471 472 See http://samtools.sourceforge.net/samtools.shtml for more details 473 474 Example: 475 476 >>> from Bio.Sequencing.Applications import SamtoolsIndexCommandline 477 >>> input = "/path/to/aln_bam" 478 >>> samtools_index_cmd = SamtoolsIndexCommandline(input_bam=input) 479 >>> print(samtools_index_cmd) 480 samtools index /path/to/aln_bam 481 482 """
483 - def __init__(self, cmd="samtools", **kwargs):
484 self.program_name = cmd 485 self.parameters = [ 486 _StaticArgument("index"), 487 _Argument(["input", "in_bam", "input_bam"], 488 "BAM file to be indexed"), 489 ] 490 AbstractCommandline.__init__(self, cmd, **kwargs)
491 492
493 -class SamtoolsIdxstatsCommandline(AbstractCommandline):
494 """Command line wrapper for samtools idxstats. 495 496 Retrieve and print stats in the index file, equivalent to:: 497 498 $ samtools idxstats <aln.bam> 499 500 See http://samtools.sourceforge.net/samtools.shtml for more details 501 502 Example: 503 504 >>> from Bio.Sequencing.Applications import SamtoolsIdxstatsCommandline 505 >>> input = "/path/to/aln_bam" 506 >>> samtools_idxstats_cmd = SamtoolsIdxstatsCommandline(input_bam=input) 507 >>> print(samtools_idxstats_cmd) 508 samtools idxstats /path/to/aln_bam 509 510 """
511 - def __init__(self, cmd="samtools", **kwargs):
512 self.program_name = cmd 513 self.parameters = [ 514 _StaticArgument("idxstats"), 515 _Argument(["input", "in_bam", "input_bam"], 516 "BAM file to be indexed") 517 ] 518 AbstractCommandline.__init__(self, cmd, **kwargs)
519 520
521 -class SamtoolsFaidxCommandline(AbstractCommandline):
522 """Command line wrapper for samtools faidx. 523 524 Retrieve and print stats in the index file, equivalent to:: 525 526 $ samtools faidx <ref.fasta> [region1 [...]] 527 528 See http://samtools.sourceforge.net/samtools.shtml for more details 529 530 Example 531 ------- 532 533 >>> from Bio.Sequencing.Applications import SamtoolsFaidxCommandline 534 >>> reference = "/path/to/reference.fasta" 535 >>> samtools_faidx_cmd = SamtoolsFaidxCommandline(reference=reference) 536 >>> print(samtools_faidx_cmd) 537 samtools faidx /path/to/reference.fasta 538 539 """
540 - def __init__(self, cmd="samtools", **kwargs):
541 self.program_name = cmd 542 self.parameters = [ 543 _StaticArgument("faidx"), 544 _Argument(["reference", "reference_fasta", "ref"], 545 "Reference FASTA to be indexed", 546 filename=True, is_required=True) 547 548 ] 549 AbstractCommandline.__init__(self, cmd, **kwargs)
550 551
552 -class SamtoolsFixmateCommandline(AbstractCommandline):
553 """Command line wrapper for samtools fixmate. 554 555 Fill in mate coordinates, ISIZE and mate related 556 flags from a name-sorted alignment, equivalent to:: 557 558 $ samtools fixmate <in.nameSrt.bam> <out.bam> 559 560 See http://samtools.sourceforge.net/samtools.shtml for more details 561 562 Example: 563 564 >>> from Bio.Sequencing.Applications import SamtoolsFixmateCommandline 565 >>> in_bam = "/path/to/in.nameSrt.bam" 566 >>> out_bam = "/path/to/out.bam" 567 >>> samtools_fixmate_cmd = SamtoolsFixmateCommandline(\ 568 input_bam=in_bam,\ 569 out_bam=out_bam) 570 >>> print(samtools_fixmate_cmd) 571 samtools fixmate /path/to/in.nameSrt.bam /path/to/out.bam 572 573 """
574 - def __init__(self, cmd="samtools", **kwargs):
575 self.program_name = cmd 576 self.parameters = [ 577 _StaticArgument("fixmate"), 578 _Argument(["in_bam", "sorted_bam", "input_bam", 579 "input", "input_file"], 580 "Name Sorted Alignment File ", 581 filename=True, is_required=True), 582 _Argument(["out_bam", "output_bam", "output", "output_file"], 583 "Output file", 584 filename=True, is_required=True) 585 ] 586 AbstractCommandline.__init__(self, cmd, **kwargs)
587 588
589 -class SamtoolsRmdupCommandline(AbstractCommandline):
590 """Command line wrapper for samtools rmdup. 591 592 Remove potential PCR duplicates, equivalent to:: 593 594 $ samtools rmdup [-sS] <input.srt.bam> <out.bam> 595 596 See http://samtools.sourceforge.net/samtools.shtml for more details 597 598 Example: 599 600 >>> from Bio.Sequencing.Applications import SamtoolsRmdupCommandline 601 >>> input_sorted_bam = "/path/to/input.srt.bam" 602 >>> out_bam = "/path/to/out.bam" 603 >>> samtools_rmdup_cmd = SamtoolsRmdupCommandline(\ 604 input_bam=input_sorted_bam,\ 605 out_bam=out_bam) 606 >>> print(samtools_rmdup_cmd) 607 samtools rmdup /path/to/input.srt.bam /path/to/out.bam 608 609 """
610 - def __init__(self, cmd="samtools", **kwargs):
611 self.program_name = cmd 612 self.parameters = [ 613 _StaticArgument("rmdup"), 614 _Switch(["-s", "s"], 615 """Remove duplicates for single-end reads. 616 617 By default, the command works for paired-end 618 reads only"""), 619 _Switch(["-S", "S"], """Treat paired-end reads 620 as single-end reads"""), 621 _Argument(["in_bam", "sorted_bam", "input_bam", 622 "input", "input_file"], 623 "Name Sorted Alignment File ", 624 filename=True, is_required=True), 625 _Argument(["out_bam", "output_bam", "output", "output_file"], 626 "Output file", filename=True, is_required=True) 627 628 ] 629 AbstractCommandline.__init__(self, cmd, **kwargs)
630 631
632 -class SamtoolsCalmdCommandline(AbstractCommandline):
633 """Command line wrapper for samtools calmd. 634 635 Generate the MD tag, equivalent to:: 636 637 $ samtools calmd [-EeubSr] [-C capQcoef] <aln.bam> <ref.fasta> 638 639 See http://samtools.sourceforge.net/samtools.shtml for more details 640 641 Example: 642 643 >>> from Bio.Sequencing.Applications import SamtoolsCalmdCommandline 644 >>> input_bam = "/path/to/aln.bam" 645 >>> reference_fasta = "/path/to/reference.fasta" 646 >>> samtools_calmd_cmd = SamtoolsCalmdCommandline(\ 647 input_bam=input_bam,\ 648 reference=reference_fasta) 649 >>> print(samtools_calmd_cmd) 650 samtools calmd /path/to/aln.bam /path/to/reference.fasta 651 652 """
653 - def __init__(self, cmd="samtools", **kwargs):
654 self.program_name = cmd 655 self.parameters = [ 656 _StaticArgument("calmd"), 657 _Switch(["-E", "E"], 658 """Extended BAQ calculation. 659 This option trades specificity for sensitivity, 660 though the effect is minor."""), 661 _Switch(["-e", "e"], 662 """Convert the read base to = if it is 663 identical to the aligned reference base. 664 665 Indel caller does not support the = bases 666 at the moment."""), 667 _Switch(["-u", "u"], "Output uncompressed BAM"), 668 _Switch(["-b", "b"], "Output compressed BAM "), 669 _Switch(["-S", "S"], "The input is SAM with header lines "), 670 _Switch(["-r", "r"], """Compute the BQ tag (without -A) 671 or cap base quality by BAQ (with -A)."""), 672 _Switch(["-A", "A"], 673 """When used jointly with -r this option overwrites 674 the original base quality"""), 675 _Option(["-C", "C"], """Coefficient to cap mapping quality 676 of poorly mapped reads. 677 678 See the pileup command for details.""", 679 equate=False, 680 checker_function=lambda x: isinstance(x, int)), 681 _Argument(["input", "input_file", "in_bam", "infile", "input_bam"], 682 "Input BAM", filename=True, is_required=True), 683 _Argument(["reference", "reference_fasta", "ref"], 684 "Reference FASTA to be indexed", 685 filename=True, is_required=True) 686 687 ] 688 AbstractCommandline.__init__(self, cmd, **kwargs)
689 690
691 -class SamtoolsTargetcutCommandline(AbstractCommandline):
692 """Command line wrapper for samtools targetcut. 693 694 This command identifies target regions by examining the continuity 695 of read depth, computes haploid consensus sequences of targets 696 and outputs a SAM with each sequence corresponding to a target, 697 equivalent to:: 698 699 $ samtools targetcut [-Q minBaseQ] [-i inPenalty] [-0 em0] 700 [-1 em1] [-2 em2] [-f ref] <in.bam> 701 702 See http://samtools.sourceforge.net/samtools.shtml for more details 703 704 Example: 705 706 >>> from Bio.Sequencing.Applications import SamtoolsTargetcutCommandline 707 >>> input_bam = "/path/to/aln.bam" 708 >>> samtools_targetcut_cmd = SamtoolsTargetcutCommandline(input_bam=input_bam) 709 >>> print(samtools_targetcut_cmd) 710 samtools targetcut /path/to/aln.bam 711 712 """
713 - def __init__(self, cmd="samtools", **kwargs):
714 self.program_name = cmd 715 self.parameters = [ 716 _StaticArgument("targetcut"), 717 _Option(["-Q", "Q"], "Minimum Base Quality ", 718 equate=False, 719 checker_function=lambda x: isinstance(x, int)), 720 _Option(["-i", "i"], "Insertion Penalty", 721 equate=False, 722 checker_function=lambda x: isinstance(x, int)), 723 _Option(["-f", "f"], "Reference Filename", 724 filename=True, equate=False, 725 checker_function=lambda x: isinstance(x, str)), 726 _Option(["-0", "em0"], "em0", equate=False, 727 checker_function=lambda x: isinstance(x, str)), 728 _Option(["-1", "em1"], "em1", equate=False, 729 checker_function=lambda x: isinstance(x, str)), 730 _Option(["-2", "em2"], "em2", equate=False, 731 checker_function=lambda x: isinstance(x, str)), 732 _Argument(["input", "input_bam", "in_bam"], 733 "Input file", 734 filename=True, is_required=True) 735 736 ] 737 AbstractCommandline.__init__(self, cmd, **kwargs)
738 739
740 -class SamtoolsPhaseCommandline(AbstractCommandline):
741 """Command line wrapper for samtools phase. 742 743 Call and phase heterozygous SNPs, equivalent to:: 744 745 $ samtools phase [-AF] [-k len] [-b prefix] 746 [-q minLOD] [-Q minBaseQ] <in.bam> 747 748 See http://samtools.sourceforge.net/samtools.shtml for more details 749 750 Example: 751 752 >>> from Bio.Sequencing.Applications import SamtoolsPhaseCommandline 753 >>> input_bam = "/path/to/in.bam" 754 >>> samtools_phase_cmd = SamtoolsPhaseCommandline(input_bam=input_bam) 755 >>> print(samtools_phase_cmd) 756 samtools phase /path/to/in.bam 757 758 """
759 - def __init__(self, cmd="samtools", **kwargs):
760 self.program_name = cmd 761 self.parameters = [ 762 _StaticArgument("phase"), 763 _Argument(["input", "input_bam", "in_bam"], "Input file", 764 filename=True, is_required=True), 765 _Switch(["-A", "A"], "Drop reads with ambiguous phase"), 766 _Option(["-b", "b"], "Prefix of BAM output", 767 filename=True, equate=False, 768 checker_function=lambda x: isinstance(x, str)), 769 _Switch(["-F", "F"], "Do not attempt to fix chimeric reads"), 770 _Option(["-k", "k"], "Maximum length for local phasing", 771 equate=False, 772 checker_function=lambda x: isinstance(x, int)), 773 _Option(["-q", "q"], """Minimum Phred-scaled LOD to 774 call a heterozygote""", 775 equate=False, 776 checker_function=lambda x: isinstance(x, int)), 777 _Option(["-Q", "Q"], """Minimum base quality to be 778 used in het calling""", 779 equate=False, 780 checker_function=lambda x: isinstance(x, int)) 781 ] 782 AbstractCommandline.__init__(self, cmd, **kwargs)
783 784 785 if __name__ == "__main__": 786 from Bio._utils import run_doctest 787 run_doctest() 788