[Bioperl-l] How to parse bl2seq report?
szhan at uoguelph.ca
szhan at uoguelph.ca
Mon Apr 19 21:59:26 EDT 2004
Hello, Bioperl experts:
I want to parse an array of pairwise protein sequences alingments with BLAST
bl2seq (standlone) in a loop. When there is no any "subject" (hit) in a bl2seq
report, the program halted. Could you please help me catch the error and
continue the loop?
Your help will be highly appreciated!
The partial code like this:
@lucsarr is array of protein sequence objects
for(my $k=0; $k<=$#lucsarr; $k++){
for(my $l=$k+1; $l<=$#lucsarr; $l++){
if($k<$l && $l<=$#lucsarr){
my $input1= $lucsarr[$k]; # take translated protein as an object
my $input2= $lucsarr[$l];
# create local factory object
my $factory = Bio::Tools::Run::StandAloneBlast->new('outfile'
=> 'bl2seq.out');
my $prgm='blastp'; # set parameter p for blastp
$factory->p($prgm);
# call executale bl2seq
my $blast_report=$factory->bl2seq($input1, $input2);
# retrieve the outputs from blast_report output file: bl2seq.out
$blast_report->sbjctName;
$blast_report->sbjctLength;
print "Score\tBits\tID_percent\tE_value\tMatch\tPositive\tLength\n";
while(my $hsp = $blast_report->next_feature) {# get all
alignments
print $hsp->score, "\t";
print $hsp->bits, "\t";
print $hsp->percent, "\t";
print $hsp->P, "\t";
print $hsp->match, "\t";
print $hsp->positive, "\t";
print $hsp->length, "\n";
}
}
}
}
The output messages as follows:
Score Bits ID_percent E_value Match Positive Length
Can't call method "nextHSP" on unblessed reference at
C:/Perl/site/lib/Bio/Tools
/BPbl2seq.pm line 236, <GEN2> line 36.
Again thank you in advance!
Joshua
More information about the Bioperl-l
mailing list