Bioperl: Flipped out on Bio::Tools::Blast

Steve A. Chervitz sac@alberich.Stanford.EDU
Tue, 25 Aug 1998 12:32:29 -0700 (PDT)


Lincoln,

This was an important oversight which has already been corrected in 
version 0.06. There is now an HSP::strand() method for retrieval of strand 
information. You can call the strand() method on HSP objects as in:

  print $hsp->strand('query');
  print $hsp->strand('sbjct');

These will return 'Minus' or 'Plus'.
Or, to get strand information for both query and sbjct (in that order):

  @strands = $hsp->strand();

As you noted, the _set_seq() method will flip start <-> end if 
start > end. This was done to facilitate HSP processing. Commenting this 
out will have consequences only for TBLAST[NX] reports and only for methods 
that rely on HSP::range() or HSP::matches(). The main effect is that the HSP 
tiling code in Sbjct.pm could break for these reports, leading to 
potentially incorrect values reported by methods of Sbjct.pm that rely on 
HSP::_tile_hsp() such as length_aln(), gaps(), frac_identical(), etc. which 
collect data across all HSPs. Tables produced by Blast::table_tiled() would 
therefore be affected as well. (Since there is now a strand() method, 
you shouldn't need to comment out any code :).

Thanks for reporting the bug in Sbjct.pm regarding sequence length 
handling. I'll include it in the next release.

Steve Chervitz
sac@genome.stanford.edu
http://genome-www.stanford.edu/perlOOP/bioperl/blast/


On Tue, 25 Aug 1998, Lincoln Stein wrote:

> When I retrieve the list of HSPs for a hit, how do I tell whether
> I'm dealing with a plus or a minus strand HSP?  The range() operator
> always sorts the endpoints of the subject and query sequences from
> lower to higher!
> ...
>
> Here's the offending code in Steve's HSP.pm.  As far as I can tell, it
> silently flips the sequence endpoints without setting any flag that
> this has happened.
> 
>     # Check for flipped sequence endpoints (from reverse complement)
>     if($self->{$seqType.'Start'} > $self->{$seqType.'Stop'}) {
>         my $tmp = $self->{$seqType.'Start'};
>         $self->{$seqType.'Start'} = $self->{$seqType.'Stop'};
>         $self->{$seqType.'Stop'} = $tmp;
>     }
> 
> If I comment out this section, what will break?
> 
> Lincoln
> -- 
> ========================================================================
> Lincoln D. Stein                           Cold Spring Harbor Laboratory
> lstein@cshl.org			                  Cold Spring Harbor, NY
> ========================================================================
> =========== Bioperl Project Mailing List Message Footer =======
> Project URL: http://www.techfak.uni-bielefeld.de/bcd/Perl/Bio/
> For info about how to (un)subscribe, where messages are archived, etc:
> http://www.techfak.uni-bielefeld.de/bcd/Perl/Bio/vsns-bcd-perl.html
> ====================================================================
> 

=========== Bioperl Project Mailing List Message Footer =======
Project URL: http://www.techfak.uni-bielefeld.de/bcd/Perl/Bio/
For info about how to (un)subscribe, where messages are archived, etc:
http://www.techfak.uni-bielefeld.de/bcd/Perl/Bio/vsns-bcd-perl.html
====================================================================