[Bioperl-l] Getting description of BLAST hits
Peter.Robinson at t-online.de
Peter.Robinson at t-online.de
Tue Apr 5 18:09:25 EDT 2005
Dear list,
I am trying to get a list of all BLAST hits to various sequences using code adapted from the bioperl website. In essence, I would like to get the descriptions
from lines such as
gb|BC021898.1| Homo sapiens adaptor-related protein complex 1, s... 2956 0.0
Using the following code:
while ( my $hit = $result->next_hit ) {
next unless ( $v > 0);
print "\thit name is ", $hit->name, "\t",
$hit->accession(), "\t", hit->description(),"\t";
leads to the error:
Can't locate object method "description" via package "hit" (perhaps you forgot to load "hit"?) at blast2SpeciesList.pl line 60.
However, the documentation for Bio::Search::Hit indicates that Hit objects should have a description() method.
What am I misunderstanding?
Thanks in advance for any tips. Just in case, I am pasting the complete code snippet at the bottom of this mail.
Peter
#!/usr/bin/perl -w
use strict;
#Remote-blast "factory object" creation and blast-parameter initialization
use Bio::Tools::Run::RemoteBlast;
my $prog = 'blastn';
my $db = 'nr';
my $e_val= '1e-10';
my @params = ( '-prog' => $prog,
'-data' => $db,
'-expect' => $e_val,
'-readmethod' => 'SearchIO' );
my $factory = Bio::Tools::Run::RemoteBlast->new(@params);
#change a paramter
#$Bio::Tools::Run::RemoteBlast::HEADER{'ENTREZ_QUERY'} = 'Homo sapiens [ORGN]';
#remove a parameter
#delete $Bio::Tools::Run::RemoteBlast::HEADER{'FILTER'};
my $v = 1;
#$v is just to turn on and off the messages
my $str = Bio::SeqIO->new(-file=>'testSeqs.fa' , '-format' => 'fasta' );
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);
#my $r = $factory->submit_blast('amino.fa');
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);
print "\nQuery Name: ", $result->query_name(), "\t",
$result->query_accession(),"\n";
while ( my $hit = $result->next_hit ) {
next unless ( $v > 0);
print "\thit name is ", $hit->name, "\t",
$hit->accession(), "\t", hit->description(),"\t";
while( my $hsp = $hit->next_hsp ) {
print "\te-val is ", $hsp->evalue, "\n";
last;
}
}
}
}
}
}
More information about the Bioperl-l
mailing list