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