[Bioperl-l] num_hits and next_hit() bug in Bio::Search::Result::BlastResult?

[Jason Stajich]

You always get a hit with BL2Seq with the Query= and Subject= in the 
report - but you'll notice for your parsed example there are no HSPs.  
The fact that the bl2seq format reports the subject means you get a hit 
no matter what.

What exactly are you trying to do this seems like a bad example. Just 
run a BLAST of a query against a database or look at the example reports 
in the t/data directory if you want to try and parse a BLAST report.

If you are using short reads to do alignments to reference sequences you 
really don't want to be using BLAST anyways.

-jason

Peng Yu wrote, On 6/6/10 11:04 PM:
> I have the following two sequences which don't have any hits. But the
> perl code gives num_hits =1 (see below). Is it a bug in parsing blast
> results?
>
> $ make
> blastn -query<(echo
> CACGACAACCCTGTATGCACACTGGTGTCTTCCAGATCGGAAGGGCGGTTCAGCAGGAATGCCGAGGCCGGTATC)
> -subject<(echo
> GGAAGACACCAGTGTGCATACAGGGGTGTCGTGAGATCGGAAGAGCGTCGTGGAGGGAAAGAGTGGAGATCTCGG)
>    
>> HWI-EAS11X_10097_4_1_1523_15064.txt
>>      
> $ cat first1.fa
>    
>> first1
>>      
> CACGACAACCCTGTATGCACACTGGTGTCTTCCAGATCGGAAGGGCGGTTCAGCAGGAATGCCGAGGCCGGTATC
> $ cat second1.fa
>    
>> second1
>>      
> GGAAGACACCAGTGTGCATACAGGGGTGTCGTGAGATCGGAAGAGCGTCGTGGAGGGAAAGAGTGGAGATCTCGG
>
>
> $ cat HWI-EAS11X_10097_4_1_1523_15064.txt
> BLASTN 2.2.23+
>
>
> Query=
> Length=75
>
> Subject=
> Length=75
>
>
> ***** No hits found *****
>
>
>
> Lambda     K      H
>      1.33    0.621     1.12
>
> Gapped
> Lambda     K      H
>      1.28    0.460    0.850
>
> Effective search space used: 4761
>
>
>
>
> Matrix: blastn matrix 1 -2
> Gap Penalties: Existence: 0, Extension: 0
>
> #####
>
> But the following code and output shows that there is one hit?
>
> $ cat main.pl
> #!/usr/bin/env perl
>
> use strict;
> use warnings;
> use Bio::Tools::Run::StandAloneBlastPlus;
> use Bio::Perl;
> use Data::Dumper;
>
> my $factory = Bio::Tools::Run::StandAloneBlastPlus->new();
>
> my $seq1 = Bio::Perl::read_sequence('first1.fa');
> my $seq2 = Bio::Perl::read_sequence('second1.fa');
>
> my $result=$factory->bl2seq(-method=>'blastn',
>    -query=>  $seq1,
>    -subject=>  $seq2
> );
>
> $factory->cleanup;
>
> print 'ref($factory): ', ref($factory), "\n";
> print 'ref($result): ', ref($result), "\n";
>
> print "num_hits: ",  $result->num_hits, "\n";
> while(my $hit = $result->next_hit()) {
>    print ref($hit), "\n";
>    print Dumper($hit), "\n";
> }
> $ ./main.pl
> ref($factory): Bio::Tools::Run::StandAloneBlastPlus
> ref($result): Bio::Search::Result::BlastResult
> num_hits: 1
> Bio::Search::Hit::BlastHit
> $VAR1 = bless( {
>                   '_hsps' =>  undef,
>                   '_iterator' =>  0,
>                   '_description' =>  '',
>                   '_query_length' =>  '75',
>                   '_accession' =>  'second1',
>                   '_length' =>  '75',
>                   '_name' =>  'second1',
>                   '_rank' =>  1,
>                   '_algorithm' =>  'BLASTN',
>                   '_root_verbose' =>  0,
>                   '_hsp_factory' =>  bless( {
>                                              'interface' =>
> 'Bio::Search::HSP::HSPI',
>                                              'type' =>
> 'Bio::Search::HSP::GenericHSP',
>                                              '_loaded_types' =>  {
>
> 'Bio::Search::HSP::GenericHSP' =>  1
>                                                                 },
>                                              '_root_verbose' =>  0
>                                            }, 'Bio::Factory::ObjectFactory' )
>                 }, 'Bio::Search::Hit::BlastHit' );
>
>
>
>