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)){