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