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