[Bioperl-guts-l] [14984] bioperl-live/trunk: fixed significant alignment parsing bug affecting hsps, and also got rid of warnings when length not known

[Senduran Balasubramaniam]

Revision: 14984
Author:   sendu
Date:     2008-11-11 13:39:20 -0500 (Tue, 11 Nov 2008)

Log Message:
-----------
fixed significant alignment parsing bug affecting hsps, and also got rid of warnings when length not known

Modified Paths:
--------------
    bioperl-live/trunk/Bio/Search/HSP/HmmpfamHSP.pm
    bioperl-live/trunk/Bio/Search/Result/HmmpfamResult.pm
    bioperl-live/trunk/t/hmmer_pull.t

Modified: bioperl-live/trunk/Bio/Search/HSP/HmmpfamHSP.pm
===================================================================
--- bioperl-live/trunk/Bio/Search/HSP/HmmpfamHSP.pm	2008-11-11 18:37:00 UTC (rev 14983)
+++ bioperl-live/trunk/Bio/Search/HSP/HmmpfamHSP.pm	2008-11-11 18:39:20 UTC (rev 14984)
@@ -137,13 +137,14 @@
     my $alignments_hash = $self->get_field('alignments');
 	
     my $identifier = $self->get_field('name').'~~~~'.$self->get_field('rank');
+	
     while (! defined $alignments_hash->{$identifier}) {
-        last unless $self->parent->parent->_next_alignment;
+		last unless $self->parent->parent->_next_alignment;
     }
     my $alignment = $alignments_hash->{$identifier};
-    
+	
     if ($alignment) {
-        # work out query, hit and homology strings, and some stats
+		# work out query, hit and homology strings, and some stats
         # (quicker to do this all at once instead of each method working on
         # $alignment string itself)
         
@@ -283,6 +284,13 @@
 sub hit {
     my $self = shift;
     unless ($self->{_created_hit}) {
+		# the full length isn't always known (given in the report), but don't
+		# warn about the missing info all the time
+		my $verbose = $self->parent->parent->parent->verbose;
+		$self->parent->parent->parent->verbose(-1);
+		my $seq_length = $self->get_field('length');
+		$self->parent->parent->parent->verbose($verbose);
+		
         $self->SUPER::hit( new  Bio::SeqFeature::Similarity
                   ('-primary'  => $self->primary_tag,
                    '-start'    => $self->get_field('hit_start'),
@@ -291,7 +299,7 @@
                    '-score'    => $self->get_field('score'),
                    '-strand'   => 1,
                    '-seq_id'   => $self->get_field('name'),
-                   '-seqlength'=> $self->get_field('length'),
+                   $seq_length ? ('-seqlength' => $seq_length) : (),
                    '-source'   => $self->get_field('algorithm'),
                    '-seqdesc'  => $self->get_field('description')
                    ) );

Modified: bioperl-live/trunk/Bio/Search/Result/HmmpfamResult.pm
===================================================================
--- bioperl-live/trunk/Bio/Search/Result/HmmpfamResult.pm	2008-11-11 18:37:00 UTC (rev 14983)
+++ bioperl-live/trunk/Bio/Search/Result/HmmpfamResult.pm	2008-11-11 18:39:20 UTC (rev 14984)
@@ -207,8 +207,7 @@
 		
 		$self->{_after_previous_alignment} = $self->_chunk_tell;
 		$self->{_next_alignment_start_text} = $chunk;
-		$self->_next_alignment;
-		return;
+		return $self->_next_alignment;
 	}
 	
 	$self->_chunk_seek($self->{_after_previous_alignment});
@@ -228,6 +227,7 @@
 	if (defined $self->{_next_alignment_start_text}) {
 		$chunk = $self->{_next_alignment_start_text}.$chunk;
 	}
+	
 	$chunk =~ s/(\S+: domain)$//;
 	$self->{_next_alignment_start_text} = $1;
 	

Modified: bioperl-live/trunk/t/hmmer_pull.t
===================================================================
--- bioperl-live/trunk/t/hmmer_pull.t	2008-11-11 18:37:00 UTC (rev 14983)
+++ bioperl-live/trunk/t/hmmer_pull.t	2008-11-11 18:39:20 UTC (rev 14984)
@@ -7,7 +7,7 @@
     use lib 't/lib';
     use BioperlTest;
     
-    test_begin(-tests => 287);
+    test_begin(-tests => 290);
 	
 	use_ok('Bio::SearchIO');
 }
@@ -194,3 +194,25 @@
 }
 
 is $searchio->result_count, 2;
+
+# bug revealed by bug 2632 - CS lines were already ignored, but we couldn't
+# parse alignments when HSPs weren't in simple order!!
+$searchio = Bio::SearchIO->new(-format => 'hmmer_pull', -file => test_input_file('hmmpfam_cs.out'), -verbose => 1);
+my $result = $searchio->next_result;
+my $hit = $result->next_hit;
+my $hsp = $hit->next_hsp;
+is $hsp->seq_str, "IPPLLAVGAVHHHLINKGLRQEASILV";
+
+# and another bug revealed: we don't always know the hit length, and
+# shouldn't complain about that with a warning
+is $hsp->hit->seqlength, 412;
+
+my $count = 0;
+while (my $hit = $result->next_hit) {
+    $count++;
+    next if $count < 6;
+    last if $count > 6;
+	my $hsp = $hit->next_hsp;
+    ok ! $hsp->hit->seqlength;
+    #*** not sure how to test for the lack of a warning though...
+}