19
19
# get a seqfeature somehow, eg, from a Sequence with Features attached
21
21
foreach $feat ( $seq->get_SeqFeatures() ) {
22
print "Feature from ", $feat->start, "to ",
23
$feat->end, " Primary tag ", $feat->primary_tag,
22
print "Feature from ", $feat->start, "to ",
23
$feat->end, " Primary tag ", $feat->primary_tag,
24
24
", produced by ", $feat->source_tag(), "\n";
26
if( $feat->strand == 0 ) {
27
print "Feature applicable to either strand\n";
29
print "Feature on strand ", $feat->strand,"\n"; # -1,1
31
print "feature location is ",$feat->start, "..",
32
$feat->end, " on strand ", $feat->strand, "\n";
33
print "easy utility to print locations in GenBank/EMBL way ",
34
$feat->location->to_FTstring(), "\n";
36
foreach $tag ( $feat->get_all_tags() ) {
37
print "Feature has tag ", $tag, " with values, ",
26
if( $feat->strand == 0 ) {
27
print "Feature applicable to either strand\n";
29
print "Feature on strand ", $feat->strand,"\n"; # -1,1
32
print "feature location is ",$feat->start, "..",
33
$feat->end, " on strand ", $feat->strand, "\n";
34
print "easy utility to print locations in GenBank/EMBL way ",
35
$feat->location->to_FTstring(), "\n";
37
foreach $tag ( $feat->get_all_tags() ) {
38
print "Feature has tag ", $tag, " with values, ",
38
39
join(' ',$feat->get_tag_values($tag)), "\n";
40
41
print "new feature\n" if $feat->has_tag('new');
41
42
# features can have sub features
42
43
my @subfeat = $feat->get_SeqFeatures();
57
58
Bioperl modules. Send your comments and suggestions preferably to one
58
59
of the Bioperl mailing lists. Your participation is much appreciated.
60
bioperl-l@bioperl.org - General discussion
61
http://bio.perl.org/MailList.html - About the mailing lists
61
bioperl-l@bioperl.org - General discussion
62
http://bioperl.org/wiki/Mailing_lists - About the mailing lists
63
64
=head2 Reporting Bugs
65
66
Report bugs to the Bioperl bug tracking system to help us keep track
66
the bugs and their resolution. Bug reports can be submitted via email
67
the bugs and their resolution. Bug reports can be submitted via the
69
bioperl-bugs@bio.perl.org
70
http://bugzilla.bioperl.org/
70
http://bugzilla.open-bio.org/
170
169
$self->throw_not_implemented();
176
Usage : $tag_exists = $self->has_tag('some_tag')
178
Returns : TRUE if the specified tag exists, and FALSE otherwise
185
my ($self,@args) = @_;
187
$self->throw_not_implemented();
191
=head2 get_tag_values
193
Title : get_tag_values
194
Usage : @values = $self->get_tag_values('some_tag')
196
Returns : An array comprising the values of the specified tag.
199
throws an exception if there is no such tag
204
shift->throw_not_implemented();
207
=head2 get_tagset_values
209
Title : get_tagset_values
210
Usage : @values = $self->get_tagset_values(qw(label transcript_id product))
212
Returns : An array comprising the values of the specified tags, in order of tags
213
Args : An array of strings
215
does NOT throw an exception if none of the tags are not present
217
this method is useful for getting a human-readable label for a
218
SeqFeatureI; not all tags can be assumed to be present, so a list of
219
possible tags in preferential order is provided
223
# interface + abstract method
224
sub get_tagset_values {
225
my ($self, @args) = @_;
227
foreach my $arg (@args) {
228
if ($self->has_tag($arg)) {
229
push(@vals, $self->get_tag_values($arg));
240
Usage : @tags = $feat->get_all_tags()
241
Function: gives all tags for this feature
242
Returns : an array of strings
249
shift->throw_not_implemented();
252
173
=head2 attach_seq
260
181
Note that it is not guaranteed that if you obtain a feature from
261
182
an object in bioperl, it will have a sequence attached. Also,
262
183
implementors of this interface can choose to provide an empty
263
implementation of this method. I.e., there is also no guarantee
184
implementation of this method. I.e., there is also no guarantee
264
185
that if you do attach a sequence, seq() or entire_seq() will not
267
188
The reason that this method is here on the interface is to enable
268
189
you to call it on every SeqFeatureI compliant object, and
269
that it will be implemented in a useful way and set to a useful
190
that it will be implemented in a useful way and set to a useful
270
191
value for the great majority of use cases. Implementors who choose
271
192
to ignore the call are encouraged to specifically state this in
272
193
their documentation.
453
340
number of N's (DNA) or X's (protein, though this is unlikely).
455
342
This function is deliberately "magical" attempting to second guess
456
what a user wants as "the" sequence for this feature
343
what a user wants as "the" sequence for this feature.
458
345
Implementing classes are free to override this method with their
459
own magic if they have a better idea what the user wants
346
own magic if they have a better idea what the user wants.
461
Args : [optional] A Bio::DB::RandomAccessI compliant object
349
-db A L<Bio::DB::RandomAccessI> compliant object if
350
one needs to retrieve remote seqs.
351
-nosort boolean if the locations should not be sorted
352
by start location. This may occur, for instance,
353
in a circular sequence where a gene span starts
354
before the end of the sequence and ends after the
355
sequence start. Example : join(15685..16260,1..207)
356
Returns : A L<Bio::PrimarySeqI> object
466
360
sub spliced_seq {
363
my ($db,$nosort) = $self->_rearrange([qw(DB NOSORT)], @args);
365
# (added 7/7/06 to allow use old API (with warnings)
366
my $old_api = (!(grep {$_ =~ /(?:nosort|db)/} @args)) ? 1 : 0;
367
if (@args && $old_api) {
368
$self->warn(q(API has changed; please use '-db' or '-nosort' ).
369
qq(for args. See POD for more details.));
370
$db = shift @args if @args;
371
$nosort = shift @args if @args;
374
if( $db && ref($db) && ! $db->isa('Bio::DB::RandomAccessI') ) {
375
$self->warn("Must pass in a valid Bio::DB::RandomAccessI object".
376
" for access to remote locations for spliced_seq");
378
} elsif( defined $db && $HasInMemory &&
379
$db->isa('Bio::DB::InMemoryCache') ) {
380
$db = new Bio::DB::InMemoryCache(-seqdb => $db);
469
383
if( ! $self->location->isa("Bio::Location::SplitLocationI") ) {
470
384
return $self->seq(); # nice and easy!
475
389
$self->throw("not atomic, not split, yikes, in trouble!");
479
393
my $seqid = $self->entire_seq->display_id;
480
394
# This is to deal with reverse strand features
481
395
# so we are really sorting features 5' -> 3' on their strand
482
396
# i.e. rev strand features will be sorted largest to smallest
483
397
# as this how revcom CDSes seem to be annotated in genbank.
484
# Might need to eventually allow this to be programable?
398
# Might need to eventually allow this to be programable?
485
399
# (can I mention how much fun this is NOT! --jason)
487
401
my ($mixed,$mixedloc, $fstrand) = (0);
488
if( $db && ref($db) && ! $db->isa('Bio::DB::RandomAccessI') ) {
489
$self->warn("Must pass in a valid Bio::DB::RandomAccessI object for access to remote locations for spliced_seq");
491
} elsif( defined $db && $HasInMemory &&
492
! $db->isa('Bio::DB::InMemoryCache') ) {
493
$db = new Bio::DB::InMemoryCache(-seqdb => $db);
495
403
if( $self->isa('Bio::Das::SegmentI') &&
496
! $self->absolute ) {
404
! $self->absolute ) {
497
405
$self->warn("Calling spliced_seq with a Bio::Das::SegmentI which does have absolute set to 1 -- be warned you may not be getting things on the correct strand");
500
my @locs = map { $_->[0] }
501
# sort so that most negative is first basically to order
502
# the features on the opposite strand 5'->3' on their strand
503
# rather than they way most are input which is on the fwd strand
505
sort { $a->[1] <=> $b->[1] } # Yes Tim, Schwartzian transformation
507
$fstrand = $_->strand unless defined $fstrand;
508
$mixed = 1 if defined $_->strand && $fstrand != $_->strand;
509
if( defined $_->seq_id ) {
510
$mixedloc = 1 if( $_->seq_id ne $seqid );
408
my @locset = $self->location->each_Location;
411
@locs = map { $_->[0] }
412
# sort so that most negative is first basically to order
413
# the features on the opposite strand 5'->3' on their strand
414
# rather than they way most are input which is on the fwd strand
416
sort { $a->[1] <=> $b->[1] } # Yes Tim, Schwartzian transformation
418
$fstrand = $_->strand unless defined $fstrand;
419
$mixed = 1 if defined $_->strand && $fstrand != $_->strand;
420
if( defined $_->seq_id ) {
421
$mixedloc = 1 if( $_->seq_id ne $seqid );
423
[ $_, $_->start * ($_->strand || 1)];
427
$self->warn("Mixed strand locations, spliced seq using the input order rather than trying to sort");
512
[ $_, $_->start* ($_->strand || 1)];
513
} $self->location->each_Location;
516
$self->warn("Mixed strand locations, spliced seq using the input order rather than trying to sort");
517
@locs = $self->location->each_Location;
518
} elsif( $mixedloc ) {
519
# we'll use the prescribed location order
520
@locs = $self->location->each_Location;
431
# use the original order instead of trying to sort
433
$fstrand = $locs[0]->strand;
523
foreach my $loc ( @locs ) {
436
foreach my $loc ( @locs ) {
524
437
if( ! $loc->isa("Bio::Location::Atomic") ) {
525
438
$self->throw("Can only deal with one level deep locations");
554
467
$called_seq = $self->entire_seq;
470
# does the called sequence make sense? Bug 1780
471
if ($called_seq->length < $loc->end) {
472
my $accession = $called_seq->accession;
474
my $length = $called_seq->length;
475
my $orig_id = $self->seq_id; # originating sequence
476
my ($locus) = $self->get_tagset_values("locus_tag");
477
$self->throw("Location end ($end) exceeds length ($length) of ".
478
"called sequence $accession.\nCheck sequence version used in ".
479
"$locus locus-tagged SeqFeature in $orig_id.");
557
482
if( $self->isa('Bio::Das::SegmentI') ) {
558
my ($s,$e) = ($loc->start,$loc->end);
483
my ($s,$e) = ($loc->start,$loc->end);
559
484
$seqstr .= $called_seq->subseq($s,$e)->seq();
561
# This is dumb subseq should work on locations...
486
# This is dumb, subseq should work on locations...
562
487
if( $loc->strand == 1 ) {
563
488
$seqstr .= $called_seq->subseq($loc->start,$loc->end);
565
$seqstr .= $called_seq->trunc($loc->start,$loc->end)->revcom->seq();
491
$seqstr = $called_seq->trunc($loc->start,$loc->end)->revcom->seq() . $seqstr;
493
$seqstr .= $called_seq->trunc($loc->start,$loc->end)->revcom->seq();
569
my $out = Bio::Seq->new( -id => $self->entire_seq->display_id . "_spliced_feat",
498
my $out = Bio::Seq->new( -id => $self->entire_seq->display_id
570
500
-seq => $seqstr);
575
=head1 RangeI methods
577
These methods are inherited from RangeI and can be used
578
directly from a SeqFeatureI interface. Remember that a
579
SeqFeature is-a RangeI, and so wherever you see RangeI you
580
can use a feature ($r in the below documentation).
585
Usage : if($feat->overlaps($r)) { do stuff }
586
if($feat->overlaps(200)) { do stuff }
587
Function: tests if $feat overlaps $r
588
Args : a RangeI to test for overlap with, or a point
589
Returns : true if the Range overlaps with the feature, false otherwise
595
Usage : if($feat->contains($r) { do stuff }
596
Function: tests whether $feat totally contains $r
597
Args : a RangeI to test for being contained
598
Returns : true if the argument is totaly contained within this range
604
Usage : if($feat->equals($r))
605
Function: test whether $feat has the same start, end, strand as $r
606
Args : a RangeI to test for equality
607
Returns : true if they are describing the same range
610
=head1 Geometrical methods
612
These methods do things to the geometry of ranges, and return
613
triplets (start, stop, strand) from which new ranges could be built.
618
Usage : ($start, $stop, $strand) = $feat->intersection($r)
619
Function: gives the range that is contained by both ranges
620
Args : a RangeI to compare this one to
621
Returns : nothing if they do not overlap, or the range that they do overlap
626
Usage : ($start, $stop, $strand) = $feat->union($r);
627
: ($start, $stop, $strand) = Bio::RangeI->union(@ranges);
628
Function: finds the minimal range that contains all of the ranges
629
Args : a range or list of ranges to find the union of
630
Returns : the range containing all of the ranges
637
508
Usage : my $location = $seqfeature->location()
638
Function: returns a location object suitable for identifying location
639
of feature on sequence or parent feature
509
Function: returns a location object suitable for identifying location
510
of feature on sequence or parent feature
640
511
Returns : Bio::LocationI object
655
526
Title : primary_id
656
527
Usage : $obj->primary_id($newval)
659
530
Returns : value of primary_id (a scalar)
660
531
Args : on set, new value (a scalar or undef, optional)
533
Primary ID is a synonym for the tag 'ID'
666
538
my $self = shift;
667
return $self->{'primary_id'} = shift if @_;
668
return $self->{'primary_id'};
539
# note from cjm@fruitfly.org:
540
# I have commented out the following 2 lines:
542
#return $self->{'primary_id'} = shift if @_;
543
#return $self->{'primary_id'};
545
#... and replaced it with the following; see
546
# http://bioperl.org/pipermail/bioperl-l/2003-December/014150.html
547
# for the discussion that lead to this change
550
if ($self->has_tag('ID')) {
551
$self->remove_tag('ID');
553
$self->add_tag_value('ID', shift);
555
my ($id) = $self->get_tagset_values('ID');
559
sub generate_unique_persistent_id {
560
# DEPRECATED - us IDHandler
562
require "Bio/SeqFeature/Tools/IDHandler.pm";
563
Bio::SeqFeature::Tools::IDHandler->new->generate_unique_persistent_id($self);
566
=head1 Bio::RangeI methods
568
These methods are inherited from RangeI and can be used
569
directly from a SeqFeatureI interface. Remember that a
570
SeqFeature is-a RangeI, and so wherever you see RangeI you
571
can use a feature ($r in the below documentation).
599
=head2 intersection()
607
=head1 Bio::AnnotatableI methods
613
B<Deprecated>. See L<Bio::AnnotatableI>
617
B<Deprecated>. See L<Bio::AnnotatableI>
619
=head2 add_tag_value()
621
B<Deprecated>. See L<Bio::AnnotatableI>
623
=head2 get_tag_values()
625
B<Deprecated>. See L<Bio::AnnotatableI>
627
=head2 get_tagset_values()
629
B<Deprecated>. See L<Bio::AnnotatableI>
631
=head2 get_all_tags()
633
B<Deprecated>. See L<Bio::AnnotatableI>