J'ai une donnée comme celle-ci
>sp|Q96A73|P33MX_HUMAN Putative monooxygenase p33MONOX OS=Homo sapiens OX=9606 GN=KIAA1191 PE=1 SV=1 RNDDDDTSVCLGTRQCSWFAGCTNRTWNSSAVPLIGLPNTQDYKWVDRNSGLTWSGNDTCLYSCQNQTKGLLYQLFRNLFCSYGLTEAHGKWRCADASITNDKGHDGHRTPTWWLTGSNLTLSVNNSGLFFLCGNGVYKGFPPKWSGRCGLGYLVPSLTRYLTLNASQITNLRSFIHKVTPHR cat input.txt | grep -v ">" | fold -w1 | sort | uniq -c 6 A 9 C 10 D 1 E 7 F 18 G 5 H 4 I 7 K 21 L 15 N 7 P 6 Q 11 R 16 S 18 T 7 V 8 W 7 Y
Je veux compter combien de chaque lettre il y a, donc si j'en ai une, je compte comme ça
>sp|Q96A73|P33MX_HUMAN Putative monooxygenase p33MONOX OS=Homo sapiens OX=9606 GN=KIAA1191 PE=1 SV=1 RNDDDDTSVCLGTRQCSWFAGCTNRTWNSSAVPLIGLPNTQDYKWVDRNSGLTWSGNDTCLYSCQNQTKGLLYQLFRNLFCSYGLTEAHGKWRCADASITNDKGHDGHRTPTWWLTGSNLTLSVNNSGLFFLCGNGVYKGFPPKWSGRCGLGYLVPSLTRYLTLNASQITNLRSFIHKVTPHR >sp|P13674|P4HA1_HUMAN Prolyl 4-hydroxylase subunit alpha-1 OS=Homo sapiens OX=9606 GN=P4HA1 PE=1 SV=2 VECCPNCRGTGMQIRIHQIGPGMVQQIQSVCMECQGHGERISPKDRCKSCNGRKIVREKKILEVHIDKGMKDGQKITFHGEGDQEPGLEPGDIIIVLDQKDHAVFTRRGEDLFMCMDIQLVEALCGFQKPISTLDNRTIVITSHPGQIVKHGDIKCVLNEGMPIYRRPYEKGRLIIEFKVNFPENGFLSPDKLSLLEKLLPERKEVEE >sp|Q7Z4N8|P4HA3_HUMAN Prolyl 4-hydroxylase subunit alpha-3 OS=Homo sapiens OX=9606 GN=P4HA3 PE=1 SV=1 MTEQMTLRGTLKGHNGWVTQIATTPQFPDMILSASRDKTIIMWKLTRDETNYGIPQRALRGHSHFVSDVVISSDGQFALSGSWDGTLRLWDLTTGTTTRRFVGHTKDVLSVAFSSDNRQIVSGSRDKTIKLWNTLGVCKYTVQDESHSEWVSCVRFSPNSSNPIIVSCGWDKLVKVWNLANCKLK >sp|P04637|P53_HUMAN Cellular tumor antigen p53 OS=Homo sapiens OX=9606 GN=TP53 PE=1 SV=4 IQVVSRCRLRHTEVLPAEEENDSLGADGTHGAGAMESAAGVLIKLFCVHTKALQDVQIRFQPQL >sp|P10144|GRAB_HUMAN Granzyme B OS=Homo sapiens OX=9606 GN=GZMB PE=1 SV=2 MQPILLLLAFLLLPRADAGEIIGGHEAKPHSRPYMAYLMIWDQKSLKRCGGFLIRDDFVLTAAHCWGSSINVTLGAHNIKEQEPTQQFIPVKRPIPHPAYNPKNFSNDIMLLQLERKAKRTRAVQPLRLPSNKAQVKPGQTCSVAGWGQTAPLGKHSHTLQEVKMTVQEDRKCES >sp|Q9UHX1|PUF60_HUMAN Poly(U)-binding-splicing factor PUF60 OS=Homo sapiens OX=9606 GN=PUF60 PE=1 SV=1 MGKDYYQTLGLARGASDEEIKRAYRRQALRYHPDKNKEPGAEEKFKEIAEAYDVLSDPRKREIFDRYGEEGLKGSGPSGGSGGGANGTSFSYTFHGDPHAMFAEFFGGRNPFDTFFGQRNGEEGMDIDDPFSGFPMGMGGFTNVNFGRSRSAQEPARKKQDPPVTHDLRVSLEEIYSGCTKKMKISHK >sp|Q06416|P5F1B_HUMAN Putative POU domain, class 5, transcription factor 1B OS=Homo sapiens OX=9606 GN=POU5F1B PE=5 SV=2 IVVKGHSTCLSEGALSPDGTVLATASHDGYVKFWQIYIEGQDEPRCLHEWKPHDGRPLSCLLFCDNHKKQDPDVPFWRFLITGADQNRELKMWCTVSWTCLQTIRFSPDIFSSVSVPPSLKVCLDLSAEYLILSDVQRKVLYVMELLQNQEEGHACFSSISEFLLTHPVLSFGIQVVSRCRLRHTEVLPAEEENDSLGADGTHGAGAMESAAGVLIKLFCVHTKALQDVQIRFQPQLNPDVVAPLPTHTAHEDFTFGESRPELGSEGLGSAAHGSQPDLRRIVELPAPADFLSLSSETKPKLMTPDAFMTPSASLQQITASPSSSSSGSSSSSSSSSSSLTAVSAMSSTSAVDPSLTRPPEELTLSPKLQLDGSLTMSSSGSLQASPRGLLPGLLPAPADKLTPKGPGQVPTATSALSLELQEVEP >sp|O14683|P5I11_HUMAN Tumor protein p53-inducible protein 11 OS=Homo sapiens OX=9606 GN=TP53I11 PE=1 SV=2 MIHNYMEHLERTKLHQLSGSDQLESTAHSRIRKERPISLGIFPLPAGDGLLTPDAQKGGETPGSEQWKFQELSQPRSHTSLKVSNSPEPQKAVEQEDELSDVSQGGSKATTPASTANSDVATIPTDTPLKEENEGFVKVTDAPNKSEISKHIEVQVAQETRNVSTGSAENEEKSEVQAIIESTPELDMDKDLSGYKGSSTPTKGIENKAFDRNTESLFEELSSAGSGLIGDVDEGADLLGMGREVENLILENTQLLETKNALNIVKNDLIAKVDELTCEKDVLQGELEAVKQAKLKLEEKNRELEEELRKARAEAEDARQKAKDDDDSDIPTAQRKRFTRVEMARVLMERNQYKERLMELQEAVRWTEMIRASRENPAMQEKKRSSIWQFFSRLFSSSSNTTKKPEPPVNLKYNAPTSHVTPSVK
4 Réponses :
Je ne sais pas à quel point ce serait plus rapide, mais si vous essayez de publier vos horaires
$ wc filey 118098 649539 16828965 filey
ce rapport compte pour le fichier, pas ligne par ligne. Comme indiqué ci-dessous, tous les awk ne prennent pas en charge le fractionnement de chaîne vide.
Voici les horaires des trois approches:
$ time grep -v ">" filey | fold -w1 | sort | uniq -c >/dev/null
real 0m11.470s
user 0m11.746s
sys 0m0.260s
$ time awk '{n=split($0,a,""); for(i=1;i<=n;i++) c[a[i]++]} END{for(k in c) print k,c[k]}' filey >/dev/null
real 0m7.441s
user 0m7.334s
sys 0m0.060s
$ time awk '{n=length($0); for(i=1;i<=n;i++) c[substr($0,i,1)]++} END{for(k in c) print k,c[k]}' filey >/dev/null
real 0m5.055s
user 0m4.979s
sys 0m0.047s
pour le fichier de test
$ awk '!/^>/ {n=split($0,a,""); for(i=1;i<=n;i++) c[a[i]]++}
END {for(k in c) print k,c[k]}' file | sort
A 6
C 9
D 10
E 1
F 7
G 18
H 5
I 4
K 7
L 21
N 15
P 7
Q 6
R 11
S 16
T 18
V 7
W 8
Y 7
cela m'a surpris que substr soit plus rapide que split . Peut-être à cause de l'allocation des tableaux.
Vous devez mentionner que le fractionnement sur un caractère nul est un comportement non défini par POSIX, donc cela ne fonctionnera que dans certains awks, par exemple GNU awk.
@Tiw Vous avez mal compris mon commentaire. Dans split ($ 0, a, "") le problème n'est pas avec $ 0 étant une chaîne nulle, c'est avec "" étant nul chaîne. Utiliser une chaîne nulle comme «expression régulière» sur laquelle se scinder est un comportement non défini par POSIX. C'est vrai que ce soit la valeur de FS ou le 3ème argument de split ().
Oh, tu veux dire que ... Je suis désolé d'avoir mal compris. Je l'ai essayé sur awk de Freebsd et cela a bien fonctionné. @EdMorton
@Tiw Oui, certains awks feront ce que vous voulez mais d'autres non. Consultez la spécification POSIX awk - lors de la discussion de FS Si FS est une chaîne nulle, le comportement n'est pas spécifié. et plus tard dans le synopsis split (..., fs) L'effet d'une chaîne nulle en tant que valeur de fs n'est pas spécifié. .
J'ai également ajouté quelques horaires à ma réponse.
Avec n'importe quel awk dans n'importe quel shell sur n'importe quelle boîte UNIX:
$ awk -v ORS= '!/^>/{gsub(/./,"&\n"); print}' file | sort | uniq -c
107 A
41 C
102 D
132 E
65 F
140 G
52 H
84 I
114 K
174 L
39 M
67 N
107 P
88 Q
101 R
168 S
115 T
101 V
27 W
30 Y
ou si vous préférez:
$ cat tst.awk
!/^>/ {
for (i=1; i<=length($0); i++) {
cnt[substr($0,i,1)]++
}
}
END {
for (char in cnt) {
print char, cnt[char]
}
}
$ awk -f tst.awk file
A 107
N 67
P 107
C 41
Q 88
D 102
E 132
R 101
F 65
S 168
G 140
T 115
H 52
I 84
V 101
W 27
K 114
Y 30
L 174
M 39
Ahhhh, tu m'as battu en quelques minutes. J'allais dormir et j'ai pensé publier cette idée, mais il semble que le monsieur Ed l'ait déjà fait ++ ve pour une belle solution.
@Ed Morton awk: erreur de syntaxe au niveau du fichier source de la ligne 1, le contexte tst.awk est cat >>> tst. <<< awk awk: renflouement à la ligne source 11
@Learner Regardez à nouveau ma réponse, vous ne faites pas ce qu'elle montre. On dirait que vous avez peut-être copié la ligne cat dans le script awk.
Je pensais à cette solution quand j'ai vu cette question, mais je pensais que plusieurs substr seraient beaucoup plus inefficaces que cette division une seule fois. Cependant, je ne suis pas sûr de cela, est-ce que je pense bien cette fois?
Le fractionnement @Tiw sur un caractère nul n'est pas portable - il se divisera en caractères dans certains awks, laissera toute la ligne telle quelle dans d'autres, et peut faire n'importe quoi d'autre et toujours être compatible POSIX puisque ce comportement n'est pas défini. Je ne sais vraiment pas si un seul split () serait plus efficace que plusieurs substr () s ou non, mais à mon humble avis, la portabilité est plus importante dans ce cas.
Le comptage des caractères dans les chaînes peut facilement être fait avec awk. Pour ce faire, vous utilisez la fonction gsub :
gsub (ere, repl [ in]) Se comporter comme
sub(voir ci-dessous), sauf qu'il remplacera toutes les occurrences de l'expression régulière (comme le substitut global de l'utilitaire ed) dans$ 0ou dans ledansargument lorsqu'il est spécifié.sub (ere, repl [ in]) Remplacez la chaîne
replà la place de la première instance de l'expression régulière étendueEREdans la chaîne dans et renvoyez le nombre de substitutions .Si dansest omis, awk utilisera l'enregistrement courant ($ 0) à sa place.source: Awk Posix Standard
Les deux fonctions suivantes effectuent le comptage de cette manière:
OP : grep,sort,uniq : 47.548 s EdMorton 1 : awk : 39.992 s EdMorton 2 : awk,sort,uniq : 53.965 s kvantour 1 : awk : 18.661 s kvantour 2 : awk : 9.309 s kvantour 3 : bioawk : 1.838 s kvantour 4 : bioawk : 1.838 s karafka : awk : 38.139 s stack0114106 1: perl : 22.754 s stack0114106 2: perl : 13.648 s stack0114106 3: perl (zdim) : 7.759 s
ou s'il peut apparaître beaucoup de caractères consécutifs égaux, la solution suivante pourrait raser quelques secondes.
bioawk -c fastx '{s=s $seq}
END { while(s!="") { c=substr(s,1,1); a[c]+=gsub(c,"",s) }
for(c in a) print c,a[c]}' file
Vous trouverez ci-dessous 4 implémentations basées sur la première fonction. Les deux premiers fonctionnent sur un awk standard, les deux derniers sur une version optimisée pour les fichiers fasta:
1. Lire la séquence et la traiter ligne par ligne:
bioawk -c fastx '{while ($seq!=""){ c=substr($seq,1,1);a[c]+=gsub(c,"",$seq) } }
END{ for(c in a) print c,a[c] }' file
2. concaténer toutes les séquences et les traiter à la fin:
awk '!/^>/{s=s $0 }
END {while(s!="") { c=substr(s,1,1); a[c]+=gsub(c,"",s) }
for(c in a) print c,a[c]}' file
3. Identique à 1 mais utilisez bioawk:
awk '!/^>/{s=$0; while(s!="") { c=substr(s,1,1); a[c]+=gsub(c,"",s) } }
END {for(c in a) print c,a[c]}' file
4. Identique à 2 mais utilisez bioawk:
function countCharacters2(str) {
n=length(str)
while(str != "") { c=substr(str,1,1); gsub(c"+","",str);
m=length(str); a[toupper[c]]+=n-m; n=m
}
}
Voici quelques résultats de synchronisation basés sur ce fichier fasta
function countCharacters(str) {
while(str != "") { c=substr(str,1,1); a[toupper[c]]+=gsub(c,"",str) }
}
est basé sur Brian Kernighan's awk qui est documenté dans " Le langage de programmation AWK ", par Al Aho, Brian Kernighan et Peter Weinberger (Addison-Wesley, 1988, ISBN 0-201-07981-X) . Je ne sais pas si cette version est compatible avec POSIX .
Merci. Malheureusement, ce n'est pas vraiment optimal C'est plus une solution O (n ^ 2) que O (n) (à cause de gsub )
Le comptage de gsub est une idée intelligente.
Wow, ce bioawk est vraiment efficace. Merci pour tout cela :)
Il m'est juste apparu que pour 1 économiser 0 $ dans une chaîne, puis match () / gsub () - utiliser CELA pourrait être plus efficace que d'opérer directement sur 0 $ puisque chaque fois que vous modifiez 0 $, vous pourriez invoquez à nouveau le fractionnement de champ (en fonction de votre implémentation awk - certains ne font que le fractionnement de champ si un champ est explicitement mentionné dans le script).
@EdMorton Je viens de lancer le test. Je n'ai pas vu de décalage horaire, mais vous avez raison. Avec une implémentation awk, ce sera différent. (mise à jour de la mise en œuvre)
Essayez cette solution Perl pour de meilleures performances.
time perl -lne'
if (not /^>/) { $l=$_; ++$kv{ substr($l,$_,1) } for 0..length($l) }
END { for ( sort keys %kv) { print "$_ $kv{$_}" }}
' chrY.fa
Encore une solution utilisant Perl, optimisée pour les performances.
if (not /^>/) { ++$kv{$_} for split //; }
Modifiez avec d'autres optimisations des performances
Toutes les durées indiquées ci-dessous sont des moyennes sur 3 à 5 exécutions sur un ordinateur, effectuées à peu près au même moment, mais permutées pour éviter des effets de cache prononcés.
Changer la boucle de style C for en for my $ i (0..length ($ _)) accélère la deuxième solution de 9,2 secondes à 6,8 secondes.
Ensuite, supprimer également un scalaire ( $ x ) à chaque opération, avec
if (not /^>/) { $l=$_; ++$kv{ substr($l,$_,1) } for 0..length($l) }
accélère cela à 5,3 secondes .
Réduction supplémentaire de l'utilisation des variables, en copiant $ _ et en libérant ainsi la boucle pour utiliser $_
if (not /^>/) { for $i (0..length($_)) { ++$kv{ substr($_,$i,1) } } }
a ajouté votre solution au calendrier: 22,745 secondes
Merci @kvantour .. J'ai ajouté une solution de plus en l'optimisant .. dans ma machine il a fallu 15,92 secondes pour votre gros fichier .. bioawk c'est super!
ajouté votre deuxième solution (13sec)
@kvantour Notez d'autres améliorations.
@ stack0114106 J'ai pris la liberté d'ajouter une section à votre réponse. Si cela ne vous convient pas, je vous présente mes excuses, veuillez me le faire savoir et je la diviserai en une réponse distincte. Je pensais que cela convenait le mieux ici car cela modifiait simplement votre code.
@zdim J'ajouterai le timing à la liste demain. Vous l'avez peut-être ajouté en tant que solution distincte. Belle explication
@zdim .. super recherche et optimisations..apprendre beaucoup .. vous êtes invités à ajouter ceci .. merci!
@kvantour " aurait pu l'ajouter en tant que solution séparée " - il modifie principalement le code de cette réponse, donc je l'ai trouvé le mieux adapté ici. (Un rappel pour l'ajouter à votre liste, si vous le souhaitez.)
@zdim a ajouté votre solution
déclarez-vous les comptes ligne par ligne ou une fois pour l'ensemble du fichier?
@karakfa pour le fichier entier
Copie possible du script Bash pour trouver le fréquence de chaque lettre d'un fichier