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

Source Code for Module Bio.Blast.Applications

   1  # Copyright 2001 Brad Chapman. 
   2  # Revisions copyright 2009-2010 by Peter Cock. 
   3  # Revisions copyright 2010 by Phillip Garland. 
   4  # All rights reserved. 
   5  # This code is part of the Biopython distribution and governed by its 
   6  # license.  Please see the LICENSE file that should have been included 
   7  # as part of this package. 
   8  """Definitions for interacting with BLAST related applications. 
   9   
  10  Wrappers for the new NCBI BLAST+ tools (written in C++): 
  11   
  12   - NcbiblastpCommandline - Protein-Protein BLAST 
  13   - NcbiblastnCommandline - Nucleotide-Nucleotide BLAST 
  14   - NcbiblastxCommandline - Translated Query-Protein Subject BLAST 
  15   - NcbitblastnCommandline - Protein Query-Translated Subject BLAST 
  16   - NcbitblastxCommandline - Translated Query-Protein Subject BLAST 
  17   - NcbipsiblastCommandline - Position-Specific Initiated BLAST 
  18   - NcbirpsblastCommandline - Reverse Position Specific BLAST 
  19   - NcbirpstblastnCommandline - Translated Reverse Position Specific BLAST 
  20   - NcbideltablastCommandline - Protein-Protein domain enhanced lookup time accelerated blast 
  21   - NcbiblastformatterCommandline - Convert ASN.1 to other BLAST output formats 
  22   
  23  For further details, see: 
  24   
  25  Camacho et al. BLAST+: architecture and applications 
  26  BMC Bioinformatics 2009, 10:421 
  27  doi:10.1186/1471-2105-10-421 
  28  """ 
  29  from __future__ import print_function 
  30   
  31  from Bio.Application import _Option, AbstractCommandline, _Switch 
  32   
  33   
