regex: matching several patterns derived from a simple string

Question

I have following task: Starting with 30 character long pattern sequence (it is actually DNA sequence, lest call it P30) I need to find in a text file all lines starting (^agacatacag... )with a exact P30, then with 29 last characters of the 30, 28 and up to 10 characters. All what I need to do is simply remove the first character of the pattern and go on searching. For the simplicity, I currently require exact matching, but allowing 1 mismatch for longer (20-30 character long patterns) would be better.

My current, rather slow solution is to create a shell file with one truncated pattern per line, and grep[1] it. Which means that I am reading huge, few GB text file 20x and this can take a day+.

I can switch to python, create a list/tuple with all required patterns and then read the file just once, looping instead 20x per each sequence, speeding things using pypy.

• Question 1: is there any regex which will be faster than such loop?
• Question 2: does it make sens to speed it up by switching to faster, compiled language? (I am trying to comprehend Dlang)

[1] Because it is DNA sequence and the input to be searched is in FASTQ format I am using fqgrep: https://github.com/indraniel/fqgrep with tre library: https://github.com/laurikari/tre/

edit_1 example of a changing (shortening pattern). Just first few steps/shorter patterns shown:

^abcde
^bcde
^cde

Or if you prefer it as DNA:

^GATACCA
^ATACCA
^TACCA

edit_2 Simple grep does not really cut it. I need to post-process each 4-lines of FASTQ format from which only line #2 matches. If I am not to use fqgrep, then I have to:

read 4 lines of the input
- check if line #2 (sequence) starts with any of the 20 patterns( P30-P10)
- if I got the match, I need to cut out the first N characters of the lines #2 and #4, where N stands for the length of the matching pattern - print out on output/write to file lines #1-$4 in no match do nothing

For the in-house solution I can try use GNU parallel splitting the input file in say 4M of lies chunks and speeding things up that way. But if I want to make it usable by others each new software I am asking end-users to install ads an extra level of complication.

** edit 3 ** Simple example of thee regexes and matching lines from Vyctor:

starting P30 regex 
^agacatacagagacatacagagacatacag
matching sequence:
^agacatacagagacatacagagacatacagGAGGACCA

P29: 
^gacatacagagacatacagagacatacag
matching sequence:
^gacatacagagacatacagagacatacagGACCACCA

P28: 
^acatacagagacatacagagacatacag
matching sequence:
^acatacagagacatacagagacatacagGATTACCA

I remove the characters/DNA bases from the left (or 5-prime end in DNA speak), because this is the way these sequences are degraded by real enzymes. The regex sequence by itself is not interesting once it is being found. The desired output is the read sequence after the regex. In the above examples it is in UPERCASE, which then can be mapped in the next step to the genome. It should be stressed that apart from this toy example I am getting way longer, a priori unknown and varied sequences after the regex pattern. In the real world, I do not have to deal with upper/lower case characters for DNA (everything is uppercase), but I am likely to encounter Ns (= unknown DNA base) in the sequences I am searching for the patterns. These can be ignored in the first approximation but for a more sensitive version of the algorithm probably should be dealt at as simple mismatches. In an ideal scenario one would not count simple mismatches at a given position but calculate more complex penalties taking into account DNA sequence quality values stored in line #4 of each 4 line long sequence record stored in FASTQ format: http://en.wikipedia.org/wiki/FASTQ_format#Quality

But this is way more complex, and so far the method "take only reads with perfect match to regex" was good enough and made the subsequent steps way easier to analyze.

Answer 1

You can generate a regex programmatically to look like below.
Its just a progressive alternation of line start or next character.

What this will give you is the ability to do a single pass search.
All you have to do is get the string length of the match to tell you
where you are.

Note - use Multi-line mode.

 #  (?:^|a)(?:^|g)(?:^|a)(?:^|c)(?:^|a)(?:^|t)(?:^|a)(?:^|c)(?:^|a)(?:^|g)(?!^)0123456789

 (?: ^ | a )      # P30
 (?: ^ | g )      # P29
 (?: ^ | a )      # P28
 (?: ^ | c )      # P27
 (?: ^ | a )      # P26
 (?: ^ | t )      # P25
 (?: ^ | a )      # P24
 (?: ^ | c )      # P23
 (?: ^ | a )      # P22
 (?: ^ | g )      # P21
                  # ..
                  # P11
 (?! ^)           # Not beginning of line
 0123456789       # P10 - P1

For example, matches these:

agacatacag0123456789
cag0123456789
gacatacag0123456789
acatacag0123456789
acag0123456789
catacag0123456789
catacag0123456789
0123456789
atacag0123456789
tacag0123456789
ag0123456789
g0123456789

But not these:

agaatacag0123456789
ca0123456789
gacataca0123456789
acaacag0123456789
acg0123456789
cataca0123456789
caacag0123456789
123456789
atacg0123456789
tcag0123456789
ag012456789
g012356789

Update

