comparison variant_effect_predictor/Bio/EnsEMBL/IndividualSlice.pm @ 0:1f6dce3d34e0

Uploaded
author mahtabm
date Thu, 11 Apr 2013 02:01:53 -0400
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:1f6dce3d34e0
1 =head1 LICENSE
2
3 Copyright (c) 1999-2012 The European Bioinformatics Institute and
4 Genome Research Limited. All rights reserved.
5
6 This software is distributed under a modified Apache license.
7 For license details, please see
8
9 http://www.ensembl.org/info/about/code_licence.html
10
11 =head1 CONTACT
12
13 Please email comments or questions to the public Ensembl
14 developers list at <dev@ensembl.org>.
15
16 Questions may also be sent to the Ensembl help desk at
17 <helpdesk@ensembl.org>.
18
19 =cut
20
21 =head1 NAME
22
23 Bio::EnsEMBL::IndividualSlice - SubClass of the Slice. Represents the
24 slice of the genome for a certain individual (applying the alleles for
25 this individual)
26
27 =head1 SYNOPSIS
28
29 $sa = $db->get_SliceAdaptor;
30
31 $slice =
32 $sa->fetch_by_region( 'chromosome', 'X', 1_000_000, 2_000_000 );
33
34 $individualSlice = $slice->get_by_Individual($individual_name);
35
36 # Get the sequence from the Individual Slice: will contain IUPAC codes
37 # for SNPs and Ensembl ambiguity codes for indels
38 my $seq = $individualSlice->seq();
39 print $seq;
40
41 # Get a subSlice of the Strain
42 my $subSlice_individual =
43 $individualSlice->sub_Slice( 5_000, 8_000, 1 )
44
45 # Compare two different individuals in the same Slice
46 my $sliceIndividual2 = $slice->get_by_Individual($individual_name2);
47 my $differences =
48 $individualSlice->get_all_differences_IndividualSlice(
49 $sliceIndividual2);
50
51 foreach my $af ( @{$differences} ) {
52 print
53 "There is a difference between $individual_name "
54 . "and $individual_name2 at ",
55 $af->start, "-", $af->end,
56 " with allele ", $af->allele_string(), "\n";
57 }
58
59 =head1 DESCRIPTION
60
61 A IndividualSlice object represents a region of a genome for a certain
62 individual. It can be used to retrieve sequence or features from a
63 individual.
64
65 =head1 METHODS
66
67 =cut
68
69 package Bio::EnsEMBL::IndividualSlice;
70 use vars qw(@ISA);
71 use strict;
72
73 use Bio::EnsEMBL::Utils::Argument qw(rearrange);
74 use Bio::EnsEMBL::Utils::Sequence qw(reverse_comp);
75 use Bio::EnsEMBL::Slice;
76 use Bio::EnsEMBL::Mapper;
77 use Bio::EnsEMBL::Utils::Exception qw(throw deprecate warning);
78
79 @ISA = qw(Bio::EnsEMBL::Slice);
80
81
82 =head2 new
83
84 Arg [1..N] : List of named arguments
85 Bio::EnsEMBL::CoordSystem COORD_SYSTEM
86 string SEQ_REGION_NAME,
87 int START,
88 int END,
89 string VERSION (optional, defaults to '')
90 int STRAND, (optional, defaults to 1)
91 Bio::EnsEMBL::DBSQL::SliceAdaptor ADAPTOR (optional)
92 Arg[N+1] : string $individual_name
93 Example : $individualSlice = Bio::EnsEMBL::IndividualSlice->new(-coord_system => $cs,
94 -start => 1,
95 -end => 10000,
96 -strand => 1,
97 -seq_region_name => 'X',
98 -seq_region_length => 12e6,
99 -individual_name => $individual_name);
100 Description : Creates a new Bio::EnsEMBL::IndividualSlice object that will contain a shallow copy of the
101 Slice object, plus additional information such as the individual this Slice refers to
102 and listref of Bio::EnsEMBL::Variation::AlleleFeatures of differences with the
103 reference sequence
104 ReturnType : Bio::EnsEMBL::IndividualSlice
105 Exceptions : none
106 Caller : general
107
108 =cut
109
110 sub new{
111 my $caller = shift;
112 my $class = ref($caller) || $caller;
113
114 #create the IndividualSlice object as the Slice, plus the individual attribute
115 my ($individual_name, $sample_id) = rearrange(['INDIVIDUAL', 'SAMPLE_ID'],@_);
116
117 my $self = $class->SUPER::new(@_);
118
119 $self->{'individual_name'} = $individual_name;
120 $self->{'sample_id'} = $sample_id;
121
122 return $self;
123
124 }
125
126 =head2 individual_name
127
128 Arg [1] : (optional) string $individual_name
129 Example : my $individual_name = $individualSlice->individual_name();
130 Description : Getter/Setter for the name of the individual in the slice
131 ReturnType : string
132 Exceptions : none
133 Caller : general
134
135 =cut
136
137 sub individual_name{
138 my $self = shift;
139 if (@_){
140 $self->{'individual_name'} = shift @_;
141 }
142 return $self->{'individual_name'};
143 }
144
145 =head2 seq
146
147 Arg [1] : none
148 Example : print "SEQUENCE = ", $strainSlice->seq();
149 Description: Returns the sequence of the region represented by this
150 StrainSlice formatted as a string.
151 Returntype : string
152 Exceptions : none
153 Caller : general
154
155 =cut
156
157 sub seq {
158 my $self = shift;
159
160 # special case for in-between (insert) coordinates
161 return '' if($self->start() == $self->end() + 1);
162
163 return $self->{'seq'} if($self->{'seq'});
164
165 if($self->adaptor()) {
166 my $seqAdaptor = $self->adaptor()->db()->get_SequenceAdaptor();
167 my $reference_sequence = $seqAdaptor->fetch_by_Slice_start_end_strand($self,1,undef,1); #get the reference sequence for that slice
168 #apply all differences to the reference sequence
169
170 # sort edits in reverse order to remove complication of
171 # adjusting downstream edits
172 my @allele_features_ordered = sort {$b->start() <=> $a->start() || $b->end() <=> $a->end()} @{$self->{'alleleFeatures'}} if (defined $self->{'alleleFeatures'});
173
174 foreach my $af (@allele_features_ordered){
175 $af->apply_edit($reference_sequence); #change, in the reference sequence, the af
176 }
177 # return substr(${$reference_sequence},0,1) if ($self->length == 1);
178 return ${$reference_sequence}; #returns the reference sequence, applying the alleleFeatures
179 }
180
181 # no attached sequence, and no db, so just return Ns
182 return 'N' x $self->length();
183 }
184
185 =head2 get_all_differences_Slice
186
187 Args : none
188 Example : my $differences = $individualSlice->get_all_differences_Slice()
189 Description : Gets all differences between the IndividualSlice object and the Slice is defined
190 ReturnType : listref of Bio::EnsEMBL::Variation::AlleleFeature
191 Exceptions : none
192 Caller : general
193
194 =cut
195
196 sub get_all_differences_Slice{
197 my $self = shift;
198 my $differences; #reference to the array with the differences between Slice and StrainSlice
199 my $ref_allele;
200 foreach my $difference (@{$self->{'alleleFeatures'}}){
201 if ($difference->length_diff == 0){
202 #the difference is a SNP, check if it is the same as the reference allele
203 $ref_allele = $self->SUPER::subseq($difference->start,$difference->end,$difference->strand);
204 $ref_allele = '-' if ($ref_allele eq '');
205 if ($ref_allele ne $difference->allele_string){
206 #when the alleleFeature is different from the reference allele, add to the differences list
207 push @{$differences},$difference;
208 }
209 }
210 else{
211 push @{$differences},$difference;
212 }
213 }
214
215 return $differences;
216
217 }
218
219 =head2 get_all_differences_IndividualSlice
220
221 Arg[1] : Bio::EnsEMBL::IndividualSlice $is
222 Example : my $differences = $individualSlice->get_all_differences_IndividualSlice($individualslice)
223 Description : Gets differences between 2 IndividualSlice objects
224 ReturnType : listref of Bio::EnsEMBL::Variation::AlleleFeature
225 Exceptions : thrown on bad argument
226 Caller : general
227
228 =cut
229
230 sub get_all_differences_IndividualSlice{
231 my $self = shift;
232 my $individualSlice = shift;
233
234 if (!ref($individualSlice) || !$individualSlice->isa('Bio::EnsEMBL::IndividualSlice')){
235 throw('Bio::EnsEMBL::IndividualSlice arg expected');
236 }
237 if ( @{$self->{'alleleFeatures'}} == 0 && @{$individualSlice->{'alleleFeatures'}} == 0){
238 return undef; #there are no differences in any of the Individuals
239
240 }
241 my $differences; #differences between individuals
242 if (@{$individualSlice->{'alleleFeatures'}} == 0){
243 #need to create a copy of alleleFeature for the first Individual
244 foreach my $difference (@{$self->{'alleleFeatures'}}){
245 my %vf = %$difference;
246 push @{$differences},bless \%vf,ref($difference);
247 }
248 }
249 elsif (@{$self->{'alleleFeatures'}} == 0){
250 #need to create a copy of AlleleFeature, but changing the allele by the allele in the reference sequence
251 foreach my $difference (@{$individualSlice->{'alleleFeatures'}}){
252 push @{$differences}, $individualSlice->_convert_difference($difference);
253 }
254 }
255 else{
256 #both individuals have differences
257 #create a hash with the differences in the first slice
258 my %allele_features_self = map {$_->start.'-'.$_->end => $_} @{$self->{'alleleFeatures'}};
259 foreach my $difference (@{$individualSlice->{'alleleFeatures'}}){
260 #there is no difference in the other individual slice, convert the allele
261 if (!defined $allele_features_self{$difference->start.'-'.$difference->end}){
262 push @{$differences},$individualSlice->_convert_difference($difference);
263 }
264 else{
265 #if it is defined and have the same allele, delete from the hash since it is not a difference
266 #between the individuals
267 if ($allele_features_self{$difference->start.'-'.$difference->end}->allele_string eq $difference->allele_string){
268 delete $allele_features_self{$difference->start.'-'.$difference->end};
269 }
270 }
271 }
272 #and finally, make a shallow copy of the differences in the first individual
273 foreach my $difference (values %allele_features_self){
274 my %vf = %$difference;
275 push @{$differences},bless \%vf,ref($difference);
276 }
277
278 }
279 #need to map differences to the first individual, self, since the coordinates are in the Slice coordinate system
280 my $mapper = $self->mapper(); #now that we have the differences, map them in the IndividualSlice
281 my @results;
282 foreach my $difference (@{$differences}){
283 @results = $mapper->map_coordinates('Slice',$difference->start,$difference->end,$difference->strand,'Slice');
284 #we can have 3 possibilities:
285 #the difference is an insertion and when mapping returns the boundaries of the insertion in the IndividualSlice
286 if (@results == 2){
287 #the first position in the result is the beginning of the insertion
288 if($results[0]->start < $results[1]->start){
289 $difference->start($results[0]->end+1);
290 $difference->end($results[1]->start-1);
291 }
292 else{
293 #it is the second position the beginning of the insertion
294 $difference->start($results[1]->end+1);
295 $difference->end($results[0]->start-1);
296 }
297 $difference->strand($results[0]->strand);
298 }
299 else{
300 #it can be either a SNP or a deletion, and we have the coordinates in the result, etither a Bio::EnsEMBL::Mapper::Coordinate
301 # or a Bio::EnsEMBL::Mapper::IndelCoordinate
302 $difference->start($results[0]->start);
303 $difference->end($results[0]->end);
304 $difference->strand($results[0]->strand);
305 }
306 }
307
308 return $differences;
309 }
310
311 #for a given AlleleFeature, converts the allele into the reference allele and returns
312 #the converted AlleleFeature
313
314 sub _convert_difference{
315 my $self = shift;
316 my $difference = shift;
317 my %new_af = %$difference; #make a copy of the alleleFeature
318 #and change the allele with the one from the reference Slice
319 $new_af{'allele_string'} = $self->SUPER::subseq($difference->start,$difference->end,$difference->strand);
320 return bless \%new_af,ref($difference);
321 }
322
323 =head2 mapper
324
325 Args : none
326 Description: Getter for the mapper between the between the IndividualSlice and the Slice it refers to.
327 It is done automatically when necessary to create subSlice or to get the differences between individuals
328 Returntype : Bio::EnsEMBL::Mapper
329 Exceptions : none
330 Caller : Internal function
331
332 =cut
333
334 sub mapper{
335 my $self = shift;
336
337 if (@_) {
338 #allow to create again the mapper
339 delete $self->{'mapper'};
340 }
341 if(!defined $self->{'mapper'}){
342 #create the mapper between the Slice and StrainSlice
343 my $mapper = Bio::EnsEMBL::Mapper->new('Slice','IndividualSlice');
344 #align with Slice
345 #get all the VariationFeatures in the Individual Slice, from start to end in the Slice
346 my @allele_features_ordered = sort {$a->start() <=> $b->start() || $b->end() <=> $a->end()} @{$self->{'alleleFeatures'}} if (defined $self->{'alleleFeatures'});
347
348 my $start_slice = 1;
349 my $end_slice;
350 my $start_individual = 1;
351 my $end_individual;
352 my $length_allele;
353 my $total_length_diff = 0;
354 #we will walk from left to right in the slice object, updating the start and end individual every time
355 #there is a new alleleFeature in the Individual
356 foreach my $allele_feature (@allele_features_ordered){
357 #we have a insertion/deletion: marks the beginning of new slice move coordinates
358 if ($allele_feature->length_diff != 0){
359 $total_length_diff += $allele_feature->length_diff;
360 $length_allele = $allele_feature->length + $allele_feature->length_diff(); #length of the allele in the Individual
361 $end_slice = $allele_feature->start() - 1; #set the end of the slice before the alleleFeature
362 if ($end_slice >= $start_slice){
363 #normal cases (not with gaps)
364 $end_individual = $end_slice - $start_slice + $start_individual; #set the end of the individual from the beginning plus the offset
365 #add the sequence that maps
366 $mapper->add_map_coordinates('Slice',$start_slice,$end_slice,1,'IndividualSlice',$start_individual,$end_individual);
367 #and add the indel
368 $mapper->add_indel_coordinates('Slice',$end_slice+1,$end_slice + $allele_feature->length,1,'IndividualSlice',$end_individual+1,$end_individual + $length_allele);
369 $start_individual = $end_individual + $length_allele + 1; #set the beginning of the individual after the allele
370 }
371 else{
372 #add the indel
373 $mapper->add_indel_coordinates('Slice',$end_slice+1,$end_slice + $allele_feature->length,1,'IndividualSlice',$end_individual+1,$end_individual + $length_allele);
374 $start_individual += $length_allele;
375 }
376 $start_slice = $end_slice + $allele_feature->length+ 1; #set the beginning of the slice after the variation feature
377 }
378 }
379 if ($start_slice <= $self->length){
380 #if we haven't reached the end of the IndividualSlice, add the final map coordinates between the individual and the slice
381 $mapper->add_map_coordinates('Slice',$start_slice,$self->length,1,'IndividualSlice',$start_individual,$start_individual + $self->length - $start_slice);
382 }
383
384 $mapper->add_map_coordinates('Slice', -$self->start+1, 0,1, 'IndividualSlice', -$self->start +1,0) if ($self->start > 0); #before individualSlice
385 $mapper->add_map_coordinates('Slice', $self->length + 1,$self->seq_region_length - ($self->length +1),1, 'IndividualSlice', $self->length + 1 + $total_length_diff,$self->seq_region_length + $total_length_diff - ($self->length +1) ) if ($self->length <= $self->seq_region_length); #after strainSlice
386 $self->{'mapper'} = $mapper;
387 }
388 return $self->{'mapper'};
389 }
390
391 =head2 sub_Slice
392
393 Arg 1 : int $start
394 Arg 2 : int $end
395 Arge [3] : int $strand
396 Example : none
397 Description: Makes another IndividualSlice that covers only part of this IndividualSlice
398 with the appropriate differences to the reference Slice
399 If a slice is requested which lies outside of the boundaries
400 of this function will return undef. This means that
401 behaviour will be consistant whether or not the slice is
402 attached to the database (i.e. if there is attached sequence
403 to the slice). Alternatively the expand() method or the
404 SliceAdaptor::fetch_by_region method can be used instead.
405 Returntype : Bio::EnsEMBL::IndividualSlice or undef if arguments are wrong
406 Exceptions : thrown when trying to get the subSlice in the middle of a
407 insertion
408 Caller : general
409
410 =cut
411
412 sub sub_Slice {
413 my ( $self, $start, $end, $strand ) = @_;
414 my $mapper = $self->mapper();
415 #map from the Individual to the Slice to get the sub_Slice, and then, apply the differences in the subSlice
416 my @results = $mapper->map_coordinates('IndividualSlice',$start,$end,$strand,'IndividualSlice');
417 my $new_start;
418 my $new_end;
419 my $new_strand;
420 my $new_seq;
421 #Get need start and end for the subSlice of the IndividualSlice
422 my @results_ordered = sort {$a->start <=> $b->start} grep {ref($_) eq 'Bio::EnsEMBL::Mapper::Coordinate'} @results;
423 $new_start = $results_ordered[0]->start();
424 $new_strand = $results_ordered[0]->strand() if (ref($results_ordered[0]) eq 'Bio::EnsEMBL::Mapper::Coordinate');
425 # $new_strand = $results_ordered[-1]->strand() if (ref($results_ordered[-1]) eq 'Bio::EnsEMBL::Mapper::Coordinate');
426 $new_end = $results_ordered[-1]->end(); #get last element of the array, the end of the slice
427
428 my $subSlice = $self->SUPER::sub_Slice($new_start,$new_end,$new_strand);
429 $subSlice->{'individual_name'} = $self->{'individual_name'};
430
431 my $new_alleles; #reference to an array that will contain the variationFeatures in the new subSlice
432 #update the VariationFeatures in the sub_Slice of the Individual
433 my %af;
434 my $new_allele_feature;
435 foreach my $alleleFeature (@{$self->{'alleleFeatures'}}){
436 $new_allele_feature = $alleleFeature->transfer($subSlice);
437 #only transfer the coordinates to the SubSlice that are within the boundaries
438 if ($new_allele_feature->start >= 1 && $new_allele_feature->end <= $subSlice->length){
439 push @{$new_alleles}, $new_allele_feature;
440 }
441 }
442 $subSlice->{'alleleFeatures'} = $new_alleles;
443 return $subSlice;
444
445 }
446
447 =head2 subseq
448
449 Arg [1] : int $startBasePair
450 relative to start of slice, which is 1.
451 Arg [2] : int $endBasePair
452 relative to start of slice.
453 Arg [3] : (optional) int $strand
454 The strand of the individual slice to obtain sequence from. Default
455 value is 1.
456 Description: returns string of dna sequence
457 Returntype : txt
458 Exceptions : end should be at least as big as start
459 strand must be set
460 Caller : general
461
462 =cut
463
464 sub subseq {
465 my ( $self, $start, $end, $strand ) = @_;
466
467 if ( $end+1 < $start ) {
468 throw("End coord + 1 is less than start coord");
469 }
470
471 # handle 'between' case for insertions
472 return '' if( $start == $end + 1);
473
474 $strand = 1 unless(defined $strand);
475
476 if ( $strand != -1 && $strand != 1 ) {
477 throw("Invalid strand [$strand] in call to Slice::subseq.");
478 }
479
480 my $subseq;
481 my $seq;
482 if($self->adaptor){
483 my $seqAdaptor = $self->adaptor()->db()->get_SequenceAdaptor();
484 $subseq = ${$seqAdaptor->fetch_by_Slice_start_end_strand($self,$start,$end,$strand)}; #get the reference sequence for that slice
485 #apply all differences to the reference sequence
486 # sort edits in reverse order to remove complication of
487 # adjusting downstream edits
488 my @allele_features_ordered = sort {$b->start() <=> $a->start() || $b->end() <=> $a->end()} @{$self->{'alleleFeatures'}} if (defined $self->{'alleleFeatures'});
489 my $af_start;
490 my $af_end;
491 foreach my $af (@allele_features_ordered){
492 if (($af->start - $start +1 > 0) && ($end - $af->end > 0)){
493 #save the current start and end of the alleleFeature before changing for apply_edit
494 $af_start = $af->start;
495 $af_end = $af->end;
496 #apply the difference if the feature is in the new slice
497 $af->start($af->start - $start +1);
498 $af->end($af->end - $start +1);
499 $af->apply_edit(\$subseq); #change, in the reference sequence, the af
500 #restore the initial values of alleleFeature start and end
501 $af->start($af_start);
502 $af->end($af_end);
503
504 }
505 }
506 }
507 else {
508 ## check for gap at the beginning and pad it with Ns
509 if ($start < 1) {
510 $subseq = "N" x (1 - $start);
511 $start = 1;
512 }
513 $subseq .= substr ($self->seq(), $start-1, $end - $start + 1);
514 ## check for gap at the end and pad it with Ns
515 if ($end > $self->length()) {
516 $subseq .= "N" x ($end - $self->length());
517 }
518 reverse_comp(\$subseq) if($strand == -1);
519 }
520 return $subseq;
521
522 }
523
524 =head2 get_all_Transcripts
525
526 Args : None
527 Example : @transcripts = @{$individualslice->get_all_Transcripts)};
528 Description: Gets all transcripts which overlap this Individual Slice. If you want to
529 specify a particular analysis or type, then you are better off
530 using get_all_Genes or get_all_Genes_by_type and iterating
531 through the transcripts of each gene.
532 Returntype : reference to a list of Bio::EnsEMBL::Transcripts
533 Exceptions : none
534 Caller : general
535
536 =cut
537
538 sub get_all_Transcripts {
539 my $self = shift;
540
541 my $transcripts = $self->SUPER::get_all_Transcripts(1);
542 $self->map_to_Individual($transcripts);
543
544 return $transcripts;
545 }
546
547
548 =head2 get_all_Exons
549
550 Arg [1] : (optional) string $dbtype
551 The dbtype of exons to obtain. This assumes that the db has
552 been added to the DBAdaptor under this name (using the
553 DBConnection::add_db_adaptor method).
554 Example : @exons = @{$individualSlice->get_all_Exons};
555 Description: Gets all exons which overlap this IndividualSlice. Note that these exons
556 will not be associated with any transcripts, so this may not
557 be terribly useful.
558 Returntype : reference to a list of Bio::EnsEMBL::Exons
559 Exceptions : none
560 Caller : general
561
562 =cut
563
564 sub get_all_Exons {
565 my $self = shift;
566 my $dbtype = shift;
567
568 my $exons = $self->SUPER::get_all_Exons($dbtype);
569 $self->map_to_Individual($exons); #map the exons to the Individual
570
571 return $exons;
572 }
573
574 =head2 get_all_Genes
575
576 Arg [1] : (optional) string $logic_name
577 The name of the analysis used to generate the genes to retrieve
578 Arg [2] : (optional) string $dbtype
579 The dbtype of genes to obtain. This assumes that the db has
580 been added to the DBAdaptor under this name (using the
581 DBConnection::add_db_adaptor method).
582 Example : @genes = @{$individualSlice->get_all_Genes};
583 Description: Retrieves all genes that overlap this slice.
584 Returntype : listref of Bio::EnsEMBL::Genes
585 Exceptions : none
586 Caller : none
587
588 =cut
589
590 sub get_all_Genes{
591 my ($self, $logic_name, $dbtype) = @_;
592
593 my $genes = $self->SUPER::get_all_Genes($logic_name, $dbtype, 1);
594
595 $self->map_to_Individual($genes);
596
597 foreach my $gene (@{$genes}){
598 $self->map_to_Individual($gene->get_all_Exons); #map the Exons to the Individual
599 $self->map_to_Individual($gene->get_all_Transcripts); #map the Transcripts to the Individual
600 }
601
602 return $genes;
603 }
604
605 =head2 map_to_Individual
606
607 Arg[1] : ref $features
608 Example : $individualSlice->map_to_Individual($exons);
609 Description : Gets the features from the Slice and maps it in the IndividualSlice, using the mapper
610 between Slice and IndividualSlice
611 ReturnType : None
612 Exceptions : None
613 Caller : general
614
615 =cut
616
617 sub map_to_Individual{
618 my $self = shift;
619 my $features = shift;
620
621 my $mapper = $self->mapper();
622 my (@results, @results_ordered, $new_start, $new_end, $new_strand);
623 #foreach of the transcripts, map them to the IndividualSlice and replace the Slice with the IndividualSlice
624 foreach my $feature (@{$features}){
625 $feature->slice($self); #replace the IndividualSlice as the Slice for this feature (the Slice plus the AlleleFeatures)
626 #map from the Slice to the Individual Slice
627 my @results = $mapper->map_coordinates('Slice',$feature->start,$feature->end,$feature->strand,'Slice');
628 #from the results, order them but filter out those that are not coordinates
629 @results_ordered = sort {$a->start <=> $b->start} grep {ref($_) eq 'Bio::EnsEMBL::Mapper::Coordinate'} @results;
630 $new_start = $results_ordered[0]->start();
631 $new_strand = $results_ordered[0]->strand();
632 $new_end = $results_ordered[-1]->end(); #get last element of the array, the end of the slice
633 $feature->start($new_start); #update new coordinates
634 $feature->end($new_end);
635 $feature->strand($new_strand);
636 }
637 }
638
639 sub alleleFeatures{
640 my $self = shift;
641 return $self->{'alleleFeatures'};
642 }
643
644 sub add_AlleleFeature{
645 my $self = shift;
646
647 if (@_){
648 if(!ref($_[0]) || !$_[0]->isa('Bio::EnsEMBL::Variation::AlleleFeature')) {
649 throw("Bio::EnsEMBL::Variation::AlleleFeature argument expected");
650 }
651 #add the alleleFeature to the individualSlice
652 push @{$self->{'alleleFeatures'}},shift;
653 }
654 }
655 1;