34 -class _NcbibaseblastCommandline(AbstractCommandline):
35 """Base Commandline object for (new) NCBI BLAST+ wrappers (PRIVATE). 36 37 This is provided for subclassing, it deals with shared options 38 common to all the BLAST tools (blastn, rpsblast, rpsblast, etc 39 AND blast_formatter). 40 """
41 - def __init__(self, cmd=None, **kwargs):
42 assert cmd is not None 43 extra_parameters = [ 44 # Core: 45 _Switch(["-h", "h"], 46 "Print USAGE and DESCRIPTION; ignore other arguments."), 47 _Switch(["-help", "help"], 48 "Print USAGE, DESCRIPTION and ARGUMENTS description; " 49 "ignore other arguments."), 50 _Switch(["-version", "version"], 51 "Print version number; ignore other arguments."), 52 # Output configuration options 53 _Option(["-out", "out"], 54 "Output file for alignment.", 55 filename=True, 56 equate=False), 57 # Formatting options: 58 _Option(["-outfmt", "outfmt"], 59 "Alignment view. Integer 0-11. Use 5 for XML output " 60 "(differs from classic BLAST which used 7 for XML).", 61 equate=False), 62 # TODO - Document and test the column options 63 _Switch(["-show_gis", "show_gis"], 64 "Show NCBI GIs in deflines?"), 65 _Option(["-num_descriptions", "num_descriptions"], 66 """Number of database sequences to show one-line descriptions for. 67 68 Integer argument (at least zero). Default is 500. 69 See also num_alignments.""", 70 equate=False), 71 _Option(["-num_alignments", "num_alignments"], 72 """Number of database sequences to show num_alignments for. 73 74 Integer argument (at least zero). Default is 200. 75 See also num_alignments.""", 76 equate=False), 77 _Option(["-line_length", "line_length"], 78 """Line length for formatting alignments (integer, at least 1, default 60). 79 80 Not applicable for outfmt > 4. 81 82 Added in BLAST+ 2.2.30. 83 """, 84 equate=False), 85 _Switch(["-html", "html"], 86 "Produce HTML output? See also the outfmt option."), 87 # Miscellaneous options 88 _Switch(["-parse_deflines", "parse_deflines"], 89 "Should the query and subject defline(s) be parsed?"), 90 ] 91 try: 92 # Insert extra parameters - at the start just in case there 93 # are any arguments which must come last: 94 self.parameters = extra_parameters + self.parameters 95 except AttributeError: 96 # Should we raise an error? The subclass should have set this up! 97 self.parameters = extra_parameters 98 AbstractCommandline.__init__(self, cmd, **kwargs)
99
100 - def _validate_incompatibilities(self, incompatibles):
101 """Used by the BLAST+ _validate method (PRIVATE).""" 102 for a in incompatibles: 103 if self._get_parameter(a): 104 for b in incompatibles[a]: 105 if self._get_parameter(b): 106 raise ValueError("Options %s and %s are incompatible." 107 % (a, b))
108 109
110 -class _NcbiblastCommandline(_NcbibaseblastCommandline):
111 """Base Commandline object for (new) NCBI BLAST+ wrappers (PRIVATE). 112 113 This is provided for subclassing, it deals with shared options 114 common to all the BLAST tools (blastn, rpsblast, rpsblast, etc). 115 """
116 - def __init__(self, cmd=None, **kwargs):
117 assert cmd is not None 118 extra_parameters = [ 119 # Input query options: 120 _Option(["-query", "query"], 121 "The sequence to search with.", 122 filename=True, 123 equate=False), # Should this be required? 124 _Option(["-query_loc", "query_loc"], 125 "Location on the query sequence (Format: start-stop)", 126 equate=False), 127 # General search options: 128 _Option(["-db", "db"], 129 "The database to BLAST against.", 130 equate=False), 131 _Option(["-evalue", "evalue"], 132 "Expectation value cutoff.", 133 equate=False), 134 _Option(["-word_size", "word_size"], 135 """Word size for wordfinder algorithm. 136 137 Integer. Minimum 2.""", 138 equate=False), 139 # BLAST-2-Sequences options: 140 # - see subclass 141 # Formatting options: 142 # - see baseclass 143 # Query filtering options 144 _Option(["-soft_masking", "soft_masking"], 145 "Apply filtering locations as soft masks (Boolean, Default = true)", 146 equate=False), 147 _Switch(["-lcase_masking", "lcase_masking"], 148 "Use lower case filtering in query and subject sequence(s)?"), 149 # Restrict search or results 150 _Option(["-gilist", "gilist"], 151 """Restrict search of database to list of GI's. 152 153 Incompatible with: negative_gilist, seqidlist, remote, subject, subject_loc""", 154 filename=True, 155 equate=False), 156 _Option(["-negative_gilist", "negative_gilist"], 157 """Restrict search of database to everything except the listed GIs. 158 159 Incompatible with: gilist, seqidlist, remote, subject, subject_loc""", 160 filename=True, 161 equate=False), 162 _Option(["-seqidlist", "seqidlist"], 163 """Restrict search of database to list of SeqID's. 164 165 Incompatible with: gilist, negative_gilist, remote, subject, subject_loc""", 166 filename=True, 167 equate=False), 168 _Option(["-entrez_query", "entrez_query"], 169 "Restrict search with the given Entrez query (requires remote).", 170 equate=False), 171 _Option(["-qcov_hsp_perc", "qcov_hsp_perc"], 172 """Percent query coverage per hsp (float, 0 to 100). 173 174 Added in BLAST+ 2.2.30. 175 """, 176 equate=False), 177 _Option(["-max_target_seqs", "max_target_seqs"], 178 "Maximum number of aligned sequences to keep (integer, at least one).", 179 equate=False), 180 # Statistical options 181 _Option(["-dbsize", "dbsize"], 182 "Effective length of the database (integer).", 183 equate=False), 184 _Option(["-searchsp", "searchsp"], 185 "Effective length of the search space (integer).", 186 equate=False), 187 _Option(["-max_hsps_per_subject", "max_hsps_per_subject"], 188 "Override max number of HSPs per subject saved for ungapped searches (integer).", 189 equate=False), 190 _Option(["-max_hsps", "max_hsps"], 191 "Set max number of HSPs saved per subject sequence (default 0 means no limit).", 192 equate=False), 193 _Switch(["-sum_statistics", "sum_statistics"], 194 "Use sum statistics."), 195 # Is -sum_stats a BLAST+ bug, why not use -sum_statistics switch? 196 _Option(["-sum_stats", "sum_stats"], 197 """Use sum statistics (boolean). 198 199 Added in BLAST+ 2.2.30. 200 """, 201 equate=False), 202 # Extension options 203 _Option(["-xdrop_ungap", "xdrop_ungap"], 204 "X-dropoff value (in bits) for ungapped extensions (float).", 205 equate=False), 206 _Option(["-xdrop_gap", "xdrop_gap"], 207 "X-dropoff value (in bits) for preliminary gapped extensions (float).", 208 equate=False), 209 _Option(["-xdrop_gap_final", "xdrop_gap_final"], 210 "X-dropoff value (in bits) for final gapped alignment (float).", 211 equate=False), 212 _Option(["-window_size", "window_size"], 213 "Multiple hits window size, use 0 to specify 1-hit algorithm (integer).", 214 equate=False), 215 # Search strategy options 216 _Option(["-import_search_strategy", "import_search_strategy"], 217 """Search strategy to use. 218 219 Incompatible with: export_search_strategy""", 220 filename=True, 221 equate=False), 222 _Option(["-export_search_strategy", "export_search_strategy"], 223 """File name to record the search strategy used. 224 225 Incompatible with: import_search_strategy""", 226 filename=True, 227 equate=False), 228 # Miscellaneous options 229 _Option(["-num_threads", "num_threads"], 230 """Number of threads to use in the BLAST search (integer, at least one). 231 232 Default is one. 233 Incompatible with: remote""", 234 equate=False), 235 _Switch(["-remote", "remote"], 236 """Execute search remotely? 237 238 Incompatible with: gilist, negative_gilist, subject_loc, num_threads, ..."""), 239 ] 240 try: 241 # Insert extra parameters - at the start just in case there 242 # are any arguments which must come last: 243 self.parameters = extra_parameters + self.parameters 244 except AttributeError: 245 # Should we raise an error? The subclass should have set this up! 246 self.parameters = extra_parameters 247 _NcbibaseblastCommandline.__init__(self, cmd, **kwargs)
248
249 - def _validate(self):
250 incompatibles = {"remote": ["gilist", "negative_gilist", "num_threads"], 251 "import_search_strategy": ["export_search_strategy"], 252 "gilist": ["negative_gilist"], 253 "seqidlist": ["gilist", "negative_gilist", "remote"]} 254 self._validate_incompatibilities(incompatibles) 255 if self.entrez_query and not self.remote: 256 raise ValueError("Option entrez_query requires remote option.") 257 AbstractCommandline._validate(self)
258 259
260 -class _Ncbiblast2SeqCommandline(_NcbiblastCommandline):
261 """Base Commandline object for (new) NCBI BLAST+ wrappers (PRIVATE). 262 263 This is provided for subclassing, it deals with shared options 264 common to all the BLAST tools supporting two-sequence BLAST 265 (blastn, psiblast, etc) but not rpsblast or rpstblastn. 266 """
267 - def __init__(self, cmd=None, **kwargs):
268 assert cmd is not None 269 extra_parameters = [ 270 # General search options: 271 _Option(["-gapopen", "gapopen"], 272 "Cost to open a gap (integer).", 273 equate=False), 274 _Option(["-gapextend", "gapextend"], 275 "Cost to extend a gap (integer).", 276 equate=False), 277 # BLAST-2-Sequences options: 278 _Option(["-subject", "subject"], 279 """Subject sequence(s) to search. 280 281 Incompatible with: db, gilist, negative_gilist. 282 See also subject_loc.""", 283 filename=True, 284 equate=False), 285 _Option(["-subject_loc", "subject_loc"], 286 """Location on the subject sequence (Format: start-stop). 287 288 Incompatible with: db, gilist, seqidlist, negative_gilist, 289 db_soft_mask, db_hard_mask, remote. 290 291 See also subject.""", 292 equate=False), 293 # Restrict search or results: 294 _Option(["-culling_limit", "culling_limit"], 295 """Hit culling limit (integer). 296 297 If the query range of a hit is enveloped by that of at 298 least this many higher-scoring hits, delete the hit. 299 300 Incompatible with: best_hit_overhang, best_hit_score_edge. 301 """, 302 equate=False), 303 _Option(["-best_hit_overhang", "best_hit_overhang"], 304 """Best Hit algorithm overhang value (float, recommended value: 0.1) 305 306 Float between 0.0 and 0.5 inclusive. 307 308 Incompatible with: culling_limit.""", 309 equate=False), 310 _Option(["-best_hit_score_edge", "best_hit_score_edge"], 311 """Best Hit algorithm score edge value (float, recommended value: 0.1) 312 313 Float between 0.0 and 0.5 inclusive. 314 315 Incompatible with: culling_limit.""", 316 equate=False), 317 ] 318 try: 319 # Insert extra parameters - at the start just in case there 320 # are any arguments which must come last: 321 self.parameters = extra_parameters + self.parameters 322 except AttributeError: 323 # Should we raise an error? The subclass should have set this up! 324 self.parameters = extra_parameters 325 _NcbiblastCommandline.__init__(self, cmd, **kwargs)
326
327 - def _validate(self):
328 incompatibles = {"subject_loc": ["db", "gilist", "negative_gilist", "seqidlist", "remote"], 329 "culling_limit": ["best_hit_overhang", "best_hit_score_edge"], 330 "subject": ["db", "gilist", "negative_gilist", "seqidlist"]} 331 self._validate_incompatibilities(incompatibles) 332 _NcbiblastCommandline._validate(self)
333 334
335 -class _NcbiblastMain2SeqCommandline(_Ncbiblast2SeqCommandline):
336 """Base Commandline object for (new) NCBI BLAST+ wrappers (PRIVATE). 337 338 This is provided for subclassing, it deals with shared options 339 common to the main BLAST tools blastp, blastn, blastx, tblastx, tblastn 340 but not psiblast, rpsblast or rpstblastn. 341 """
342 - def __init__(self, cmd=None, **kwargs):
343 assert cmd is not None 344 extra_parameters = [ 345 # Restrict search or results: 346 _Option(["-db_soft_mask", "db_soft_mask"], 347 """Filtering algorithm for soft masking (integer). 348 349 Filtering algorithm ID to apply to the BLAST database as soft masking. 350 351 Incompatible with: db_hard_mask, subject, subject_loc""", 352 equate=False), 353 _Option(["-db_hard_mask", "db_hard_mask"], 354 """Filtering algorithm for hard masking (integer). 355 356 Filtering algorithm ID to apply to the BLAST database as hard masking. 357 358 Incompatible with: db_soft_mask, subject, subject_loc""", 359 equate=False), 360 ] 361 try: 362 # Insert extra parameters - at the start just in case there 363 # are any arguments which must come last: 364 self.parameters = extra_parameters + self.parameters 365 except AttributeError: 366 # Should we raise an error? The subclass should have set this up! 367 self.parameters = extra_parameters 368 _Ncbiblast2SeqCommandline.__init__(self, cmd, **kwargs)
369
370 - def _validate(self):
371 incompatibles = {"db_soft_mask": ["db_hard_mask", "subject", "subject_loc"], 372 "db_hard_mask": ["db_soft_mask", "subject", "subject_loc"]} 373 self._validate_incompatibilities(incompatibles) 374 _Ncbiblast2SeqCommandline._validate(self)
375 376
377 -class NcbiblastpCommandline(_NcbiblastMain2SeqCommandline):
378 """Create a commandline for the NCBI BLAST+ program blastp (for proteins). 379 380 With the release of BLAST+ (BLAST rewritten in C++ instead of C), the NCBI 381 replaced the old blastall tool with separate tools for each of the searches. 382 This wrapper therefore replaces BlastallCommandline with option -p blastp. 383 384 >>> from Bio.Blast.Applications import NcbiblastpCommandline 385 >>> cline = NcbiblastpCommandline(query="rosemary.pro", db="nr", 386 ... evalue=0.001, remote=True, ungapped=True) 387 >>> cline 388 NcbiblastpCommandline(cmd='blastp', query='rosemary.pro', db='nr', evalue=0.001, remote=True, ungapped=True) 389 >>> print(cline) 390 blastp -query rosemary.pro -db nr -evalue 0.001 -remote -ungapped 391 392 You would typically run the command line with cline() or via the Python 393 subprocess module, as described in the Biopython tutorial. 394 """
395 - def __init__(self, cmd="blastp", **kwargs):
396 self.parameters = [ 397 # General search options: 398 _Option(["-task", "task"], 399 "Task to execute (string, blastp (default), blastp-fast or blastp-short).", 400 checker_function=lambda value: value in ["blastp", 401 "blastp-fast", 402 "blastp-short"], 403 equate=False), 404 _Option(["-matrix", "matrix"], 405 "Scoring matrix name (default BLOSUM62)."), 406 _Option(["-threshold", "threshold"], 407 "Minimum score for words to be added to the BLAST lookup table (float).", 408 equate=False), 409 _Option(["-comp_based_stats", "comp_based_stats"], 410 """Use composition-based statistics (string, default 2, i.e. True). 411 412 0, F or f: no composition-based statistics 413 2, T or t, D or d : Composition-based score adjustment as in 414 Bioinformatics 21:902-911, 2005, conditioned on sequence properties 415 416 Note that tblastn also supports values of 1 and 3.""", 417 checker_function=lambda value: value in "0Ft2TtDd", 418 equate=False), 419 # Query filtering options: 420 _Option(["-seg", "seg"], 421 """Filter query sequence with SEG (string). 422 423 Format: "yes", "window locut hicut", or "no" to disable. 424 Default is "12 2.2 2.5""", 425 equate=False), 426 # Extension options: 427 _Switch(["-ungapped", "ungapped"], 428 "Perform ungapped alignment only?"), 429 # Miscellaneous options: 430 _Switch(["-use_sw_tback", "use_sw_tback"], 431 "Compute locally optimal Smith-Waterman alignments?"), 432 ] 433 _NcbiblastMain2SeqCommandline.__init__(self, cmd, **kwargs)
434 435
436 -class NcbiblastnCommandline(_NcbiblastMain2SeqCommandline):
437 """Wrapper for the NCBI BLAST+ program blastn (for nucleotides). 438 439 With the release of BLAST+ (BLAST rewritten in C++ instead of C), the NCBI 440 replaced the old blastall tool with separate tools for each of the searches. 441 This wrapper therefore replaces BlastallCommandline with option -p blastn. 442 443 For example, to run a search against the "nt" nucleotide database using the 444 FASTA nucleotide file "m_code.fasta" as the query, with an expectation value 445 cut off of 0.001, saving the output to a file in XML format: 446 447 >>> from Bio.Blast.Applications import NcbiblastnCommandline 448 >>> cline = NcbiblastnCommandline(query="m_cold.fasta", db="nt", strand="plus", 449 ... evalue=0.001, out="m_cold.xml", outfmt=5) 450 >>> cline 451 NcbiblastnCommandline(cmd='blastn', out='m_cold.xml', outfmt=5, query='m_cold.fasta', db='nt', evalue=0.001, strand='plus') 452 >>> print(cline) 453 blastn -out m_cold.xml -outfmt 5 -query m_cold.fasta -db nt -evalue 0.001 -strand plus 454 455 You would typically run the command line with cline() or via the Python 456 subprocess module, as described in the Biopython tutorial. 457 """
458 - def __init__(self, cmd="blastn", **kwargs):
459 self.parameters = [ 460 # Input query options: 461 _Option(["-strand", "strand"], 462 """Query strand(s) to search against database/subject. 463 464 Values allowed are "both" (default), "minus", "plus".""", 465 checker_function=lambda value: value in ["both", 466 "minus", 467 "plus"], 468 equate=False), 469 # General search options: 470 _Option(["-task", "task"], 471 """Task to execute (string, default 'megablast') 472 473 Allowed values 'blastn', 'blastn-short', 'dc-megablast', 'megablast' 474 (the default), or 'vecscreen'.""", 475 checker_function=lambda value: value in ['blastn', 476 'blastn-short', 477 'dc-megablast', 478 'megablast', 479 'vecscreen'], 480 equate=False), 481 _Option(["-penalty", "penalty"], 482 "Penalty for a nucleotide mismatch (integer, at most zero).", 483 equate=False), 484 _Option(["-reward", "reward"], 485 "Reward for a nucleotide match (integer, at least zero).", 486 equate=False), 487 _Option(["-use_index", "use_index"], 488 "Use MegaBLAST database index (Boolean, Default = False)", 489 equate=False), 490 _Option(["-index_name", "index_name"], 491 "MegaBLAST database index name.", 492 equate=False), 493 # Query filtering options: 494 _Option(["-dust", "dust"], 495 """Filter query sequence with DUST (string). 496 497 Format: 'yes', 'level window linker', or 'no' to disable. 498 Default = '20 64 1'. 499 """, 500 equate=False), 501 _Option(["-filtering_db", "filtering_db"], 502 "BLAST database containing filtering elements (i.e. repeats).", 503 equate=False), 504 _Option(["-window_masker_taxid", "window_masker_taxid"], 505 "Enable WindowMasker filtering using a Taxonomic ID (integer).", 506 equate=False), 507 _Option(["-window_masker_db", "window_masker_db"], 508 "Enable WindowMasker filtering using this repeats database (string).", 509 equate=False), 510 # Restrict search or results: 511 _Option(["-perc_identity", "perc_identity"], 512 "Percent identity (real, 0 to 100 inclusive).", 513 equate=False), 514 # Discontiguous MegaBLAST options 515 _Option(["-template_type", "template_type"], 516 """Discontiguous MegaBLAST template type (string). 517 518 Allowed values: 'coding', 'coding_and_optimal' or 'optimal' 519 Requires: template_length.""", 520 checker_function=lambda value: value in ['coding', 'coding_and_optimal', 'optimal'], 521 equate=False), 522 _Option(["-template_length", "template_length"], 523 """Discontiguous MegaBLAST template length (integer). 524 525 Allowed values: 16, 18, 21 526 527 Requires: template_type.""", 528 checker_function=lambda value: value in [16, 18, 21, '16', '18', '21'], 529 equate=False), 530 # Extension options: 531 _Switch(["-no_greedy", "no_greedy"], 532 "Use non-greedy dynamic programming extension"), 533 _Option(["-min_raw_gapped_score", "min_raw_gapped_score"], 534 "Minimum raw gapped score to keep an alignment in the " 535 "preliminary gapped and traceback stages (integer).", 536 equate=False), 537 _Switch(["-ungapped", "ungapped"], 538 "Perform ungapped alignment only?"), 539 _Option(["-off_diagonal_range", "off_diagonal_range"], 540 """Number of off-diagonals to search for the 2nd hit (integer). 541 542 Expects a positive integer, or 0 (default) to turn off. 543 544 Added in BLAST 2.2.23+ 545 """, 546 equate=False), 547 ] 548 _NcbiblastMain2SeqCommandline.__init__(self, cmd, **kwargs)
549
550 - def _validate(self):
551 if (self.template_type and not self.template_length) \ 552 or (self.template_length and not self.template_type): 553 raise ValueError("Options template_type and template_type require each other.") 554 _NcbiblastMain2SeqCommandline._validate(self)
555 556
557 -class NcbiblastxCommandline(_NcbiblastMain2SeqCommandline):
558 """Wrapper for the NCBI BLAST+ program blastx (nucleotide query, protein database). 559 560 With the release of BLAST+ (BLAST rewritten in C++ instead of C), the NCBI 561 replaced the old blastall tool with separate tools for each of the searches. 562 This wrapper therefore replaces BlastallCommandline with option -p blastx. 563 564 >>> from Bio.Blast.Applications import NcbiblastxCommandline 565 >>> cline = NcbiblastxCommandline(query="m_cold.fasta", db="nr", evalue=0.001) 566 >>> cline 567 NcbiblastxCommandline(cmd='blastx', query='m_cold.fasta', db='nr', evalue=0.001) 568 >>> print(cline) 569 blastx -query m_cold.fasta -db nr -evalue 0.001 570 571 You would typically run the command line with cline() or via the Python 572 subprocess module, as described in the Biopython tutorial. 573 """
574 - def __init__(self, cmd="blastx", **kwargs):
575 self.parameters = [ 576 # Input query options: 577 _Option(["-task", "task"], 578 "Task to execute (string, blastx (default) or blastx-fast).", 579 checker_function=lambda value: value in ["blastx", 580 "blastx-fast"], 581 equate=False), 582 _Option(["-strand", "strand"], 583 """Query strand(s) to search against database/subject. 584 585 Values allowed are "both" (default), "minus", "plus".""", 586 checker_function=lambda value: value in ["both", "minus", "plus"], 587 equate=False), 588 # Input query options: 589 _Option(["-query_gencode", "query_gencode"], 590 "Genetic code to use to translate query (integer, default 1).", 591 equate=False), 592 # General search options: 593 _Option(["-frame_shift_penalty", "frame_shift_penalty"], 594 """Frame shift penalty (integer, at least 1, default ignored) (OBSOLETE). 595 596 This was removed in BLAST 2.2.27+""", 597 equate=False), 598 _Option(["-max_intron_length", "max_intron_length"], 599 """Maximum intron length (integer). 600 601 Length of the largest intron allowed in a translated nucleotide 602 sequence when linking multiple distinct alignments (a negative 603 value disables linking). Default zero.""", 604 equate=False), 605 _Option(["-matrix", "matrix"], 606 "Scoring matrix name (default BLOSUM62).", 607 equate=False), 608 _Option(["-threshold", "threshold"], 609 "Minimum score for words to be added to the BLAST lookup table (float).", 610 equate=False), 611 _Option(["-comp_based_stats", "comp_based_stats"], 612 """Use composition-based statistics for blastp, blastx, or tblastn: 613 614 D or d: default (equivalent to 2 ) 615 0 or F or f: no composition-based statistics 616 1: Composition-based statistics as in NAR 29:2994-3005, 2001 617 2 or T or t : Composition-based score adjustment as in 618 Bioinformatics 21:902-911, 2005, conditioned on sequence properties 619 3: Composition-based score adjustment as in 620 Bioinformatics 21:902-911, 2005, unconditionally 621 622 For programs other than tblastn, must either be absent or be D, F or 0 623 Default = `2' 624 """, 625 equate=False), 626 # Query filtering options: 627 _Option(["-seg", "seg"], 628 """Filter query sequence with SEG (string). 629 630 Format: "yes", "window locut hicut", or "no" to disable. 631 Default is "12 2.2 2.5""", 632 equate=False), 633 # Extension options: 634 _Switch(["-ungapped", "ungapped"], 635 "Perform ungapped alignment only?"), 636 _Switch(["-use_sw_tback", "use_sw_tback"], 637 "Compute locally optimal Smith-Waterman alignments?"), 638 ] 639 _NcbiblastMain2SeqCommandline.__init__(self, cmd, **kwargs)
640 641
642 -class NcbitblastnCommandline(_NcbiblastMain2SeqCommandline):
643 """Wrapper for the NCBI BLAST+ program tblastn. 644 645 With the release of BLAST+ (BLAST rewritten in C++ instead of C), the NCBI 646 replaced the old blastall tool with separate tools for each of the searches. 647 This wrapper therefore replaces BlastallCommandline with option -p tblastn. 648 649 >>> from Bio.Blast.Applications import NcbitblastnCommandline 650 >>> cline = NcbitblastnCommandline(help=True) 651 >>> cline 652 NcbitblastnCommandline(cmd='tblastn', help=True) 653 >>> print(cline) 654 tblastn -help 655 656 You would typically run the command line with cline() or via the Python 657 subprocess module, as described in the Biopython tutorial. 658 """
659 - def __init__(self, cmd="tblastn", **kwargs):
660 self.parameters = [ 661 # General search options: 662 _Option(["-task", "task"], 663 "Task to execute (string, tblastn (default) or tblastn-fast).", 664 checker_function=lambda value: value in ["tblastn", 665 "tblastn-fast"], 666 equate=False), 667 _Option(["-db_gencode", "db_gencode"], 668 "Genetic code to use to translate query (integer, default 1).", 669 equate=False), 670 _Option(["-frame_shift_penalty", "frame_shift_penalty"], 671 """Frame shift penalty (integer, at least 1, default ignored) (OBSOLETE). 672 673 This was removed in BLAST 2.2.27+""", 674 equate=False), 675 _Option(["-max_intron_length", "max_intron_length"], 676 """Maximum intron length (integer). 677 678 Length of the largest intron allowed in a translated nucleotide 679 sequence when linking multiple distinct alignments (a negative 680 value disables linking). Default zero.""", 681 equate=False), 682 _Option(["-matrix", "matrix"], 683 "Scoring matrix name (default BLOSUM62).", 684 equate=False), 685 _Option(["-threshold", "threshold"], 686 "Minimum score for words to be added to the BLAST lookup table (float).", 687 equate=False), 688 _Option(["-comp_based_stats", "comp_based_stats"], 689 """Use composition-based statistics (string, default 2, i.e. True). 690 691 0, F or f: no composition-based statistics 692 1: Composition-based statistics as in NAR 29:2994-3005, 2001 693 2, T or t, D or d : Composition-based score adjustment as in 694 Bioinformatics 21:902-911, 2005, conditioned on sequence properties 695 3: Composition-based score adjustment as in Bioinformatics 21:902-911, 696 2005, unconditionally 697 698 Note that only tblastn supports values of 1 and 3.""", 699 checker_function=lambda value: value in "0Ft12TtDd3", 700 equate=False), 701 # Query filtering options: 702 _Option(["-seg", "seg"], 703 """Filter query sequence with SEG (string). 704 705 Format: "yes", "window locut hicut", or "no" to disable. 706 Default is "12 2.2 2.5""", 707 equate=False), 708 # Extension options: 709 _Switch(["-ungapped", "ungapped"], 710 "Perform ungapped alignment only?"), 711 # Miscellaneous options: 712 _Switch(["-use_sw_tback", "use_sw_tback"], 713 "Compute locally optimal Smith-Waterman alignments?"), 714 # PSI-TBLASTN options: 715 _Option(["-in_pssm", "in_pssm"], 716 """PSI-BLAST checkpoint file 717 718 Incompatible with: remote, query""", 719 filename=True, 720 equate=False), 721 ] 722 _NcbiblastMain2SeqCommandline.__init__(self, cmd, **kwargs)
723 724
725 -class NcbitblastxCommandline(_NcbiblastMain2SeqCommandline):
726 """Wrapper for the NCBI BLAST+ program tblastx. 727 728 With the release of BLAST+ (BLAST rewritten in C++ instead of C), the NCBI 729 replaced the old blastall tool with separate tools for each of the searches. 730 This wrapper therefore replaces BlastallCommandline with option -p tblastx. 731 732 >>> from Bio.Blast.Applications import NcbitblastxCommandline 733 >>> cline = NcbitblastxCommandline(help=True) 734 >>> cline 735 NcbitblastxCommandline(cmd='tblastx', help=True) 736 >>> print(cline) 737 tblastx -help 738 739 You would typically run the command line with cline() or via the Python 740 subprocess module, as described in the Biopython tutorial. 741 """
742 - def __init__(self, cmd="tblastx", **kwargs):
743 self.parameters = [ 744 # Input query options: 745 _Option(["-strand", "strand"], 746 """Query strand(s) to search against database/subject. 747 748 Values allowed are "both" (default), "minus", "plus".""", 749 checker_function=lambda value: value in ["both", "minus", "plus"], 750 equate=False), 751 # Input query options: 752 _Option(["-query_gencode", "query_gencode"], 753 "Genetic code to use to translate query (integer, default 1).", 754 equate=False), 755 # General search options: 756 _Option(["-db_gencode", "db_gencode"], 757 "Genetic code to use to translate query (integer, default 1).", 758 equate=False), 759 _Option(["-max_intron_length", "max_intron_length"], 760 """Maximum intron length (integer). 761 762 Length of the largest intron allowed in a translated nucleotide 763 sequence when linking multiple distinct alignments (a negative 764 value disables linking). Default zero.""", 765 equate=False), 766 _Option(["-matrix", "matrix"], 767 "Scoring matrix name (default BLOSUM62).", 768 equate=False), 769 _Option(["-threshold", "threshold"], 770 "Minimum score for words to be added to the BLAST lookup table (float).", 771 equate=False), 772 # Query filtering options: 773 _Option(["-seg", "seg"], 774 """Filter query sequence with SEG (string). 775 776 Format: "yes", "window locut hicut", or "no" to disable. 777 Default is "12 2.2 2.5""", 778 equate=False), 779 ] 780 _NcbiblastMain2SeqCommandline.__init__(self, cmd, **kwargs)
781 782
783 -class NcbipsiblastCommandline(_Ncbiblast2SeqCommandline):
784 """Wrapper for the NCBI BLAST+ program psiblast. 785 786 With the release of BLAST+ (BLAST rewritten in C++ instead of C), the NCBI 787 replaced the old blastpgp tool with a similar tool psiblast. This wrapper 788 therefore replaces BlastpgpCommandline, the wrapper for blastpgp. 789 790 >>> from Bio.Blast.Applications import NcbipsiblastCommandline 791 >>> cline = NcbipsiblastCommandline(help=True) 792 >>> cline 793 NcbipsiblastCommandline(cmd='psiblast', help=True) 794 >>> print(cline) 795 psiblast -help 796 797 You would typically run the command line with cline() or via the Python 798 subprocess module, as described in the Biopython tutorial. 799 """
800 - def __init__(self, cmd="psiblast", **kwargs):
801 self.parameters = [ 802 # General search options: 803 _Option(["-matrix", "matrix"], 804 "Scoring matrix name (default BLOSUM62).", 805 equate=False), 806 _Option(["-threshold", "threshold"], 807 "Minimum score for words to be added to the BLAST lookup table (float).", 808 equate=False), 809 _Option(["-comp_based_stats", "comp_based_stats"], 810 """Use composition-based statistics (string, default 2, i.e. True). 811 812 0, F or f: no composition-based statistics 813 2, T or t, D or d : Composition-based score adjustment 814 as in Bioinformatics 21:902-911, 2005, conditioned on 815 sequence properties 816 817 Note that tblastn also supports values of 1 and 3.""", 818 checker_function=lambda value: value in "0Ft2TtDd", 819 equate=False), 820 # Query filtering options: 821 _Option(["-seg", "seg"], 822 """Filter query sequence with SEG (string). 823 824 Format: "yes", "window locut hicut", or "no" to disable. 825 Default is "12 2.2 2.5""", 826 equate=False), 827 # Extension options: 828 _Option(["-gap_trigger", "gap_trigger"], 829 "Number of bits to trigger gapping (float, default 22).", 830 equate=False), 831 # Miscellaneous options: 832 _Switch(["-use_sw_tback", "use_sw_tback"], 833 "Compute locally optimal Smith-Waterman alignments?"), 834 # PSI-BLAST options: 835 _Option(["-num_iterations", "num_iterations"], 836 """Number of iterations to perform (integer, at least one). 837 838 Default is one. 839 Incompatible with: remote""", 840 equate=False), 841 _Option(["-out_pssm", "out_pssm"], 842 "File name to store checkpoint file.", 843 filename=True, 844 equate=False), 845 _Option(["-out_ascii_pssm", "out_ascii_pssm"], 846 "File name to store ASCII version of PSSM.", 847 filename=True, 848 equate=False), 849 _Switch(["-save_pssm_after_last_round", "save_pssm_after_last_round"], 850 "Save PSSM after the last database search."), 851 _Switch(["-save_each_pssm", "save_each_pssm"], 852 """Save PSSM after each iteration 853 854 File name is given in -save_pssm or -save_ascii_pssm options. 855 """), 856 _Option(["-in_msa", "in_msa"], 857 """File name of multiple sequence alignment to restart PSI-BLAST. 858 859 Incompatible with: in_pssm, query""", 860 filename=True, 861 equate=False), 862 _Option(["-msa_master_idx", "msa_master_idx"], 863 """Index of sequence to use as master in MSA. 864 865 Index (1-based) of sequence to use as the master in the 866 multiple sequence alignment. If not specified, the first 867 sequence is used.""", 868 equate=False), 869 _Option(["-in_pssm", "in_pssm"], 870 """PSI-BLAST checkpoint file. 871 872 Incompatible with: in_msa, query, phi_pattern""", 873 filename=True, 874 equate=False), 875 # PSSM engine options: 876 _Option(["-pseudocount", "pseudocount"], 877 """Pseudo-count value used when constructing PSSM. 878 879 Integer. Default is zero.""", 880 equate=False), 881 _Option(["-inclusion_ethresh", "inclusion_ethresh"], 882 "E-value inclusion threshold for pairwise alignments (float, default 0.002).", 883 equate=False), 884 _Switch(["-ignore_msa_master", "ignore_msa_master"], 885 """Ignore the master sequence when creating PSSM 886 887 Requires: in_msa 888 Incompatible with: msa_master_idx, in_pssm, query, query_loc, phi_pattern 889 """), 890 # PHI-BLAST options: 891 _Option(["-phi_pattern", "phi_pattern"], 892 """File name containing pattern to search. 893 894 Incompatible with: in_pssm""", 895 filename=True, 896 equate=False), 897 ] 898 _Ncbiblast2SeqCommandline.__init__(self, cmd, **kwargs)
899
900 - def _validate(self):
901 incompatibles = {"num_iterations": ["remote"], 902 "in_msa": ["in_pssm", "query"], 903 "in_pssm": ["in_msa", "query", "phi_pattern"], 904 "ignore_msa_master": ["msa_master_idx", "in_pssm", 905 "query", "query_loc", "phi_pattern"], 906 } 907 self._validate_incompatibilities(incompatibles) 908 _Ncbiblast2SeqCommandline._validate(self)
909 910
911 -class NcbirpsblastCommandline(_NcbiblastCommandline):
912 """Wrapper for the NCBI BLAST+ program rpsblast. 913 914 With the release of BLAST+ (BLAST rewritten in C++ instead of C), the NCBI 915 replaced the old rpsblast tool with a similar tool of the same name. This 916 wrapper replaces RpsBlastCommandline, the wrapper for the old rpsblast. 917 918 >>> from Bio.Blast.Applications import NcbirpsblastCommandline 919 >>> cline = NcbirpsblastCommandline(help=True) 920 >>> cline 921 NcbirpsblastCommandline(cmd='rpsblast', help=True) 922 >>> print(cline) 923 rpsblast -help 924 925 You would typically run the command line with cline() or via the Python 926 subprocess module, as described in the Biopython tutorial. 927 """
928 - def __init__(self, cmd="rpsblast", **kwargs):
929 # TODO - remove the -word_size argument as per BLAST+ 2.2.30 930 # (BLAST team say it should never have been included, since 931 # the word size is set when building the domain database.) 932 # This likely means reviewing the class hierarchy again. 933 self.parameters = [ 934 # Query filtering options: 935 _Option(["-seg", "seg"], 936 """Filter query sequence with SEG (string). 937 938 Format: "yes", "window locut hicut", or "no" to disable. 939 Default is "12 2.2 2.5""", 940 equate=False), 941 # Restrict search or results: 942 _Option(["-culling_limit", "culling_limit"], 943 """Hit culling limit (integer). 944 945 If the query range of a hit is enveloped by that of at 946 least this many higher-scoring hits, delete the hit. 947 948 Incompatible with: best_hit_overhang, best_hit_score_edge. 949 """, 950 equate=False), 951 _Option(["-best_hit_overhang", "best_hit_overhang"], 952 """Best Hit algorithm overhang value (recommended value: 0.1) 953 954 Float between 0.0 and 0.5 inclusive. 955 956 Incompatible with: culling_limit.""", 957 equate=False), 958 _Option(["-best_hit_score_edge", "best_hit_score_edge"], 959 """Best Hit algorithm score edge value (recommended value: 0.1) 960 961 Float between 0.0 and 0.5 inclusive. 962 963 Incompatible with: culling_limit.""", 964 equate=False), 965 # General search options: 966 _Option(["-comp_based_stats", "comp_based_stats"], 967 """Use composition-based statistics. 968 969 D or d: default (equivalent to 0 ) 970 0 or F or f: Simplified Composition-based statistics as in 971 Bioinformatics 15:1000-1011, 1999 972 1 or T or t: Composition-based statistics as in NAR 29:2994-3005, 2001 973 974 Default = `0' 975 """, 976 checker_function=lambda value: value in "Dd0Ff1Tt", 977 equate=False), 978 # Misc options: 979 _Switch(["-use_sw_tback", "use_sw_tback"], 980 "Compute locally optimal Smith-Waterman alignments?"), 981 ] 982 _NcbiblastCommandline.__init__(self, cmd, **kwargs)
983
984 - def _validate(self):
985 incompatibles = {"culling_limit": ["best_hit_overhang", "best_hit_score_edge"]} 986 self._validate_incompatibilities(incompatibles) 987 _NcbiblastCommandline._validate(self)
988 989
990 -class NcbirpstblastnCommandline(_NcbiblastCommandline):
991 """Wrapper for the NCBI BLAST+ program rpstblastn. 992 993 With the release of BLAST+ (BLAST rewritten in C++ instead of C), the NCBI 994 replaced the old rpsblast tool with a similar tool of the same name, and a 995 separate tool rpstblastn for Translated Reverse Position Specific BLAST. 996 997 >>> from Bio.Blast.Applications import NcbirpstblastnCommandline 998 >>> cline = NcbirpstblastnCommandline(help=True) 999 >>> cline 1000 NcbirpstblastnCommandline(cmd='rpstblastn', help=True) 1001 >>> print(cline) 1002 rpstblastn -help 1003 1004 You would typically run the command line with cline() or via the Python 1005 subprocess module, as described in the Biopython tutorial. 1006 """
1007 - def __init__(self, cmd="rpstblastn", **kwargs):
1008 # TODO - remove the -word_size argument as per BLAST+ 2.2.30 1009 # (BLAST team say it should never have been included, since 1010 # the word size is set when building the domain database.) 1011 # This likely means reviewing the class hierarchy again. 1012 self.parameters = [ 1013 # Input query options: 1014 _Option(["-strand", "strand"], 1015 """Query strand(s) to search against database/subject. 1016 1017 Values allowed are "both" (default), "minus", "plus".""", 1018 checker_function=lambda value: value in ["both", 1019 "minus", 1020 "plus"], 1021 equate=False), 1022 # Input query options: 1023 _Option(["-query_gencode", "query_gencode"], 1024 "Genetic code to use to translate query (integer, default 1).", 1025 equate=False), 1026 # Query filtering options: 1027 _Option(["-seg", "seg"], 1028 """Filter query sequence with SEG (string). 1029 1030 Format: "yes", "window locut hicut", or "no" to disable. 1031 Default is "12 2.2 2.5""", 1032 equate=False), 1033 # Extension options: 1034 _Switch(["-ungapped", "ungapped"], 1035 "Perform ungapped alignment only?"), 1036 ] 1037 _NcbiblastCommandline.__init__(self, cmd, **kwargs)
1038 1039
1040 -class NcbiblastformatterCommandline(_NcbibaseblastCommandline):
1041 """Wrapper for the NCBI BLAST+ program blast_formatter. 1042 1043 With the release of BLAST 2.2.24+ (i.e. the BLAST suite rewritten in C++ 1044 instead of C), the NCBI added the ASN.1 output format option to all the 1045 search tools, and extended the blast_formatter to support this as input. 1046 1047 The blast_formatter command allows you to convert the ASN.1 output into 1048 the other output formats (XML, tabular, plain text, HTML). 1049 1050 >>> from Bio.Blast.Applications import NcbiblastformatterCommandline 1051 >>> cline = NcbiblastformatterCommandline(archive="example.asn", outfmt=5, out="example.xml") 1052 >>> cline 1053 NcbiblastformatterCommandline(cmd='blast_formatter', out='example.xml', outfmt=5, archive='example.asn') 1054 >>> print(cline) 1055 blast_formatter -out example.xml -outfmt 5 -archive example.asn 1056 1057 You would typically run the command line with cline() or via the Python 1058 subprocess module, as described in the Biopython tutorial. 1059 1060 Note that this wrapper is for the version of blast_formatter from BLAST 1061 2.2.24+ (or later) which is when the NCBI first announced the inclusion 1062 this tool. There was actually an early version in BLAST 2.2.23+ (and 1063 possibly in older releases) but this did not have the -archive option 1064 (instead -rid is a mandatory argument), and is not supported by this 1065 wrapper. 1066 """
1067 - def __init__(self, cmd="blast_formatter", **kwargs):
1068 self.parameters = [ 1069 # Input options 1070 _Option(["-rid", "rid"], 1071 "BLAST Request ID (RID), not compatible with archive arg", 1072 equate=False), 1073 _Option(["-archive", "archive"], 1074 "Archive file of results, not compatible with rid arg.", 1075 filename=True, 1076 equate=False), 1077 # Restrict search or results 1078 _Option(["-max_target_seqs", "max_target_seqs"], 1079 "Maximum number of aligned sequences to keep", 1080 checker_function=lambda value: value >= 1, 1081 equate=False), 1082 ] 1083 _NcbibaseblastCommandline.__init__(self, cmd, **kwargs)
1084
1085 - def _validate(self):
1086 incompatibles = {"rid": ["archive"]} 1087 self._validate_incompatibilities(incompatibles) 1088 _NcbibaseblastCommandline._validate(self)
1089 1090
1091 -class NcbideltablastCommandline(_Ncbiblast2SeqCommandline):
1092 """Create a commandline for the NCBI BLAST+ program deltablast (for proteins). 1093 1094 This is a wrapper for the deltablast command line command included in 1095 the NCBI BLAST+ software (not present in the original BLAST). 1096 1097 >>> from Bio.Blast.Applications import NcbideltablastCommandline 1098 >>> cline = NcbideltablastCommandline(query="rosemary.pro", db="nr", 1099 ... evalue=0.001, remote=True) 1100 >>> cline 1101 NcbideltablastCommandline(cmd='deltablast', query='rosemary.pro', db='nr', evalue=0.001, remote=True) 1102 >>> print(cline) 1103 deltablast -query rosemary.pro -db nr -evalue 0.001 -remote 1104 1105 You would typically run the command line with cline() or via the Python 1106 subprocess module, as described in the Biopython tutorial. 1107 """
1108 - def __init__(self, cmd="deltablast", **kwargs):
1109 self.parameters = [ 1110 # General search options: 1111 _Option(["-matrix", "matrix"], 1112 "Scoring matrix name (default BLOSUM62)."), 1113 _Option(["-threshold", "threshold"], 1114 "Minimum score for words to be added to the BLAST lookup table (float).", 1115 equate=False), 1116 _Option(["-comp_based_stats", "comp_based_stats"], 1117 """Use composition-based statistics (string, default 2, i.e. True). 1118 1119 0, F or f: no composition-based statistics. 1120 2, T or t, D or d : Composition-based score adjustment as in 1121 Bioinformatics 21:902-911, 2005, conditioned on sequence properties 1122 1123 Note that tblastn also supports values of 1 and 3.""", 1124 checker_function=lambda value: value in "0Ft2TtDd", 1125 equate=False), 1126 # Query filtering options: 1127 _Option(["-seg", "seg"], 1128 """Filter query sequence with SEG (string). 1129 1130 Format: "yes", "window locut hicut", or "no" to disable. 1131 Default is "12 2.2 2.5""", 1132 equate=False), 1133 # Extension options: 1134 _Option(["-gap_trigger", "gap_trigger"], 1135 "Number of bits to trigger gapping Default = 22", 1136 equate=False), 1137 # Miscellaneous options: 1138 _Switch(["-use_sw_tback", "use_sw_tback"], 1139 "Compute locally optimal Smith-Waterman alignments?"), 1140 # PSI-BLAST options 1141 _Option(["-num_iterations", "num_iterations"], 1142 """Number of iterations to perform. (integer >=1, Default is 1) 1143 1144 Incompatible with: remote""", 1145 equate=False), 1146 _Option(["-out_pssm", "out_pssm"], 1147 "File name to store checkpoint file.", 1148 filename=True, 1149 equate=False), 1150 _Option(["-out_ascii_pssm", "out_ascii_pssm"], 1151 "File name to store ASCII version of PSSM.", 1152 filename=True, 1153 equate=False), 1154 _Switch(["-save_pssm_after_last_round", "save_pssm_after_last_round"], 1155 "Save PSSM after the last database search."), 1156 _Switch(["-save_each_pssm", "save_each_pssm"], 1157 """Save PSSM after each iteration 1158 1159 File name is given in -save_pssm or -save_ascii_pssm options. 1160 """), 1161 # PSSM engine options 1162 _Option(["-pseudocount", "pseudocount"], 1163 "Pseudo-count value used when constructing PSSM (integer, default 0).", 1164 equate=False), 1165 _Option(["-domain_inclusion_ethresh", "domain_inclusion_ethresh"], 1166 """E-value inclusion threshold for alignments with conserved domains. 1167 1168 (float, Default is 0.05)""", 1169 equate=False), 1170 _Option(["-inclusion_ethresh", "inclusion_ethresh"], 1171 "Pairwise alignment e-value inclusion threshold (float, default 0.002).", 1172 equate=False), 1173 # DELTA-BLAST options 1174 _Option(["-rpsdb", "rpsdb"], 1175 "BLAST domain database name (dtring, Default = 'cdd_delta').", 1176 equate=False), 1177 _Switch(["-show_domain_hits", "show_domain_hits"], 1178 """Show domain hits? 1179 1180 Incompatible with: remote, subject""") 1181 ] 1182 _Ncbiblast2SeqCommandline.__init__(self, cmd, **kwargs)
1183 1184
1185 -def _test():
1186 """Run the Bio.Blast.Applications module's doctests.""" 1187 import doctest 1188 doctest.testmod(verbose=1)
1189 1190 if __name__ == "__main__": 1191 # Run the doctests 1192 _test() 1193