LynHsiong | 6 Jun 05:22 2013
Picon

Errors in a bioperl script, pls help me to have a check, thanks in advance

Dear friends:

I have a blast report, and want to extract following information for each
result(query): the Query_name, hit_number, name and description of the
hit(HSP) with the highest identity.

my bioperl script is:

#!/usr/bin/perl -w
use strict;
use warnings;
use Bio::SearchIO;

if ( <at> ARGV != 2){die "Usage: $0 <BLAST-report-file> <output-file>\n"};

my ($infile,$outfile) =  <at> ARGV;
print "Parsing the BLAST result ...";
my $blast = new Bio::SearchIO(
-format => 'blast',
-file   => $infile);
open (OUT,">$outfile") or die "Cannot open $outfile: $!";

print OUT "query\tNumber of hits\tGreatest identity %\tAccession (identity
%)\tDescription (identity %)\n";

while (my $result = $blast->next_result){
   print OUT $result->query_name . "\t";
   print OUT $result->num_hits. "\t";
        if (my $hit = &sort_hit($result)){
                         if (my $hsp=$hit->hsp){
                                        print OUT $hsp->percent_identity.
"\t"; 
                                        print OUT $hit->accession. "\t";
                                        print OUT $hit->description. "\n";
                            }                   
           }
    }
 close OUT;
 print " DONE!!!\n";

 sub sort_hit{
                  my $result = shift;
                  my  <at> hits;
                  while (my $hit = $result->next_hit) {
                        push  <at> hits, $hit;
                        }
                  my  <at> hit_sort = sort { $b-> hsp ->percent_identity * 100
<=> $a-> hsp ->percent_identity * 100 }  <at> hits;
                  $result->rewind;
                  my $flag=0;
                  my $temp_hit;
                  while (my $hit=$result->next_hit){
                        if($flag==0){
                           $temp_hit=$hit;
                     $flag++;
                           }
                  return $temp_hit
                        }
         }

error information:

-------------EXCEPTION: Bio::Root::Exception -------------
MSG: Can't get HSPs: data not collected.
STACK: Error::throw
STACK: Bio::Root::Root::throw
/home/xiongtl/Bioperl/dag/lib/perl5/Bio/Root/Root.pm:368
STACK: Bio::Search::Hit::GenericHit::hsp
/home/xiongtl/Bioperl/dag/lib/perl5/Bio/Search/Hit/GenericHit.pm:688
STACK: main::sort_hit xtl_new.pl:37
STACK: xtl_new.pl:20
-----------------------------------------------------------
Parsing the BLAST result ...

please give me some clues, all your words are welcome!

--
View this message in context: http://bioperl.996286.n3.nabble.com/Errors-in-a-bioperl-script-pls-help-me-to-have-a-check-thanks-in-advance-tp16986.html
Sent from the Bioperl-L mailing list archive at Nabble.com.

Gmane