comparison variant_effect_predictor/Bio/Variation/IO/flat.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 # $Id: flat.pm,v 1.12 2002/10/22 07:38:50 lapp Exp $
2 # BioPerl module for Bio::Variation::IO::flat
3 #
4 # Cared for by Heikki Lehvaslaiho <Heikki@ebi.ac.uk>
5 #
6 # Copyright Heikki Lehvaslaiho
7 #
8 # You may distribute this module under the same terms as perl itself
9 #
10
11 # POD documentation - main docs before the code
12
13 =head1 NAME
14
15 Bio::Variation::IO::flat - flat file sequence variation input/output stream
16
17 =head1 SYNOPSIS
18
19 Do not use this module directly. Use it via the Bio::Variation::IO class.
20
21 =head1 DESCRIPTION
22
23 This object can transform Bio::Variation::SeqDiff objects to and from
24 flat file databases. The format used is EMBL like extension of what is
25 used by the "EBI Mutation Checker" at
26 http://www.ebi.ac.uk/cgi-bin/mutations/check.cgi and will eventually
27 replace it.
28
29 More information of the attributes and values use can be found at
30 http://www.ebi.ac.uk/mutations/recommendations/.
31
32 =head1 FEEDBACK
33
34 =head2 Mailing Lists
35
36 User feedback is an integral part of the evolution of this and other
37 Bioperl modules. Send your comments and suggestions preferably to the
38 Bioperl mailing lists Your participation is much appreciated.
39
40 bioperl-l@bioperl.org - General discussion
41 http://bio.perl.org/MailList.html - About the mailing lists
42
43 =head2 Reporting Bugs
44
45 report bugs to the Bioperl bug tracking system to help us keep track
46 the bugs and their resolution. Bug reports can be submitted via
47 email or the web:
48
49 bioperl-bugs@bio.perl.org
50 http://bugzilla.bioperl.org/
51
52 =head1 AUTHOR - Heikki Lehvaslaiho
53
54 Email: heikki@ebi.ac.uk
55 Address:
56
57 EMBL Outstation, European Bioinformatics Institute
58 Wellcome Trust Genome Campus, Hinxton
59 Cambs. CB10 1SD, United Kingdom
60
61
62
63 =head1 APPENDIX
64
65 The rest of the documentation details each of the object
66 methods. Internal methods are usually preceded with a _
67
68 =cut
69
70 # Let the code begin...
71
72 package Bio::Variation::IO::flat;
73 my $VERSION=1.0;
74 use vars qw(@ISA);
75 use strict;
76
77 use Text::Wrap;
78 use Bio::Variation::IO;
79 use Bio::Variation::SeqDiff;
80 use Bio::Variation::DNAMutation;
81 use Bio::Variation::RNAChange;
82 use Bio::Variation::AAChange;
83 use Bio::Variation::Allele;
84
85
86 @ISA = qw(Bio::Variation::IO);
87
88 sub new {
89 my($class, @args) = @_;
90 my $self = bless {}, $class;
91 $self->_initialize(@args);
92 return $self;
93 }
94
95 sub _initialize {
96 my($self,@args) = @_;
97 return unless $self->SUPER::_initialize(@args);
98 }
99
100 =head2 next
101
102
103 Title : next
104 Usage : $haplo = $stream->next()
105 Function: returns the next seqDiff in the stream
106 Returns : Bio::Variation::SeqDiff object
107 Args : NONE
108
109 =cut
110
111 sub next {
112 my( $self ) = @_;
113 local $/ = '//';
114 return unless my $entry = $self->_readline;
115
116 return if $entry =~ /^\s+$/;
117
118 $entry =~ /\s*ID\s+\S+/ || $self->throw("We do need an ID!");
119
120 my ($id, $offset, $alphabet) = $entry =~ /\s*ID +([^:]+)..(\d+)[^\)]*.\[?([cg])?/
121 or $self->throw("Can't parse ID line");
122 # $self->throw("$1|$2|$3");
123 my $h =Bio::Variation::SeqDiff->new(-id => $id,
124 -offset => $offset,
125 );
126 if ($alphabet) {
127 if ($alphabet eq 'g') {
128 $alphabet = 'dna';
129 }
130 elsif ($alphabet eq 'c') {
131 $alphabet = 'rna';
132 }
133 $h->alphabet($alphabet);
134 }
135 #
136 # DNA
137 #
138 my @dna = split ( / DNA;/, $entry );
139 shift @dna;
140 my $prevdnaobj;
141 foreach my $dna (@dna) {
142 $dna =~ s/Feature[ \t]+//g;
143 ($dna) = split "RNA; ", $dna;
144 #$self->warn("|$dna|") ;
145 #exit;
146 my ($mut_number, $proof, $location, $upflank, $change, $dnflank) =
147 $dna =~ m|\W+([\d\.]+).+/proof: (\w+).+/location: ([^ \n]+).+/upflank: ([ \n\w]+).+/change: ([^ /]+).+/dnflank: ([ \n\w]+)|s;
148 $change =~ s/[ \n]//g;
149 my ($ori, $mut) = split /[>\|]/, $change;
150 my ($variation_number, $change_number) = split /\./, $mut_number;
151 #$self->warn("|$mut_number|>|$variation_number|$change_number|");
152 my $dnamut;
153 if ($change_number and $change_number > 1 ) {
154 my $a3 = Bio::Variation::Allele->new;
155 $a3->seq($mut) if $mut;
156 #$dnamut->add_Allele($a3);
157 $prevdnaobj->add_Allele($a3);
158 } else {
159 $upflank =~ s/[ \n]//g;
160 $dnflank =~ s/[ \n]//g;
161 my ($region, $junk, $region_value, $junk2, $region_dist) =
162 $dna =~ m|.+/region: ([\w\']+)(; )?(\w+)?( ?\(\+?)?(-?\d+)?|s;
163 #my $s = join ("|", $mut_number, $proof, $location, $upflank,
164 # $change, $dnflank, $region, $region_value, $region_dist, $1,$2,$3,$4,$5);
165 #$self->warn($s);
166 #exit;
167 my ($start, $sep, $end) = $location =~ /(-?\d+)(.)?\D?(-?\d+)?/;
168 $end = $start if not $end ;
169 my ($len) = $end - $start +1;
170 $len = 0, $start = $end if defined $sep and $sep eq '^';
171 my $ismut = 0;
172 $ismut = 1 if $change =~ m/>/;
173
174 $dnamut = Bio::Variation::DNAMutation->new
175 ('-start' => $start,
176 '-end' => $end,
177 '-length' => $len,
178 '-upStreamSeq' => $upflank,
179 '-dnStreamSeq' => $dnflank,
180 '-proof' => $proof,
181 '-mut_number' => $mut_number
182 );
183 $prevdnaobj = $dnamut;
184 my $a1 = Bio::Variation::Allele->new;
185 $a1->seq($ori) if $ori;
186 $dnamut->allele_ori($a1);
187 my $a2 = Bio::Variation::Allele->new;
188 $a2->seq($mut) if $mut;
189 $dnamut->add_Allele($a2);
190 if ($ismut) {
191 $dnamut->isMutation(1);
192 $dnamut->allele_mut($a2);
193 }
194 $dnamut->region($region) if defined $region;
195 $dnamut->region_value($region_value) if defined $region_value;
196 $dnamut->region_dist($region_dist) if defined $region_dist;
197
198 $h->add_Variant($dnamut);
199 $dnamut->SeqDiff($h);
200 }
201 }
202
203 #
204 # RNA
205 #
206 my @rna = split ( / RNA;/, $entry );
207 shift @rna;
208 my $prevrnaobj;
209 foreach my $rna (@rna) {
210 $rna = substr ($rna, 0, index($rna, 'Feature AA'));
211 $rna =~ s/Feature[ \t]+//g;
212 ($rna) = split "DNA; ", $rna;
213 #$self->warn("|$rna|") ;
214 my ($mut_number, $proof, $location, $upflank, $change, $dnflank) =
215 $rna =~ m|\W+([\d\.]+).+/proof: (\w+).+/location: ([^ \n]+).+/upflank: (\w+).+/change: ([^/]+).+/dnflank: (\w+)|s ;#'
216 my ($region, $junk, $region_value, $junk2, $region_dist) =
217 $rna =~ m|.+/region: ([\w\']+)(; )?(\w+)?( ?\(\+?)?(-?\d+)?|s;
218 #my $s = join ("|", $mut_number, $proof, $location, $upflank,
219 # $change, $dnflank, $region, $region_value, $region_dist, $1,$2,$3,$4,$5);
220 #$self->warn($s);
221 #exit;
222 $change =~ s/[ \n]//g;
223 my ($ori, $mut) = split /[>\|]/, $change;
224 my $rnamut;
225 my ($variation_number, $change_number) = split /\./, $mut_number;
226 if ($change_number and $change_number > 1 ) {
227 my $a3 = Bio::Variation::Allele->new;
228 $a3->seq($mut) if $mut;
229 #$rnamut->add_Allele($a3);
230 $prevrnaobj->add_Allele($a3);
231 } else {
232 my ($start, $sep, $end) = $location =~ /(-?\d+)(.)?\D?(-?\d+)?/;
233 $end = $start if not $end ;
234 my ($len) = $end - $start + 1;
235 $len = 0, $start = $end if defined $sep and $sep eq '^';
236 my $ismut;
237 $ismut = 1 if $change =~ m/>/;
238 my ($codon_table) = $rna =~ m|.+/codon_table: (\d+)|s;
239 my ($codon_pos) = $rna =~ m|.+/codon:[^;]+; ([123])|s;
240
241 $rnamut = Bio::Variation::RNAChange->new
242 ('-start' => $start,
243 '-end' => $end,
244 '-length' => $len,
245 '-upStreamSeq' => $upflank,
246 '-dnStreamSeq' => $dnflank,
247 '-proof' => $proof,
248 '-mut_number' => $mut_number
249
250 );
251 $prevrnaobj = $rnamut;
252 my $a1 = Bio::Variation::Allele->new;
253 $a1->seq($ori) if $ori;
254 $rnamut->allele_ori($a1);
255 my $a2 = Bio::Variation::Allele->new;
256 $a2->seq($mut) if $mut;
257 $rnamut->add_Allele($a2);
258 if ($ismut) {
259 $rnamut->isMutation(1);
260 $rnamut->allele_mut($a2);
261 }
262 $rnamut->region($region) if defined $region;
263 $rnamut->region_value($region_value) if defined $region_value;
264 $rnamut->region_dist($region_dist) if defined $region_dist;
265
266 $rnamut->codon_table($codon_table) if $codon_table;
267 $rnamut->codon_pos($codon_pos) if $codon_pos;
268 $h->add_Variant($rnamut);
269 foreach my $mut ($h->each_Variant) {
270 if ($mut->isa('Bio::Variation::DNAMutation') ) {
271 if ($mut->mut_number == $rnamut->mut_number) {
272 $rnamut->DNAMutation($mut);
273 $mut->RNAChange($rnamut);
274 }
275 }
276 }
277 }
278 }
279 #
280 # AA
281 #
282 my @aa = split ( / AA;/, $entry );
283 shift @aa;
284 my $prevaaobj;
285 foreach my $aa (@aa) {
286 $aa = substr ($aa, 0, index($aa, 'Feature AA'));
287 $aa =~ s/Feature[ \t]+//g;
288 ($aa) = split "DNA; ", $aa;
289 #$self->warn("|$aa|") ;
290 my ($mut_number, $proof, $location, $change) =
291 $aa =~ m|\W+([\d\.]+).+/proof: (\w+).+/location: ([^ \n]+)./change: ([^/;]+)|s;
292 $change =~ s/[ \n]//g;
293 #my $s = join ("|", $mut_number, $proof, $location, $change);
294 #$self->warn($s);
295 #exit;
296 $change =~ s/[ \n]//g;
297 $change =~ s/DNA$//;
298 my ($ori, $mut) = split /[>\|]/, $change;
299 #print "------$location----$ori-$mut-------------\n";
300 my ($variation_number, $change_number) = split /\./, $mut_number;
301 my $aamut;
302 if ($change_number and $change_number > 1 ) {
303 my $a3 = Bio::Variation::Allele->new;
304 $a3->seq($mut) if $mut;
305 $prevaaobj->add_Allele($a3);
306 } else {
307 my ($start, $sep, $end) = $location =~ /(-?\d+)(.)?\D?(-?\d+)?/;
308 $end = $start if not $end ;
309 my ($len) = $end - $start + 1;
310 $len = 0, $start = $end if defined $sep and $sep eq '^';
311 my $ismut;
312 $ismut = 1 if $change =~ m/>/;
313 my ($region) = $aa =~ m|.+/region: (\w+)|s ;
314 $aamut = Bio::Variation::AAChange->new
315 ('-start' => $start,
316 '-end' => $end,
317 '-length' => $len,
318 '-proof' => $proof,
319 '-mut_number' => $mut_number
320 );
321 $prevaaobj = $aamut;
322 my $a1 = Bio::Variation::Allele->new;
323 $a1->seq($ori) if $ori;
324 $aamut->allele_ori($a1);
325 my $a2 = Bio::Variation::Allele->new;
326 $a2->seq($mut) if $mut;
327 $aamut->add_Allele($a2);
328 if ($ismut) {
329 $aamut->isMutation(1);
330 $aamut->allele_mut($a2);
331 }
332 $region && $aamut->region($region);
333 $h->add_Variant($aamut);
334 foreach my $mut ($h->each_Variant) {
335 if ($mut->isa('Bio::Variation::RNAChange') ) {
336 if ($mut->mut_number == $aamut->mut_number) {
337 $aamut->RNAChange($mut);
338 $mut->AAChange($aamut);
339 }
340 }
341 }
342
343 }
344 }
345 return $h;
346 }
347
348 =head2 write
349
350 Title : write
351 Usage : $stream->write(@seqDiffs)
352 Function: writes the $seqDiff object into the stream
353 Returns : 1 for success and 0 for error
354 Args : Bio::Variation::SeqDiff object
355
356
357 =cut
358
359 sub write {
360 my ($self,@h) = @_;
361
362 #$columns = 75; #default for Text::Wrap
363 my %tag =
364 (
365 'ID' => 'ID ',
366 'Description' => 'Description ',
367 'FeatureKey' => 'Feature ',
368 'FeatureQual' => "Feature ",
369 'FeatureWrap' => "Feature ",
370 'ErrorComment' => 'Comment '
371 #'Comment' => 'Comment -!-',
372 #'CommentLine' => 'Comment ',
373 );
374
375 if( !defined $h[0] ) {
376 $self->throw("Attempting to write with no information!");
377 }
378
379 foreach my $h (@h) {
380
381 my @entry =();
382
383 my ($text, $tmp, $tmp2, $sep);
384 my ($count) = 0;
385
386
387 $text = $tag{ID};
388
389 $text .= $h->id;
390 $text .= ":(". $h->offset;
391 $text .= "+1" if $h->sysname =~ /-/;
392 $text .= ")". $h->sysname;
393 $text .= "; ". $h->trivname if $h->trivname;
394 push (@entry, $text);
395
396 #Variants need to be ordered accoding to mutation_number attribute
397 #put them into a hash of arrays holding the Variant objects
398 #This is necessary for cases like several distict mutations present
399 # in the same sequence.
400 my @allvariants = $h->each_Variant;
401 my %variants = ();
402 foreach my $mut ($h->each_Variant) {
403 push @{$variants{$mut->mut_number} }, $mut;
404 }
405 #my ($variation_number, $change_number) = split /\./, $mut_number;
406 foreach my $var (sort keys %variants) {
407 #print $var, ": ", join (" ", @{$variants{$var}}), "\n";
408
409 foreach my $mut (@{$variants{$var}}) {
410 #
411 # DNA
412 #
413 if ( $mut->isa('Bio::Variation::DNAMutation') ) {
414 #collect all non-reference alleles
415 $self->throw("allele_ori needs to be defined in [$mut]")
416 if not $mut->allele_ori;
417 if ($mut->isMutation) {
418 $sep = '>';
419 } else {
420 $sep = '|';
421 }
422 my @alleles = $mut->each_Allele;
423 #push @alleles, $mut->allele_mut if $mut->allele_mut;
424 my $count = 0; # two alleles
425 foreach my $allele (@alleles) {
426 $count++;
427 my ($variation_number, $change_number) = split /\./, $mut->mut_number;
428 if ($change_number and $change_number != $count){
429 $mut->mut_number("$change_number.$count");
430 }
431 $mut->allele_mut($allele);
432 push (@entry,
433 $tag{FeatureKey}. 'DNA'. "; ". $mut->mut_number
434 );
435 #label
436 $text=$tag{FeatureQual}. '/label: '. $mut->label;
437 push (@entry, $text);
438
439 #proof
440 if ($mut->proof) {
441 $text = $tag{FeatureQual}. '/proof: '. $mut->proof;
442 push (@entry, $text) ;
443 }
444 #location
445 $text = $tag{FeatureQual}. '/location: ';
446 #$mut->id. '; '. $mut->start;
447 if ($mut->length > 1 ) {# if ($mut->end - $mut->start ) {
448 my $l = $mut->start + $mut->length -1;
449 $text .= $mut->start. '..'. $l;
450 }
451 elsif ($mut->length == 0) {
452 my $tmp_start = $mut->start - 1;
453 $tmp_start-- if $tmp_start == 0;
454 $text .= $tmp_start. '^'. $mut->end;
455 } else {
456 $text .= $mut->start;
457 }
458
459 if ($h->alphabet && $h->alphabet eq 'dna') {
460 $tmp = $mut->start + $h->offset;
461 $tmp-- if $tmp <= 0;
462 $mut->start < 1 && $tmp++;
463 #$text.= ' ('. $h->id. '::'. $tmp;
464 $tmp2 = $mut->end + $h->offset;
465 if ( $mut->length > 1 ) {
466 $mut->end < 1 && $tmp2++;
467 $text.= ' ('. $h->id. '::'. $tmp. "..". $tmp2;
468 }
469 elsif ($mut->length == 0) {
470 $tmp--;
471 $tmp-- if $tmp == 0;
472 $text .= ' ('. $h->id. '::'. $tmp. '^'. $tmp2;
473 } else {
474 $text.= ' ('. $h->id. '::'. $tmp;
475 }
476 $text .= ')';
477 }
478 push (@entry, $text);
479 #sequence
480 push (@entry,
481 $tag{FeatureQual}. '/upflank: '. $mut->upStreamSeq
482 );
483 $text = '';
484 $text = $mut->allele_ori->seq if $mut->allele_ori->seq;
485 $text .= $sep;
486 $text .= $mut->allele_mut->seq if $mut->allele_mut->seq;
487 push (@entry,
488 wrap($tag{FeatureQual}. '/change: ', $tag{FeatureWrap},
489 $text)
490 );
491
492 push (@entry,
493 $tag{FeatureQual}. '/dnflank: '. $mut->dnStreamSeq
494 );
495 #restriction enzyme
496 if ($mut->restriction_changes ne '') {
497 $text = $mut->restriction_changes;
498 $text = wrap($tag{FeatureQual}. '/re_site: ', $tag{FeatureWrap}, $text);
499 push (@entry,
500 $text
501 );
502 }
503 #region
504 if ($mut->region ) {
505 $text = $tag{FeatureQual}. '/region: '. $mut->region;
506 $text .= ';' if $mut->region_value or $mut->region_dist;
507 $text .= ' '. $mut->region_value if $mut->region_value;
508 if ($mut->region_dist ) {
509 $tmp = '';
510 $tmp = '+' if $mut->region_dist > 1;
511 $text .= " (". $tmp. $mut->region_dist. ')';
512 }
513 push (@entry, $text);
514 }
515 #CpG
516 if ($mut->CpG) {
517 push (@entry,
518 $tag{FeatureQual}. "/CpG"
519 );
520 }
521 }
522 }
523 #
524 # RNA
525 #
526 elsif ($mut->isa('Bio::Variation::RNAChange') ) {
527 #collect all non-reference alleles
528 $self->throw("allele_ori needs to be defined in [$mut]")
529 if not $mut->allele_ori;
530 my @alleles = $mut->each_Allele;
531 #push @alleles, $mut->allele_mut if $mut->allele_mut;
532 if ($mut->isMutation) {
533 $sep = '>';
534 } else {
535 $sep = '|';
536 }
537
538 my $count = 0; # two alleles
539 foreach my $allele (@alleles) {
540 $count++;
541 my ($variation_number, $change_number) = split /\./, $mut->mut_number;
542 if ($change_number and $change_number != $count){
543 $mut->mut_number("$change_number.$count");
544 }
545 $mut->allele_mut($allele);
546 push (@entry,
547 $tag{FeatureKey}. 'RNA'. "; ". $mut->mut_number
548 );
549 #label
550 $text=$tag{FeatureQual}. '/label: '. $mut->label;
551 push (@entry, $text);
552 #proof
553 if ($mut->proof) {
554 $text = $tag{FeatureQual}. '/proof: '. $mut->proof;
555 push (@entry, $text) ;
556 }
557 #location
558 $text = $tag{FeatureQual}. '/location: ' ;
559 if ($mut->length > 1 ) {
560 $text .= $mut->start. '..'. $mut->end;
561 $tmp2 = $mut->end + $h->offset;
562 }
563 elsif ($mut->length == 0) {
564 my $tmp_start = $mut->start;
565 $tmp_start--;
566 $tmp_start-- if $tmp_start == 0;
567 $text .= $tmp_start. '^'. $mut->end;
568 } else {
569 $text .= $mut->start;
570 }
571
572 if ($h->alphabet && $h->alphabet eq 'rna') {
573 $tmp = $mut->start + $h->offset;
574 $tmp-- if $tmp <= 0;
575 #$mut->start < 1 && $tmp++;
576 #$text.= ' ('. $h->id. '::'. $tmp;
577 $tmp2 = $mut->end + $h->offset;
578 #$mut->end < 1 && $tmp2++;
579 if ( $mut->length > 1 ) {
580 $text.= ' ('. $h->id. '::'. $tmp. "..". $tmp2;
581 }
582 elsif ($mut->length == 0) {
583 $tmp--;
584 $text .= ' ('. $h->id. '::'. $tmp. '^'. $tmp2;
585 } else {
586 $text.= ' ('. $h->id. '::'. $tmp;
587 }
588
589 $text .= ')';
590 }
591 push (@entry, $text);
592
593 #sequence
594 push (@entry,
595 $tag{FeatureQual}. '/upflank: '. $mut->upStreamSeq
596 );
597 $text = '';
598 $text = $mut->allele_ori->seq if $mut->allele_ori->seq;
599 $text .= $sep;
600 $text .= $mut->allele_mut->seq if $mut->allele_mut->seq;
601 push (@entry,
602 wrap($tag{FeatureQual}. '/change: ', $tag{FeatureWrap},
603 $text)
604 );
605 push (@entry,
606 $tag{FeatureQual}. '/dnflank: '. $mut->dnStreamSeq
607 );
608 #restriction
609 if ($mut->restriction_changes ne '') {
610 $text = $mut->restriction_changes;
611 $text = wrap($tag{FeatureQual}. '/re_site: ', $tag{FeatureWrap}, $text);
612 push (@entry,
613 $text
614 );
615 }
616 #coding
617 if ($mut->region eq 'coding') {
618 #codon table
619 $text = $tag{FeatureQual}. '/codon_table: ';
620 $text .= $mut->codon_table;
621 push (@entry, $text);
622 #codon
623
624 $text = $tag{FeatureQual}. '/codon: '. $mut->codon_ori. $sep;
625 if ($mut->DNAMutation->label =~ /.*point/) {
626 $text .= $mut->codon_mut;
627 }
628 else {
629 $text .= '-';
630 }
631 $text .= "; ". $mut->codon_pos;
632 push (@entry, $text);
633 }
634 #region
635 if ($mut->region ) {
636 $text = $tag{FeatureQual}. '/region: '. $mut->region;
637 $text .= ';' if $mut->region_value or $mut->region_dist;
638 $text .= ' '. $mut->region_value if $mut->region_value;
639 if ($mut->region_dist ) {
640 $tmp = '';
641 $tmp = '+' if $mut->region_dist > 1;
642 $text .= " (". $tmp. $mut->region_dist. ')';
643 }
644 push (@entry, $text);
645 }
646 }
647 }
648 #
649 # AA
650 #
651 elsif ($mut->isa('Bio::Variation::AAChange')) {
652 #collect all non-reference alleles
653 $self->throw("allele_ori needs to be defined in [$mut]")
654 if not $mut->allele_ori;
655 if ($mut->isMutation) {
656 $sep = '>';
657 } else {
658 $sep = '|';
659 }
660 my @alleles = $mut->each_Allele;
661 #push @alleles, $mut->allele_mut if $mut->allele_mut;
662 my $count = 0; # two alleles
663 foreach my $allele (@alleles) {
664 $count++;
665 my ($variation_number, $change_number) = split /\./, $mut->mut_number;
666 if ($change_number and $change_number != $count){
667 $mut->mut_number("$change_number.$count");
668 }
669 $mut->allele_mut($allele);
670 push (@entry,
671 $tag{FeatureKey}. 'AA'. "; ". $mut->mut_number
672 );
673 #label
674 $text=$tag{FeatureQual}. '/label: '. $mut->label;
675 push (@entry, $text) ;
676 #proof
677 if ($mut->proof) {
678 $text = $tag{FeatureQual}. '/proof: '. $mut->proof;
679 push (@entry, $text) ;
680 }
681 #location
682 $text = $tag{FeatureQual}. '/location: '.
683 #$mut->id. '; '. $mut->start;
684 $mut->start;
685 if ($mut->length > 1 ) {
686 $tmp = $mut->start + $mut->length -1;
687 $text .= '..'. $tmp;
688 }
689 push (@entry, $text);
690 #sequence
691 $text = '';
692 $text = $mut->allele_ori->seq if $mut->allele_ori->seq;
693 $text .= $sep;
694 $text .= $mut->allele_mut->seq if $mut->allele_mut->seq;
695 push (@entry,
696 wrap($tag{FeatureQual}. '/change: ', $tag{FeatureWrap},
697 $text)
698 );
699 #region
700 if ($mut->region ) {
701 $text = $tag{FeatureQual}. '/region: '. $mut->region;
702 $text .= ';' if $mut->region_value or $mut->region_dist;
703 $text .= ' '. $mut->region_value if $mut->region_value;
704 if ($mut->region_dist ) {
705 $tmp = '';
706 $tmp = '+' if $mut->region_dist > 1;
707 $text .= " (". $tmp. $mut->region_dist. ')';
708 }
709 push (@entry, $text);
710 }
711 }
712 }
713 }
714 }
715 push (@entry,
716 "//"
717 );
718 my $str = join ("\n", @entry). "\n";
719 $str =~ s/\t/ /g;
720 $self->_print($str);
721 }
722 return 1;
723 }
724
725 1;