source: lib/ProccesAllicVariants.pm @ 12

Last change on this file since 12 was 12, checked in by maarten, 10 years ago

feature ok version

File size: 15.1 KB
Line 
1package ProccesAllicVariants;
2use strict;
3use warnings;
4use Carp;
5use Data::Dumper;
6
7sub new {
8        my $pubmedlocation = $_;
9        my $self           = {};
10        $self->{Omim2Pubmed} = "";
11        $self->{Omim2DBsnp}  = "";
12        $self->{av}          = "";
13
14        bless($self);    # but see below
15        return $self;
16}
17sub saveAllelicVariants{
18        my ($self,$avref)=@_;
19        #create a hash
20        my %avhash=();
21
22        foreach my $av (@$avref){
23                $av->number()=~/(\d{4})/i;
24                $avhash{$1}=$av;
25        }
26        #save the hash
27        $self->{av}=\%avhash;
28
29       
30        }
31
32
33sub trim($) {
34        my $string = shift;
35        if ( length($string) > 0 ) {
36                $string =~ s/^\s+//;
37                $string =~ s/\s+$//;
38        }
39        return $string;
40}
41
42sub getfilename {
43        my ( $title, $alt ) = @_;
44        my $protein = "";
45        #print($title);
46        #print("\n\n".$alt);
47        if ( $title =~ m/\;/ ) {
48                if ( $title =~ /([\w\-]{2,10})\s{0,1}$/ ) {
49                        $protein = $1;
50
51                }
52        }
53
54        if ( $protein eq "" ) {
55
56                my $titlelong = $title . ";" . $alt . ";";
57
58                #remove breaks
59                $titlelong =~ s/\r\n|\n|\r/;/g;
60
61                if ( $protein eq "" ) {
62
63                        if ( $alt =~ /^(\w{2,10})$/ ) {
64                                $protein = $alt;
65
66                        }
67
68                        else {
69                                if ( $titlelong =~ /[\;\:][ ]*?([\w\-]{2,10})[ ]*?\;/ ) {
70                                        $protein = $1;
71                                }
72                        }
73                }
74        }
75        if (length(trim( $protein))==0){
76                if($alt =~ /\;\s([\d\w]{2,8}),\sINCLUDED/){
77                      $protein = $1;
78                     
79                }
80        }
81
82        return ( trim $protein);
83}
84
85sub getAAMutation {
86        my $av = shift;
87
88        my $mutation = $av->aa_ori() . $av->position() . $av->aa_mut();
89#       if ( length($mutation) == 0 ) {
90#               $mutation =
91#                 formatInDelMutation( $av->additional_mutations(),
92#                       $av->description() );
93#
94#       }
95        if ( length($mutation) > 0 ) {
96                $mutation="p.".$mutation;
97        }
98        return ($mutation);
99
100}
101sub getDNAMutation {
102        my $av = shift;
103        my $mutation ="";
104        my $aamutation = $av->aa_ori() . $av->position() . $av->aa_mut();
105        if ( length($mutation) == 0 ) {
106                $mutation =$av->additional_mutations();
107#                 formatInDelMutation( $av->additional_mutations(),
108#                       $av->description() );
109
110        }
111        return ($mutation);
112
113}
114
115sub addOmim2pubmed {
116        my ( $self, $o2p ) = @_;
117        $self->{Omim2Pubmed} = $o2p;
118
119}
120
121sub addOmim2dbsnp {
122        my ( $self, $o2dbsnp ) = @_;
123        $self->{Omim2DBsnp} = $o2dbsnp;
124
125}
126
127sub extractwithposition {
128        my ( $desc, $aapos, $returnDesc ) = @_;
129        my $result = "";
130
131        #110G-T transversion
132        #110G-A transition
133        #
134        my $regex = "(($aapos)[ACTG]-[ACTG] trans[a-z]{4,9})";
135        if ( $desc =~ /$regex/i ) {
136                $result = $1;
137                $desc =~ s/$regex/ /i;
138
139                #print($result);
140        }
141        if ( $returnDesc == 1 ) {
142                $result = $desc;
143        }
144
145        return ($result);
146}
147
148sub extractnuclitdemutation {
149        my ( $desc, $aapos, $returnDesc ) = @_;
150        my $result = "";
151
152        #check corresponding nucl position from aa pos
153        if ( $aapos =~ /[0-9]{1,6}/i ) {
154
155                my $nuclposition = "";
156                $nuclposition = ( $aapos * 3 ) . "|";
157                $nuclposition .= ( ( $aapos * 3 ) - 1 ) . "|" . ( ( $aapos * 3 ) - 2 );
158                if ( $desc =~ /($nuclposition)/ ) {
159
160                        $result = extractwithposition( $desc, $nuclposition, $returnDesc );
161                }
162        }
163
164        if ( length($result) == 0 ) {
165
166                #heterozygous A-to-G transition in exon 3
167                #heterozygous G-to-C transversion in exon 2
168                # T-to-C transition in exon 3
169                my $regex =
170"(([a-z]{3,9}gous)? [ATCG]-to-[ATCG] trans[a-z]{5,8} in exon [0-9]{1,2})";
171
172                #n A-to-G transition at nucleotide 1730 in exon 13, resulting in
173                my $regex2 =
174"([ACTG]-to-[ACTG] trans[a-z]{4,9})[\\w\\s\\-\\.]{0,10}(nucleotide \\d{1,5})[\\w\\s\\-\\.]{0,90}((exon|intron)\\s\\d{1,3})";
175
176                #codon 163 was changed from GTG (val) to CTG (leu)
177                #find 1110A-C transversion with wrong number (!= AA*3)
178                #ACTA1 .0009
179                #codon 163 was changed from GTG (val) to CTG (leu)
180
181                my $regex3 = "(codon \\d{1,4})";
182
183                #find 1110A-C transversion with wrong number (!= AA*3)
184                #ACTA1 .0009
185                my $regex4 = "([0-9]{1,5}[ACTG]-[ACTG] trans[a-z]{4,9})";
186
187                my $regex5 = "([ACTG]{3}-to-[ACTG]{3})";
188
189                my $regex6 =
190"([ACTG]-to-[ACTG] trans[a-z]{4,9})[\\w\\s\\-\\.]{0,90}((exon|intron)\\s\\d{1,3})";
191
192                my $regex7 = "([ACTG]{3} to [ACTG]{3})";
193
194                #heterozygous A-to-G transition in exon 3
195                #heterozygous G-to-C transversion in exon 2
196                # T-to-C transition in exon 3
197                if ( $desc =~ /$regex/i ) {
198                        $result = $1;
199                        $desc =~ s/$regex/ /i;
200
201                }
202
203                #n A-to-G transition at nucleotide 1730 in exon 13, resulting in
204
205                elsif ( $desc =~ /$regex2/i ) {
206                        $result = $1 . " at " . $2 . " in " . $3;
207                        $desc =~ s/$regex2/ /i;
208
209                        #print($result);
210
211                }
212
213                #find 1110A-C transversion with wrong number (!= AA*3)
214                #ACTA1 .0009
215                elsif ( $desc =~ /$regex4/i ) {
216                        $result = $1;
217                        $desc =~ s/$regex4/ /i;
218
219                        #print($result);
220
221                        #T-to-C transition BLABLA in exon 30
222                }
223                elsif ( $desc =~ /$regex5/i ) {
224                        $result = $1;
225                        $desc =~ s/$regex5/ /i;
226
227                        #print($result);
228
229                        #T-to-C transition BLABLA in exon 30
230                }
231
232                elsif ( $desc =~ /$regex6/i ) {
233                        $result = $1 . " in " . $2;
234                        $desc =~ s/$regex6/ /i;
235                }
236                elsif ( $desc =~ /$regex7/ ) {
237                        $result = $1;
238                        $desc =~ s/$regex7/ /i;
239                }
240               
241        }
242        if ($returnDesc==1){
243                $result=$desc;
244        }
245       
246        return ($result);
247}
248
249sub getAuthor {
250
251        my $desc   = shift;
252        my $authAV = "";
253        my $yearAV = "";
254        $desc =~ tr/\n/ /;
255        $desc = trim($desc);
256   my $descOriginal=$desc;
257        my @regexp = ();
258
259        push( @regexp, "([\\w-]*.)\\set al.[,\\s]*?\\(([\\d\\s,]{4,17})\\)" );
260
261        # add regexp for authors with  And between names
262
263        push( @regexp,
264                "\\b([A-Za-z\\-]*.[^\\)\\,]\\sand\\s[\\w-]*.)\\s*.\\(([\\d]{4})\\)" );
265
266        #TODO: add solution for mulitple authors (author 2009; auth3 2002 ,etc )
267
268        #Bradley et al., 1975;
269        push( @regexp, "([\\w-]*.)\\set al.[,\\s]*?([\\d]{4})" );
270
271        #See Brennan (1985)
272
273        push( @regexp, "([\\w-]*.)\\s*.\\(([\\d]{4})\\)" );
274
275        #Tiller et al. (1993, 1995)
276        push( @regexp, "([\\w-]*.)\\set al.\\s*?([\\d\\s,]{4,17})\\)" );
277
278        #de Vries and de Wet (1986, 1987)
279        push( @regexp,
280"\\b([A-Za-z\\-]{3,11}\\sand\\s[\\w\\s-]{3,11})\\s*.\\(([\\d\\s,]{4,17})\\)"
281        );
282
283        #See 123456.1234
284        push( @regexp, "See\\s(\\d{6}).(\\d{4})" );
285
286        my $reg         = "";
287        my @returnarray = ();
288        my @srcs        = ();
289        foreach $reg (@regexp) {
290
291                if ( @srcs = ( $desc =~ /$reg/g ) ) {
292
293                        #print($reg);
294                        push( @returnarray, @srcs );
295                        $desc =~ s/$reg//g;
296                }
297        }
298
299        #check if return has a author and year
300        if ( scalar(@returnarray) % 2 ) {
301                carp( " length of authors is not a multiple from 2 possible error\n"
302                          . Dumper(@returnarray) );
303        }
304        if ( scalar(@returnarray) == 0 ) {
305                carp( "No authors found in:\n\n" . $descOriginal . "\n" );
306        }
307
308        return ( \@returnarray );
309}
310
311sub getRefenceAndClean {
312        my ( $self, $id ) = @_;
313
314 my $cleanDesc="";
315 my $desc="";
316        #carp(Dumper($self->{av}));
317        if (exists($self->{av}->{$id})){
318                my $av=$self->{av}->{$id};
319                $desc= $av->description();
320                $desc =~ tr/\n/ /;
321                $cleanDesc = extractnuclitdemutation( $desc, $av->number(),1 );
322        }else{
323                carp("/n/nReferecne $id not found/n/n");
324        }
325       
326
327        return ($cleanDesc);
328
329}
330
331sub getPubmedID {
332
333        my ( $self, $refsref, $authAV, $yearAV, $numb ) = @_;
334        my $result = "";
335        if ( length($authAV) != 0 ) {
336
337                my @refs    = @$refsref;
338                my $counter = 1;
339                foreach my $ref (@refs) {
340                        my $auth = $ref->authors();
341                        if ( $auth =~ m/^$authAV/ ) {
342                                my $authY      = $ref->location();
343                                my $tempresult = "";
344                                if ( $authY =~ m/$yearAV/ ) {
345
346                                        $tempresult =
347                                          $self->{Omim2Pubmed}->getpubmed( $numb, $counter );
348                                        if ( $tempresult ne "" ) {
349                                                $result .= "{PMID"
350                                                  . $tempresult . ":"
351                                                  . $authAV . " ("
352                                                  . $yearAV . ")}";
353                                        }
354                                        else {
355                                                $result .=
356                                                  $authAV . "(" . $yearAV . ") " . $ref->location();
357
358                                        }
359                                }
360                        }
361                        $counter = $counter + 1;
362                }
363        }
364
365        return ($result);
366}
367sub formatdnamutation{
368        my ($dnaUnformated) = shift;
369        if($dnaUnformated=~/c\.(.*)/){
370                $dnaUnformated=$1;
371        }
372        my $formated="";
373                if ( $dnaUnformated =~ /^\s?(IVS\d{1,5})\w{0,5},\s?([\w\-]{3})\s?,\s?([+\-\d]{2,9})\s?$/ ) {
374                        my $substition=$2;
375                        $substition =~ tr/-/>/;
376                        $formated=$1.$3.$substition;
377        }elsif( $dnaUnformated =~ /^\s?(\d{1,4})-BP (DEL|INS), [A-z]*(\d{1,5})([ACTG]{0,10})\s?$/ ){ 
378                my $amount=$1;
379                my $event=lc($2);
380                my $position=$3;
381                my $sequence=$4;
382                my $postion2="";
383                if ($event eq "del"){
384                         $postion2=($position+1)-$amount;
385                }else{
386                        $postion2=($position+1);
387                }
388                $formated=$position."_".$postion2.$event.$sequence;
389       
390       
391        }else{
392                $formated=$dnaUnformated;
393        }
394        if (length($formated)>0){
395                $formated="c.".$formated;
396        }
397        return($formated);
398}
399
400
401sub fomatYears {
402        my ($yearsUnformated) = shift;
403        my @years = ();
404        if ( $yearsUnformated =~ /^(\d{4})$/ ) {
405                push( @years, $1 );
406        }
407        else {
408                my @catch = ();
409                if ( @catch = ( $yearsUnformated =~ /\d{4}/g ) ) {
410                        push( @years, @catch );
411
412                        #$yearsUnformated =~ s/\d{4}//g;
413                }
414        }
415        if ( scalar(@years) == 0 ) {
416                carp( " No years found in formatYears in string: \n"
417                          . Dumper($yearsUnformated) );
418        }
419
420        return ( \@years );
421
422}
423
424sub formatAuthors {
425        my ($authref)       = shift;
426        my @authorsAndYears = @$authref;
427        my @formatedAuth    = ();
428
429        for ( my $i = 0 ; $i < scalar(@authorsAndYears) ; $i = $i + 2 ) {
430                my $years = fomatYears( $authorsAndYears[ $i + 1 ] );
431                foreach my $year (@$years) {
432                        push( @formatedAuth, $authorsAndYears[$i], $year );
433                }
434        }
435        return ( \@formatedAuth );
436}
437
438sub formatNucleotideMutation{
439        my ($mutation,$postion)       = @_;
440        my $formated="";
441        #6088c-t --> 6088C>T
442        if (length($mutation)>0){
443          if ( $mutation =~ /([\d]+)([actg])-([actg])/i ) {
444                  $formated=$1.$2.">".$3;
445         
446                 
447          }#T-to-C transition at nucleotide 1622 in exon 12
448          elsif($mutation =~ /([actgT])\-to\-([aCctg]).*nucleotide\s([\d]+)/i ){
449                  $formated=$3.$1.">".$2 ;
450                 
451          }elsif($postion=~/^\d+$/){
452                $formated=calculateNuclitideMutation($mutation,$postion); 
453          }
454        }
455        if (length($formated)>0){
456                $formated="c.".$formated;
457        }
458        return ($formated);
459}
460
461sub calculateNuclitideMutation{
462        my ($mutation,$postion)       = @_;
463        my $formated="";
464         if ( $mutation =~ /[actg]{3}[-\s]to[-\s]([actg]{3})/i ) {
465                  $formated=(($postion*3)-2)."_".($postion*3)."delins".$1;
466         
467          }
468        return ($formated);
469}
470
471
472sub getDeceaseID {
473        my ($text) = shift;
474        my $deceaseID = "";
475
476        if ( $text =~ /\([\W\;\S]{0,15}(\d{6})\)/g ) {
477                $deceaseID = $1;
478
479                                $deceaseID =~ s/^see //i;
480        }
481       
482        return ($deceaseID);
483}
484
485#sub formatAAMutation {
486#       my $mutation = @_;
487#       #additional_mutations
488#       ##TODO determ mutation is deletion, substitution of insertion.
489#       if  (! $mutation =~ /\w{3}\d{1,4}\w{w}/ ) {
490#               carp("non mutation\t".$mutation);
491#       }
492#       return ($mutation);
493#}
494
495sub formatInDelMutation {
496        my ( $mutation, $desc ) = @_;
497
498        #create return vraiable
499        my $formatMutation = "";
500
501        #check a postion is availble
502        if ( !( $mutation =~ /,/ ) ) {
503
504                #carp("no postion\n $mutation \n");
505                if ( $desc =~ /loss of nucleotide(s | )(\d{1,5})/i ) {
506                        $mutation = $mutation . "," . $2;
507                }
508        }
509
510        #codon check
511        if ( $mutation =~ /,[\w ]*?(\d+)/ ) {
512                my $codonnumber = $1;
513
514                #               carp("test 1 \n");
515                if ( $desc =~ /codon $codonnumber/ ) {
516
517                        #                       carp("test 2 \n");
518                        my $nuclitdenumber = ( $codonnumber * 3 ) - 2;
519                        $mutation =~ s/$codonnumber/$nuclitdenumber/;
520
521                }
522        }
523
524        #check indels
525        if ( $mutation =~ /(DEL|INS).*(DEL|INS)/ ) {
526
527                #               carp("INDEL 2 \n");
528                my $indelregex =
529                  "(\\d{1,5})-BP\\s*?(DEL|INS)\\/(\\d{1,5})-BP\\s*?(DEL|INS)";
530                if ( $mutation =~ /$indelregex/ ) {
531
532                        #                               carp("INDEL 3 \n");
533
534                        my $del = 0;
535                        my $ins = 0;
536
537                        #derm length of del en ins and asign to right variable
538                        if ( $2 eq "DEL" ) {
539                                $del = $1;
540                                $ins = $3;
541                        }
542                        else {
543                                $del = $3;
544                                $ins = $1;
545                        }
546                        my $mutationclean = $mutation;
547                        $mutationclean =~ s/$indelregex/ /;
548
549                        #                       carp("\n\n $mutationclean \n\n");
550                        if ( my @results =
551                                $mutationclean =~ /([ACTGN]*)(\d{1,5})([ACTGN]*)/i )
552                        {
553
554                                #                                               carp("INDEL 4 \n");
555                                $formatMutation =
556                                    $results[1] . "_"
557                                  . ( ( $results[1] ) + ( $del - 1 ) ) . "del"
558                                  . $results[0] . "ins"
559                                  . $results[2];
560                        }
561                }
562        }
563
564        #               print("'".$mutation."'");
565
566        if ( $mutation =~ /^1-BP\s*?(DEL)\,\s*?(\d{1,5})([A-Z])/i ) {
567
568                $formatMutation = $2 . lc($1) . $3;
569        }
570        elsif ( $mutation =~
571                /^(\d{1,5})-BP\s*?(DEL|INS)\,[\sNT]*?(\d{1,5})([A-Z]*.?)/i )
572        {
573
574                #correction for calculating deletions
575                my $endposition = $1 + $3 - 1;
576                if ( lc($2) eq "ins" ) {
577                        $endposition = $3 + 1;
578                }
579
580                $formatMutation = $3 . "_" . $endposition . lc($2) . $4;
581        }
582        elsif ( $mutation =~ /^(\d{1,5})([A-Z]*?)\s*?(DEL)/i ) {
583                $formatMutation = $1 . lc($3) . $2;
584        }
585        elsif ( $mutation =~ /(INT(\d{0,3})-(\d{0,3}) (DEL))/i ) {
586
587                $formatMutation = "IVS" . $2 . "+?_IVS" . $3 . "+?" . lc($4);
588                $mutation =~ s/INT\d{0,3}-\d{0,3} DEL//;
589                if ( length($mutation) > 5 ) {
590                        $formatMutation .= "," . formatInDelMutation( $mutation, $desc );
591                }
592#               IVS9, G-T,-1
593        }elsif( $mutation =~ /(IVS\d{1,2})[AS]{0,2},\s([AZTG])-([ACTG]),\s*?([\+\-\d]{1,3})/i ){
594               
595                        my ($ivs,$from,$to,$length)=($1,$2,$3,$4);
596                        #267-1G-C
597                        my $regexp="(\\d{1,5})".$length.$from."-".$to;
598       
599                        if($desc =~/$regexp/i){
600               
601                                $ivs=$1;
602                        }
603                        $formatMutation=$ivs.$length.$from.">".$to;
604                       
605        }
606        elsif ( $formatMutation eq "" ) {
607                $formatMutation = $mutation;
608        }
609
610        return ($formatMutation);
611}
612
613sub createAuthorsOutput {
614        my ( $self, $refsref, $authref, $numb ) = @_;
615        my $result      = "";
616        my @authAndYear = @$authref;
617        my ( $author, $year, $pubmedidref ) = "";
618        for ( my $i = 0 ; $i < scalar(@authAndYear) ; $i = $i + 2 ) {
619                $author = $authAndYear[$i];
620                $year   = $authAndYear[ $i + 1 ];
621                $result .= $self->getPubmedID( $refsref, $author, $year, $numb );
622
623        }
624       
625
626        return ($result);
627}
628
629sub proccesAllicVariants {
630        my ( $self, $avsref, $refsref, $numb ) = @_;
631        my @avs    = @$avsref;
632        my $result =   "omimrecord" . "\t"
633                      . "av number" . "\t"
634                  . "title" . "\t"
635                  . "decease id" . "\t"
636                  . "aa mutation". "\t"
637                  . "dna mutation"."\t"
638                 . "formated dna mutation"."\t"
639                  . "DBsnp" . "\t"
640                  . "references" . "\t"
641                  . "nucleotide mutation" . "\t"
642                  . "formated nucleotide mutation"."\n";
643
644
645        #print(Dumper(@avs ));
646        #todo save $avsref so it  can be accesed lateron
647
648        $self->saveAllelicVariants($avsref);
649       
650        foreach my $av (@avs) {
651
652                #print($av->number() . "\t");
653                if ( $av->title() ne "REMOVED FROM DATABASE" ) {
654
655        # check if a record reference with see 123456.1234 reference to other record
656                        my $desc = $av->description();
657
658                        $desc =~ tr/\n/ /;
659
660                        if ( $desc =~ /see (\d{6}).(\d{4})/i ) {
661                                #carp($av->number());
662                                #todo clean  referenced record and append to record
663                                if ( $1 eq $numb ) {
664                                        $desc.=$self->getRefenceAndClean($2);
665                                }
666                                else {
667                                        carp("reference found to $1 . $2 from but is not loaded when in $numb " .$av->number());
668                                }
669                        }
670
671                        #my $mutation = getMutation($av);
672                        my $aamutation=getAAMutation($av);
673                        my $dnamutation=getDNAMutation($av);
674                        my $formateddnamutation=formatdnamutation($dnamutation);
675                        my $authorsAndYears = getAuthor( $av->description() );
676                        $authorsAndYears = formatAuthors($authorsAndYears);
677                        my $articles =
678                          $self->createAuthorsOutput( $refsref, $authorsAndYears, $numb );
679
680                        #lool up allelic variant in dbsnp and return accesionnumber
681                        my $dbsnp = $self->{Omim2DBsnp}->getdbsnp( $numb, $av->number() );
682
683                        my $deceaseID = getDeceaseID($desc);
684
685                        my $nucleotidemutation =
686                          extractnuclitdemutation( $desc, $av->position(),0);
687                        my $formatedNucleotidemutation =formatNucleotideMutation($nucleotidemutation,$av->position());
688
689                        $result .=
690                            $numb . "\t"
691                          . $av->number() . "\t"
692                          . $av->title() . "\t"
693                          . $deceaseID . "\t"
694                          . $aamutation . "\t"
695                          . $dnamutation."\t"
696                        . $formateddnamutation."\t"
697                          . $dbsnp . "\t"
698                          . $articles . "\t"
699                          . $nucleotidemutation . "\t"
700                          . $formatedNucleotidemutation."\n";
701
702                        #                       }
703                }
704        }
705
706        return ($result);
707}
7081;
Note: See TracBrowser for help on using the repository browser.