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