source: lib/ProccesAllicVariants.pm @ 11

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

fixed latest requests

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