Package Bio :: Package Align :: Package Applications :: Module _Probcons
[hide private]
[frames] | no frames]

Source Code for Module Bio.Align.Applications._Probcons

  1  # Copyright 2009 by Cymon J. Cox.  All rights reserved. 
  2  # This code is part of the Biopython distribution and governed by its 
  3  # license.  Please see the LICENSE file that should have been included 
  4  # as part of this package. 
  5  """Command line wrapper for the multiple alignment program PROBCONS. 
  6  """ 
  7   
  8  from __future__ import print_function 
  9   
 10  __docformat__ = "restructuredtext en"  # Don't just use plain text in epydoc API pages! 
 11   
 12  from Bio.Application import _Option, _Switch, _Argument, AbstractCommandline 
 13   
 14   
15 -class ProbconsCommandline(AbstractCommandline):
16 """Command line wrapper for the multiple alignment program PROBCONS. 17 18 http://probcons.stanford.edu/ 19 20 Example: 21 -------- 22 23 To align a FASTA file (unaligned.fasta) with the output in ClustalW 24 format, and otherwise default settings, use: 25 26 >>> from Bio.Align.Applications import ProbconsCommandline 27 >>> probcons_cline = ProbconsCommandline(input="unaligned.fasta", 28 ... clustalw=True) 29 >>> print(probcons_cline) 30 probcons -clustalw unaligned.fasta 31 32 You would typically run the command line with probcons_cline() or via 33 the Python subprocess module, as described in the Biopython tutorial. 34 Note that PROBCONS will write the alignment to stdout, which you may 35 want to save to a file and then parse, e.g.:: 36 37 stdout, stderr = probcons_cline() 38 with open("aligned.aln", "w") as handle: 39 handle.write(stdout) 40 from Bio import AlignIO 41 align = AlignIO.read("aligned.fasta", "clustalw") 42 43 Alternatively, to parse the output with AlignIO directly you can 44 use StringIO to turn the string into a handle:: 45 46 stdout, stderr = probcons_cline() 47 from StringIO import StringIO 48 from Bio import AlignIO 49 align = AlignIO.read(StringIO(stdout), "clustalw") 50 51 Citations: 52 ---------- 53 54 Do, C.B., Mahabhashyam, M.S.P., Brudno, M., and Batzoglou, S. 2005. 55 PROBCONS: Probabilistic Consistency-based Multiple Sequence Alignment. 56 Genome Research 15: 330-340. 57 58 Last checked against version: 1.12 59 """
60 - def __init__(self, cmd="probcons", **kwargs):
61 self.parameters = \ 62 [ 63 # Note that some options cannot be assigned via properties using the 64 # original documented option (because hyphens are not valid for names in 65 # python), e.g cmdline.pre-training = 3 will not work 66 # In these cases the shortened option name should be used 67 # cmdline.pre = 3 68 _Switch(["-clustalw", "clustalw"], 69 "Use CLUSTALW output format instead of MFA"), 70 _Option(["-c", "c", "--consistency", "consistency"], 71 "Use 0 <= REPS <= 5 (default: 2) passes of consistency transformation", 72 checker_function=lambda x: x in range(0, 6), 73 equate=False), 74 _Option(["-ir", "--iterative-refinement", "iterative-refinement", "ir"], 75 "Use 0 <= REPS <= 1000 (default: 100) passes of " 76 "iterative-refinement", 77 checker_function=lambda x: x in range(0, 1001), 78 equate=False), 79 _Option(["-pre", "--pre-training", "pre-training", "pre"], 80 "Use 0 <= REPS <= 20 (default: 0) rounds of pretraining", 81 checker_function=lambda x: x in range(0, 21), 82 equate=False), 83 _Switch(["-pairs", "pairs"], 84 "Generate all-pairs pairwise alignments"), 85 _Switch(["-viterbi", "viterbi"], 86 "Use Viterbi algorithm to generate all pairs " 87 "(automatically enables -pairs)"), 88 _Switch(["-verbose", "verbose"], 89 "Report progress while aligning (default: off)"), 90 _Option(["-annot", "annot"], 91 "Write annotation for multiple alignment to FILENAME", 92 equate=False), 93 _Option(["-t", "t", "--train", "train"], 94 "Compute EM transition probabilities, store in FILENAME " 95 "(default: no training)", 96 equate=False), 97 _Switch(["-e", "e", "--emissions", "emissions"], 98 "Also reestimate emission probabilities (default: off)"), 99 _Option(["-p", "p", "--paramfile", "paramfile"], 100 "Read parameters from FILENAME", 101 equate=False), 102 _Switch(["-a", "--alignment-order", "alignment-order", "a"], 103 "Print sequences in alignment order rather than input " 104 "order (default: off)"), 105 # Input file name 106 _Argument(["input"], 107 "Input file name. Must be multiple FASTA alignment "+ 108 "(MFA) format", 109 filename=True, 110 is_required=True), 111 ] 112 AbstractCommandline.__init__(self, cmd, **kwargs)
113 114
115 -def _test():
116 """Run the module's doctests (PRIVATE).""" 117 print("Running modules doctests...") 118 import doctest 119 doctest.testmod() 120 print("Done")
121 122 if __name__ == "__main__": 123 _test() 124