[Bioperl-l] Parsing Blast Output

Brian Osborne brian_osborne@cognia.com
Tue, 11 Jun 2002 09:52:22 -0400


Mat and Anthony,

I'd recently modified the example code in bptutorial.pl so that it runs,
I've inserted it below. Note there's no call to next_result(). I'd also
modified the example code in RemoteBlast.pm, it should run as well but I'm
guessing these fixes are only in version bioperl-live and version 1.0.1.
Perhaps you're using 1.0?

Brian O.


$run_remoteblast = sub {
    print "\nBeginning run_remoteblast example... \n";
    eval { require Bio::Tools::Run::RemoteBlast; };

    if ( $@ ) {
      print STDERR "Cannot load Bio::Tools::Run::RemoteBlast\n";
      print STDERR "Cannot run run_remoteblast example:\n$@\n";
    } else {
      my (@params, $remote_blast_object, $blast_file, $r, $rc,
         $database);

      $database =  'ecoli';
      @params = ('-prog'   => 'blastp',
                 '-data'   => $database,
                 '-expect' => '1e-10',
                 '-readmethod' => 'BPlite' );

      $remote_blast_object = Bio::Tools::Run::RemoteBlast->new(@params);
      $blast_file = Bio::Root::IO->catfile("t","data","ecolitst.fa");
      $r = $remote_blast_object->submit_blast( $blast_file);
      print "submitted Blast job\n";
      while ( my @rids = $remote_blast_object->each_rid ) {
        foreach my $rid ( @rids ) {
          $rc = $remote_blast_object->retrieve_blast($rid);
          "retrieving results...\n";
          if( !ref($rc) ) {   # $rc not a reference => either error
                              # or job not yet finished
            if( $rc < 0 ) {
              $remote_blast_object->remove_rid($rid);
              print "Error return code for BlastID code $rid ... \n";
            }
            sleep 5;
          } else {
            $remote_blast_object->remove_rid($rid);
            while ( my $sbjct = $rc->nextSbjct ) {
              print "sbjct name is ", $sbjct->name, "\n";
              while ( my $hsp = $sbjct->nextHSP ) {
                print "score is ", $hsp->score, "\n";
              }
            }
          }
        }
      }
    }
  return 1;
} ;




-----Original Message-----
From: bioperl-l-admin@bioperl.org [mailto:bioperl-l-admin@bioperl.org]On
Behalf Of Wiepert, Mathieu
Sent: Tuesday, June 11, 2002 9:12 AM
To: 'bioperl-l@bioperl.org'
Subject: RE: [Bioperl-l] Parsing Blast Output

I got a bounce message, if this shows up twice, I apologize.


Hi,

I see that you copied the example from the RemoteBlast.pm, and it isn't
working.  I modified your code, and got an error from BPLite myself:

Can't locate object method "next_result" via package "Bio::Tools::BPlite" at
testremoteblast.pl line 32

Which is true (according to the documentation, I didn't look at the code for
the module).  But, the docs say basically say that your line of code:

my $rc = $factory->retrieve_blast($rid);

is already returning the BPLite object, so you don't need the line:

my $result = $rc->next_result (you already have it in $rc)

I modified the code to work (there were other errors in the script, the POD
should be fixed perhaps, to have working code).  My command line was
(because I don't have 1.0 bioperl installed):

perl -I  ~/bioperl_latest/lib/site_perl/5.6.0 testremoteblast.pl


#!/usr/bin/perl -w
use Bio::Tools::Run::RemoteBlast;
use strict;
my $v = 1;

my $prog = 'blastn';
my $db   = 'nr';
my $e_val= '1e-10';

my @params = ( '-prog' => $prog,
                 '-data' => $db,
                 '-expect' => $e_val );

  my $factory = Bio::Tools::Run::RemoteBlast->new(@params);
  $v = 1;
  my $str = Bio::SeqIO->new(-file=>'comt.cdna.fa' , '-format' => 'fasta' );
  my $input = $str->next_seq();
  #  Blast a sequence against a database:
  my $r = $factory->submit_blast($input);
  print STDERR "waiting..." if( $v > 0 );
  while ( my @rids = $factory->each_rid ) {
      foreach my $rid ( @rids ) {
          my $rc = $factory->retrieve_blast($rid);
          if( !ref($rc) ) {
              if( $rc < 0 )

                      $factory->remove_rid($rid);
                  }
              print STDERR "." if ( $v > 0 );
              sleep 5;
          } else {
              $factory->remove_rid($rid);
#             my $result = $rc->nextSbjct;
              print "db is ", $rc->database(), "\n";
              print "query is ", $rc->query(), "\n";
              my $count = 0;
              while( my $hit = $rc->nextSbjct )

                          $count++;
                          next unless ( $v > 0);
                          print "hit name is ", $hit->name, "\n";
                        while( my $hsp = $hit->nextHSP ) {
                          print "score is ", $hsp->score, "\n";
                        }
              }
          }
      }
  }


If there are errors in my code, I apologize, I am not the best at this
either :-)

-Mat
_______________________________________________
Bioperl-l mailing list
Bioperl-l@bioperl.org
http://bioperl.org/mailman/listinfo/bioperl-l