This is a graphic illustration that a single regex can replace all 30.
There is really no need for 30 separate regex's, all you need is 1 constant regex.
In this example the cluster group is replaced with capture groups so you can see what it is doing.

 # (^|G)(^|A)(^|T)(^|A)(^|C)(^|C)(?!^)A

 ( ^ | G )        # (1)
 ( ^ | A )        # (2)
 ( ^ | T )        # (3)
 ( ^ | A )        # (4)
 ( ^ | C )        # (5)
 ( ^ | C )        # (6)
 (?! ^ )          # Not beginning of line
 A

Input, 6 lines:

GATACCA
ATACCA
TACCA
ACCA
CCA
CA

Output:

 **  Grp 0 -  ( pos 0 , len 7 ) 
GATACCA  
 **  Grp 1 -  ( pos 0 , len 1 ) 
G  
 **  Grp 2 -  ( pos 1 , len 1 ) 
A  
 **  Grp 3 -  ( pos 2 , len 1 ) 
T  
 **  Grp 4 -  ( pos 3 , len 1 ) 
A  
 **  Grp 5 -  ( pos 4 , len 1 ) 
C  
 **  Grp 6 -  ( pos 5 , len 1 ) 
C  

------------------------------

 **  Grp 0 -  ( pos 9 , len 6 ) 
ATACCA  
 **  Grp 1 -  ( pos 9 , len 0 )  EMPTY 
 **  Grp 2 -  ( pos 9 , len 1 ) 
A  
 **  Grp 3 -  ( pos 10 , len 1 ) 
T  
 **  Grp 4 -  ( pos 11 , len 1 ) 
A  
 **  Grp 5 -  ( pos 12 , len 1 ) 
C  
 **  Grp 6 -  ( pos 13 , len 1 ) 
C  

------------------------------

 **  Grp 0 -  ( pos 17 , len 5 ) 
TACCA  
 **  Grp 1 -  ( pos 17 , len 0 )  EMPTY 
 **  Grp 2 -  ( pos 17 , len 0 )  EMPTY 
 **  Grp 3 -  ( pos 17 , len 1 ) 
T  
 **  Grp 4 -  ( pos 18 , len 1 ) 
A  
 **  Grp 5 -  ( pos 19 , len 1 ) 
C  
 **  Grp 6 -  ( pos 20 , len 1 ) 
C  

------------------------------

 **  Grp 0 -  ( pos 24 , len 4 ) 
ACCA  
 **  Grp 1 -  ( pos 24 , len 0 )  EMPTY 
 **  Grp 2 -  ( pos 24 , len 0 )  EMPTY 
 **  Grp 3 -  ( pos 24 , len 0 )  EMPTY 
 **  Grp 4 -  ( pos 24 , len 1 ) 
A  
 **  Grp 5 -  ( pos 25 , len 1 ) 
C  
 **  Grp 6 -  ( pos 26 , len 1 ) 
C  

------------------------------

 **  Grp 0 -  ( pos 30 , len 3 ) 
CCA  
 **  Grp 1 -  ( pos 30 , len 0 )  EMPTY 
 **  Grp 2 -  ( pos 30 , len 0 )  EMPTY 
 **  Grp 3 -  ( pos 30 , len 0 )  EMPTY 
 **  Grp 4 -  ( pos 30 , len 0 )  EMPTY 
 **  Grp 5 -  ( pos 30 , len 1 ) 
C  
 **  Grp 6 -  ( pos 31 , len 1 ) 
C  

------------------------------

 **  Grp 0 -  ( pos 35 , len 2 ) 
CA  
 **  Grp 1 -  ( pos 35 , len 0 )  EMPTY 
 **  Grp 2 -  ( pos 35 , len 0 )  EMPTY 
 **  Grp 3 -  ( pos 35 , len 0 )  EMPTY 
 **  Grp 4 -  ( pos 35 , len 0 )  EMPTY 
 **  Grp 5 -  ( pos 35 , len 0 )  EMPTY 
 **  Grp 6 -  ( pos 35 , len 1 ) 
C  

Answer 2

If I understand you correctly you have preset you want to match:

agacatacagagacatacagagacatacag

And then match lines that match:

re: agacatacagagacatacagagacatacag
30: agacatacagagacatacagagacatacag
29: agacatacagagacatacagagacatacac
28: agacatacagagacatacagagacataccc

You don't really need regexps for this, you just need to find the difference between two rows and since it's DNA I assume that string abcde and aebcd has difference of 4 because all sequences need to be at the right places.

If the order doesn't matter and you want to just search for rows that match at least at 28 characters you can just count differences in the strings.

reg = 'agacatacagagacatacagagacatacag'
for row in file:
    letters = diff_letters(reg, row.strip())
    if letters == 30:   # complete match
    elif letters == 29: # one different character
                        # so on

If you need match that actually starts with proper sequence, you can just get points of difference between strings and if first difference is at point >=28

reg = 'agacatacagagacatacagagacatacag'
for row in file:
    diffs = list(i for i,(a1,a2) in enumerate(zip(s1,s2)) if a1!=a2)
    if not len(diffs):
        difference = len(reg)
    else:
        difference = diffs[0]

    if difference == 30: # First difference is at last offset