[Bioperl-l] BioPerl newbie question about exporting sub_SeqFe atures

[Jason Stajich]

okay I fixed this - it was a bug in that Bio::RangeI::union($feature) was
returning $self->new(-start => $start, -end => $end, -strand => $strand)
rather than ( $start, $end, $strand) as what used to be the case I assume.
Bio::SeqFeature::Generic now understands this and your code should work.

BTW: If you are trying to represent a single feature that has multiple
locations try the Bio::Location::Split class to set your location
points.  This will output the following in the GenBank Feature table
if that is the behaviour you want.

testfeature 	join(2..4, 6..8)

Otherwise your code below should work if you sync off of CVS (have to
wait an hour or so for the CVS code to propigate to anonymous.  The
diff of what I changed is included below if you don't want to wait.


Index: Bio/SeqFeature/Generic.pm
===================================================================
RCS file: /home/repository/bioperl/bioperl-live/Bio/SeqFeature/Generic.pm,v
retrieving revision 1.46
diff -r1.46 Generic.pm
818d817
<
828,831c827,830
<       my ($start, $end, $strand) = $self->union($feat);
<       $self->start($start);
<       $self->end($end);
<       $self->strand($strand);
---
>       my $range = $self->union($feat);
>       $self->start($range->start);
>       $self->end($range->end);
>       $self->strand($range->strand);

-Jason

On Tue, 17 Apr 2001, Gosink, Mark (SEA) wrote:

> Thanks for the suggestion, unfortunately it didn't seem to work for more
> than 1 sub_feature. Here is the script with multiple sub_features:
>
> use Bio::SeqFeature::Generic;
> use Bio::Seq;
> use Bio::SeqIO;
>
> my $feat = Bio::SeqFeature::Generic->new(-primary => "test");
> my $seq = Bio::Seq->new(-seq => "agctacgaccttgatcgcta", -id =>"test01");
>
> # add sub-features to the feature
> $feat->add_sub_SeqFeature(Bio::SeqFeature::Generic->new(-start => 2, -end =>
> 4, -primary => "sub1"), EXPAND);
> $feat->add_sub_SeqFeature(Bio::SeqFeature::Generic->new(-start => 6, -end =>
> 8, -primary => "sub2"), EXPAND);
> # could add more
>
> # attach to sequence
> $seq->add_SeqFeature($feat);
>
> # print in GenBank format (note that there is lots of stuff missing, like
> species etc)
> my $seqout = Bio::SeqIO->new(-fh => \*STDOUT, -format => "genbank");
>
> $seqout->write_seq($seq);
>
> $seqout->close();
>
> And here's the results:
>
> LOCUS       test01            20 bp     DNA  UNK
> DEFINITION
> ACCESSION   unknown
> FEATURES             Location/Qualifiers
>      sub1            2..4
>      sub2            6..8
>      test            Bio::SeqFeature::Generic=HASH(0xd58d0)..4
> BASE COUNT        5 a      6 c      4 g      5 t
> ORIGIN
>         1 agctacgacc ttgatcgcta
> //
>
> > -----Original Message-----
> > From: Hilmar Lapp [mailto:hilmarl@yahoo.com]
> > Sent: Tuesday, April 17, 2001 11:06 AM
> > To: Gosink, Mark (SEA)
> > Cc: 'bioperl-l@bioperl.org'
> > Subject: Re: [Bioperl-l] BioPerl newbie question about exporting
> > sub_SeqFeatures
> >
> >
> > "Gosink, Mark (SEA)" wrote:
> > >
> > > Hello all,
> > >         Does anyone have any example of adding a feature
> > with sub_features
> > > to a sequence and then saving that as a GenBank/EMBL/other?
> > When I try, I
> > > get a HASH feature.
> > > Thanks for any tips,
> >
> > use Bio::SeqFeature::Generic;
> > use Bio::Seq;
> > use Bio::SeqIO;
> >
> > my $feat = Bio::SeqFeature::Generic->new(-start => 1, -end => 5,
> > -primary => "test");
> > my $seq = Bio::Seq->new(-seq => "agctacgaccttgatcgcta", -id =>
> > "test01");
> >
> > # add sub-features to the feature
> > $feat->add_SeqFeature(Bio::SeqFeature::Generic->new(-start => 2,
> > -end => 4, -primary => "sub1"));
> > # could add more
> >
> > # attach to sequence
> > $seq->add_SeqFeature($feat);
> >
> > # print in GenBank format (note that there is lots of stuff
> > missing, like species etc)
> > my $seqout = Bio::SeqIO->new(-fh => \*STDOUT, -format =>
> > "genbank");
> > $seqout->write_seq($seq);
> >
> > $seqout->close();
> >
> >
> > This works for me. The two features come out as top-level features
> > though, because Genbank does not have subfeatures of features (in
> > Genbank features can only have annotation in the form of tag/value
> > pairs).
> >
> > Does this answer your question? If not, please provide your code
> > and the output it produces on which input.
> >
> > 	Hilmar
> > --
> > -----------------------------------------------------------------
> > Hilmar Lapp                              email: hilmarl@yahoo.com
> > GNF, San Diego, Ca. 92122                phone: +1 858 812 1757
> > -----------------------------------------------------------------
> >
>
>
> The information contained in this email is intended for the
> personal and confidential use of the addressee only. It may
> also be privileged information. If you are not the intended
> recipient then you are hereby notified that you have received
> this document in error and that any review, distribution or
> copying of this document is strictly prohibited. If you have
> received  this communication in error, please notify Celltech
> Group immediately on:
>
> +44 (0)1753 534655, or email 'is@celltech.co.uk'
>
> Celltech Group plc
> 216 Bath Road, Slough, SL1 4EN, Berkshire, UK
>
> Registered Office as above. Registered in England No. 2159282
> _______________________________________________
> Bioperl-l mailing list
> Bioperl-l@bioperl.org
> http://bioperl.org/mailman/listinfo/bioperl-l
>

Jason Stajich
jason@chg.mc.duke.edu
Center for Human Genetics
Duke University Medical Center
http://www.chg.duke.edu/