#!/usr/bin/perl

  # Authors: DEEPAK K GUPTA
  # Company: SBiK-F, Senckenberg.

## This perl script finds effect of nucleotide polymorphism (SNP) on Amino Acid Sequences: Synonymous and NonSynonymous SNP Prediction
## This script takes TAB DELIMITED fasta sequence file and SNP list as given below
use strict;
#use warnings;
my ($in_seq, $snp) = @ARGV; # Input Tab-delimited fasta sequence file and list of SNPs 
my $type;my $type1;my $type2; my $codon1;my $pos_1;my $pos1; my $pos2; my $AA1; my $AA2; my $codon2; my $len1; my $len2;my $cod2; my $var;my $codon3;my $codon4;my $cod3; my $cod4;my $AA3; my $AA4;
open F, "$in_seq" or die $!; # Opening Seqeunce File
do{

  chomp(my $fasta = <F>);
  my @f=split(/\t/, $fasta);
  $f[1] =~ s/\s//g;
  $f[1] =~ s/\n//g;
   open N, "$snp" or die $!;
     while(<N>)
	{
	 chomp;
	 my @s = split(/\t/, $_);
#--------------------------------------
#EXTRACT CODON IN WHICH SNP IS PRESENT
#--------------------------------------
	 if($s[0] =~ /\b$f[0]\b/)
           {
        	$pos1 = int($s[4]/3);
        	$pos2 = $s[4]%3;
        	$pos_1 = $pos1*3;
		$len1 = length($f[1]);
        	if($pos2 != 0)
         	{
       	  	 $codon1 = substr($f[1],$pos_1,3);
         	}
         	
       		elsif($pos2 == 0)
         	{
         	 $codon1 = substr($f[1],($pos_1-3),3);
         	} 
	    $AA1 = codon2aa($codon1);# Translation of Codon into Amino Acid
	    $codon4 = $codon3 = $codon2 = $codon1;
	    my @cod2 = split(//, $codon2);
	    my @cod3 = split(//, $codon3);
 	    my @cod4 = split(//, $codon4);	
	    $len2 = length($s[3]);
	    my @var = split(//, $s[3]);
#-----------------------------------------------------
# CHECK OF EFFECT OF SNP (1ST BASE CALL) ON AMINO ACID
#-----------------------------------------------------
            if ($len2 == 2)
	    {	    
	     if($pos2 == 1)
	 	{
		  if($cod2[0] =~  $var[1])
		    {$cod2[0] =~ s/$cod2[0]/$var[0]/g;}
		  elsif($cod2[0] ==  $var[0])
		    {$cod2[0] =~ s/$cod2[0]/$var[1]/g;}
		}
             elsif($pos2 == 2)
                {
                  if($cod2[1] =~  $var[1])
                    {$cod2[1] =~ s/$cod2[1]/$var[0]/g;}
                  elsif($cod2[1] =~  $var[0])
                    {$cod2[1] =~ s/$cod2[1]/$var[1]/g;}
                }
	     elsif($pos2 == 0)
		{
		  if($cod2[2] =~ $var[1])
                    {$cod2[2] =~ s/$cod2[2]/$var[0]/g;}
                  elsif($cod2[2] =~  $var[0])
                    {$cod2[2] =~ s/$cod2[2]/$var[1]/g;}
		}
           }
#-------------------------------------------------------------
# CHECK OF EFFECT OF SNP (1ST AND 2ND BASE CALL) ON AMINO ACID
#-------------------------------------------------------------
	   if($len2 == 3)
	   { 	
	    if($pos2 == 1)
                {
                  if($cod2[0] =~  $var[0])
                    {
			$cod2[0] =~ s/$cod2[0]/$var[1]/g;
		        $cod3[0] =~ s/$cod3[0]/$var[2]/g;
                    }
                  elsif($cod2[0] =~ $var[1])
                    {
			$cod2[0] =~ s/$cod2[0]/$var[0]/g;
		    	$cod3[0] =~ s/$cod3[0]/$var[2]/g;
		    }
        	  elsif($cod2[0] =~  $var[2])
                    {
                        $cod2[0] =~ s/$cod2[0]/$var[1]/g;
                        $cod3[0] =~ s/$cod3[0]/$var[0]/g;
                    }
	        }
            elsif($pos2 == 2)
                {
  		  if($cod2[1] =~  $var[0])
                    {
                        $cod2[1] =~ s/$cod2[1]/$var[1]/g;
                        $cod3[1] =~ s/$cod3[1]/$var[2]/g;
                    }
                  elsif($cod2[1] =~  $var[1])
                    {
                        $cod2[1] =~ s/$cod2[1]/$var[0]/g;
                        $cod3[1] =~ s/$cod3[1]/$var[2]/g;
                    }
                  elsif($cod2[1] =~ $var[2])
                    {
                        $cod2[1] =~ s/$cod2[1]/$var[1]/g;
                        $cod3[1] =~ s/$cod3[1]/$var[0]/g;
                    }
                }
            elsif($pos2 == 0)
	        {
                  if($cod2[2] =~  $var[0])
                    {
                        $cod2[2] =~ s/$cod2[2]/$var[1]/g;
                        $cod3[2] =~ s/$cod3[2]/$var[2]/g;
                    }
                  elsif($cod2[2] =~  $var[1])
                    {
                        $cod2[2] =~ s/$cod2[2]/$var[0]/g;
                        $cod3[2] =~ s/$cod3[2]/$var[2]/g;
                    }
                  elsif($cod2[2] =~  $var[2])
                    {
                        $cod2[2] =~ s/$cod2[2]/$var[1]/g;
                        $cod3[2] =~ s/$cod3[2]/$var[0]/g;
                    }
                }
             }
#------------------------------------------------------------------
# CHECK OF EFFECT OF SNP (1ST,2ND AND 3RD BASE CALL) ON AMINO ACID
#------------------------------------------------------------------

	    if($len2 == 4)
	    { 
              if($pos2 == 1)
                {
                  if($cod2[0] =~  $var[0])
                    {
                        $cod2[0] =~ s/$cod2[0]/$var[1]/g;
                        $cod3[0] =~ s/$cod3[0]/$var[2]/g;
                        $cod4[0] =~ s/$cod4[0]/$var[3]/g;
                    }
                  elsif($cod2[0] =~ $var[1])
                    {
                        $cod2[0] =~ s/$cod2[0]/$var[0]/g;
                        $cod3[0] =~ s/$cod3[0]/$var[2]/g;
                        $cod4[0] =~ s/$cod4[0]/$var[3]/g;
                    }
                  elsif($cod2[0] =~  $var[2])
                    {
                        $cod2[0] =~ s/$cod2[0]/$var[1]/g;
                        $cod3[0] =~ s/$cod3[0]/$var[0]/g;
                        $cod4[0] =~ s/$cod4[0]/$var[3]/g;
                    }
                  elsif($cod2[0] =~  $var[3])
                    {
                        $cod2[0] =~ s/$cod2[0]/$var[1]/g;
                        $cod3[0] =~ s/$cod3[0]/$var[0]/g;
                        $cod4[0] =~ s/$cod4[0]/$var[2]/g;
                    }
                }
            elsif($pos2 == 2)
                {
                  if($cod2[1] =~  $var[0])
                    {
                        $cod2[1] =~ s/$cod2[1]/$var[1]/g;
                        $cod3[1] =~ s/$cod3[1]/$var[2]/g;
                        $cod4[1] =~ s/$cod4[1]/$var[3]/g;
                    }
                  elsif($cod2[1] =~  $var[1])
                    {
                        $cod2[1] =~ s/$cod2[1]/$var[0]/g;
                        $cod3[1] =~ s/$cod3[1]/$var[2]/g;
                        $cod4[1] =~ s/$cod4[1]/$var[3]/g;
                    }
                  elsif($cod2[1] =~ $var[2])
                    {
                        $cod2[1] =~ s/$cod2[1]/$var[1]/g;
                        $cod3[1] =~ s/$cod3[1]/$var[0]/g;
                        $cod4[1] =~ s/$cod4[1]/$var[3]/g;
                    }
                  elsif($cod2[1] =~ $var[3])
                    {
                        $cod2[1] =~ s/$cod2[1]/$var[1]/g;
                        $cod3[1] =~ s/$cod3[1]/$var[0]/g;
                        $cod4[1] =~ s/$cod4[1]/$var[2]/g;
                    }
	        }
            elsif($pos2 == 0)
                {
                  if($cod2[2] =~  $var[0])
                    {
                        $cod2[2] =~ s/$cod2[2]/$var[1]/g;
                        $cod3[2] =~ s/$cod3[2]/$var[2]/g;
                        $cod4[2] =~ s/$cod4[2]/$var[3]/g;
                    }
                  elsif($cod2[2] =~  $var[1])
                    {
                        $cod2[2] =~ s/$cod2[2]/$var[0]/g;
                        $cod3[2] =~ s/$cod3[2]/$var[2]/g;
                        $cod4[2] =~ s/$cod4[2]/$var[3]/g;
                    }
                  elsif($cod2[2] =~  $var[2])
                    {
                        $cod2[2] =~ s/$cod2[2]/$var[1]/g;
                        $cod3[2] =~ s/$cod3[2]/$var[0]/g;
                        $cod4[2] =~ s/$cod4[2]/$var[3]/g;
                    }
	          elsif($cod2[2] =~  $var[3])
                    {
                        $cod2[2] =~ s/$cod2[2]/$var[1]/g;
                        $cod3[2] =~ s/$cod3[2]/$var[0]/g;
                        $cod4[2] =~ s/$cod4[2]/$var[2]/g;
                    }
                }
                }
	
		$codon2 = join(undef,@cod2);
		$codon3 = join(undef,@cod3);
		$codon4 = join(undef,@cod4);
		$AA2 = codon2aa($codon2);
		$AA3 = codon2aa($codon3);
		$AA4 = codon2aa($codon4);
		if($AA1 =~ $AA2)
		{$type = 'SYNONYMOUS';}
		else
		{$type = 'NON-SYNONYMOUS';}
		if($AA1 =~ $AA3)
                {$type1 = 'SYNONYMOUS';}
                else
                {$type1 = 'NON-SYNONYMOUS';}
		if($AA1 =~ $AA4)
                {$type2 = 'SYNONYMOUS';}
                else
                {$type2 = 'NON-SYNONYMOUS';}

#----------------------
# PRINTING THE RESULTS
#----------------------
	   if($len2 == 2){print "$_\t$pos_1,$pos2\t$codon1\t$AA1\t$codon2\t$AA2\t$type\t--\t--\t--\t--\t--\t--\n";}
	elsif($len2 == 3){print "$_\t$pos_1,$pos2\t$codon1\t$AA1\t$codon2\t$AA2\t$type\t$codon3\t$AA3\t$type1\t--\t--\t--\n";}
	elsif($len2 == 4){print "$_\t$pos_1,$pos2\t$codon1\t$AA1\t$codon2\t$AA2\t$type\t$codon3\t$AA3\t$type1\t$codon4\t$AA4\t$type2\n";}
	   }
	} close N;

}until eof(F);
 close F;
#-----------------------------
# SUBROUTINE FOR GENETIC CODES
#-----------------------------
sub codon2aa {
    my($codon) = @_;

    $codon = uc $codon;
 
    my(%genetic_code) = (
    
    'TCA' => 'S',    # Serine
    'TCC' => 'S',    # Serine
    'TCG' => 'S',    # Serine
    'TCT' => 'S',    # Serine
    'TTC' => 'F',    # Phenylalanine
    'TTT' => 'F',    # Phenylalanine
    'TTA' => 'L',    # Leucine
    'TTG' => 'L',    # Leucine
    'TAC' => 'Y',    # Tyrosine
    'TAT' => 'Y',    # Tyrosine
    'TAA' => '_',    # Stop
    'TAG' => '_',    # Stop
    'TGC' => 'C',    # Cysteine
    'TGT' => 'C',    # Cysteine
    'TGA' => '_',    # Stop
    'TGG' => 'W',    # Tryptophan
    'CTA' => 'L',    # Leucine
    'CTC' => 'L',    # Leucine
    'CTG' => 'L',    # Leucine
    'CTT' => 'L',    # Leucine
    'CCA' => 'P',    # Proline
    'CCC' => 'P',    # Proline
    'CCG' => 'P',    # Proline
    'CCT' => 'P',    # Proline
    'CAC' => 'H',    # Histidine
    'CAT' => 'H',    # Histidine
    'CAA' => 'Q',    # Glutamine
    'CAG' => 'Q',    # Glutamine
    'CGA' => 'R',    # Arginine
    'CGC' => 'R',    # Arginine
    'CGG' => 'R',    # Arginine
    'CGT' => 'R',    # Arginine
    'ATA' => 'I',    # Isoleucine
    'ATC' => 'I',    # Isoleucine
    'ATT' => 'I',    # Isoleucine
    'ATG' => 'M',    # Methionine
    'ACA' => 'T',    # Threonine
    'ACC' => 'T',    # Threonine
    'ACG' => 'T',    # Threonine
    'ACT' => 'T',    # Threonine
    'AAC' => 'N',    # Asparagine
    'AAT' => 'N',    # Asparagine
    'AAA' => 'K',    # Lysine
    'AAG' => 'K',    # Lysine
    'AGC' => 'S',    # Serine
    'AGT' => 'S',    # Serine
    'AGA' => 'R',    # Arginine
    'AGG' => 'R',    # Arginine
    'GTA' => 'V',    # Valine
    'GTC' => 'V',    # Valine
    'GTG' => 'V',    # Valine
    'GTT' => 'V',    # Valine
    'GCA' => 'A',    # Alanine
    'GCC' => 'A',    # Alanine
    'GCG' => 'A',    # Alanine
    'GCT' => 'A',    # Alanine
    'GAC' => 'D',    # Aspartic Acid
    'GAT' => 'D',    # Aspartic Acid
    'GAA' => 'E',    # Glutamic Acid
    'GAG' => 'E',    # Glutamic Acid
    'GGA' => 'G',    # Glycine
    'GGC' => 'G',    # Glycine
    'GGG' => 'G',    # Glycine
    'GGT' => 'G',    # Glycine
    );

    if(exists $genetic_code{$codon}) {
        return $genetic_code{$codon};
    }else{

            print STDERR "Bad codon \"$codon\"!!\n";
            
    }
}



# Example SNP file:
#FUN_00578       #       #       AG      30      +
#FUN_00578       #       #       C       32      +
#FUN_00578       #       #       C       58      +
