Recherche:Répétition des bases dans l'ADN des procaryotes/Annexe/Perl
Compilation des codons
modifier# ********** Codon.pl *********** # #!/usr/bin/env perl #use strict; #use warnings; #use v5.22; ## Décompte des codons des gènes de protéines. Documentation, Paris le 9.9.16 # # Ce programme décompte les 64 codons d'un gène d'une protéine d'un organisme donné. Les décomptes de plusieurs gènes du même organisme # peuvent être cumulés etle résultat affiche une ligne pour cet organisme. # − cumul: Il a éte écrit au début pour afficher les cumuls de plusieurs organismes pour étudier l'évolution de la fréquence d'un codon # donné en fonction du contenu en GC (%GC) de l'ADN des procaryotes. # − décompte: Il est surtout utilisé pour étudier l'évolution de la fréquence d'un codon donné d'un seul gène en fonction du contenu en # GC (%GC) de l'ADN des procaryotes ou de celui du gène lui-même. # fichiers en entrée: xxxppp d'après la liste de fichiers @x sous forme d'un tableau de tableaux pour pouvoir faire des cumuls sur les # protéines si nécessaire. # - En lecture, longueur des lignes fixe de 60 caractères minuscules,format KEGG. Si majuscules pensez à basuler le tableau @cod en # majuscules sous writer ou word. Copier la séquence de nucléotides "NT seq". # - xxx, nom codé KEGG du procaryote. # - ppp, nom de la protéine sous format EC (Enzyme Commission numbers) sans les séparations par point. Dans le cas de non-enzyme le nom # est empreinté à KEGG dans son entrée "KEGG orthology" que je signale dans les références. # - pour récupérer le fichier dans KEGG il suffit de rechercher le code EC ( avec points ) dans un moteur de recherche et prendre une # adresse KEGG. # + rechercher chaque procaryote par son code KEGG sous la forme 'xxx:'. Il peut y avoir plusieurs gènes pour le même enzyme. # Voir les références suivantes (le tableau des références est dans l'article wikipédia): # + fr.wikiversity.org/wiki/Recherche:Répétition_des_bases_dans_l'ADN_des_procaryotes. # + la recherche par code "orthology" dans un moteur de recherche est plus simple parce qu'on accède directement aux gènes # après avoir choisi le procaryote. # fichier en sortie, résultat du comptage des codons: codonsc.txt (création et mise à jour) contient les tableaux "1 gène/n procaryotes" # et "1 procaryote/ n gènes". # fichier en sortie, sommation des codons en acides aminés: codonsp.txt (création et mise à jour), tableaux "1 gène/n procaryotes" # et "1 procaryote/ n gènes". # # Note: pour la programmatioon Perl #$x0=@{$x[0]}; @z0=@{$x[0]}; # $x0=@z0 équivalent à $x0=@{$x[0]} #$x1=@{$x[1]}; @z1=@{$x[1]}; # @z0=@{$x[0]} est la 1ère sous-liste #$x2=@{$x[2]}; @z2=@{$x[2]}; # # liste avec cumuls: #my @x=([a1653, a1653m, a1653n], [a2777A, a27771, a27761, a27762, a63551, a63552], [a6114, a6115, a27140, aseca, a4213]); # liste décompte, sans cumuls: #my @x=([a1653],[a1653m],[a1653n],[a2777A],[a27771],[a27761],[a27762],[a63551],[a63552],[a6114],[a6115],[a27140],[aseca],[a4213]); # my @x=([a1653],[a1653m],[a1653n],[a2777A],[a27771],[a27761],[a27762],[a63551],[a63552],[a6114],[a6115],[a27140],[aseca],[a4213]); $t=@x; $t1=$t-1; @noms=(age, amd, roa, saci, sma, sgr, sho, sus ); # liste des organismes $n=@noms; # nombre d'organismes @u = ([[(0)x64]x$t]x$n); # **************************************************** my $filenamer = "codonsc.txt"; # Sous Linux open ($fhr, '>>', $filenamer) or die "Impossible d'ouvrir le fichier $filenamer en ajout"; my $filenamera = "codonsa.txt"; # Sous Linux open ($fhra, '>>', $filenamera) or die "Impossible d'ouvrir le fichier $filenamera en ajout"; # **************************************************** # foreach $r (0..$t1) { @y=@{$x[$r]}; # liste de plusieurs gène préfixés par a, à traiter un à un. protein(@y);} # traiter un seul gène pour plusieurs organismes à partir d'une liste de plusieurs gènes ptprot(@u, @x); # 1 seul organisme pour plusieurs gènes à partir d'une liste d'organisme close $fhr; close $fhra; # sub protein { # *** sub protein *** @z=@y; # # *** Tableau des codons ****** # liste des organismes @cod = (tta,ttg,ctt,cta,ctg,ctc, # Leu 6 codons Leu 0..5 aaa,aag,ttt,ttc, # Lys Phe 2 codons Lys et 2 codons Phe 6..9 tct,tca,agt,tcg,tcc,agc, # Ser 6 codons Ser 10..15 tat,tac,aat,aac, # Tyr Asn 2 codons Tyr et 2 Codons Asn 16..19 aga,cga,cgt,agg,cgc,cgg, # Arg 6 codons Arg 20..25 cat,cac,caa,cag, # His Gln 2 codons His et 2 codons Gln 26..29 gga,ggt,ggc,ggg, # Gly 4 codons Gly 30..33 gaa,gag,gat,gac, # Glu Asp 2 codons Glu et 2 codons Asp 34..37 gca,gct,gcg,gcc, # Ala 4 codons Ala 38..41 tgt,tgc,atg, # Cys Met 2 codons Cys et 1 codon Met 42..44 cca,cct,ccg,ccc, # Pro 4 codons Pro 45..48 tgg, # Trp 1 codon Trp 49 gtt,gta,gtg,gtc, # Val 4 codons Val 50..53 act,aca,acg,acc, # Thr 4 codons Thr 54..57 att,ata,atc, # Ile 3 codons Ile 58..60 taa,tag,tga); # Stp 3 codons Stp 61..64 # $m=@z; # nombre de gènes $pt=$m-1; $nt=$n-1; @q = ([(0)x64]x$n); # tableau des compteurs de 64 codons par organisme @fic = (a)x$m; # tableau intermédiaire des noms d'organismes $i=0; $j=0; @p = (0)x64; # tableau de 64 compteurs intermédiaires # tableau des compteurs de 64 codons par organisme foreach $l (0..$nt) { # traitement de chaque organisme (age, amd, roa, saci, sma, sgr, sho, sus) @p = (0)x64; # foreach $k (0..$pt) { # constitution des noms de ses fichiers de gènes présents ou non pour les cumuls $lz=length($z[$k])-1; $zi=substr($z[$k],1,$lz); # sauter le préfixe du nom $fic[$k]=$noms[$l].$zi; # concaténation organisme..gène: les fichiers KEGG de 60 cractères du gène de cet organisme, # exemple ade2777A, organisme ade gène 2777A. } # la liste des noms des fichiers de l'organisme est constituée. foreach $k (0..$pt) { comptfich($fic[$k]);} # traitement de la liste des fichiers de l'organisme sub comptfich { # routine de décompte d'un fichier comme age1653 $filename=$fic[$k]; # récupération du nom du fichier open ($fh, '<', $filename) or die "Impossible d'ouvrir le fichier $filename en lecture"; while( defined( $lis = <$fh> ) ) { # traitement ligne par ligne jusqu'à la fin chomp $lis; while ($i <= 60) { $lu = substr($lis,$i,3 ); # un codon foreach $j (0..63) { # traitement du codon lu if( $lu eq $cod[$j] ) { $p[$j]+=1; } # +1 dans son compteur } $i=$i+3; } # avancer de 3 caractères et relire $i = 0; } # fin de la ligne, remise à zéro du compteur de cractères close $fh; # # fin de la routine } # fin d'un organisme @{$q[$l]}=@p; } # sauvegarde des compteurs intermédiaires dans les compteurs affectés à l'organisme dans @{$u[$r]}=@q; # le tableau général @q. Organisme suivant. # ******* écriture de l'en-tête codons ******* # foreach $k (0..$pt) { # écriture des noms des gènes traités print $fhr ($z[$k],";"); } print $fhr ("\n"); print $fhr ( "Génome;"); # organisme foreach $i (0..9) { print $fhr ($cod[$i],";");} # Leu, Lys, Phe. # foreach $i (10..19) { print $fhr ($cod[$i],";");} # Ser, Tyr, Asn # foreach $i (20..29) { print $fhr ($cod[$i],";");} # Arg, His, Gln # foreach $i (30..37) { print $fhr ($cod[$i],";");} # Gly, Glu, Asp # foreach $i (38..44) { print $fhr ($cod[$i],";");} # Ala, Cys, Met print $fhr ($cod[49],";"); # Trp # foreach $i (45..48) { print $fhr ($cod[$i],";");} # Pro # foreach $i (50..53) { print $fhr ($cod[$i],";");} # Val # foreach $i (54..57) { print $fhr ($cod[$i],";");} # Thr # foreach $i (58..60) { print $fhr ($cod[$i],";");} # Ile # foreach $i (61..63) { print $fhr ($cod[$i],";"); } # Stp # print $fhr ("\n"); # # # ******* écriture des compteurs codons de tous les organismes ******* # foreach $i (0..$nt) { print $fhr ( $noms[$i],";"); # nom de l'organisme @w=@{$q[$i]}; # récupération des compteurs de l'organisme foreach $k (0..9) { # compteurs des codons de Leu, Lys, Phe print $fhr ( $w[$k],";");} # foreach $k (10..19) { # compteurs des codons de Ser, Tyr, Asn print $fhr ( $w[$k],";");} # foreach $k (20..29) { # compteurs des codons de Arg, His, Gln print $fhr ( $w[$k],";");} # foreach $k (30..37) { # compteurs des codons de Gly, Glu, Asp print $fhr ( $w[$k],";");} # foreach $k (38..44) { # compteurs des codons de Ala, Cys, Met, Trp print $fhr ( $w[$k],";");} print $fhr ( $w[49],";"); # foreach $k (45..48) { # compteurs des codons de Pro print $fhr ( $w[$k],";");} # foreach $k (50..53) { # compteurs des codons de Val print $fhr ( $w[$k],";");} # foreach $k (54..57) { # compteurs des codons de Thr print $fhr ( $w[$k],";");} # foreach $k (58..60) { # compteurs des codons de Ile print $fhr ( $w[$k],";");} # foreach $k (61..63) { # compteurs des codons de Stp print $fhr ( $w[$k],";");} # print $fhr ("\n"); } # # ******* écriture de l'en-tête aas ******* # foreach $k (0..$pt) { # écriture des noms des gènes traités print $fhra ($z[$k],";"); } print $fhra ("\n"); print $fhra ( "Génome;"); # organisme print $fhra ( "Leu;"); # Leu print $fhra ( "Lys;Phe;Ser;"); # Lys, Phe, Ser print $fhra ( "Tyr;Asn;Arg;"); # Tyr, Asn, Arg print $fhra ( "His;Gln;Gly;"); # His, Gln, Gly print $fhra ( "Glu;Asp;Ala;"); # Glu, Asp, Ala print $fhra ( "Cys;Met;Trp;Pro;"); # Cys, Met, Trp, Pro print $fhra ( "Val;"); # Val print $fhra ( "Thr;"); # Thr print $fhra ( "Ile;"); # Ile print $fhra ( "Stp;"); # Stp print $fhra ("\n"); # # # ******* écriture des compteurs des aas de tous les organismes ******* # foreach $i (0..$nt) { print $fhra ( $noms[$i],";"); # nom de l'organisme @w=@{$q[$i]}; # récupération des compteurs de l'organisme #foreach $k (0..9) $leu=$w[0]+$w[1]+$w[2]+$w[3]+$w[4]+$w[5]; print $fhra ($leu,";"); # Leu total $lys=$w[6]+$w[7]; print $fhra ($lys,";"); # Lys total $phe=$w[8]+$w[9]; print $fhra ($phe,";"); # Phe total #foreach $k (10..19) $ser=$w[10]+$w[11]+$w[12]+$w[13]+$w[14]+$w[15]; print $fhra ($ser,";"); # Ser total $tyr=$w[16]+$w[17]; print $fhra ($tyr,";"); # Tyr total $asn=$w[18]+$w[19]; print $fhra ($asn,";"); # Asn total #foreach $k (20..29) $arg=$w[20]+$w[21]+$w[22]+$w[23]+$w[24]+$w[25]; print $fhra ($arg,";"); # Arg total $his=$w[26]+$w[27]; print $fhra ($his,";"); # His total $gln=$w[28]+$w[29]; print $fhra ($gln,";"); # Gln total #foreach $k (30..37) $gly=$w[30]+$w[31]+$w[32]+$w[33]; print $fhra ($gly,";"); # Gly total $glu=$w[34]+$w[35]; print $fhra ($glu,";"); # Glu total $asp=$w[36]+$w[37]; print $fhra ($asp,";"); # Asp total #foreach $k (38..44) $ala=$w[38]+$w[39]+$w[40]+$w[41]; print $fhra ($ala,";"); # Ala total $cys=$w[42]+$w[43]; print $fhra ($cys,";"); # Cys total $met=$w[44]; print $fhra ($met,";"); # Met $trp=$w[49]; print $fhra ($trp,";"); # Trp #foreach $k (45..48) $pro=$w[45]+$w[46]+$w[47]+$w[48]; print $fhra ($pro,";"); # Pro total #foreach $k (50..53) $val=$w[50]+$w[51]+$w[52]+$w[53]; print $fhra ($val,";"); # Val total #foreach $k (54..57) $thr=$w[54]+$w[55]+$w[56]+$w[57]; print $fhra ($thr,";"); # Thr total #foreach $k (58..60) $ile=$w[58]+$w[59]+$w[60]; print $fhra ($ile,";"); # Ile total #foreach $k (61..63) $stp=$w[61]+$w[62]+$w[63]; print $fhra ($stp,";"); # Stp total print $fhra ("\n"); } # } # *** sub protein *** # # ***sub ptprot**** écriture 1 tableau par organisme pour n organismes ********sub ptprot*********** # sub ptprot { $n1=$n-1; $t=@x; $t1=$t-1; # *** sub ptprot *** # foreach $h (0..$n1) { # 1 organisme codons ptprot print $fhr ($noms[$h],";"); # nom print $fhr ("\n"); # ******* écriture de l'en-tête codons ******* print $fhr ( "Gène;"); #Gène de protéine foreach $i (0..9) { print $fhr ($cod[$i],";");} # Leu, Lys, Phe. foreach $i (10..19) { print $fhr ($cod[$i],";");} # Ser, Tyr, Asn foreach $i (20..29) { print $fhr ($cod[$i],";");} # Arg, His, Gln foreach $i (30..37) { print $fhr ($cod[$i],";");} # Gly, Glu, Asp foreach $i (38..44) { print $fhr ($cod[$i],";");} # Ala, Cys, Met print $fhr ($cod[49],";"); # Trp foreach $i (45..48) { print $fhr ($cod[$i],";");} # Pro foreach $i (50..53) { print $fhr ($cod[$i],";");} # Val foreach $i (54..57) { print $fhr ($cod[$i],";");} # Thr foreach $i (58..60) { print $fhr ($cod[$i],";");} # Ile foreach $i (61..63) { print $fhr ($cod[$i],";"); } # Stp # print $fhr ("\n"); # ******* écriture des compteurs de toutes les protéines codons******* # foreach $v (0..$t1) { $xl=@{$x[$v]}; $lot=$x[$v][0]; $xl1=$xl-1; foreach $j (1..$xl1) { $lot=$lot.",".$x[$v][$j];} print $fhr ( $lot,";"); # nom de l'organisme sous l'en-tête Leu (1 place pour le total) #@u = ([[(0)x64]x$n]x$t); # récupération des compteurs de l'organisme foreach $k (0..9) { # compteurs des codons de Leu, Lys, Phe print $fhr ( $u[$v][$h][$k],";");} # foreach $k (10..19) { # compteurs des codons Ser, Tyr, Asn print $fhr ( $u[$v][$h][$k],";");} # foreach $k (20..29) { # compteurs des codons Arg, His, Gln print $fhr ( $u[$v][$h][$k],";");} # foreach $k (30..37) { # compteurs des codons Gly, Glu, Asp print $fhr ( $u[$v][$h][$k],";");} # foreach $k (38..44) { # compteurs des codons Ala, Cys, Met, Trp print $fhr ( $u[$v][$h][$k],";");} print $fhr ( $u[$v][$h][49],";"); # foreach $k (45..48) { # compteurs des codons Pro print $fhr ( $u[$v][$h][$k],";");} # foreach $k (50..53) { # compteurs des codons Val print $fhr ( $u[$v][$h][$k],";");} # foreach $k (54..57) { # compteurs des codons Thr print $fhr ( $u[$v][$h][$k],";");} # foreach $k (58..60) { # compteurs des codons Ile print $fhr ( $u[$v][$h][$k],";");} # foreach $k (61..63) { # compteurs des codons Stp print $fhr ( $u[$v][$h][$k],";");} # print $fhr ("\n"); } } # 1 organisme codons ptprot # # ********************* acides aminés ********************* # foreach $h (0..$n1) { # écriture de l'organisme # 1 organisme aas ptprot print $fhra ($noms[$h],";"); # nom print $fhra ("\n"); # ******* écriture de l'en-tête aas ******* print $fhra ( "Gène;"); # Gène de protéine print $fhra ( "Leu;"); # Leu print $fhra ( "Lys;Phe;Ser;"); # Lys, Phe, Ser print $fhra ( "Tyr;Asn;Arg;"); # Tyr, Asn, Arg print $fhra ( "His;Gln;Gly;"); # His, Gln, Gly print $fhra ( "Glu;Asp;Ala;"); # Glu, Asp, Ala print $fhra ( "Cys;Met;Trp;Pro;"); # Cys, Met, Trp, Pro print $fhra ( "Val;"); # Val print $fhra ( "Thr;"); # Thr print $fhra ( "Ile;"); # Ile print $fhra ( "Stp;"); # Stp print $fhra ("\n"); # foreach $v (0..$t1) { $xl=@{$x[$v]}; $lot=$x[$v][0]; $xl1=$xl-1; # 1 gène aas ptprot foreach $j (1..$xl1) { $lot=$lot.",".$x[$v][$j];} print $fhra ( $lot,";"); # nom de l'organisme sous l'en-tête Leu (1 place pour le total) #@u = ([[(0)x64]x$n]x$t); # récupération des compteurs de l'organisme $leu=0; $ser=0; $arg=0; $gly=0; $ala=0; $pro=0; $val=0; $thr=0; foreach $k (0..5) { $leu=$leu+$u[$v][$h][$k];} print $fhra ($leu,";"); # Leu total $lys=$u[$v][$h][6]+$u[$v][$h][7]; print $fhra ($lys,";"); # Lys total $phe=$u[$v][$h][8]+$u[$v][$h][9]; print $fhra ($phe,";"); # Phe total foreach $k (10..15) { $ser=$ser+$u[$v][$h][$k];} print $fhra ($ser,";"); # Ser total $tyr=$u[$v][$h][16]+$u[$v][$h][17]; print $fhra ($tyr,";"); # Tyr total $asn=$u[$v][$h][18]+$u[$v][$h][19]; print $fhra ($asn,";"); # Asn total foreach $k (20..25) { $Arg=$Arg+$u[$v][$h][$k];} print $fhra ($arg,";"); # Arg total $his=$u[$v][$h][26]+$u[$v][$h][27]; print $fhra ($his,";"); # His total $gln=$u[$v][$h][28]+$u[$v][$h][29]; print $fhra ($gln,";"); # Gln total foreach $k (30..33) { $gly=$gly+$u[$v][$h][$k];} print $fhra ($gly,";"); # Gly total $glu=$u[$v][$h][34]+$u[$v][$h][35]; print $fhra ($glu,";"); # Glu total $asp=$u[$v][$h][36]+$u[$v][$h][37]; print $fhra ($asp,";"); # Asp total foreach $k (38..41) { $ala=$ala+$u[$v][$h][$k];} print $fhra ($ala,";"); # Ala total $cys=$u[$v][$h][42]+$u[$v][$h][43]; print $fhra ($cys,";"); # Cys total $met=$u[$v][$h][44]; print $fhra ($met,";"); # Met $trp=$u[$v][$h][49]; print $fhra ($trp,";"); # Trp foreach $k (45..48) { $pro=$pro+$u[$v][$h][$k];} print $fhra ($pro,";"); # Pro total foreach $k (50..53) { $val=$val+$u[$v][$h][$k];} print $fhra ($val,";"); # Val total foreach $k (54..57) { $thr=$thr+$u[$v][$h][$k];} print $fhra ($thr,";"); # Thr total $ile=$u[$v][$h][58]+$u[$v][$h][59]+$u[$v][$h][60]; print $fhra ($ile,";"); # Ile total $stp=$u[$v][$h][61]+$u[$v][$h][62]+$u[$v][$h][63]; print $fhra ($stp,";"); # Stp total print $fhra ("\n"); } # 1 gène aas ptprot # } # 1 organisme aas ptprot } # **** écriture 1 tableau par organisme pour n organismes **** #
Décompte des répétitions
modifier- repete.pl
# Décompte des répétitions des bases dans l'ADN des procaryotes. Documentation, Paris le 18.9.16 # # fichiers en entrée: xxx99 d'après la liste de fichiers @noms. En lecture, longueur des lignes fixe de 70 caractères (format FASTA de NCBI). # fichier en sortie, résultat du comptage: repete.txt (création et mise à jour) # fichier en sortie, repete1.txt (création et mise à jour), pour exploitation sous Calc avec l'en-tête : # KEGG;effectif;%GC;>4T;>4A;>4C;>4G;<4T;<4A;<4C;<4G;4T;4A;4C;4G;Ctrl # KEGG: les 1ères lettres de l'organisme d'après KEGG. C'est aussi le prefixe des fichiers xxx99 créés avec le logiciel gedit ( sans le .txt); # effectif: taille du fichier xxx99 calculé avec le paragraphe "calcul de la taille du fichier"; # Le total calculé à partir du comptage des répétitions sert à contrôler le nombre de caractères différents de T,A,C,G; # Ctrl: est la différence "effectif-total". Ctrl négatif veut dire que la dernière ligne n'a pas été comptée parce qu'il y a des lignes à blanc à la # fin du fichier. # Nettoyer et recommencer. Ctrl positif veut dire que le programme n'a pas compté certains caractères: soit ils sont différents de AGTC et à ce moment # retirer le fichier de la liste ou le conserver si la différence n'a pas d'impact sur les interprétations. Par contre une autre # raison de Ctrl>0 c'est quele comptage n' a pas consigné les répétitions > 20. Il faut conserver le fichier et rechercher manuellement ces répétitions. # %GC: 100x(total des GC)/total. # >4T,A,C,G: total des répétitions supérieures à 4 de T,A,C,G. # <4T,A,C,G: total des répétitions 2 et 3 de T,A,C,G. # 4T,A,C,G: total des répétitions 4 de T,A,C,G. # # Fabrication des fichiers xxx99 (99 est modfiable, penser à le changer dans "my @z=") sous gedit: # - se connecter à KEGG http://www.genome.jp/kegg/pathway.html # - allez à genome puis organisms # - rechercher les 3 lettres KEGG: ctrl F pour les noms de 3 lettres à faible occurrences ou les noms à plus de 3 lettres. # Pour les noms à forte occurrences mais seulement les noms à 3 lettres # * retour à http://www.genome.jp/kegg/pathway.html # * mettre le nom dans le champs "Select prefix" puis entrée # * choisir le 1er schéma affiché et y sélectionner un gène dont le rectangle est coloré # * on retrouve le nom en 3 lettres # - cliquer sur le nom en 3 lettres:retrouver la ligne contenant la référence à la séquence du chromosome commençant par "GB:", # * "Sequence RS: NC_000913 (GB: U00096)" # * cliquer sur "GB: U00096" on obtient une entrée de NCBI qui est dans ce cas u00096: http://www.ncbi.nlm.nih.gov/nuccore/U00096 # * cliquer alors sur FASTA puis copier la séquence du gène dans xxx99. # * En fin de fichier supprimer les lignes à blanc s'il y a lieu. Puis sauvegarder. # # Exécution du programme Perl (sous distribution Linusx Ubuntu 15.4 ou plus avec Perl 2.7 ou plus installé): # * mettre repete.pl et les fichiers xxx99 dans un même dossier # * dans un terminal accéder à ce dossier ( attention aux autorisations pour créer et écrire les .txt dans ce répertoire) et lancer # * perl repete.pl # # Résultat fichier repete.txt: KEGG donne 4 641 652 paires de bases. Chercher dans organismes de KEGG eco (Escherichia Coli). # # eco # n t a c g # 1 572076 575694 684256 684448 # 2 155479 155226 185070 183945 # 3 50279 50147 31599 31526 −−−−> Il y a 31 526 séquences de ggg. # 4 15523 15380 6067 6027 # 5 5929 5810 1086 1056 # 6 1929 1891 153 154 # 7 475 472 31 25 # 8 97 109 7 3 # 9 11 7 1 0 # 10 0 0 0 1 #.....jusqu'à 20 répétitions # # resultat fichier repete1.txt: # KEGG;effectif;%GC;>4T;>4A;>4C;>4G;<4T;<4A;<4C;<4G;4T;4A;4C;4G;Ctrl # eco; 4641652; 50.7907098593346; 8441; 8289; 1278; 1239; 205758; 205373; 216669; 215471; 15523; 15380; 6067; 6027; 0 # my @noms=(aae); # **************************************************** my $filenamer = "repete.txt"; # Sous Linux my $filenamer1 = "repete1.txt"; # Sous Linux my $fhr; $fhr1; open ($fhr, '>>', $filenamer) or die "Impossible d'ouvrir le fichier $filenamer en ajout"; open ($fhr1, '>>', $filenamer1) or die "Impossible d'ouvrir le fichier $filenamer1 en ajout"; print $fhr1 ( "KEGG;effectif;%GC;>4T;>4A;>4C;>4G;<4T;<4A;<4C;<4G;4T;4A;4C;4G;Ctrl \n"); # **************************************************** my @z=( a99 ); my $m=@z; my $pt=$m-1; my $o=@noms; my $nt=$o-1; my @fic = (a)x$m; foreach $l (0..$nt) { # un fichier foreach $kt (0..$pt) { $lz=length($z[$kt])-1; $zi=substr($z[$kt],1,$lz); $fic[$kt]=$noms[$l].$zi; #format du nom du fichier; @noms.99 ou xxx99 } foreach $kt (0..$pt) { comptfich($fic[$kt]);} # boucle pour tous les fichiers de @noms sub comptfich { # début de sub comptfich: comptage dans 1 fichier $filename=$fic[$kt]; @Valeurs = (); # Initialisation du tableau [4 bases x 20 compteurs de répétition] foreach $b (0..3) {foreach $n (0..19) { $Valeurs [$b] [$n] = 0;}} # augmenter les bornes supérieures si nécessaire: séquences spécifiques $NbElem = @Valeurs; #********************************** calcul de la taille du fichier ************* (enlever les lignes à blancs de fin de fichier) open (my $fh, '<', $filename) or die "Impossible d'ouvrir le fichier '$filename' en lecture"; my $compteur = 0; while (my $ligne = <$fh>) { $lig=$ligne; if (length($lig)>1) { ++$compteur;} } $lu = length($lig); $eff=($compteur-1)*70+$lu-1; print ("$noms[$l] $eff octets \n"); close $fh; # fin calcul taille du fichier # #*************** Initialisation des compteurs de travail ********** my $i = 0, $k=0; my @e= (1,2,3,4); my $j=0; my @p = (0)x4; # ************* Comptage de la 1ère base des chromosomes circulaires ************ open ($fh, '<', $filename); $lis1 = <$fh>; $car = substr($lis1,0,1 ); # 1ère base sauvegardée dans $car $j=1; for ($i=1; $i<=70; $i+=1) { # 70 à changer si on doute que la 1ère base est répétée sur une séquence plus grande que 70 $lu = substr($lis1,$i,1 ); if ($lu eq $car) { $j+=1;} # ajout au compteur si répétition de la 1ère base else { $i=70;}} # fin de la répétition de la 1ère base $d=$j; # sauvegarde du nombre de répétition de la 1ère base dans $d close $fh; # fin 1ère base # ************ Comptage de toutes les séquences répétitives. ************ $i=0; $j=0; open ($fh, '<', $filename); while( defined( $lis = <$fh> ) ) { # lecture d'une ligne chomp $lis; while ($i <= 70) { # extrait d'une base; 70, longueur des lignes du fichier en entrée. À changer en 60 pour les fichiers KEGG et non NCBI. $lu = substr($lis,$i,1 ); if( $lu eq "T" ) { # Traitement de la base "T" if($p[1] == 0) { # on passe d'une autre base à la base "T" foreach $j (@e, $e!=1) { # ajout au compteur de l'ancienne base if( $p[$j] !=0 ) { $n=$p[$j]-1; $b=$j-1; $Valeurs [$b] [$n] +=1; $p[$j]=0; } } $p[1]+=1;} # ajout pour la 1ère fois au compteur de "T" else { $p[1]+=1; }} # ajout au compteur de "T" pour les fois suivantes if( $lu eq "A" ) { # Traitement de la base "A" if($p[2] == 0) { foreach $j (@e, $e!=2) { if( $p[$j] !=0 ) { $n=$p[$j]-1; $b=$j-1; $Valeurs [$b] [$n] +=1; $p[$j]=0; }} $p[2]+=1;} else {$p[2]+=1;}} if( $lu eq "C" ) { # Traitement de la base "C" if($p[3] == 0) { foreach $j (@e, $e!=3) { if( $p[$j] !=0 ) { $n=$p[$j]-1; $b=$j-1; $Valeurs [$b] [$n] +=1; $p[$j]=0; }} $p[3]+=1;} else {$p[3]+=1;}} if( $lu eq "G" ) { # Traitement de la base "G" if($p[4] == 0) { foreach $j (@e, $e!=4) { if( $p[$j] !=0 ) { $n=$p[$j]-1; $b=$j-1; $Valeurs [$b] [$n] +=1; $p[$j]=0; }} $p[4]+=1;} else {$p[4]+=1;}} $i+=1;}; # Base suivante $i = 0;}; # Ligne suivante foreach $j (@e) { # ajout au compteur de la dernière base du fichier if( $p[$j] !=0 ) { $n=$p[$j]-1; $b=$j-1; $Valeurs [$b] [$n] +=1; $p[$j]=0; } }; # # ********* Le chromosome est circulaire. Joindre les 2 bouts **************** @c=("T","A","C","G"); foreach $k (0..3) { if ($c[$k] eq $car) { $x=$k;}}; #adaptation au tableau if ( $x == $b ) { # $b et $n sont ceux de la dernière base du fichier $Valeurs [$b] [$n] -=1; # retrait de 1 du compteur de la dernière répétition de cette base si la dernière base est la même que la 1ère, $car. $n = $n + $d; # ajout de la répétition de la 1ère base ($d) au compteur de la répétition de la dernière base $Valeurs [$b] [$n] +=1; # ajout de 1 au compteur de la nouvelle répétition de la dernière base $n = $d-1; # compteur de la 1ère répétition de la 1ère base, adaptation $Valeurs [$b] [$n] -=1;} # retrait de 1 du compteur de la 1ère répétition de la base # # ********* Fin, le chromosome est circulaire. ****************************** # my $v=0, $w=0, $t1=0, $t2=0, $t3=0, $a1=0, $a2=0, $a3=0, $c1=0, $c2=0, $c3=0, $g1=0, $g2=0, $g3=0; # # ********* écriture des répétitions pour le procaryote traité dans repete.txt ********** print $fhr ($noms[$l], "\n"); print $fhr ("n t a c g \n"); $t2=$Valeurs [0] [1]+$Valeurs [0] [2]; $a2=$Valeurs [1] [1]+$Valeurs [1] [2]; # total <4T,A $c2=$Valeurs [2] [1]+$Valeurs [2] [2]; $g2=$Valeurs [3] [1]+$Valeurs [3] [2]; # total <4C,G $t3=$Valeurs [0] [3]; $a3=$Valeurs [1] [3]; # total =4T,A $c3=$Valeurs [2] [3]; $g3=$Valeurs [3] [3]; # total =4C,G foreach $k (0..19) { $j=$k+1; print $fhr ( $j," ", $Valeurs [0] [$k]," ", $Valeurs [1] [$k]," ", $Valeurs [2] [$k]," ", $Valeurs [3] [$k], "\n"); # # ********* Mise en forme de la ligne correspondant à l'en-tête du fichier repete1.txt ********** Sommations avec le foreach: $v=$v+($Valeurs [0] [$k]+$Valeurs [1] [$k]+$Valeurs [2] [$k]+ $Valeurs [3] [$k])*$j; # total $w=$w+($Valeurs [2] [$k]+ $Valeurs [3] [$k])*$j; # total G+C if ($j > 4) { $t1=$t1+$Valeurs [0] [$k]; $a1=$a1+$Valeurs [1] [$k]; # total >4T,A $c1=$c1+$Valeurs [2] [$k]; $g1=$g1+$Valeurs [3] [$k];} # total >4C,G }; $diff=$eff-$v; $w=100*$w/$v; # fin des sommations et de l'écriture dans repete.txt # **** écriture de la ligne dans repete1.txt ****** print $fhr1 ($noms[$l],"; ",$eff,"; ",$w,"; ",$t1,"; ",$a1,"; ",$c1,"; ",$g1,"; ", $t2,"; ",$a2,"; ",$c2,"; ",$g2,"; ", $t3,"; ",$a3,"; ",$c3,"; ",$g3,"; ",$diff, "\n"); } # ********** fin de sub comptfich ************ close $fh; # fermeture du fichier du procaryote traité *********** } # ********** un fichier procaryote a été traité********* close $fhr; close $fhr1;
Répétitions aléatoires
modifierrepete-alea.pl
# Comptage des répétitions des bases dans l'ADN des procaryotes. Aléas. Documentation, Paris le 4.9.16 # # fichier en sortie, résultat du comptage: repete.txt (création et mise à jour) # fichier en sortie, repete1.txt (création et mise à jour), pour exploitation sous Calc avec l'en-tête : d-n;effectif;%GC;>4GC;>4AT;<5GC;<5AT;diff GC. # d-n: pour 15% souhaité on aura avec un dénominateur 100, "100,15". # effectif: un chromosome fictif et constant pour calculer l'équation de l'aléa; un chromosome d'un procaryote pour comparaison dans repete.txt. # %GC: 100x(total des GC)/effectif. # >4GC et >4AT: 10000x(total des répétitions supérieures à 4 de G et de C)/total; de même pour A et T. # <5GC et <5AT: 1000x(total des répétitions 2, 3 et 4 de G et de C)/total; de même pour A et T. # diff GC: différence entre le %GC souhaité (15%) et le %GC obtenu par l'aléa. # Ce n'est pas nécessaire de refaire les courbes de l'aléa (sauf pour les courbes en base unique A ou T ou G ou C): # elle a pour équation en AT -0.00083x^3+0.16x^2-12.16x+302.4 # et en GC -0.00083x^3+0.16x^2-12.16x+302.4.. # Les paramètres pour comparaison ave un procaryote donné, résultant de repete.pl: # my @suit=(15,18,20,23,25,28,30,33,35,38,40,43,45,48,50,53,55,58,60,63,65,68,70,73,75,78,80,83,85); # # pour courbes de l'aléa, mettre à jour la ligne "foreach $pk (0..28) {". # my $tet=100; $efc=5000000; # pour courbes de l'aléa. # # my @suit=(21); # pour comparaison d'un seul chromosome avec l'aléa, mettre à jour la ligne "foreach $pk (0..0) {". %GC=21/47. # my $tet=47; $efc=2123000; # pour courbes de l'aléa. La taille du chromosomen est celle du procaryote à comparer issue de repete.pl. # # Exécution du programme Perl (sous distribution Linusx Ubuntu 15.4 ou plus avec Perl 2.7 ou plus installé): # * dans un terminal accéder au dossier home et lancer # * perl repete.pl # my @suit=(25); # pour comparaison d'un seul chromosome avec l'aléa, mettre à jour la ligne "foreach $pk (0..0) {". %GC=21/47. my $tet=98; $efc=751719; # pour courbes de l'aléa. La taille du chromosomen est celle du procaryote à comparer issue de repete.pl. # **************************************************** my $filenamer = "repete.txt"; # Sous Linux my $filenamer1 = "repete1.txt"; # Sous Linux my $fhr; $fhr1; open ($fhr, '>>', $filenamer) or die "Impossible d'ouvrir le fichier $filenamer en ajout"; open ($fhr1, '>>', $filenamer1) or die "Impossible d'ouvrir le fichier $filenamer1 en ajout"; print $fhr1 ( "d-n;effectif;%GC;>4T;>4A;>4C;>4G;<4T;<4A;<4C;<4G;4T;4A;4C;4G;diff GC \n"); # **************************************************** foreach $pk (0..0) { my $ll=$tet; $gc=$suit[$pk]; my $at=$ll-$gc; my $gc1=$gc-1; $at1=$at-1; my $gc2=$gc+3; $gc3=$gc2+1; $gc4=$gc3+$gc-2; my $at2=$gc4+1; $at3=$at2+$at-2; $at4=$at3+1; $at5=$at4+$at-2; my $l=$ll*2; # le max de l'aléa *** my $i=0; my $r=0; my $k=0; my $s=0; my $v=0; my $a=0; my $e=0; # Initialisation du tableau pour les fréquences @Valeurs = (); foreach $b (0..3) {foreach $n (0..19) { $Valeurs [$b] [$n] = 0;}} # ********* augmenter les bornes supérieures si nécessaire: séquences spécifiques ************ my @m= (1,2,3,4); my @p = (0)x4; # my @u=(1,$at2..$at3); @d=(2,$at4..$at5); @t=(3,5..$gc2); @q=(4,$gc3..$gc4); for ($k=1; $k<=$efc; $k+=1) { $r = int rand($l) + 1; if ($s==0) { foreach $i (0..$at1) { if ($r == $u[$i]) { $s=1; }}} if ($s==0) { foreach $i (0..$at1) { if ($r == $d[$i]) { $s=2; }}} if ($s==0) { foreach $i (0..$gc1) { if ($r == $t[$i]) { $s=3; }}} if ($s==0) { foreach $i (0..$gc1) { if ($r == $q[$i]) { $s=4; }}} Freq($s); if($v == 0 ) { if($a==0) {$a=$s;} if($s==$a) { $e +=1;} else {$v=1; }} $s=0;} $s=0; Freq($s); # sub Freq { # my ($s) = @_; my $j=0; # ********* Comptage de toutes les séquences répétitives. À compléter pour les séquences spécifiques. ************ # if( $s == 1 ) { if($p[1] == 0) { foreach $j (@m, $m[$j]!=1) { if ( $p[$j] !=0 ) { $n=$p[$j]-1; $b=$j-1; $Valeurs [$b] [$n] +=1; $p[$j]=0; }} $p[1]+=1;} else { $p[1]+=1;}} if( $s == 2 ) { if($p[2] == 0) { foreach $j (@m, $m[$j]!=2) { if( $p[$j] !=0 ) { $n=$p[$j]-1; $b=$j-1; $Valeurs [$b] [$n] +=1; $p[$j]=0; }} $p[2]+=1;} else {$p[2]+=1;}} if( $s == 3 ) { if($p[3] == 0) { foreach $j (@m, $m[$j]!=3) { if( $p[$j] !=0 ) { $n=$p[$j]-1; $b=$j-1; $Valeurs [$b] [$n] +=1; $p[$j]=0; }} $p[3]+=1;} else {$p[3]+=1;}} if($s == 4 ) { if($p[4] == 0) { foreach $j (@m, $m[$j]!=4) { if( $p[$j] !=0 ) { $n=$p[$j]-1; $b=$j-1; $Valeurs [$b] [$n] +=1; $p[$j]=0; }} $p[4]+=1;} else {$p[4]+=1;}} if($s != 0 ) { return} if($s == 0 ) { foreach $j (@m) { if ( $p[$j] !=0 ) { $n=$p[$j]-1; $b=$j-1; $Valeurs [$b] [$n] +=1; $p[$j]=0; } }} @c=(1..4); foreach $k (0..3) { if ($c[$k] == $a) { $x=$k;}}; if ( $x == $b ) { $Valeurs [$b] [$n] -=1; $n = $n + $e; $Valeurs [$b] [$n] +=1; $n = $e-1; $Valeurs [$b] [$n] -=1;} # #**************************************************************** my $v=0, $w=0, $t1=0, $t2=0, $a1=0, $a2=0, $c1=0, $c2=0, $g1=0, $g2=0; # print $fhr ($tet, " ",$suit[$pk]," ", $efc, "\n"); print $fhr ("n t a c g \n"); $t2=$Valeurs [0] [1]+$Valeurs [0] [2]; $a2=$Valeurs [1] [1]+$Valeurs [1] [2]; # total <4T,A $c2=$Valeurs [2] [1]+$Valeurs [2] [2]; $g2=$Valeurs [3] [1]+$Valeurs [3] [2]; # total <4C,G $t3=$Valeurs [0] [3]; $a3=$Valeurs [1] [3]; # total =4T,A $c3=$Valeurs [2] [3]; $g3=$Valeurs [3] [3]; # total =4C,G foreach $k (0..19) { $j=$k+1; print $fhr ( $j," ", $Valeurs [0] [$k]," ", $Valeurs [1] [$k]," ", $Valeurs [2] [$k]," ", $Valeurs [3] [$k], "\n"); $v=$v+($Valeurs [0] [$k]+$Valeurs [1] [$k]+$Valeurs [2] [$k]+ $Valeurs [3] [$k])*$j; $w=$w+($Valeurs [2] [$k]+ $Valeurs [3] [$k])*$j; if ($j > 4) { $t1=$t1+$Valeurs [0] [$k]; $a1=$a1+$Valeurs [1] [$k]; # total >4T,A $c1=$c1+$Valeurs [2] [$k]; $g1=$g1+$Valeurs [3] [$k];} # total >4C,G }; $w=100*$w/$v; $diff=$w-100*$suit[$pk]/$tet; # controle %GC print $fhr1 ( $tet,", ", $suit[$pk],"; ",$v,"; ",$w,"; ",$t1,"; ",$a1,"; ",$c1,"; ",$g1,"; ",$t2,"; ",$a2,"; ",$c2,"; ",$g2,"; ",$t3,"; ",$a3,"; ",$c3,"; ",$g3,"; ",$diff, "\n"); return;} # } close $fhr; close $fhr1;