Mercurial > repos > mgarnier > pangenome_cog_analysis
changeset 12:b36506d26a43 draft
Uploaded
author | mgarnier |
---|---|
date | Tue, 06 Jul 2021 15:51:16 +0000 |
parents | 9ac1a58fda89 |
children | 45cc191a3290 |
files | pangenomeCogAnalysis_V1.pl |
diffstat | 1 files changed, 66 insertions(+), 11 deletions(-) [+] |
line wrap: on
line diff
--- a/pangenomeCogAnalysis_V1.pl Tue Jul 06 15:51:02 2021 +0000 +++ b/pangenomeCogAnalysis_V1.pl Tue Jul 06 15:51:16 2021 +0000 @@ -93,6 +93,7 @@ my $coregene_line; my %coregenes2 = (); # HASH -> key1: colonne i ; key2: gène ; val: orthogroupe +my %specificgenes2 = (); # HASH -> key1: colonne i ; key2: gène ; val: orthogroupe my %Genes_of_OG = (); # HASH -> key1: orthogroupe ; key2: colonne i ; val: gène @@ -109,6 +110,7 @@ my $combi_abs = ""; my $val; my $gene_random; + my $unique_col_detected; for (my $i=1; $i <= $#infos; $i++){ # on travaille par ligne puis dans chaque ligne (while(<M>)), cellule par cellule (cette boucle for) @@ -119,6 +121,7 @@ $combi_prs .= "_".$i; # ...on va concaténer notre chaine $combi_prs pour que cela forme une combinaison $nb_found++; # on incrémente le compteur qui permet de savoir cb de fois notre orthogroupe est présent (le but sera de l'utiliser quand nb_found == 9) $gene_random=$val; # on récupère la valeur de la case (les gènes) + $unique_col_detected = $i; my @table_genes = split (',', $val); my $premier_gene = $table_genes[0]; @@ -145,19 +148,26 @@ my $first_gene = $list_of_genes[0]; # prend la valeur du premier gène uniquement ! $coregenes{$first_gene}= $orthogroup; # on va récupérer ce premier gène qu'on met dans un hash (pour y avoir accès facilement, d'où val = 1, ici ça n'a pas d'importance) $coregenes2{$i}{$first_gene}= $orthogroup; - - - - + } if (!$coregene_line){ $coregene_line = $line; } } elsif ($nb_found == 1) { # si on a un gène spé + + # print "$gene_random\n"; + # print "$line\n"; + # print "$unique_col_detected\n"; + + + my @list_of_genes = split (',', $gene_random); # idem, on ne veut qu'un seul gène donc on crée la liste my $first_gene = $list_of_genes[0]; # on ne prend que le premier + # print "$first_gene\n"; + # exit; $specificgenes{$first_gene}= $orthogroup; # et pareil on crée la table de hash + $specificgenes2{$unique_col_detected}{$first_gene}= $orthogroup; } else { # là c'est le génome accessoire, i.e tout le reste ! @@ -189,12 +199,12 @@ } } # exit; -# foreach my $i (sort keys (%coregenes2)){ # parcours de la table %hCount2 au niveau des catégories -# foreach my $gene (keys %{$coregenes2{$i} }){ # parcours de la table %hCount2 au niveau des espèces -# print "$i\t$gene\t".$coregenes2{$i}{$gene}."\n"; +# foreach my $i (sort keys (%specificgenes2)){ # parcours de la table %hCount2 au niveau des catégories +# foreach my $gene (keys %{$specificgenes2{$i} }){ # parcours de la table %hCount2 au niveau des espèces +# print "$i\t$gene\t".$specificgenes2{$i}{$gene}."\n"; # } # } - +# exit; # while (my ($k,$v) = each(%strain_specie)) { # print "i=$k strain=$v\n"; # } @@ -478,6 +488,7 @@ my %hgff_order = (); #HASH -> key: un fichier GFF ; val: un nom de souche (ces 2 données sont entrées en input = $annotation_GFF et $order_GFF) my %Gene_position = (); my %Cat_genes = (); +my %Cat_genes2 = (); my %hash_of_genes = (); @@ -495,8 +506,8 @@ my $type = $table_gff[2]; - - if ($type && ($type eq "mRNA" or $type eq "CDS") && $gene_name =~ /ID=([^;]+);/){ +#or $type eq "CDS" + if ($type && $type eq "mRNA" && $gene_name =~ /ID=([^;]+);/){ my $gene = $1; # print $gene."\n"; # exit; @@ -507,7 +518,12 @@ $Cat_genes{$gene}=$cog; } } - + foreach my $cog_bis (keys (%hSpecific_Cat)){ + if ($hSpecific_Cat{$cog_bis} eq $gene){ + $Cat_genes2{$gene}=$cog_bis; + } + } + $Gene_position{$gene}="$chr\t$start\t$end"; } @@ -561,6 +577,45 @@ close (OUT5); } +mkdir("StrainSpecific"); +foreach my $i (keys (%specificgenes2)){ + + if (!$hCol_Annotated{$i}) { # si le fichier GFF n'existe pas + next; + } + + + my $strain_name = $hCol_Annotated{$i}; + + my $specie_name = $hSpecies{$strain_name}; + + + + open (OUT7, "> StrainSpecific/$strain_name.$specie_name.txt") or die "Cannot create file $!\n"; + print OUT7 "Orthogroup\tGene\tChromosome\tStart\tEnd\tCOG categories\tNumber assigned\n"; + + my $refspecificgenes2 = $specificgenes2{$i}; + my %subhash = %$refspecificgenes2; + foreach my $gene (keys (%subhash)){ + # print "$gene\n"; exit; + my $cat = "unknown"; + if ($Cog_of_gene{$gene}){ + $cat = $Cog_of_gene{$gene}; + } + # if (!$Gene_position{$gene}){ + # print "$gene\n coucou"; exit; + # } + + # if (!$subhash{$gene}){ + # print "$gene\n"; + # } + print OUT7 $subhash{$gene}."\t"."$gene\t".$Gene_position{$gene}."\t".$cat."\t".$Hash_Convert{$cat}."\n"; + + } + + close (OUT7); +} + mkdir("GroupSpecific"); foreach my $i (keys (%Genes_of_OG)){