Apologies!!!
I want to find a "pattern" in a sequence and count how many of those patterns occur in it. I know fairly that 'grep' can serve the purpose but I have a single file with multiple sequences in .fasta format. I want to have my output (number of patterns in each of the sequence) for individual sequences in a list. How can we do this using 'awk'?
Thanks a lot in advance for help!
3 answers
Perhaps not the most efficient, but this perl 1 liner should do it:
cat INPUT_FILE | perl -ne 'if($_ =~ />/){print "\n$_"; next;} else {chomp $_; print $_;}' | sed 1d | perl -ne 'if($_ =~ />/){print "$_"; next;} else {$count = () = $_ =~ /PATTERN/gi; print "$count\n"; $count=0;}'
Description:
Print sequence header followed by the corresponding sequence on next line (for each sequence):
perl -ne 'if($_ =~ />/){print "\n$_"; next;} else {chomp $_; print $_;}'
Remove the erroneous new line introduced above:
sed 1d
If header is encountered print else print the number of times the "pattern" is seen
perl -ne 'if($_ =~ />/){print "$_"; next;} else {$count = () = $_ =~ /PATTERN/gi; print "$count\n"; $count=0;}'
•
written
by
zlskidmore
Edit: I'm not Pierre either, so apologies for that!
You could probably use bioawk and a substitution count sub(/<pattern>/,"") to count the number of occurrences per sequence, assuming no two matches overlap.
If you would like to deal with overlaps, use a Perl script instead. my @matches = $string ~ m/pattern/g with a scalar(@matches) should do it.
Edit: Je suis pas Pierre (or something like that, I speak English and German).
Just linearize the file (google "fasta linearize" for awk and biopieces methods...or just see a couple here) and then pipe that to grep.
•
written
by
Devon Ryan
I think it'll just count once in that case, in which case your bioawk solution would be more appropriate.
•
written
by
Devon Ryan
Yeah, I suspect the best method is biopython/bioperl based, but I guess we should wait for Pierre to reply.
•
written
by
Devon Ryan
•
written
by
Pierre Lindenbaum
Haha.. This should be the official meme of Biostars for an accepted answer :)
•
written
by
poisonAlien
Log in to answer this question.
-
Conservation scores at each location of few motifs in MSA
written by wannabio •
-
Convert phased fasta sequences to IUPAC sequence
written by quadras •
-
Extract sequences from a fasta file with specific nucleotide repetition
written by shabbas12 •
-
How to extract a sequence below a fasta header (">") from a fasta file ?
written by v.berriosfarias •
-
How to remove just sequence ID's from a fasta file with multiple sequences
written by Nyksubuz •
-
Counting base and nucleotide frequency of multifasta file
written by saadleeshehreen
-
Pattren counter in fasta file
written by adityagt15 •
-
Insert a pattern in a sequence
written by jaqy •
-
filter out fasta File by pattern
written by Assa Yeroslaviz
-
Is there a way to split single .txt file with multiple fasta sequences into individual sequences in…
written by Alpha •
I wouldn't have thought grep would be capable of this unless your fasta has one line per sequence
I found seqinR does pattern searching very nicely for one single fasta sequence. But, I have multiple sequences (say 100) in a single file, and wants the results (no#of patterns found) of individual sequence separately in a tab delimited file.
Thanks you in advance for help!
pipe split output to grep
Why is this addressed to Pierre? I like the direct approach though, cut out the middle men :)
Don't worry about the comments, people were just having a laugh because you addressed Pierre specifically and it's funny because he probably is the best person to ask about awk. Nothing other than good-hearted fun was intended.