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  __docformat__ = "restructuredtext en" 
  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 SamtoolsSortCommandline(AbstractCommandline):
319 """Command line wrapper for samtools 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 SamtoolsSortCommandline 331 >>> input_bam = "/path/to/input_bam" 332 >>> out_prefix = "/path/to/out_prefix" 333 >>> samtools_sort_cmd = SamtoolsSortCommandline(\ 334 input_bam=input_bam,\ 335 out_prefix=out_prefix) 336 >>> print(samtools_sort_cmd) 337 samtools sort /path/to/input_bam /path/to/out_prefix 338 339 """
340 - def __init__(self, cmd="samtools", **kwargs):
341 self.program_name = cmd 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_bam"], "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 SamtoolsMergeCommandline(AbstractCommandline):
360 """Command line wrapper for samtools merge. 361 362 Merge multiple sorted alignments, equivalent to:: 363 364 $ samtools merge [-nur1f] [-h inh.sam] [-R reg] 365 <out.bam> <in1.bam> <in2.bam> [...] 366 367 See http://samtools.sourceforge.net/samtools.shtml for more details 368 369 Example 370 ------- 371 372 >>> from Bio.Sequencing.Applications import SamtoolsMergeCommandline 373 >>> out_bam = "/path/to/out_bam" 374 >>> in_bam = ["/path/to/input_bam1", "/path/to/input_bam2"] 375 >>> samtools_merge_cmd = SamtoolsMergeCommandline(\ 376 out_bam=out_bam,\ 377 input_bam=in_bam) 378 >>> print(samtools_merge_cmd) 379 samtools merge /path/to/out_bam /path/to/input_bam1 /path/to/input_bam2 380 381 """
382 - def __init__(self, cmd="samtools", **kwargs):
383 self.program_name = cmd 384 self.parameters = [ 385 _StaticArgument("merge"), 386 _Switch(["-n", "n"], 387 """The input alignments are sorted by read names 388 rather than by chromosomal coordinates"""), 389 _Switch(["-r", "r"], """Attach an RG tag to each alignment. 390 The tag value is inferred from file names"""), 391 _Switch(["-u", "u"], "Uncompressed BAM output"), 392 _Switch(["-1", "fast_bam"], """Use zlib compression level 1 393 to compress the output"""), 394 _Switch(["-f", "f"], """Force to overwrite the 395 output file if present"""), 396 _Option(["-h", "h"], """Use the lines of FILE as '@' 397 headers to be copied to out.bam""", 398 filename=True, equate=False, 399 checker_function=lambda x: isinstance(x, str)), 400 _Option(["-R", "R"], 401 "Merge files in the specified region indicated by STR", 402 equate=False, 403 checker_function=lambda x: isinstance(x, str)), 404 _Argument(["output_bam", "out_bam", "out", "output"], 405 "Output BAM file", 406 filename=True, is_required=True), 407 _ArgumentList(["input_bam", "in_bam", "input", "bam"], 408 "Input BAM", 409 filename=True, is_required=True), 410 ] 411 AbstractCommandline.__init__(self, cmd, **kwargs)
412 413
414 -class SamtoolsIndexCommandline(AbstractCommandline):
415 """Command line wrapper for samtools index. 416 417 Index sorted alignment for fast random access, equivalent to:: 418 419 $ samtools index <aln.bam> 420 421 See http://samtools.sourceforge.net/samtools.shtml for more details 422 423 Example: 424 425 >>> from Bio.Sequencing.Applications import SamtoolsIndexCommandline 426 >>> input = "/path/to/aln_bam" 427 >>> samtools_index_cmd = SamtoolsIndexCommandline(input_bam=input) 428 >>> print(samtools_index_cmd) 429 samtools index /path/to/aln_bam 430 431 """
432 - def __init__(self, cmd="samtools", **kwargs):
433 self.program_name = cmd 434 self.parameters = [ 435 _StaticArgument("index"), 436 _Argument(["input", "in_bam", "input_bam"], 437 "BAM file to be indexed"), 438 ] 439 AbstractCommandline.__init__(self, cmd, **kwargs)
440 441
442 -class SamtoolsIdxstatsCommandline(AbstractCommandline):
443 """Command line wrapper for samtools idxstats. 444 445 Retrieve and print stats in the index file, equivalent to:: 446 447 $ samtools idxstats <aln.bam> 448 449 See http://samtools.sourceforge.net/samtools.shtml for more details 450 451 Example: 452 453 >>> from Bio.Sequencing.Applications import SamtoolsIdxstatsCommandline 454 >>> input = "/path/to/aln_bam" 455 >>> samtools_idxstats_cmd = SamtoolsIdxstatsCommandline(input_bam=input) 456 >>> print(samtools_idxstats_cmd) 457 samtools idxstats /path/to/aln_bam 458 459 """
460 - def __init__(self, cmd="samtools", **kwargs):
461 self.program_name = cmd 462 self.parameters = [ 463 _StaticArgument("idxstats"), 464 _Argument(["input", "in_bam", "input_bam"], 465 "BAM file to be indexed") 466 ] 467 AbstractCommandline.__init__(self, cmd, **kwargs)
468 469
470 -class SamtoolsFaidxCommandline(AbstractCommandline):
471 """Command line wrapper for samtools faidx. 472 473 Retrieve and print stats in the index file, equivalent to:: 474 475 $ samtools faidx <ref.fasta> [region1 [...]] 476 477 See http://samtools.sourceforge.net/samtools.shtml for more details 478 479 Example 480 ------- 481 482 >>> from Bio.Sequencing.Applications import SamtoolsFaidxCommandline 483 >>> reference = "/path/to/reference.fasta" 484 >>> samtools_faidx_cmd = SamtoolsFaidxCommandline(reference=reference) 485 >>> print(samtools_faidx_cmd) 486 samtools faidx /path/to/reference.fasta 487 488 """
489 - def __init__(self, cmd="samtools", **kwargs):
490 self.program_name = cmd 491 self.parameters = [ 492 _StaticArgument("faidx"), 493 _Argument(["reference", "reference_fasta", "ref"], 494 "Reference FASTA to be indexed", 495 filename=True, is_required=True) 496 497 ] 498 AbstractCommandline.__init__(self, cmd, **kwargs)
499 500
501 -class SamtoolsFixmateCommandline(AbstractCommandline):
502 """Command line wrapper for samtools fixmate. 503 504 Fill in mate coordinates, ISIZE and mate related 505 flags from a name-sorted alignment, equivalent to:: 506 507 $ samtools fixmate <in.nameSrt.bam> <out.bam> 508 509 See http://samtools.sourceforge.net/samtools.shtml for more details 510 511 Example: 512 513 >>> from Bio.Sequencing.Applications import SamtoolsFixmateCommandline 514 >>> in_bam = "/path/to/in.nameSrt.bam" 515 >>> out_bam = "/path/to/out.bam" 516 >>> samtools_fixmate_cmd = SamtoolsFixmateCommandline(\ 517 input_bam=in_bam,\ 518 out_bam=out_bam) 519 >>> print(samtools_fixmate_cmd) 520 samtools fixmate /path/to/in.nameSrt.bam /path/to/out.bam 521 522 """
523 - def __init__(self, cmd="samtools", **kwargs):
524 self.program_name = cmd 525 self.parameters = [ 526 _StaticArgument("fixmate"), 527 _Argument(["in_bam", "sorted_bam", "input_bam", 528 "input", "input_file"], 529 "Name Sorted Alignment File ", 530 filename=True, is_required=True), 531 _Argument(["out_bam", "output_bam", "output", "output_file"], 532 "Output file", 533 filename=True, is_required=True) 534 ] 535 AbstractCommandline.__init__(self, cmd, **kwargs)
536 537
538 -class SamtoolsRmdupCommandline(AbstractCommandline):
539 """Command line wrapper for samtools rmdup. 540 541 Remove potential PCR duplicates, equivalent to:: 542 543 $ samtools rmdup [-sS] <input.srt.bam> <out.bam> 544 545 See http://samtools.sourceforge.net/samtools.shtml for more details 546 547 Example: 548 549 >>> from Bio.Sequencing.Applications import SamtoolsRmdupCommandline 550 >>> input_sorted_bam = "/path/to/input.srt.bam" 551 >>> out_bam = "/path/to/out.bam" 552 >>> samtools_rmdup_cmd = SamtoolsRmdupCommandline(\ 553 input_bam=input_sorted_bam,\ 554 out_bam=out_bam) 555 >>> print(samtools_rmdup_cmd) 556 samtools rmdup /path/to/input.srt.bam /path/to/out.bam 557 558 """
559 - def __init__(self, cmd="samtools", **kwargs):
560 self.program_name = cmd 561 self.parameters = [ 562 _StaticArgument("rmdup"), 563 _Switch(["-s", "s"], 564 """Remove duplicates for single-end reads. 565 566 By default, the command works for paired-end 567 reads only"""), 568 _Switch(["-S", "S"], """Treat paired-end reads 569 as single-end reads"""), 570 _Argument(["in_bam", "sorted_bam", "input_bam", 571 "input", "input_file"], 572 "Name Sorted Alignment File ", 573 filename=True, is_required=True), 574 _Argument(["out_bam", "output_bam", "output", "output_file"], 575 "Output file", filename=True, is_required=True) 576 577 ] 578 AbstractCommandline.__init__(self, cmd, **kwargs)
579 580
581 -class SamtoolsCalmdCommandline(AbstractCommandline):
582 """Command line wrapper for samtools calmd. 583 584 Generate the MD tag, equivalent to:: 585 586 $ samtools calmd [-EeubSr] [-C capQcoef] <aln.bam> <ref.fasta> 587 588 See http://samtools.sourceforge.net/samtools.shtml for more details 589 590 Example: 591 592 >>> from Bio.Sequencing.Applications import SamtoolsCalmdCommandline 593 >>> input_bam = "/path/to/aln.bam" 594 >>> reference_fasta = "/path/to/reference.fasta" 595 >>> samtools_calmd_cmd = SamtoolsCalmdCommandline(\ 596 input_bam=input_bam,\ 597 reference=reference_fasta) 598 >>> print(samtools_calmd_cmd) 599 samtools calmd /path/to/aln.bam /path/to/reference.fasta 600 601 """
602 - def __init__(self, cmd="samtools", **kwargs):
603 self.program_name = cmd 604 self.parameters = [ 605 _StaticArgument("calmd"), 606 _Switch(["-E", "E"], 607 """Extended BAQ calculation. 608 This option trades specificity for sensitivity, 609 though the effect is minor."""), 610 _Switch(["-e", "e"], 611 """Convert the read base to = if it is 612 identical to the aligned reference base. 613 614 Indel caller does not support the = bases 615 at the moment."""), 616 _Switch(["-u", "u"], "Output uncompressed BAM"), 617 _Switch(["-b", "b"], "Output compressed BAM "), 618 _Switch(["-S", "S"], "The input is SAM with header lines "), 619 _Switch(["-r", "r"], """Compute the BQ tag (without -A) 620 or cap base quality by BAQ (with -A)."""), 621 _Switch(["-A", "A"], 622 """When used jointly with -r this option overwrites 623 the original base quality"""), 624 _Option(["-C", "C"], """Coefficient to cap mapping quality 625 of poorly mapped reads. 626 627 See the pileup command for details.""", 628 equate=False, 629 checker_function=lambda x: isinstance(x, int)), 630 _Argument(["input", "input_file", "in_bam", "infile", "input_bam"], 631 "Input BAM", filename=True, is_required=True), 632 _Argument(["reference", "reference_fasta", "ref"], 633 "Reference FASTA to be indexed", 634 filename=True, is_required=True) 635 636 ] 637 AbstractCommandline.__init__(self, cmd, **kwargs)
638 639
640 -class SamtoolsTargetcutCommandline(AbstractCommandline):
641 """Command line wrapper for samtools targetcut. 642 643 This command identifies target regions by examining the continuity 644 of read depth, computes haploid consensus sequences of targets 645 and outputs a SAM with each sequence corresponding to a target, 646 equivalent to:: 647 648 $ samtools targetcut [-Q minBaseQ] [-i inPenalty] [-0 em0] 649 [-1 em1] [-2 em2] [-f ref] <in.bam> 650 651 See http://samtools.sourceforge.net/samtools.shtml for more details 652 653 Example: 654 655 >>> from Bio.Sequencing.Applications import SamtoolsTargetcutCommandline 656 >>> input_bam = "/path/to/aln.bam" 657 >>> samtools_targetcut_cmd = SamtoolsTargetcutCommandline(input_bam=input_bam) 658 >>> print(samtools_targetcut_cmd) 659 samtools targetcut /path/to/aln.bam 660 661 """
662 - def __init__(self, cmd="samtools", **kwargs):
663 self.program_name = cmd 664 self.parameters = [ 665 _StaticArgument("targetcut"), 666 _Option(["-Q", "Q"], "Minimum Base Quality ", 667 equate=False, 668 checker_function=lambda x: isinstance(x, int)), 669 _Option(["-i", "i"], "Insertion Penalty", 670 equate=False, 671 checker_function=lambda x: isinstance(x, int)), 672 _Option(["-f", "f"], "Reference Filename", 673 filename=True, equate=False, 674 checker_function=lambda x: isinstance(x, str)), 675 _Option(["-0", "em0"], "em0", equate=False, 676 checker_function=lambda x: isinstance(x, str)), 677 _Option(["-1", "em1"], "em1", equate=False, 678 checker_function=lambda x: isinstance(x, str)), 679 _Option(["-2", "em2"], "em2", equate=False, 680 checker_function=lambda x: isinstance(x, str)), 681 _Argument(["input", "input_bam", "in_bam"], 682 "Input file", 683 filename=True, is_required=True) 684 685 ] 686 AbstractCommandline.__init__(self, cmd, **kwargs)
687 688
689 -class SamtoolsPhaseCommandline(AbstractCommandline):
690 """Command line wrapper for samtools phase. 691 692 Call and phase heterozygous SNPs, equivalent to:: 693 694 $ samtools phase [-AF] [-k len] [-b prefix] 695 [-q minLOD] [-Q minBaseQ] <in.bam> 696 697 See http://samtools.sourceforge.net/samtools.shtml for more details 698 699 Example: 700 701 >>> from Bio.Sequencing.Applications import SamtoolsPhaseCommandline 702 >>> input_bam = "/path/to/in.bam" 703 >>> samtools_phase_cmd = SamtoolsPhaseCommandline(input_bam=input_bam) 704 >>> print(samtools_phase_cmd) 705 samtools phase /path/to/in.bam 706 707 """
708 - def __init__(self, cmd="samtools", **kwargs):
709 self.program_name = cmd 710 self.parameters = [ 711 _StaticArgument("phase"), 712 _Argument(["input", "input_bam", "in_bam"], "Input file", 713 filename=True, is_required=True), 714 _Switch(["-A", "A"], "Drop reads with ambiguous phase"), 715 _Option(["-b", "b"], "Prefix of BAM output", 716 filename=True, equate=False, 717 checker_function=lambda x: isinstance(x, str)), 718 _Switch(["-F", "F"], "Do not attempt to fix chimeric reads"), 719 _Option(["-k", "k"], "Maximum length for local phasing", 720 equate=False, 721 checker_function=lambda x: isinstance(x, int)), 722 _Option(["-q", "q"], """Minimum Phred-scaled LOD to 723 call a heterozygote""", 724 equate=False, 725 checker_function=lambda x: isinstance(x, int)), 726 _Option(["-Q", "Q"], """Minimum base quality to be 727 used in het calling""", 728 equate=False, 729 checker_function=lambda x: isinstance(x, int)) 730 ] 731 AbstractCommandline.__init__(self, cmd, **kwargs)
732 733
734 -def _test():
735 """Run the module's doctests (PRIVATE).""" 736 print("Running modules doctests...") 737 import doctest 738 doctest.testmod() 739 print("Done")
740 741 if __name__ == "__main__": 742 _test() 743