Skip to content

HTTPS clone URL

Subversion checkout URL

You can clone with HTTPS or Subversion.

Download ZIP
Browse files

Bio/SeqFeature/Generic.pm: Fix for issue #76, adjusted "seq"

subroutine to correctly handle circular features cut by origin,
which were making it crash, and added tests
  • Loading branch information...
commit 409a7a5b20ddce0bfd295e78e68167b18bf1c40d 1 parent f3e8778
@fjossandon fjossandon authored
Showing with 44 additions and 7 deletions.
  1. +21 −6 Bio/SeqFeature/Generic.pm
  2. +23 −1 t/SeqIO/Splicedseq.t
View
27 Bio/SeqFeature/Generic.pm
@@ -703,17 +703,32 @@ sub seq {
# assumming our seq object is sensible, it should not have to yank
# the entire sequence out here.
+ my $seq;
+ my $start = $self->start;
+ my $end = $self->end;
+ # Check circular sequences cut by origin (e.g. join(2006035..2007700,1..257))
+ if ( $self->{'_gsf_seq'}->is_circular
+ and $self->location->isa('Bio::Location::SplitLocationI')
+ and $start > $end
+ ) {
+ my $primary_seq_length = $self->{'_gsf_seq'}->length;
- my $seq = $self->{'_gsf_seq'}->trunc($self->start(), $self->end());
+ # Get duplicate object with the first sequence piece using trunc()
+ $seq = $self->{'_gsf_seq'}->trunc($start, $primary_seq_length);
+ # Get post-origin sequence and build the complete sequence
+ my $post_origin = $self->{'_gsf_seq'}->subseq(1, $end);
+ my $complete_seq = $seq->seq() . $post_origin;
- if ( defined $self->strand &&
- $self->strand == -1 ) {
+ # Add complete sequence to object
+ $seq->seq($complete_seq);
+ }
+ else {
+ $seq = $self->{'_gsf_seq'}->trunc($start, $end);
+ }
- # ok. this does not work well (?)
- #print STDERR "Before revcom", $seq->str, "\n";
+ if ( defined $self->strand && $self->strand == -1 ) {
$seq = $seq->revcom;
- #print STDERR "After revcom", $seq->str, "\n";
}
return $seq;
View
24 t/SeqIO/Splicedseq.t
@@ -7,8 +7,9 @@ BEGIN {
use lib '.';
use Bio::Root::Test;
- test_begin(-tests => 19);
+ test_begin(-tests => 25);
+ use_ok('Bio::Seq');
use_ok('Bio::SeqIO');
}
@@ -67,6 +68,27 @@ warnings_like { $len_nodb = length($feats[1]->spliced_seq()->seq); }
"appropriate warning if db not provided for remote sequence";
ok($len_nodb == 374, "correct number of Ns added if remote sequence not provided");
+# Test for cut by origin features
+my $seq_obj = Bio::Seq->new(-display_id => 'NC_008309',
+ -seq => 'AAAAACCCCCGGGGGTTTTT');
+$seq_obj->is_circular(1);
+my $loc_obj = Bio::Factory::FTLocationFactory->from_string('join(16..20,1..2)');
+my $cut_feat = Bio::SeqFeature::Generic->new(-primary_tag => 'CDS',
+ -location => $loc_obj,
+ -tag => { locus_tag => 'HS_1792',
+ product => 'hypothetical protein',
+ protein_id => 'YP_718205.1',
+ } );
+$seq_obj->add_SeqFeature($cut_feat);
+is $cut_feat->seq->seq, 'TTTTTAA', 'cut by origin sequence';
+is $cut_feat->start, 16, 'cut by origin start using $feat->start';
+is $cut_feat->end, 2, 'cut by origin end using $feat->end';
+TODO: {
+ local $TODO = 'wrong coords for $feat->location->start/end';
+ is $cut_feat->location->start, 16, 'cut by origin start using $feat->location->start';
+ is $cut_feat->location->end, 2, 'cut by origin end using $feat->location->end';
+}
+
SKIP: {
test_skip(-tests => 3, -requires_module => 'LWP::UserAgent', -requires_networking => 1);
my $db_in;
Please sign in to comment.
Something went wrong with that request. Please try again.