File Coverage

Bio/SeqFeature/Gene/Transcript.pm
Criterion Covered Total %
statement 164 222 73.8
branch 48 86 55.8
condition 25 44 56.8
subroutine 24 32 75.0
pod 20 23 86.9
total 281 407 69.0


line stmt bran cond sub pod time code
1             #
2             # BioPerl module for Bio::SeqFeature::Gene::Transcript
3             #
4             # Please direct questions and support issues to
5             #
6             # Cared for by Hilmar Lapp
7             #
8             # Copyright Hilmar Lapp
9             #
10             # You may distribute this module under the same terms as perl itself
11              
12             # POD documentation - main docs before the code
13              
14             =head1 NAME
15              
16             Bio::SeqFeature::Gene::Transcript - A feature representing a transcript
17              
18             =head1 SYNOPSIS
19              
20             # See documentation of methods.
21              
22             =head1 DESCRIPTION
23              
24             A feature representing a transcript.
25              
26             =head1 FEEDBACK
27              
28             =head2 Mailing Lists
29              
30             User feedback is an integral part of the evolution of this and other
31             Bioperl modules. Send your comments and suggestions preferably to one
32             of the Bioperl mailing lists. Your participation is much appreciated.
33              
34             bioperl-l@bioperl.org - General discussion
35             http://bioperl.org/wiki/Mailing_lists - About the mailing lists
36              
37             =head2 Support
38              
39             Please direct usage questions or support issues to the mailing list:
40              
41             I
42              
43             rather than to the module maintainer directly. Many experienced and
44             reponsive experts will be able look at the problem and quickly
45             address it. Please include a thorough description of the problem
46             with code and data examples if at all possible.
47              
48             =head2 Reporting Bugs
49              
50             Report bugs to the Bioperl bug tracking system to help us keep track
51             the bugs and their resolution. Bug reports can be submitted via the
52             web:
53              
54             https://github.com/bioperl/bioperl-live/issues
55              
56             =head1 AUTHOR - Hilmar Lapp
57              
58             Email hlapp@gmx.net
59              
60             =head1 APPENDIX
61              
62             The rest of the documentation details each of the object methods. Internal methods are usually preceded with a _
63              
64             =cut
65              
66              
67             # Let the code begin...
68              
69             package Bio::SeqFeature::Gene::Transcript;
70 9     9   927 use strict;
  9         11  
  9         218  
71              
72              
73 9     9   625 use Bio::PrimarySeq;
  9         9  
  9         175  
74              
75 9     9   28 use base qw(Bio::SeqFeature::Generic Bio::SeqFeature::Gene::TranscriptI);
  9         11  
  9         2900  
