[Bioperl-l] Bio::SeqFeature::Generic unexpected behavior

Jonathan Flowers jonathanmflowers at gmail.com
Wed Oct 7 22:13:46 UTC 2009


Hi,

I am building Bio::SeqFeature::Gene::Transcript objects to facilitate easy
parsing of cds and intron sequences.  I ran across an unexpected problem
when creating two or more Transcript objects and attaching independent
Bio::Seq objects to them. In the following code, 2 transcript objects are
constructed independently (not cloned).  However, extracting the cds for one
transcript gives an errant (and unexpected result).  Curiously, extracting
the introns from the same transcript objects gives the correct result.  I
have highlighted the unexpected outcome below.  Finally, I noticed that
reversing the order of the attach_seq commands determines which cds is
printed.

The following example illustrates the problem.

Can someone explain this?

Jonathan

#! usr/bin/env perl
use strict;
use warnings;
use Bio::FeatureIO::gff;
use Bio::SeqFeature::Annotated;
use Bio::SeqFeature::Gene::Transcript;
use Bio::SeqFeature::Gene::Exon;
#use Bio::SeqFeature::Gene::UTR;
#use Bio::SeqFeature::Gene::Intron;

#Build the objects...
my $seq1 = Bio::Seq->new(-seq => 'AAAAA');
my $seq2 = Bio::Seq->new(-seq => 'TTTTT');

my $tscript1 = Bio::SeqFeature::Gene::Transcript->new(-start => 1, -end =>
5, -strand => +1);
my $tscript2 = Bio::SeqFeature::Gene::Transcript->new(-start => 1, -end =>
5, -strand => +1);

my $exon1 = Bio::SeqFeature::Gene::Exon->new(-start => 1 , -end => 2,
-strand => +1);
my $exon2 = Bio::SeqFeature::Gene::Exon->new(-start => 4 , -end => 5,
-strand => +1);

$tscript1->add_exon($exon1);
$tscript1->add_exon($exon2);

$tscript2->add_exon($exon1);
$tscript2->add_exon($exon2);

$tscript1->attach_seq($seq1);
$tscript2->attach_seq($seq2);

#now extract the desired features and print
print $tscript1->cds->seq,"\n"; #prints 'TTTT' <- unexpected !!!
print $tscript2->cds->seq,"\n"; #prints 'TTTT' <- as expected

foreach($tscript1->introns){
    print $_->seq->seq,"\n";     #prints 'A'   <- as expected
}

foreach($tscript2->introns){
    print $_->seq->seq,"\n";    #prints 'T'  <- as expected
}



More information about the Bioperl-l mailing list