1
$\begingroup$

Note: this question was also asked on biostars

are you aware of any tool that is able to perform error-tolerant pattern-matching search on protein FASTA files?

For example, I want to know, which proteins in my fasta file to matchADNG..C.G regexp (which represents ADNGCG pattern). However, I want to be tolerant to matching errors, meaning that I'm good with any protein that differs in any 1 letter in the motif: ADNG..C.D, ADNG..D.G, MDNG..C.G, etc are all good. Running all possible variants through grep is possible, but exponentially long for longer patterns (I usually have 15-20 letters in pattern, and scan ~10e6-7 sequences).

I am aware of the tool agrep (docs. As far as I know, it's not supported anymore, and also does not distinguish letters in which it's allowed to make error. Also, it does not support long enough patterns (with more than 9 errors -- and yes, I tried to recompile it with myself from here and could not get it work). Also, it's not designed for proteins specifically.

I am also aware of ScanProsite and Protein Pattern Find services. The first one is not error-prone, the second one is seems to be a web interface for grep.

$\endgroup$
  • $\begingroup$ As correctly mentioned on biostars, where I cross-posted it, I should add cross reference to that topic: biostars.org/p/409111/#409142 $\endgroup$ – marinegor Nov 22 '19 at 2:22
  • $\begingroup$ Thanks for the link. I personally have no problem with cross-posting, especially if you're linking to other sites and acknowledging the source of any gained knowledge that's added here. $\endgroup$ – gringer Nov 22 '19 at 5:24
0
$\begingroup$

Python solution:
e<=2 will allow up to 2 insertions, 2 deletions, or 2 SNPs in any 2 combinations

import regex
seq_list=["GGGAGGG", "GGGAAGGG", "GGGAAAGGG", "TTTATTT", "TTTAATTT", "TTTAAATTT", "GGGG", "GGG"]
for sequence in seq_list:
    if regex.search("(GGGGGG|TTTTTT){e<=2}", sequence):
        print(sequence)

GGGAGGG
GGGAAGGG
TTTATTT
TTTAATTT
GGGG

Can also include a wildcard

seq_list=["GGGAAAAAGGG", "GGGAAAAATGG", "GGGAAAAATTG", "GGGAAAAATTT"]
for sequence in seq_list:
    if regex.search("(GGG.*GGG|TTTTTT){e<=2}", sequence):
        print(sequence)

GGGAAAAAGGG
GGGAAAAATGG
GGGAAAAATTG

| improve this answer | |
$\endgroup$
  • $\begingroup$ Agreed on the approach, it can be solved with an outer loop forcing s floating wildcard [AGTC] or [A-Z] across the regex search term $\endgroup$ – Michael Dec 7 '19 at 5:03
  • $\begingroup$ Thanks for the solution. The problem is, it is unadequately slow -- e.g. matching whole uniprot with that is very slow. Also, I am not aware of it's computational complexity for python regex fuzzy matching -- are you aware of its description? $\endgroup$ – marinegor Jan 8 at 14:51

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.