76              
77             sub new {
78 96     96 1 199 my ($caller, @args) = @_;
79 96         242 my $self = $caller->SUPER::new(@args);
80 96         243 $self->_register_for_cleanup(\&transcript_destroy);
81 96         240 my ($primary) = $self->_rearrange([qw(PRIMARY)],@args);
82              
83 96 100       197 $primary = 'transcript' unless $primary;
84 96         193 $self->primary_tag($primary);
85 96 50       165 $self->strand(0) if(! defined($self->strand()));
86 96         211 return $self;
87             }
88              
89              
90             =head2 promoters
91              
92             Title : promoters()
93             Usage : @proms = $transcript->promoters();
94             Function: Get the promoter features/sites of this transcript.
95              
96             Note that OO-modeling of regulatory elements is not stable yet.
97             This means that this method might change or even disappear in a
98             future release. Be aware of this if you use it.
99              
100             Returns : An array of Bio::SeqFeatureI implementing objects representing the
101             promoter regions or sites.
102             Args :
103              
104              
105             =cut
106              
107             sub promoters {
108 9     9 1 11 my ($self) = @_;
109 9         14 return $self->get_feature_type('Bio::SeqFeature::Gene::Promoter');
110             }
111              
112             =head2 add_promoter
113              
114             Title : add_promoter()
115             Usage : $transcript->add_promoter($feature);
116             Function: Add a promoter feature/site to this transcript.
117              
118              
119             Note that OO-modeling of regulatory elements is not stable yet.
120             This means that this method might change or even disappear in a
121             future release. Be aware of this if you use it.
122              
123             Returns :
124             Args : A Bio::SeqFeatureI implementing object.
125              
126              
127             =cut
128              
129             sub add_promoter {
130 9     9 1 13 my ($self, $fea) = @_;
131 9         19 $self->_add($fea,'Bio::SeqFeature::Gene::Promoter');
132             }
133              
134             =head2 flush_promoters
135              
136             Title : flush_promoters()
137             Usage : $transcript->flush_promoters();
138             Function: Remove all promoter features/sites from this transcript.
139              
140             Note that OO-modeling of regulatory elements is not stable yet.
141             This means that this method might change or even disappear in a
142             future release. Be aware of this if you use it.
143              
144             Returns : the removed features as a list
145             Args : none
146              
147              
148             =cut
149              
150             sub flush_promoters {
151 0     0 1 0 my ($self) = @_;
152 0         0 return $self->_flush('Bio::SeqFeature::Gene::Promoter');
153             }
154              
155             =head2 exons
156              
157             Title : exons()
158             Usage : @exons = $gene->exons();
159             ($inital_exon) = $gene->exons('Initial');
160             Function: Get all exon features or all exons of specified type of this
161             transcript.
162              
163             Exon type is treated as a case-insensitive regular expression and
164             is optional. For consistency, use only the following types:
165             initial, internal, terminal.
166              
167             Returns : An array of Bio::SeqFeature::Gene::ExonI implementing objects.
168             Args : An optional string specifying the primary_tag of the feature.
169              
170              
171             =cut
172              
173             sub exons {
174 111     111 1 630 my ($self, $type) = @_;
175 111         203 return $self->get_unordered_feature_type('Bio::SeqFeature::Gene::ExonI',
176             $type);
177             }
178              
179             =head2 exons_ordered
180              
181             Title : exons_ordered
182             Usage : @exons = $gene->exons_ordered();
183             @exons = $gene->exons_ordered("Internal");
184             Function: Get an ordered list of all exon features or all exons of specified
185             type of this transcript.
186              
187             Exon type is treated as a case-insensitive regular expression and
188             is optional. For consistency, use only the following types:
189              
190             Returns : An array of Bio::SeqFeature::Gene::ExonI implementing objects.
191             Args : An optional string specifying the primary_tag of the feature.
192              
193             =cut
194              
195             sub exons_ordered {
196 4     4 1 5 my ($self,$type) = @_;
197 4         12 return $self->get_feature_type('Bio::SeqFeature::Gene::ExonI', $type);
198             }
199              
200             =head2 add_exon
201              
202             Title : add_exon()
203             Usage : $transcript->add_exon($exon,'initial');
204             Function: Add a exon feature to this transcript.
205              
206             The second argument denotes the type of exon. Mixing exons with and
207             without a type is likely to cause trouble in exons(). Either
208             leave out the type for all exons or for none.
209              
210             Presently, the following types are known: initial, internal,
211             terminal, utr, utr5prime, and utr3prime (all case-insensitive).
212             UTR should better be added through utrs()/add_utr().
213              
214             If you wish to use other or additional types, you will almost
215             certainly have to call exon_type_sortorder() in order to replace
216             the default sort order, or mrna(), cds(), protein(), and exons()
217             may yield unexpected results.
218              
219             Returns :
220             Args : A Bio::SeqFeature::Gene::ExonI implementing object.
221             A string indicating the type of the exon (optional).
222              
223              
224             =cut
225              
226             sub add_exon {
227 509     509 1 533 my ($self, $fea, $type) = @_;
228 509 50       1255 if(! $fea->isa('Bio::SeqFeature::Gene::ExonI') ) {
229 0         0 $self->throw("$fea does not implement Bio::SeqFeature::Gene::ExonI");
230             }
231 509         871 $self->_add($fea,'Bio::SeqFeature::Gene::Exon', $type);
232             }
233              
234             =head2 flush_exons
235              
236             Title : flush_exons()
237             Usage : $transcript->flush_exons();
238             $transcript->flush_exons('terminal');
239             Function: Remove all or a certain type of exon features from this transcript.
240              
241             See add_exon() for documentation about types.
242              
243             Calling without a type will not flush UTRs. Call flush_utrs() for
244             this purpose.
245             Returns : the deleted features as a list
246             Args : A string indicating the type of the exon (optional).
247              
248              
249             =cut
250              
251             sub flush_exons {
252 0     0 1 0 my ($self, $type) = @_;
253 0         0 return $self->_flush('Bio::SeqFeature::Gene::Exon',$type);
254             }
255              
256             =head2 introns
257              
258             Title : introns()
259             Usage : @introns = $gene->introns();
260             Function: Get all intron features this gene structure.
261              
262             Note that this implementation generates these features
263             on-the-fly, that is, it simply treats all regions between
264             exons as introns, assuming that exons do not overlap. A
265             consequence is that a consistent correspondence between the
266             elements in the returned array and the array that exons()
267             returns will exist only if the exons are properly sorted
268             within their types (forward for plus- strand and reverse
269             for minus-strand transcripts). To ensure correctness the
270             elements in the array returned will always be sorted.
271              
272             Returns : An array of Bio::SeqFeature::Gene::Intron objects representing
273             the intron regions.
274             Args :
275              
276              
277             =cut
278              
279             sub introns {
280 1     1 1 298 my ($self) = @_;
281 1         2 my @introns = ();
282 1         3 my @exons = $self->exons();
283 1         2 my ($strand, $rev_order);
284              
285             # if there's 1 or less exons we're done
286 1 50       4 return () unless($#exons > 0);
287             # record strand and order (a minus-strand transcript is likely to have
288             # the exons stacked in reverse order)
289 1         2 foreach my $exon (@exons) {
290 1         4 $strand = $exon->strand();
291 1 50       3 last if $strand; # we're done if we've got 1 or -1
292             }
293 1 50       3 $rev_order = ($exons[0]->end() < $exons[1]->start() ? 0 : 1);
294              
295             # Make sure exons are sorted. Because we assume they don't overlap, we
296             # simply sort by start position.
297 1 50 33     5 if((! defined($strand)) || ($strand != -1) || (! $rev_order)) {
      33        
298             # always sort forward for plus-strand transcripts, and for negative-
299             # strand transcripts that appear to be unsorted or forward sorted
300 4         6 @exons = map { $_->[0] } sort { $a->[1] <=> $b->[1] }
  4         6  
301 1   50     2 map { [ $_, $_->start * ($_->strand || 1)] } @exons;
  4         5  
302             } else {
303             # sort in reverse order for transcripts on the negative strand and
304             # found to be in reverse order
305 0         0 @exons = map { $_->[0] } sort { $b->[1] <=> $a->[1] } map { [ $_, $_->start()] } @exons;
  0         0  
  0         0  
  0         0  
306             }
307             # loop over all intervening gaps
308 1   66     7 while ((my $exonA = shift (@exons)) &&(my $exonB = shift(@exons))){
309 3         18 my $intron = Bio::SeqFeature::Gene::Intron->new(-primary=>'intron');
310 3         6 $intron->upstream_Exon($exonA);
311 3         5 $intron->downstream_Exon($exonB);
312 3 50       5 $intron->attach_seq($self->entire_seq) if $self->entire_seq;
313 3         4 unshift(@exons,$exonB);
314 3         12 push @introns,$intron;
315             }
316 1         3 return @introns;
317             }
318              
319             =head2 poly_A_site
320              
321             Title : poly_A_site()
322             Usage : $polyAsite = $transcript->poly_A_site();
323             Function: Get/set the poly-adenylation feature/site of this transcript.
324             Returns : A Bio::SeqFeatureI implementing object representing the
325             poly-adenylation region.
326             Args : A Bio::SeqFeatureI implementing object on set, or FALSE to flush
327             a previously set object.
328              
329              
330             =cut
331              
332             sub poly_A_site {
333 31     31 1 29 my ($self, $fea) = @_;
334 31 100       54 if ($fea) {
335 11         30 $self->_add($fea,'Bio::SeqFeature::Gene::Poly_A_site');
336             }
337 31         56 return ($self->get_feature_type('Bio::SeqFeature::Gene::Poly_A_site'))[0];
338             }
339              
340             =head2 utrs
341              
342             Title : utrs()
343             Usage : @utr_sites = $transcript->utrs('utr3prime');
344             @utr_sites = $transcript->utrs('utr5prime');
345             @utr_sites = $transcript->utrs();
346             Function: Get the features representing untranslated regions (UTR) of this
347             transcript.
348              
349             You may provide an argument specifying the type of UTR. Currently
350             the following types are recognized: utr5prime utr3prime for UTR on the
351             5' and 3' end of the CDS, respectively.
352              
353             Returns : An array of Bio::SeqFeature::Gene::UTR objects
354             representing the UTR regions or sites.
355             Args : Optionally, either utr3prime, or utr5prime for the the type of UTR
356             feature.
357              
358              
359             =cut
360              
361             sub utrs {
362 3     3 1 4 my ($self, $type) = @_;
363 3         6 return $self->get_feature_type('Bio::SeqFeature::Gene::UTR',$type);
364              
365             }
366              
367             =head2 add_utr
368              
369             Title : add_utr()
370             Usage : $transcript->add_utr($utrobj, 'utr3prime');
371             $transcript->add_utr($utrobj);
372             Function: Add a UTR feature/site to this transcript.
373              
374             The second parameter is optional and denotes the type of the UTR
375             feature. Presently recognized types include 'utr5prime' and 'utr3prime'
376             for UTR on the 5' and 3' end of a gene, respectively.
377              
378             Calling this method is the same as calling
379             add_exon($utrobj, 'utr'.$type). In this sense a UTR object is a
380             special exon object, which is transcribed, not spliced out, but
381             not translated.
382              
383             Note that the object supplied should return FALSE for is_coding().
384             Otherwise cds() and friends will become confused.
385              
386             Returns :
387             Args : A Bio::SeqFeature::Gene::UTR implementing object.
388              
389              
390             =cut
391              
392             sub add_utr {
393 2     2 1 3 my ($self, $fea, $type) = @_;
394 2         5 $self->_add($fea,'Bio::SeqFeature::Gene::UTR',$type);
395             }
396              
397             =head2 flush_utrs
398              
399             Title : flush_utrs()
400             Usage : $transcript->flush_utrs();
401             $transcript->flush_utrs('utr3prime');
402             Function: Remove all or a specific type of UTR features/sites from this
403             transcript.
404              
405             Cf. add_utr() for documentation about recognized types.
406             Returns : a list of the removed features
407             Args : Optionally a string denoting the type of UTR feature.
408              
409              
410             =cut
411              
412             sub flush_utrs {
413 0     0 1 0 my ($self, $type) = @_;
414 0         0 return $self->_flush('Bio::SeqFeature::Gene::UTR',$type);
415             }
416              
417             =head2 sub_SeqFeature
418              
419             Title : sub_SeqFeature
420             Usage : @feats = $transcript->sub_SeqFeature();
421             Function: Returns an array of all subfeatures.
422              
423             This method is defined in Bio::SeqFeatureI. We override this here
424             to include the exon etc features.
425              
426             Returns : An array Bio::SeqFeatureI implementing objects.
427             Args : none
428              
429              
430             =cut
431              
432             sub sub_SeqFeature {
433 9     9 1 12 my ($self) = @_;
434 9         9 my @feas;
435            
436             # get what the parent already has
437 9         41 @feas = $self->SUPER::sub_SeqFeature();
438             # add the features we have in addition
439 9         25 push(@feas, $self->exons()); # this includes UTR features
440 9         30 push(@feas, $self->promoters());
441 9 100       23 push(@feas, $self->poly_A_site()) if($self->poly_A_site());
442 9         28 return @feas;
443             }
444              
445             =head2 flush_sub_SeqFeature
446              
447             Title : flush_sub_SeqFeature
448             Usage : $transcript->flush_sub_SeqFeature();
449             $transcript->flush_sub_SeqFeature(1);
450             Function: Removes all subfeatures.
451              
452             This method is overridden from Bio::SeqFeature::Generic to flush
453             all additional subfeatures like exons, promoters, etc., which is
454             almost certainly not what you want. To remove only features added
455             through $transcript->add_sub_SeqFeature($feature) pass any
456             argument evaluating to TRUE.
457              
458             Example :
459             Returns : none
460             Args : Optionally, an argument evaluating to TRUE will suppress flushing
461             of all transcript-specific subfeatures (exons etc.).
462              
463              
464             =cut
465              
466             sub flush_sub_SeqFeature {
467 0     0 1 0 my ($self,$fea_only) = @_;
468              
469 0         0 $self->SUPER::flush_sub_SeqFeature();
470 0 0       0 if(! $fea_only) {
471 0         0 $self->flush_promoters();
472 0         0 $self->flush_exons();
473 0         0 $self->flush_utrs();
474 0         0 $self->poly_A_site(0);
475             }
476             }
477              
478              
479             =head2 cds
480              
481             Title : cds
482             Usage : $seq = $transcript->cds();
483             Function: Returns the CDS (coding sequence) as defined by the exons
484             of this transcript and the attached sequence.
485              
486             If no sequence is attached this method will return false.
487              
488             Note that the implementation provided here returns a
489             concatenation of all coding exons, thereby assuming that
490             exons do not overlap.
491              
492             Note also that you cannot set the CDS via this method. Set
493             a single CDS feature as a single exon, or derive your own
494             class if you want to store a predicted CDS.
495              
496             Example :
497             Returns : A Bio::PrimarySeqI implementing object.
498             Args :
499              
500             =cut
501              
502             sub cds {
503 4     4 1 19 my ($self) = @_;
504 4         53 my @exons = $self->exons_ordered(); #this is always sorted properly according to strand
505 4         10 my $strand;
506              
507 4 50       16 return unless(@exons);
508             # record strand (a minus-strand transcript must have the exons sorted in
509             # reverse order)
510 4         12 foreach my $exon (@exons) {
511 40 100 100     53 if(defined($exon->strand()) && (! $strand)) {
512 5         9 $strand = $exon->strand();
513             }
514 40 50 66     51 if($exon->strand() && (($exon->strand() * $strand) < 0)) {
515 0         0 $self->throw("Transcript mixes coding exons on plus and minus ".
516             "strand. This makes no sense.");
517             }
518             }
519 4         17 my $cds = $self->_make_cds(@exons);
520 4 50       11 return unless $cds;
521 4         18 return Bio::PrimarySeq->new('-id' => $self->seq_id(),
522             '-seq' => $cds,
523             '-alphabet' => "dna");
524             }
525              
526             =head2 protein
527              
528             Title : protein()
529             Usage : $protein = $transcript->protein();
530             Function: Get the protein encoded by the transcript as a sequence object.
531              
532             The implementation provided here simply calls translate() on the
533             object returned by cds().
534              
535             Returns : A Bio::PrimarySeqI implementing object.
536             Args :
537              
538              
539             =cut
540              
541             sub protein {
542 0     0 1 0 my ($self) = @_;
543 0         0 my $seq;
544              
545 0         0 $seq = $self->cds();
546 0 0       0 return $seq->translate() if $seq;
547 0         0 return;
548             }
549              
550             =head2 mrna
551              
552             Title : mrna()
553             Usage : $mrna = $transcript->mrna();
554             Function: Get the mRNA of the transcript as a sequence object.
555              
556             The difference to cds() is that the sequence object returned by
557             this methods will also include UTR and the poly-adenylation site,
558             but not promoter sequence (TBD).
559              
560             HL: do we really need this method?
561              
562             Returns : A Bio::PrimarySeqI implementing object.
563             Args :
564              
565              
566             =cut
567              
568             sub mrna {
569 1     1 1 336 my ($self) = @_;
570 1         2 my ($seq, $mrna, $elem);
571              
572             # get the coding part
573 1         3 $seq = $self->cds();
574 1 50       4 if(! $seq) {
575 0         0 $seq = Bio::PrimarySeq->new('-id' => $self->seq_id(),
576             '-alphabet' => "rna",
577             '-seq' => "");
578             }
579             # get and add UTR sequences
580 1         1 $mrna = "";
581 1         3 foreach $elem ($self->utrs('utr5prime')) {
582 1         3 $mrna .= $elem->seq()->seq();
583             }
584 1         3 $seq->seq($mrna . $seq->seq());
585 1         2 $mrna = "";
586 1         2 foreach $elem ($self->utrs('utr3prime')) {
587 1         3 $mrna .= $elem->seq()->seq();
588             }
589 1         3 $seq->seq($seq->seq() . $mrna);
590 1 50       2 if($self->poly_A_site()) {
591 1         4 $seq->seq($seq->seq() . $self->poly_A_site()->seq()->seq());
592             }
593 1 50       3 return if($seq->length() == 0);
594 1         3 return $seq;
595             }
596              
597             sub _get_typed_keys {
598 0     0   0 my ($self, $keyprefix, $type) = @_;
599 0         0 my @keys = ();
600 0         0 my @feas = ();
601              
602             # make case-insensitive
603 0 0       0 $type = ($type ? lc($type) : "");
604             # pull out all feature types that exist and match
605 0         0 @keys = grep { /^_$keyprefix$type/i; } (keys(%{$self}));
  0         0  
  0         0  
606 0         0 return @keys;
607             }
608              
609             sub _make_cds {
610 4     4   8 my ($self,@exons) = @_;
611 4         5 my $cds = "";
612              
613 4         7 foreach my $exon (@exons) {
614 40 100 33     85 next if((! defined($exon->seq())) || (! $exon->is_coding()));
615 38         57 my $phase = length($cds) % 3;
616             # let's check the simple case
617 38 100 100     76 if((! defined($exon->frame())) || ($phase == $exon->frame())) {
618             # this one fits exactly, or frame of the exon is undefined (should
619             # we warn about that?); we bypass the $exon->cds() here (hmm,
620             # not very clean style, but I don't see where this screws up)
621 37         59 $cds .= $exon->seq()->seq();
622             } else {
623             # this one is probably from exon shuffling and needs some work
624 1         4 my $seq = $exon->cds(); # now $seq is guaranteed to be in frame 0
625 1 50       3 next if(! $seq);
626 1         3 $seq = $seq->seq();
627             # adjustment needed?
628 1 50       3 if($phase > 0) {
629             # how many Ns can we chop off the piece to be added?
630 0         0 my $n_crop = 0;
631 0 0       0 if($seq =~ /^(n+)/i) {
632 0         0 $n_crop = length($1);
633             }
634 0 0       0 if($n_crop >= $phase) {
635             # chop off to match the phase
636 0         0 $seq = substr($seq, $phase);
637             } else {
638             # fill in Ns
639 0         0 $seq = ("n" x (3-$phase)) . $seq;
640             }
641             }
642 1         3 $cds .= $seq;
643             }
644             }
645 4         20 return $cds;
646             }
647              
648             =head2 features
649              
650             Title : features
651             Usage : my @features=$transcript->features;
652             Function: returns all the features associated with this transcript
653             Returns : a list of SeqFeatureI implementing objects
654             Args : none
655              
656              
657             =cut
658              
659              
660             sub features {
661 254     254 1 229 my $self = shift;
662 254 50       221 return grep { defined } @{$self->{'_features'} || []};
  1674         1767  
  254         597  
663             }
664              
665             =head2 features_ordered
666              
667             Title : features_ordered
668             Usage : my @features=$transcript->features_ordered;
669             Function: returns all the features associated with this transcript,
670             in order by feature start, according to strand
671             Returns : a list of SeqFeatureI implementing objects
672             Args : none
673              
674              
675             =cut
676              
677             sub features_ordered{
678 0     0 1 0 my ($self) = @_;
679 0 0       0 return $self->_stranded_sort(@{$self->{'_features'} || []});
  0         0  
680             }
681              
682              
683             sub get_unordered_feature_type{
684 158     158 0 157 my ($self, $type, $pri)=@_;
685 158         121 my @list;
686 158         234 foreach ( $self->features) {
687 1143 100       2221 if ($_->isa($type)) {
688 773 100 100     1024 if ($pri && $_->primary_tag !~ /$pri/i) {
689 35         30 next;
690             }
691 738         632 push @list,$_;
692             }
693             }
694 158         422 return @list;
695              
696             }
697              
698             sub get_feature_type {
699 47     47 0 63 my ($self)=shift;
700 47         73 return $self->_stranded_sort($self->get_unordered_feature_type(@_));
701             }
702              
703             #This was fixed by Gene Cutler - the indexing on the list being reversed
704             #fixed a bad bug. Thanks Gene!
705             sub _flush {
706 0     0   0 my ($self, $type, $pri)=@_;
707 0         0 my @list=$self->features;
708 0         0 my @cut;
709 0         0 for (reverse (0..$#list)) {
710 0 0 0     0 if (defined $list[$_] &&
711             $list[$_]->isa($type)) {
712 0 0 0     0 if ($pri && $list[$_]->primary_tag !~ /$pri/i) {
713 0         0 next;
714             }
715 0         0 push @cut, splice @list, $_, 1; #remove the element of $type from @list
716             #and return each of them in @cut
717             }
718             }
719 0         0 $self->{'_features'}=\@list;
720 0         0 return reverse @cut;
721             }
722              
723             sub _add {
724 531     531   510 my ($self, $fea, $type, $pri)=@_;
725 531         5210 require Bio::SeqFeature::Gene::Promoter;
726 531         3121 require Bio::SeqFeature::Gene::UTR;
727 531         1101 require Bio::SeqFeature::Gene::Exon;
728 531         3176 require Bio::SeqFeature::Gene::Intron;
729 531         3340 require Bio::SeqFeature::Gene::Poly_A_site;
730              
731 531 50       1229 if(! $fea->isa('Bio::SeqFeatureI') ) {
732 0         0 $self->throw("$fea does not implement Bio::SeqFeatureI");
733             }
734 531 100 100     2089 if(! $fea->isa($type) || $pri) {
735 382         554 $fea=$self->_new_of_type($fea,$type,$pri);
736             }
737 531 100       847 if (! $self->strand) {
738 87         148 $self->strand($fea->strand);
739             } else {
740 444 50       666 if ($self->strand * $fea->strand == -1) {
741 0         0 $self->throw("$fea is on opposite strand from $self");
742             }
743             }
744              
745 531         1232 $self->_expand_region($fea);
746 531 0 33     1049 if(defined($self->entire_seq()) && (! defined($fea->entire_seq())) &&
      33        
747             $fea->can('attach_seq')) {
748 0         0 $fea->attach_seq($self->entire_seq());
749             }
750 531 50       1170 if (defined $self->parent) {
751 0         0 $self->parent->_expand_region($fea);
752             }
753 531         394 push(@{$self->{'_features'}}, $fea);
  531         934  
754 531         1026 1;
755             }
756              
757             sub _stranded_sort {
758 47     47   54 my ($self,@list)=@_;
759 47         37 my $strand;
760 47         49 foreach my $fea (@list) {
761 78 100       113 if($fea->strand()) {
762             # defined and != 0
763 64 100       99 $strand = $fea->strand() if(! $strand);
764 64 50       83 if(($fea->strand() * $strand) < 0) {
765 0         0 $strand = undef;
766 0         0 last;
767             }
768             }
769             }
770 47 100 100     127 if (defined $strand && $strand == - 1) { #reverse strand
771 20         27 return map { $_->[0] } sort {$b->[1] <=> $a->[1]} map { [$_, $_->start] } @list;
  52         101  
  32         27  
  52         78  
772             } else { #undef or forward strand
773 27         38 return map { $_->[0] } sort {$a->[1] <=> $b->[1]} map { [$_, $_->start] } @list;
  26         68  
  7         14  
  26         42  
774             }
775             }
776              
777             sub _new_of_type {
778 382     382   374 my ($self, $fea, $type, $pri)= @_;
779 382         272 my $primary;
780 382 100       464 if ($pri) {
781 363         323 $primary = $pri; #can set new primary tag if desired
782             } else {
783 19         88 ($primary) = $type =~ /.*::(.+)/; #or else primary is just end of type string
784             }
785 382         406 bless $fea,$type;
786 382         656 $fea->primary_tag($primary);
787 382         486 return $fea;
788             }
789              
790             sub transcript_destroy {
791 96     96 0 83 my $self = shift;
792             # We're going to be really explicit to insure memory leaks
793             # don't occur
794 96         151 foreach my $f ( $self->features ) {
795 531         362 $f = undef;
796             }
797 96         210 $self->parent(undef);
798             }
799              
800             1;