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 maximum number of HSPs per subject to save for ungapped searches (integer).", 189 equate=False), 190 _Option(["-max_hsps", "max_hsps"], 191 "Set maximum number of HSPs per subject sequence to save (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 Bioinformatics 21:902-911, 2005, conditioned on sequence properties 618 3: Composition-based score adjustment as in Bioinformatics 21:902-911, 2005, unconditionally 619 620 For programs other than tblastn, must either be absent or be D, F or 0 621 Default = `2' 622 """, 623 equate=False), 624 # Query filtering options: 625 _Option(["-seg", "seg"], 626 """Filter query sequence with SEG (string). 627 628 Format: "yes", "window locut hicut", or "no" to disable. 629 Default is "12 2.2 2.5""", 630 equate=False), 631 # Extension options: 632 _Switch(["-ungapped", "ungapped"], 633 "Perform ungapped alignment only?"), 634 _Switch(["-use_sw_tback", "use_sw_tback"], 635 "Compute locally optimal Smith-Waterman alignments?"), 636 ] 637 _NcbiblastMain2SeqCommandline.__init__(self, cmd, **kwargs)
638 639
640 -class NcbitblastnCommandline(_NcbiblastMain2SeqCommandline):
641 """Wrapper for the NCBI BLAST+ program tblastn. 642 643 With the release of BLAST+ (BLAST rewritten in C++ instead of C), the NCBI 644 replaced the old blastall tool with separate tools for each of the searches. 645 This wrapper therefore replaces BlastallCommandline with option -p tblastn. 646 647 >>> from Bio.Blast.Applications import NcbitblastnCommandline 648 >>> cline = NcbitblastnCommandline(help=True) 649 >>> cline 650 NcbitblastnCommandline(cmd='tblastn', help=True) 651 >>> print(cline) 652 tblastn -help 653 654 You would typically run the command line with cline() or via the Python 655 subprocess module, as described in the Biopython tutorial. 656 """
657 - def __init__(self, cmd="tblastn", **kwargs):
658 self.parameters = [ 659 # General search options: 660 _Option(["-task", "task"], 661 "Task to execute (string, tblastn (default) or tblastn-fast).", 662 checker_function=lambda value: value in ["tblastn", 663 "tblastn-fast"], 664 equate=False), 665 _Option(["-db_gencode", "db_gencode"], 666 "Genetic code to use to translate query (integer, default 1).", 667 equate=False), 668 _Option(["-frame_shift_penalty", "frame_shift_penalty"], 669 """Frame shift penalty (integer, at least 1, default ignored) (OBSOLETE). 670 671 This was removed in BLAST 2.2.27+""", 672 equate=False), 673 _Option(["-max_intron_length", "max_intron_length"], 674 """Maximum intron length (integer). 675 676 Length of the largest intron allowed in a translated nucleotide 677 sequence when linking multiple distinct alignments (a negative 678 value disables linking). Default zero.""", 679 equate=False), 680 _Option(["-matrix", "matrix"], 681 "Scoring matrix name (default BLOSUM62).", 682 equate=False), 683 _Option(["-threshold", "threshold"], 684 "Minimum score for words to be added to the BLAST lookup table (float).", 685 equate=False), 686 _Option(["-comp_based_stats", "comp_based_stats"], 687 """Use composition-based statistics (string, default 2, i.e. True). 688 689 0, F or f: no composition-based statistics 690 1: Composition-based statistics as in NAR 29:2994-3005, 2001 691 2, T or t, D or d : Composition-based score adjustment as in 692 Bioinformatics 21:902-911, 2005, conditioned on sequence properties 693 3: Composition-based score adjustment as in Bioinformatics 21:902-911, 694 2005, unconditionally 695 696 Note that only tblastn supports values of 1 and 3.""", 697 checker_function=lambda value: value in "0Ft12TtDd3", 698 equate=False), 699 # Query filtering options: 700 _Option(["-seg", "seg"], 701 """Filter query sequence with SEG (string). 702 703 Format: "yes", "window locut hicut", or "no" to disable. 704 Default is "12 2.2 2.5""", 705 equate=False), 706 # Extension options: 707 _Switch(["-ungapped", "ungapped"], 708 "Perform ungapped alignment only?"), 709 # Miscellaneous options: 710 _Switch(["-use_sw_tback", "use_sw_tback"], 711 "Compute locally optimal Smith-Waterman alignments?"), 712 # PSI-TBLASTN options: 713 _Option(["-in_pssm", "in_pssm"], 714 """PSI-BLAST checkpoint file 715 716 Incompatible with: remote, query""", 717 filename=True, 718 equate=False), 719 ] 720 _NcbiblastMain2SeqCommandline.__init__(self, cmd, **kwargs)
721 722
723 -class NcbitblastxCommandline(_NcbiblastMain2SeqCommandline):
724 """Wrapper for the NCBI BLAST+ program tblastx. 725 726 With the release of BLAST+ (BLAST rewritten in C++ instead of C), the NCBI 727 replaced the old blastall tool with separate tools for each of the searches. 728 This wrapper therefore replaces BlastallCommandline with option -p tblastx. 729 730 >>> from Bio.Blast.Applications import NcbitblastxCommandline 731 >>> cline = NcbitblastxCommandline(help=True) 732 >>> cline 733 NcbitblastxCommandline(cmd='tblastx', help=True) 734 >>> print(cline) 735 tblastx -help 736 737 You would typically run the command line with cline() or via the Python 738 subprocess module, as described in the Biopython tutorial. 739 """
740 - def __init__(self, cmd="tblastx", **kwargs):
741 self.parameters = [ 742 # Input query options: 743 _Option(["-strand", "strand"], 744 """Query strand(s) to search against database/subject. 745 746 Values allowed are "both" (default), "minus", "plus".""", 747 checker_function=lambda value: value in ["both", "minus", "plus"], 748 equate=False), 749 # Input query options: 750 _Option(["-query_gencode", "query_gencode"], 751 "Genetic code to use to translate query (integer, default 1).", 752 equate=False), 753 # General search options: 754 _Option(["-db_gencode", "db_gencode"], 755 "Genetic code to use to translate query (integer, default 1).", 756 equate=False), 757 _Option(["-max_intron_length", "max_intron_length"], 758 """Maximum intron length (integer). 759 760 Length of the largest intron allowed in a translated nucleotide 761 sequence when linking multiple distinct alignments (a negative 762 value disables linking). Default zero.""", 763 equate=False), 764 _Option(["-matrix", "matrix"], 765 "Scoring matrix name (default BLOSUM62).", 766 equate=False), 767 _Option(["-threshold", "threshold"], 768 "Minimum score for words to be added to the BLAST lookup table (float).", 769 equate=False), 770 # Query filtering options: 771 _Option(["-seg", "seg"], 772 """Filter query sequence with SEG (string). 773 774 Format: "yes", "window locut hicut", or "no" to disable. 775 Default is "12 2.2 2.5""", 776 equate=False), 777 ] 778 _NcbiblastMain2SeqCommandline.__init__(self, cmd, **kwargs)
779 780
781 -class NcbipsiblastCommandline(_Ncbiblast2SeqCommandline):
782 """Wrapper for the NCBI BLAST+ program psiblast. 783 784 With the release of BLAST+ (BLAST rewritten in C++ instead of C), the NCBI 785 replaced the old blastpgp tool with a similar tool psiblast. This wrapper 786 therefore replaces BlastpgpCommandline, the wrapper for blastpgp. 787 788 >>> from Bio.Blast.Applications import NcbipsiblastCommandline 789 >>> cline = NcbipsiblastCommandline(help=True) 790 >>> cline 791 NcbipsiblastCommandline(cmd='psiblast', help=True) 792 >>> print(cline) 793 psiblast -help 794 795 You would typically run the command line with cline() or via the Python 796 subprocess module, as described in the Biopython tutorial. 797 """
798 - def __init__(self, cmd="psiblast", **kwargs):
799 self.parameters = [ 800 # General search options: 801 _Option(["-matrix", "matrix"], 802 "Scoring matrix name (default BLOSUM62).", 803 equate=False), 804 _Option(["-threshold", "threshold"], 805 "Minimum score for words to be added to the BLAST lookup table (float).", 806 equate=False), 807 _Option(["-comp_based_stats", "comp_based_stats"], 808 """Use composition-based statistics (string, default 2, i.e. True). 809 810 0, F or f: no composition-based statistics 811 2, T or t, D or d : Composition-based score adjustment 812 as in Bioinformatics 21:902-911, 2005, conditioned on 813 sequence properties 814 815 Note that tblastn also supports values of 1 and 3.""", 816 checker_function=lambda value: value in "0Ft2TtDd", 817 equate=False), 818 # Query filtering options: 819 _Option(["-seg", "seg"], 820 """Filter query sequence with SEG (string). 821 822 Format: "yes", "window locut hicut", or "no" to disable. 823 Default is "12 2.2 2.5""", 824 equate=False), 825 # Extension options: 826 _Option(["-gap_trigger", "gap_trigger"], 827 "Number of bits to trigger gapping (float, default 22).", 828 equate=False), 829 # Miscellaneous options: 830 _Switch(["-use_sw_tback", "use_sw_tback"], 831 "Compute locally optimal Smith-Waterman alignments?"), 832 # PSI-BLAST options: 833 _Option(["-num_iterations", "num_iterations"], 834 """Number of iterations to perform (integer, at least one). 835 836 Default is one. 837 Incompatible with: remote""", 838 equate=False), 839 _Option(["-out_pssm", "out_pssm"], 840 "File name to store checkpoint file.", 841 filename=True, 842 equate=False), 843 _Option(["-out_ascii_pssm", "out_ascii_pssm"], 844 "File name to store ASCII version of PSSM.", 845 filename=True, 846 equate=False), 847 _Switch(["-save_pssm_after_last_round", "save_pssm_after_last_round"], 848 "Save PSSM after the last database search."), 849 _Switch(["-save_each_pssm", "save_each_pssm"], 850 """Save PSSM after each iteration 851 852 File name is given in -save_pssm or -save_ascii_pssm options. 853 """), 854 _Option(["-in_msa", "in_msa"], 855 """File name of multiple sequence alignment to restart PSI-BLAST. 856 857 Incompatible with: in_pssm, query""", 858 filename=True, 859 equate=False), 860 _Option(["-msa_master_idx", "msa_master_idx"], 861 """Index of sequence to use as master in MSA. 862 863 Index (1-based) of sequence to use as the master in the 864 multiple sequence alignment. If not specified, the first 865 sequence is used.""", 866 equate=False), 867 _Option(["-in_pssm", "in_pssm"], 868 """PSI-BLAST checkpoint file. 869 870 Incompatible with: in_msa, query, phi_pattern""", 871 filename=True, 872 equate=False), 873 # PSSM engine options: 874 _Option(["-pseudocount", "pseudocount"], 875 """Pseudo-count value used when constructing PSSM. 876 877 Integer. Default is zero.""", 878 equate=False), 879 _Option(["-inclusion_ethresh", "inclusion_ethresh"], 880 "E-value inclusion threshold for pairwise alignments (float, default 0.002).", 881 equate=False), 882 _Switch(["-ignore_msa_master", "ignore_msa_master"], 883 """Ignore the master sequence when creating PSSM 884 885 Requires: in_msa 886 Incompatible with: msa_master_idx, in_pssm, query, query_loc, phi_pattern 887 """), 888 # PHI-BLAST options: 889 _Option(["-phi_pattern", "phi_pattern"], 890 """File name containing pattern to search. 891 892 Incompatible with: in_pssm""", 893 filename=True, 894 equate=False), 895 ] 896 _Ncbiblast2SeqCommandline.__init__(self, cmd, **kwargs)
897
898 - def _validate(self):
899 incompatibles = {"num_iterations": ["remote"], 900 "in_msa": ["in_pssm", "query"], 901 "in_pssm": ["in_msa", "query", "phi_pattern"], 902 "ignore_msa_master": ["msa_master_idx", "in_pssm", 903 "query", "query_loc", "phi_pattern"], 904 } 905 self._validate_incompatibilities(incompatibles) 906 _Ncbiblast2SeqCommandline._validate(self)
907 908
909 -class NcbirpsblastCommandline(_NcbiblastCommandline):
910 """Wrapper for the NCBI BLAST+ program rpsblast. 911 912 With the release of BLAST+ (BLAST rewritten in C++ instead of C), the NCBI 913 replaced the old rpsblast tool with a similar tool of the same name. This 914 wrapper replaces RpsBlastCommandline, the wrapper for the old rpsblast. 915 916 >>> from Bio.Blast.Applications import NcbirpsblastCommandline 917 >>> cline = NcbirpsblastCommandline(help=True) 918 >>> cline 919 NcbirpsblastCommandline(cmd='rpsblast', help=True) 920 >>> print(cline) 921 rpsblast -help 922 923 You would typically run the command line with cline() or via the Python 924 subprocess module, as described in the Biopython tutorial. 925 """
926 - def __init__(self, cmd="rpsblast", **kwargs):
927 # TODO - remove the -word_size argument as per BLAST+ 2.2.30 928 # (BLAST team say it should never have been included, since 929 # the word size is set when building the domain database.) 930 # This likely means reviewing the class hierarchy again. 931 self.parameters = [ 932 # Query filtering options: 933 _Option(["-seg", "seg"], 934 """Filter query sequence with SEG (string). 935 936 Format: "yes", "window locut hicut", or "no" to disable. 937 Default is "12 2.2 2.5""", 938 equate=False), 939 # Restrict search or results: 940 _Option(["-culling_limit", "culling_limit"], 941 """Hit culling limit (integer). 942 943 If the query range of a hit is enveloped by that of at 944 least this many higher-scoring hits, delete the hit. 945 946 Incompatible with: best_hit_overhang, best_hit_score_edge. 947 """, 948 equate=False), 949 _Option(["-best_hit_overhang", "best_hit_overhang"], 950 """Best Hit algorithm overhang value (recommended value: 0.1) 951 952 Float between 0.0 and 0.5 inclusive. 953 954 Incompatible with: culling_limit.""", 955 equate=False), 956 _Option(["-best_hit_score_edge", "best_hit_score_edge"], 957 """Best Hit algorithm score edge value (recommended value: 0.1) 958 959 Float between 0.0 and 0.5 inclusive. 960 961 Incompatible with: culling_limit.""", 962 equate=False), 963 # General search options: 964 _Option(["-comp_based_stats", "comp_based_stats"], 965 """Use composition-based statistics. 966 967 D or d: default (equivalent to 0 ) 968 0 or F or f: Simplified Composition-based statistics as in 969 Bioinformatics 15:1000-1011, 1999 970 1 or T or t: Composition-based statistics as in NAR 29:2994-3005, 2001 971 972 Default = `0' 973 """, 974 checker_function=lambda value: value in "Dd0Ff1Tt", 975 equate=False), 976 # Misc options: 977 _Switch(["-use_sw_tback", "use_sw_tback"], 978 "Compute locally optimal Smith-Waterman alignments?"), 979 ] 980 _NcbiblastCommandline.__init__(self, cmd, **kwargs)
981
982 - def _validate(self):
983 incompatibles = {"culling_limit": ["best_hit_overhang", "best_hit_score_edge"]} 984 self._validate_incompatibilities(incompatibles) 985 _NcbiblastCommandline._validate(self)
986 987
988 -class NcbirpstblastnCommandline(_NcbiblastCommandline):
989 """Wrapper for the NCBI BLAST+ program rpstblastn. 990 991 With the release of BLAST+ (BLAST rewritten in C++ instead of C), the NCBI 992 replaced the old rpsblast tool with a similar tool of the same name, and a 993 separate tool rpstblastn for Translated Reverse Position Specific BLAST. 994 995 >>> from Bio.Blast.Applications import NcbirpstblastnCommandline 996 >>> cline = NcbirpstblastnCommandline(help=True) 997 >>> cline 998 NcbirpstblastnCommandline(cmd='rpstblastn', help=True) 999 >>> print(cline) 1000 rpstblastn -help 1001 1002 You would typically run the command line with cline() or via the Python 1003 subprocess module, as described in the Biopython tutorial. 1004 """
1005 - def __init__(self, cmd="rpstblastn", **kwargs):
1006 # TODO - remove the -word_size argument as per BLAST+ 2.2.30 1007 # (BLAST team say it should never have been included, since 1008 # the word size is set when building the domain database.) 1009 # This likely means reviewing the class hierarchy again. 1010 self.parameters = [ 1011 # Input query options: 1012 _Option(["-strand", "strand"], 1013 """Query strand(s) to search against database/subject. 1014 1015 Values allowed are "both" (default), "minus", "plus".""", 1016 checker_function=lambda value: value in ["both", 1017 "minus", 1018 "plus"], 1019 equate=False), 1020 # Input query options: 1021 _Option(["-query_gencode", "query_gencode"], 1022 "Genetic code to use to translate query (integer, default 1).", 1023 equate=False), 1024 # Query filtering options: 1025 _Option(["-seg", "seg"], 1026 """Filter query sequence with SEG (string). 1027 1028 Format: "yes", "window locut hicut", or "no" to disable. 1029 Default is "12 2.2 2.5""", 1030 equate=False), 1031 # Extension options: 1032 _Switch(["-ungapped", "ungapped"], 1033 "Perform ungapped alignment only?"), 1034 ] 1035 _NcbiblastCommandline.__init__(self, cmd, **kwargs)
1036 1037
1038 -class NcbiblastformatterCommandline(_NcbibaseblastCommandline):
1039 """Wrapper for the NCBI BLAST+ program blast_formatter. 1040 1041 With the release of BLAST 2.2.24+ (i.e. the BLAST suite rewritten in C++ 1042 instead of C), the NCBI added the ASN.1 output format option to all the 1043 search tools, and extended the blast_formatter to support this as input. 1044 1045 The blast_formatter command allows you to convert the ASN.1 output into 1046 the other output formats (XML, tabular, plain text, HTML). 1047 1048 >>> from Bio.Blast.Applications import NcbiblastformatterCommandline 1049 >>> cline = NcbiblastformatterCommandline(archive="example.asn", outfmt=5, out="example.xml") 1050 >>> cline 1051 NcbiblastformatterCommandline(cmd='blast_formatter', out='example.xml', outfmt=5, archive='example.asn') 1052 >>> print(cline) 1053 blast_formatter -out example.xml -outfmt 5 -archive example.asn 1054 1055 You would typically run the command line with cline() or via the Python 1056 subprocess module, as described in the Biopython tutorial. 1057 1058 Note that this wrapper is for the version of blast_formatter from BLAST 1059 2.2.24+ (or later) which is when the NCBI first announced the inclusion 1060 this tool. There was actually an early version in BLAST 2.2.23+ (and 1061 possibly in older releases) but this did not have the -archive option 1062 (instead -rid is a mandatory argument), and is not supported by this 1063 wrapper. 1064 """
1065 - def __init__(self, cmd="blast_formatter", **kwargs):
1066 self.parameters = [ 1067 # Input options 1068 _Option(["-rid", "rid"], 1069 "BLAST Request ID (RID), not compatible with archive arg", 1070 equate=False), 1071 _Option(["-archive", "archive"], 1072 "Archive file of results, not compatible with rid arg.", 1073 filename=True, 1074 equate=False), 1075 # Restrict search or results 1076 _Option(["-max_target_seqs", "max_target_seqs"], 1077 "Maximum number of aligned sequences to keep", 1078 checker_function=lambda value: value >= 1, 1079 equate=False), 1080 ] 1081 _NcbibaseblastCommandline.__init__(self, cmd, **kwargs)
1082
1083 - def _validate(self):
1084 incompatibles = {"rid": ["archive"]} 1085 self._validate_incompatibilities(incompatibles) 1086 _NcbibaseblastCommandline._validate(self)
1087 1088
1089 -class NcbideltablastCommandline(_Ncbiblast2SeqCommandline):
1090 """Create a commandline for the NCBI BLAST+ program deltablast (for proteins). 1091 1092 This is a wrapper for the deltablast command line command included in 1093 the NCBI BLAST+ software (not present in the original BLAST). 1094 1095 >>> from Bio.Blast.Applications import NcbideltablastCommandline 1096 >>> cline = NcbideltablastCommandline(query="rosemary.pro", db="nr", 1097 ... evalue=0.001, remote=True) 1098 >>> cline 1099 NcbideltablastCommandline(cmd='deltablast', query='rosemary.pro', db='nr', evalue=0.001, remote=True) 1100 >>> print(cline) 1101 deltablast -query rosemary.pro -db nr -evalue 0.001 -remote 1102 1103 You would typically run the command line with cline() or via the Python 1104 subprocess module, as described in the Biopython tutorial. 1105 """
1106 - def __init__(self, cmd="deltablast", **kwargs):
1107 self.parameters = [ 1108 # General search options: 1109 _Option(["-matrix", "matrix"], 1110 "Scoring matrix name (default BLOSUM62)."), 1111 _Option(["-threshold", "threshold"], 1112 "Minimum score for words to be added to the BLAST lookup table (float).", 1113 equate=False), 1114 _Option(["-comp_based_stats", "comp_based_stats"], 1115 """Use composition-based statistics (string, default 2, i.e. True). 1116 1117 0, F or f: no composition-based statistics. 1118 2, T or t, D or d : Composition-based score adjustment as in 1119 Bioinformatics 21:902-911, 2005, conditioned on sequence properties 1120 1121 Note that tblastn also supports values of 1 and 3.""", 1122 checker_function=lambda value: value in "0Ft2TtDd", 1123 equate=False), 1124 # Query filtering options: 1125 _Option(["-seg", "seg"], 1126 """Filter query sequence with SEG (string). 1127 1128 Format: "yes", "window locut hicut", or "no" to disable. 1129 Default is "12 2.2 2.5""", 1130 equate=False), 1131 # Extension options: 1132 _Option(["-gap_trigger", "gap_trigger"], 1133 "Number of bits to trigger gapping Default = 22", 1134 equate=False), 1135 # Miscellaneous options: 1136 _Switch(["-use_sw_tback", "use_sw_tback"], 1137 "Compute locally optimal Smith-Waterman alignments?"), 1138 # PSI-BLAST options 1139 _Option(["-num_iterations", "num_iterations"], 1140 """Number of iterations to perform. (integer >=1, Default is 1) 1141 1142 Incompatible with: remote""", 1143 equate=False), 1144 _Option(["-out_pssm", "out_pssm"], 1145 "File name to store checkpoint file.", 1146 filename=True, 1147 equate=False), 1148 _Option(["-out_ascii_pssm", "out_ascii_pssm"], 1149 "File name to store ASCII version of PSSM.", 1150 filename=True, 1151 equate=False), 1152 _Switch(["-save_pssm_after_last_round", "save_pssm_after_last_round"], 1153 "Save PSSM after the last database search."), 1154 _Switch(["-save_each_pssm", "save_each_pssm"], 1155 """Save PSSM after each iteration 1156 1157 File name is given in -save_pssm or -save_ascii_pssm options. 1158 """), 1159 # PSSM engine options 1160 _Option(["-pseudocount", "pseudocount"], 1161 "Pseudo-count value used when constructing PSSM (integer, default 0).", 1162 equate=False), 1163 _Option(["-domain_inclusion_ethresh", "domain_inclusion_ethresh"], 1164 """E-value inclusion threshold for alignments with conserved domains. 1165 1166 (float, Default is 0.05)""", 1167 equate=False), 1168 _Option(["-inclusion_ethresh", "inclusion_ethresh"], 1169 "E-value inclusion threshold for pairwise alignments (float, Default is 0.002).", 1170 equate=False), 1171 # DELTA-BLAST options 1172 _Option(["-rpsdb", "rpsdb"], 1173 "BLAST domain database name (dtring, Default = 'cdd_delta').", 1174 equate=False), 1175 _Switch(["-show_domain_hits", "show_domain_hits"], 1176 """Show domain hits? 1177 1178 Incompatible with: remote, subject""") 1179 ] 1180 _Ncbiblast2SeqCommandline.__init__(self, cmd, **kwargs)
1181 1182
1183 -def _test():
1184 """Run the Bio.Blast.Applications module's doctests.""" 1185 import doctest 1186 doctest.testmod(verbose=1)
1187 1188 if __name__ == "__main__": 1189 # Run the doctests 1190 _test() 1191