Mercurial > repos > mgarnier > pangenome_cog_analysis
changeset 14:ef012cfc41e8 draft
Deleted selected files
| author | mgarnier |
|---|---|
| date | Thu, 19 Aug 2021 13:39:24 +0000 |
| parents | 45cc191a3290 |
| children | ae68b2a86229 |
| files | pangenomeCogAnalysis_V1.pl |
| diffstat | 1 files changed, 0 insertions(+), 649 deletions(-) [+] |
line wrap: on
line diff
--- a/pangenomeCogAnalysis_V1.pl Thu Aug 19 13:39:14 2021 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,649 +0,0 @@ -#!/usr/bin/perl - -use strict; -use warnings; - -my $num_args = $#ARGV + 1; -if ($num_args != 9) { - print "Il n'y a pas le bon nombre d'arguments !\n"; - exit; -} - -# INPUT_ -my $matrix_file = $ARGV[0]; # fichier tabulé : une liste d'orthogroupes qui se retrouvent ou non dans les différentes souches -my $species_file = $ARGV[1]; # association de chaque souche à son espèce (fichier tabulé également) -my $annotation = $ARGV[2]; # collection de fichiers tabulés qui contiennent pour chaque gène la ou les catégories de COG associée(s) -my $order = $ARGV[3]; # cette entrée correspond simplement au nom des souches qui sont rentrées dans le même ordre que les fichiers d'annotation : cela permet de savoir pour un fichier COG à quelle souche et donc plus tard à quelle espèce il correspond -my $annotation_GFF = $ARGV[4]; # fichiers avec les GFF -# my $order_GFF = $ARGV[5]; - -# OUTPUT_ -my $output = $ARGV[5]; # liste des espèces avec leurs orthogroupes (présence-absence) -my $output2 = $ARGV[6]; # fichier des moyennes -my $output3 = $ARGV[7]; # fichier de la liste des valeurs pour chaque catégorie de COG et pour chaque espèce -my $output4 = $ARGV[8]; # fichier avec les catégories de COG pour core-génome / génome accessoire / gènes spé - - -# print "ok\n"; -# exit; - -my @list_gff = split(',', $annotation_GFF); # liste des différents fichiers GFF (qui se retrouvent dans le dossier Annotation Maker) -my %hSpecies = (); # HASH -> key: N_Id (ex NF_AR12) ; val: nom de l'esp (ex Naegleria Fowleri) - -######################## LE SPECIES_FILE ########################### -open (S, $species_file); -while (my $line = <S>){ - - $line =~s/\n//g; $line =~s/\r//g; - my @sp = split('\t', $line); - # print "$line\n"; - # exit; - $hSpecies{$sp[0]} = $sp[1]; # HASH -> key: N_Id ; val: name - -} - my $nbr = keys (%hSpecies); #compter le nombre de souches max - # = taille de la table de hash -# print "J'ai $nbr clés\n"; -# exit; - -close (S); - -#/////////////////////////////////////////////////////////////////////////////////////////////////// - -############################################ LA MATRICE ############################################ - -open(M, $matrix_file); - -my $first_line = <M>; -$first_line =~s/\n//g; $first_line =~s/\r//g; # ne garder que la première ligne du tableau -my @samples = split(/\t/,$first_line); # mettre dans une liste (@samples) chaque intitulé de colonne = N_Id -# print "$first_line\n"; -# exit; - -# Le but ici est de récupérer les combinaisons associées à chaque espèce : NF, NG et NL -my %hCombination =(); # HASH -> key: N_Id ; val: combinaison - -for (my $i=1; $i <= $#samples; $i++){ # on parcourt chaque colonne ($i) mais on ne regarde que le N_Id - my $header = $samples[$i]; # on récupère le N_Id dans $header (soit le nom de la colonne i) - my $species = $hSpecies{$header}; # on regarde dans la table avec N_Id => Nom esp et on attribue à chaque header (qui est ici une clé) sa valeur donc son nom d'esp correspondant - $hCombination{$species} .= "_".$i; # à chaque tour de boucle, pour une $species spé va ajouter le n° de colonne $i pour avoir la combinaison spé à chaque esp - # print "$header\n"; - # exit; -} - - -# foreach my $species (keys (%hCombination)){ -# my $combination = $hCombination{$species}; -# # print "$species $combination\n"; -# } - - -# exit; - -# orthogrp présents : -my %hCombination_prs = (); # HASH -> key: combinaison ; val: liste des orthogroupes -# orthogrp absents : -my %hCombination_abs = (); # idem - - - -my %coregenes = (); # HASH -> key: gene ; val: orthogroupe (pour core-genome) -my %specificgenes = (); # HASH -> key: gene ; val: orthogroupe (pour gènes spécifiques) -my %accessorygenes = (); # HASH -> key: gene ; val: orthogroupe (pour génome accessoire) - -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 - - -while(<M>) { - - my $line = $_; - $line =~s/\n//g; $line =~s/\r//g; - my $nb_found = 0; - my @infos = split(/\t/,$line); - my $orthogroup = $infos[0]; # on récupère le nom de l'orthogroupe dans $orthogroup - my $first_column = $infos[1]; # ici on récupère les gènes de la première colonne qui vont nous servir pour le core-génome - my $combi_prs = ""; - 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) - - $val = $infos[$i]; # on récupère l'information contenue dans la case $i - - if ($val =~/\w/){ # s'il cette cellule contient qq chose... - $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]; - $Genes_of_OG{$i}{$orthogroup} = $premier_gene; # pour chaque orthorgoupe de chaque colonne, on récupère le premier gène - } - - else { # si jamais il n'y a rien dans la cellule... - $combi_abs .= "_".$i; # ... on fait la même chose mais avec $combi_abs - } - - } - - # $hCount{$combi}++; - $hCombination_prs{$combi_prs}.=$orthogroup."\n"; # à la fin de chaque ligne, on va ajouter notre orthogroupe à la combinaison qui lui correspond - $hCombination_abs{$combi_abs}.=$orthogroup."\n"; - - - - if ($nb_found == $#infos){ # si nb_found = au nombre de souche, c'est qu'on a à faire à un core-génome - # print "$orthogroup\n"; - # print "$nb_found\n=================\n"; - for (my $i=1; $i <= $#infos; $i++){ - my @list_of_genes = split (',', $infos[$i]); # ici va séparer tous les gènes (qui se présentent comme une liste, séparés par des ',') - 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 ! - my @list_of_genes = split (',', $gene_random); - my $first_gene = $list_of_genes[0]; - $accessorygenes{$first_gene}= $orthogroup; - } - -} - -my %hCol_Annotated = (); # HASH -> key: colonne ; val: 1 (colonnes pour lesquelles les GFF sont présents) - -# Le but ici est de ne garder que les colonnes (donc les souches) qui ont un fichier GFF associé -my @list_column = split ('\t', $coregene_line); -for (my $i=1; $i <= $#list_column; $i++){ - my @list_genes = split (', ', $list_column[$i]); - my $premier_gene = $list_genes[0]; - my $strain = $samples[$i]; # récupérer le nom de la souche - - - foreach my $gff (@list_gff){ - my $result_grep = `grep $premier_gene $gff`; - - if ($result_grep){ - $hCol_Annotated{$i}=$strain; - - } - # print "$result_grep\n"; - } -} -# exit; -# 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"; -# } -# exit; -# foreach my $oups (keys (%coregenes)) { -# print "$oups\n"; -# } - # exit; - -close (M); - -my %Hash_Specific = (); - -open (OUT, '>', $output) or die $!; -print OUT "$annotation\n"; -foreach my $species (keys (%hCombination)){ # parcours de la table de hash %hCombination (key: nom esp ; val: combi) - my $combination = $hCombination{$species}; # on récupère dans la variable $combination la valeur de chaque clé {species} (= nom esp) de la table de hash %hCombination - my $ortho_presents = $hCombination_prs{$combination}; # $ortho_presents prend la valeur de chaque clé {combination} (récupérée juste au-dessus) de la table de hash %hCombination - my $ortho_absents = $hCombination_abs{$combination}; # en somme on a 3 combi possibles (_1_2_3_4_5 | _6 | _7_8_9) donc pour ces 3 combi-là, qui sont les clés de %hCombination_prs ou_abs, on va retrouver la liste des orthogroupes qui correspondent - - # open (OUT,">results.list.txt"); - - if ($ortho_presents){ - print OUT "> $species - present\n"; - print OUT "$ortho_presents\n"; - my @orthogroups_name = split ('\n', $ortho_presents); - foreach my $ortho (@orthogroups_name){ - $Hash_Specific{$ortho} = $species; - } - } - - if ($ortho_absents){ - # open (OUT2,">$species.$combination.absents.list.txt"); - print OUT "> $species - absent\n"; - print OUT "$ortho_absents\n"; - } - -# close(OUT2); -} - -close(OUT); - -#////////////////////////////////////////////////////////////////////////////////////////////////// - -############################################### COG ############################################### - -# STEP 1 : CORRESPONDANCE ENTRE LES DIFFERENTS FICHIERS DE COG ET L'ORDRE -------------------------------------------- -my @files = split(',', $annotation); # liste des différents fichiers COG (qui se retrouvent dans le dossier Naegleria) -my @list = split(',', $order); # liste de l'ordre des souches -#my ($f,$l); - -my %hCorrespondance = (); #HASH -> key: un fichier COG ; val: un nom de souche (ces 2 données sont entrées en input = $annotation et $order) - -# ++++++++++++ parcours de 2 listes en même temps ++++++++++++ # -my $l = 1; -foreach my $f (@files){ - $hCorrespondance{$f} = $list[$l]; # on fait correspondre pour chaque fichier de COG, un nom de souche - $l++; -} - -# #Affichage du hash -# foreach my $f (keys %hCorrespondance){ -# print $f."=>".$hCorrespondance{$f}."\n" -# } -# exit; - -# STEP 2 : POUR CHAQUE FICHIER DE COG, FAIRE CORRESPONDRE L'ESPECE (ET NON LA SOUCHE) ------------------------------------- -my %hCorresp_file_species = (); # HASH -> key: un fichier de COG ; val: une espèce -my %species_names; # HASH -> key: nom d'espèce ; val: 1 - -foreach my $h (keys (%hCorrespondance)){ # parcours de la table de hash {fichier COG => nom souche} - my $smpl = $hCorrespondance{$h}; # $smpl prend la valeur de la clé (donc d'un nom de souche) - my $espece = $hSpecies{$smpl}; # on regarde la correspondance entre ce $smpl et les nom qu'on a dans notre table de hash %hSpecies (fichier "species.txt") pour avoir le nom de l'espèce dans $espece - $species_names{$espece} = 1; # on garde sous le coude nos nom d'espèce dans cette nouvelle table de hash - $hCorresp_file_species{$h} = $espece; # BUT ATTEINT : on donne pour chaque fichier de COG le nom de l'espèce qui lui correspond -} -# while (my ($k,$v) = each(%hCorresp_file_species)) { -# print "file=$k sp=$v\n"; -# } -# exit; - -# STEP 3 : COMPTAGE DES CATEGORIES DE COG ------------------------------------------------------------------------------ -my %hCount2 = (); # HASH -> key1: catégorie de COG ; key2: espèce associée ; val: comptage - -# comptage du core-genome / des gènes spé / du génome accessoire -my %hCore_Count = (); # HASH -> key: catégorie de COG ; val: comptage (ce hash ne sera utilisé que pour le core-genome) -my %hSpecific_Count = (); # HASH -> key: catégorie de COG ; val: comptage -my %hAccessory_Count = (); # HASH -> key: catégorie de COG ; val: comptage - -# hash pour récupérer le gène -my %hCore_Cat = (); # HASH -> key: catégorie de COG ; val: gène -my %hAccessory_Cat = (); # HASH -> key: catégorie de COG ; val: gène -my %hSpecific_Cat = (); # HASH -> key: catégorie de COG ; val: gène - -# hash pour récupérer le gène -my %hCore_Cat_Esp = (); # HASH -> key1: catégorie de COG ; key2: espèce ; val: gène -my %hAccessory_Cat_Esp = (); # HASH -> key1: catégorie de COG ; key2: espèce ; val: gène -my %hSpecific_Cat_Esp = (); # HASH -> key1: catégorie de COG ; key2: espèce ; val: gène - -my %Cog_of_gene = (); # HASH -> key: gène ; val: cat de COG -my %Specie_of_gene = (); # HASH -> key: gène ; val: souche - -foreach my $file(@files){ # parcours de la liste des fichiers - my $esp = $hCorresp_file_species{$file}; # on récupère l'espèce pour chaque fichier de COG dans $esp - # print $esp."\n"; - # exit; - - my %hCount = (); # HASH -> key: catégorie de COG ; val: comptage - - - open (A, $file); # on va parcourir maintenant chaque fichier un à un - - while (my $line2 = <A>){ - - $line2 =~s/\n//g; $line2 =~s/\r//g; # on procède ligne par ligne - my @Genes = split('\t', $line2); - my $gene = $Genes[0]; - my $first_cat = $Genes[2]; - $Cog_of_gene{$gene} = $first_cat; - - for (my $j=2; $j <= $#Genes; $j++) { - my $cat = $Genes[$j]; # on récupère la ou les catégorie(s) de COG - $hCount{$cat}++; # pour la catégorie donnée, on incrémente son nb d'occurences - - if ($coregenes{$gene}){ # si le $gene fait bien partie du core-genome (donc de notre table de hash %coregenes) - $hCore_Count{$cat}++; # on incrémente le hash - $hCore_Cat{$cat}=$gene; # on récupère le nom du gène - } - if ($accessorygenes{$gene}){ # s'il fait partie des gènes accessoires - $hAccessory_Count{$cat}++; - $hAccessory_Cat{$cat}=$gene; - } - if ($specificgenes{$gene}){ # s'il fait partie des gènes spécifiques - $hSpecific_Count{$cat}++; - $hSpecific_Cat{$cat}=$gene; - } - # $hCount2{$cat}{$esp}++; # TABLE DE HASH AVEC CLES=CAT DE COG + ESPECE VAL=COMPTAGE - } - - } - close (A); - - # print "$file $esp\n=============\n"; - while (my ($k,$v) = each(%hCount)) { # parcours de la table de hash de comptage - # print "cat=$k nb=$v\n"; - $hCount2{$k}{$esp}.= "$v,"; # pour un $k (= une catégorie de COG) on lui associe son espèce et on donne la valeur du comptage qui vient de %hCount - # le but ici est en fait pour une espèce et une catégorie données on veut le nombre d'occurences par souche (pour NF par ex on aura 5 valeurs car il y a 5 souches) - } - - # Récupérer les gènes du core-génome - while (my ($cat_core,$gene_core) = each(%hCore_Cat)) { - $hCore_Cat_Esp{$cat_core}{$esp}=$gene_core; - } - # Récupérer les gènes du génome-accessoire - while (my ($cat_acc,$gene_acc) = each(%hAccessory_Cat)) { - $hAccessory_Cat_Esp{$cat_acc}{$esp}=$gene_acc; - } - # Récupérer les gènes spécifique - while (my ($cat_spe,$gene_spe) = each(%hSpecific_Cat)) { - $hSpecific_Cat_Esp{$cat_spe}{$esp}=$gene_spe; - } - -} -# foreach my $category (sort keys (%hSpecific_Cat_Esp)) { # parcours au niveau de la 1ere clé - -# foreach my $especeee (keys %{$hSpecific_Cat_Esp{$category} }) { # parcours au niveau de la 2e clé pour la $category donnée - -# print "$category\t$especeee\t$hSpecific_Cat_Esp{$category}{$especeee}\n"; # on crée une sortie qui affiche en somme notre hash %hCount2 -# } -# } -# exit; - -# STEP 4 : AFFICHAGE DANS LE FICHIER DE SORTIE ------------------------------------------------------------------------------ -open (OUT4, ">$output4") or die $!; - -print OUT4 "Species"."\t"."COG categories"."\t"."Core-genome"."\t"."Accessory genome"."\t"."Specific genes"."\n"; - -foreach my $category (sort keys (%hCount2)){ # parcours de la table %hCount2 au niveau des catégories - foreach my $especeee (keys %{$hCount2{$category} }){ # parcours de la table %hCount2 au niveau des espèces - print OUT4 "$especeee\t$category\t"; # affichage des esp puis des cat - - # if ($hCore_Cat_Esp{$category}{$especeee}) { - # print OUT4 "$hCore_Cat_Esp{$category}{$especeee}\t"; - # } - my $c = 0; - if ($hCore_Count{$category}){ # si cette catégorie existe dans le core-génome - $c = ($hCore_Count{$category}/scalar keys (%coregenes))*100; # calcul du % du comptage - } - print OUT4 "$c\t"; # affichage du % - - # if ($hAccessory_Cat_Esp{$category}{$especeee}) { - # print OUT4 "$hAccessory_Cat_Esp{$category}{$especeee}\t"; - # } - my $acc = 0; - if ($hAccessory_Count{$category}){ # si cette catégorie existe dans le génome accessoire - $acc = ($hAccessory_Count{$category}/scalar keys (%accessorygenes))*100; # calcul du % du comptage - } - print OUT4 "$acc\t"; # affichage du % - - # if ($hSpecific_Cat_Esp{$category}{$especeee}) { - # print OUT4 "$hSpecific_Cat_Esp{$category}{$especeee}\t"; - # } - my $s = 0; - if ($hSpecific_Count{$category}){ # si cette catégorie existe dans les gènes spécifiques - $s = ($hSpecific_Count{$category}/scalar keys (%specificgenes))*100; # calcul du % du comptage - } - print OUT4 "$s\n"; # affichage du % - } -} -close (OUT4); - -open (OUT3, ">$output3") or die $!; -foreach my $category (sort keys (%hCount2)) { # parcours au niveau de la 1ere clé - - foreach my $especeee (keys %{$hCount2{$category} }) { # parcours au niveau de la 2e clé pour la $category donnée - - print OUT3 "$category\t$especeee\t$hCount2{$category}{$especeee}\n"; # on crée une sortie qui affiche en somme notre hash %hCount2 - } -} - -close (OUT3); - - -open (OUT2, ">$output2") or die $!; - -print OUT2 "category"; -foreach my $e (sort keys (%species_names)){ # on parcours le hash d'espèces... - print OUT2 "\t".$e; #... où on récupère le nom de celles-ci -} -print OUT2 "\n"; - -foreach my $category (sort keys (%hCount2)) { # on parcourt de nouveau les catégories de notre hash à 2 clés - print OUT2 $category; - - foreach my $especes (sort keys (%species_names)) { # on parcourt également le hash d'espèces - - my $nbr = 0; - if ($hCount2{$category}{$especes}) { # si pour une catégorie et une espèce données, on a un nombre : $nbr prend la valeur de ce dernier - $nbr = $hCount2{$category}{$especes}; - } - # $nbr =~s/\n//g; $nbr =~s/\r//g; - - - my @liste = split(',', $nbr); # vu qu'il peut y avoir plusieurs nombres on les dissocie - - my $somme=0; - my $n=0; - my $moyenne=0; - #print "\nma liste de $nbr: ".join("%",@liste)."\n"; - foreach my $x (@liste) { # on parcourt nos nombres - $somme=$somme+$x; - $n=$n+1; - } - - if ($n>0){ - $moyenne = $somme/$n; # on fait le calcul de la moyenne - } - # print "$category, $especes: $hCount2{$category}{$especes}\t"; - # print "moyenne = $moyenne\n=============\n"; - - print OUT2 "\t".$moyenne; # fichier de sortie - } -print OUT2 "\n"; -} - -close (OUT2); - -# foreach my $cat (keys (%hCore_Cat)){ -# print OUT4 $c_gene."\t"; -# } - - -#////////////////////////////////////////////////////////////////////////////////////////////////// - -############################################### GFF ############################################### - - -# my @order_gff = split(',', $order_GFF); # liste de l'ordre des souches -my ($g,$o); - -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 = (); - - -foreach $g (@list_gff){ - # print "$g\n"; - # $hgff_order{$g} = $order_gff[$o++]; # on fait correspondre pour chaque fichier GFF, un nom de souche - open (G, $g); - while (<G>) { - my @table_gff = split (/\t/, $_); - my $chr = $table_gff[0]; - my $start = $table_gff[3]; - my $end = $table_gff[4]; - my $gene_name = $table_gff[8]; - my $type = $table_gff[2]; - - -#or $type eq "CDS" - if ($type && $type eq "mRNA" && $gene_name =~ /ID=([^;]+);/){ - my $gene = $1; - # print $gene."\n"; - # exit; - $hash_of_genes{$gene}=1; - - foreach my $cog (keys (%hCore_Cat)){ - if ($hCore_Cat{$cog} eq $gene){ - $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"; - } - - # foreach my $gene (keys (%hash_of_genes)){ - # my $orthogrp = $hGene_OG{$gene}; - # print "$orthogrp\n"; - # } - } - - close (G); -} - -my %Hash_Convert = ( "A"=>1, "B"=>2, "C"=>3, "D"=>4, "E"=>5, "F"=>6, "G"=>7, "H"=>8, "I"=>9, "J"=>10, "K"=>11, "L"=>12, "M"=>13, "N"=>14, "O"=>15, "P"=>16, "Q"=>17, "R"=>18,"S"=>19, "T"=>20, "U"=>21, "V"=>22, "W"=>23, "X"=>24, "Y"=>25, "Z"=>26, "unknown"=>27); - -mkdir("Core"); -foreach my $i (keys (%coregenes2)){ - - 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 (OUT5, "> Core/$strain_name.$specie_name.txt") or die "Cannot create file $!\n"; - print OUT5 "Orthogroup\tGene\tChromosome\tStart\tEnd\tCOG categories\tNumber assigned\n"; - - my $refcoregenes2 = $coregenes2{$i}; - my %subhash = %$refcoregenes2; - foreach my $gene (keys (%subhash)){ - # print "$gene\n"; - 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 OUT5 $subhash{$gene}."\t"."$gene\t".$Gene_position{$gene}."\t".$cat."\t".$Hash_Convert{$cat}."\n"; - - } - - 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)){ - 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 (OUT6, "> GroupSpecific/$strain_name.$specie_name.txt") or die "Cannot create file $!\n"; - print OUT6 "Orthogroup\tGene\tChromosome\tStart\tEnd\tCOG categories\tNumber assigned\n"; - - my $refGenes_of_OG = $Genes_of_OG{$i}; - my %subhash = %$refGenes_of_OG; - - foreach my $orthogroup (keys (%subhash)){ - if ($Hash_Specific{$orthogroup} && $Hash_Specific{$orthogroup} eq $specie_name){ - my $gene = $subhash{$orthogroup}; - - my $cat = "unknown"; - if ($Cog_of_gene{$gene}){ - $cat = $Cog_of_gene{$gene}; - } - print OUT6 $orthogroup."\t".$subhash{$orthogroup}."\t".$Gene_position{$gene}."\t".$cat."\t".$Hash_Convert{$cat}."\n"; - } - - } - close (OUT6); -}