comparison FunAnn/annotate.pl @ 1:40b81942976a draft default tip

Uploaded

0:7e7168ebc150

$fss=13468;

unless (-e "CorGAT")

{
   system("wget -i https://raw.githubusercontent.com/matteo14c/CorGAT_galaxy/dev/ann.txt");
   system("gzip -d GCF_009858895.2_ASM985889v3_genomic.fna.gz");
}

$gen_code="genetic_code";
die("need genetic code file in the current folder\n") unless -e "genetic_code";
open(IN,$gen_code);

{
	($triplet,$oneL)=(split());
	$code{$triplet}=$oneL;
}

$genome="GCF_009858895.2_ASM985889v3_genomic.fna";
die("need reference genome file in the current folder\n") unless -e "GCF_009858895.2_ASM985889v3_genomic.fna";
open(IN,$genome);
while(<IN>)
{
	next if $_=~/^>/;

	{
		$pos=$i+1;
		$res=$seq_res[$i];
	}
}
%AF_data=%{read_simple_table("af_data_new.csv")};
%MFE_data=%{read_simple_table("MFE_annot.csv")};
%epi_data=%{read_epitopes("epitopes_annot.csv")};
%hyphy_data=%{read_hyphy("hyphy_novel.csv")};



$var_File=shift;#""cl7.csv";#"phenetic_indels_sars_cov2.csv";
$out_File=shift;

			die("2\n $triplet b:$Bs[1] r:$ref")unless ($Bs[1] eq $ref);
			$Bs[1]=$alt;
		}elsif ($mod==0){
			$triplet=substr($seq,$pos-3,3);
			@Bs=split('',$triplet);
			die("3\n $triplet b:$Bs[2] r:$ref")unless ($Bs[2] eq $ref);
			$Bs[2]=$alt;
		}
		#print "$pos_inG $relpos $mod @Bs\n";
		$Atriplet=join("",@Bs);
		$A1=$code{$triplet};

				$eff="inframe";
				if ($ref=~/\./)
				{
					$eff.="Ins";
					$Talt=(translate($Calt,\%code))[1];
				}
                                if ($alt=~/\./)
				{
					$eff.="Del";
					$Tref=(translate($Cref,\%code))[1];
				}
				if ($eff=~/Del/)
				{
					@Tref=split('',$Ttref);
					for ($i=0;$i<=$#Tref;$i++)

	{
		chomp();
		($gene,$pos,@annot)=(split(/\,/));
		foreach $a (@annot)
		{
			($ref,$key)=(split(/\:/,$a));
			$ref=~s/{//g;
			$ref=~s/}//g;
			$ref=~s/"//g;
			$key=~s/"//g;
			#print "$gene $pos $ref $key\n";
			if ($kw{$ref})
			{

1:40b81942976a

$fss=13468;

unless (-e "GCF_009858895.2_ASM985889v3_genomic.fna")
{
	
	system("wget -i https://raw.githubusercontent.com/matteo14c/CorGAT_galaxy/master/ann.txt");
	system("gzip -d GCF_009858895.2_ASM985889v3_genomic.fna.gz");
}

$gen_code="genetic_code";
die("need genetic code file in the current folder\n") unless -e "genetic_code";
open(IN,$gen_code);

{
	($triplet,$oneL)=(split());
	$code{$triplet}=$oneL;
}

$genome="GCF_009858895.2_ASM985889v3_genomic.fna ";
die("need reference genome file in the current folder\n") unless -e "GCF_009858895.2_ASM985889v3_genomic.fna";
open(IN,$genome);
while(<IN>)
{
	next if $_=~/^>/;

	{
		$pos=$i+1;
		$res=$seq_res[$i];
	}
}
%AF_data=%{read_simple_table("AF_current.csv")};
%MFE_data=%{read_simple_table("MFE_annot.csv")};
%epi_data=%{read_epitopes("epitopes_annot.csv")};
%hyphy_data=%{read_hyphy("hyphy_current.csv")};



$var_File=shift;#""cl7.csv";#"phenetic_indels_sars_cov2.csv";
$out_File=shift;

			die("2\n $triplet b:$Bs[1] r:$ref")unless ($Bs[1] eq $ref);
			$Bs[1]=$alt;
		}elsif ($mod==0){
			$triplet=substr($seq,$pos-3,3);
			@Bs=split('',$triplet);
			die("$_ 3\n $triplet b:$Bs[2] r:$ref")unless ($Bs[2] eq $ref);
			$Bs[2]=$alt;
		}
		#print "$pos_inG $relpos $mod @Bs\n";
		$Atriplet=join("",@Bs);
		$A1=$code{$triplet};

				$eff="inframe";
				if ($ref=~/\./)
				{
					$eff.="Ins";
					$Talt=(translate($Calt,\%code))[1];
					if ($Cref=~/[ACTG]{1,}/)
                                        {
                                                $Cref=~s/\.//g;
                                                $Tref=(translate($Cref,\%code))[1];
                                        }

				}
                                if ($alt=~/\./)
				{
					$eff.="Del";
					$Tref=(translate($Cref,\%code))[1];
					if ($Calt=~/[ACTG]{1,}/)
                                        {
                                                $Calt=~s/\.//g;
                                                $Talt=(translate($Calt,\%code))[1];
                                        }
				}
				if ($eff=~/Del/)
				{
					@Tref=split('',$Ttref);
					for ($i=0;$i<=$#Tref;$i++)

	{
		chomp();
		($gene,$pos,@annot)=(split(/\,/));
		foreach $a (@annot)
		{
			$a=~s/{//g;
			$a=~s/}//g;
			($ref,$key)=(split(/\:/,$a));
			$ref=~s/"//g;
			$key=~s/"//g;
			#print "$gene $pos $ref $key\n";
			if ($kw{$ref})
			{