[Bioperl-l] Remote blast

Roopa Raghuveer rtbio.2009 at gmail.com
Wed Dec 2 12:07:08 UTC 2009


Hello everyone,

I have a problem. I am new to Bioperl. I am working on RNAi tool wherein a
cgi script was written which connects to NCBI blast using remote blast
program,i.e.,

The input sequence given in the html page is taken as input and Remote blast
is performed on this based on the code for Remote blast.But,I have a problem
in the Remote blast code.

My code goes like this

@compseqs=blastcode($in{'Inputseq'});

sub blastcode
{
$input1= $_[0];

open(NUC,'>',$nuc);
print NUC $input1;
close(NUC);

 my $prog = 'blastn';
 my $db   = 'refseq_rna';
 my $e_val= '1e-10';
 my $organism= 'Trypanosoma Brucei';

$gb = new Bio::DB::GenBank;

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

my $factory = Bio::Tools::Run::RemoteBlast->new(@params);

  #change a paramter
 $Bio::Tools::Run::RemoteBlast::HEADER{'ENTREZ_QUERY'} = 'Trypanosoma
brucei[ORGN]';

  my $v = 1;
  #$v is just to turn on and off the messages

 my $str = Bio::SeqIO->new(-file => $nuc , '-format' => 'fasta' ,
'-organism' => 'Trypanosoma Brucei' );


 while (my $input = $str->next_seq())

{
   #Blast a sequence against a database:
    #Alternatively, you could  pass in a file with many
    #sequences rather than loop through sequence one at a time
    #Remove the loop starting 'while (my $input = $str->next_seq())'
    #and swap the two lines below for an example of that.

   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 {
          my $result = $rc->next_result();
         #save the output
          my $filename = $result->query_name()."\.out";
           $factory->save_output($filename);
          $factory->remove_rid($rid);
         #       open(BLASTDEBUGFILE,'>',$blastdebugfile);
  #     print BLASTDEBUGFILE "Test1  $result";
   #     close(BLASTDEBUGFILE);

     open(OUTFILE,'>',$outfile);
     print OUTFILE "Test2 $result->database_name()";
     close(OUTFILE);

    while ( my $hit = $result->next_hit ) {
            next unless ( $v > 0);

              # open(OUTFILE,'>',$outfile);
              # print OUTFILE "in while hits";
              #close(OUTFILE);

       my $sequ = $gb->get_Seq_by_version($hit->name);
           my $dna = $sequ->seq();        # get the sequence as a string
                  push(@seqs,$dna);
          }
        }
      }
    }
}
# open(OUTFILE,'>',$outfile);
  #print OUTFILE $seqs[0];
 # close(OUTFILE);

return(@seqs);
}

Here in the above code,my program is able to go till the 'else' part and
writing the output file i.e.,this step.
my $filename = $result->query_name()."\.out";

But when I tried to enter in to the next while loop where I can get the
hits,the program is not entering into the while loop i.e.,

Not entering into this
while ( my $hit = $result->next_hit ) {
            next unless ( $v > 0);


Hence I am unable to get any hits for my query.
Ex:-If the query's accession number is Tb11.02.2210, I could just get a file
Tb11.02.2210.out file,it is just displaying the file name on the browser.

Please help me in solving this problem and mail me regarding any confusions.

Regards,
Roopa.



More information about the Bioperl-l mailing list