examples/searchio/hspwriter.pl

   1 #!/usr/bin/perl
   2
   3 # Demonstrates the use of a SearchIO Blast parser and a SearchWriterI object
   4 # for producing tab-delimited output of HSP data from a Blast report
   5 # input stream.
   6 #
   7 # Each row in the output represents data for a single HSP.
   8 #
   9 # This parser represents a new and improved version of Bio::Tools::Blast.
  10 #
  11 # Usage:
  12 #   STDIN:  stream containing one or more BLAST or PSI-BLAST reports.
  13 #   STDOUT: none, but generates an output file "hspwriter.out"
  14 #           containing tab-delimited data on a per-HSP basis.
  15 #   STDERR: Progress info and any errors.
  16 #
  17 # In this example, we create a SearchIO parser that screens out hits
  18 # based on expect (or P) scores and a default HSPTableWriter. This writer
  19 # provides the same functionality as the original Bio::Tools::Blast::table2()
  20 # function (i.e., a tab-delimited summary of each hit per row).
  21 # HSPTableWriter, however, is customizable so you can specify just the columns
  22 # you want to have in the output table.
  23 #
  24 # For more documentation about the writer, including
  25 # a complete list of columns, execute:
  26 #   perldoc Bio::SearchIO::Writer::HSPTableWriter.
  27 #
  28 # For more documentation about working with Blast result objects,
  29 # see docs for these modules:
  30 #   Bio::Search::Result::BlastResult
  31 #   Bio::Search::Iteration::IterationI
  32 #   Bio::Search::Hit::BlastHit
  33 #   Bio::Search::HSP::BlastHSP
  34 #
  35 # For more documentation about the Blast parser, see docs for
  36 #   Bio::SearchIO
  37 #
  38 # Author: Steve Chervitz <sac@bioperl.org>
  39
  40 use strict;
  41
  42 use Bio::SearchIO;
  43 use Bio::SearchIO::Writer::HSPTableWriter;
  44
  45 # These are the columns that will be in the output table of BLAST results.
  46 my @columns = qw(
  47                  query_name
  48                  query_length
  49                  hit_name
  50                  hit_length
  51                  rank
  52                  expect
  53                  frac_identical_query
  54                  length_aln_query
  55                  gaps_total
  56                  strand_query
  57                  strand_hit
  58                 );
  59
  60
  61 print STDERR "\nUsing SearchIO->new()\n";
  62
  63 # Note that all parameters for the $in, $out, and $writer objects are optional.
  64 # Default in = STDIN; Default out = STDOUT; Default writer = all columns
  65 # In this example, we're reading from STDIN and  writing to a STDOUT
  66 my $in     = Bio::SearchIO->new( -format => 'blast',
  67                                  -fh => \*ARGV
  68                                );
  69 my $writer = Bio::SearchIO::Writer::HSPTableWriter->new( -columns => \@columns );
  70 my $out    = Bio::SearchIO->new( -format => 'blast',
  71                                  -writer => $writer,
  72                                  -file   => ">hspwriter.out" );
  73
  74 while ( my $result = $in->next_result() ) {
  75   printf STDERR "\nReport %d: $result\n", $in->result_count;
  76
  77   if( $result->hits ) {
  78     $out->write_result($result, ($in->result_count - 1 ? 0 : 1) );
  79   }
  80   else {
  81     print STDERR "Hitless Blast Report: $result ";
  82     print STDERR ($result->no_hits_found ? "\n" : "(filtered)\n");
  83   }
  84
  85   ## For a simple progress monitor, uncomment this line:
  86   #print STDERR "."; print STDERR "\n" if $in->result_count % 50 == 0;
  87 }
  88
  89 printf STDERR "\n%d Blast report(s) processed.\n", $in->result_count;
  90 printf STDERR "Output sent to file: %s\n",  $out->file if $out->file;
  91