18
$\begingroup$

I have an array of 100,000 strings, all of length $k$. I want to compare each string to every other string to see if any two strings differ by 1 character. Right now, as I add each string to the array, I'm checking it against every string already in the array, which has a time complexity of $\frac{n(n-1)}{2} k$.

Is there a data structure or algorithm that can compare strings to each other faster than what I'm already doing?

Some additional information:

  • Order matters: abcde and xbcde differ by 1 character, while abcde and edcba differ by 4 characters.

  • For each pair of strings that differ by one character, I will be removing one of those strings from the array.

  • Right now, I'm looking for strings that differ by only 1 character, but it would be nice if that 1 character difference could be increased to, say, 2, 3, or 4 characters. However, in this case, I think efficiency is more important than the ability to increase the character-difference limit.

  • $k$ is usually in the range of 20-40.

$\endgroup$
  • 4
    $\begingroup$ Searching a string dictionary with 1 error is a fairly well-known problem, eg cs.nyu.edu/~adi/CGL04.pdf $\endgroup$ – KWillets Jun 25 '18 at 20:48
  • 1
    $\begingroup$ 20-40mers can use a fair bit of space. You might look at a Bloom filter (en.wikipedia.org/wiki/Bloom_filter) to test if degenerate strings — the set of all mers from one, two or more substitutions on a test mer — are "maybe-in" or "definitely-not-in" a set of kmers. If you get a "maybe-in", then further compare the two strings to determine whether or not it is a false positive. The "definitely-not-in" cases are true negatives that will reduce the overall number of letter-by-letter comparisons you have to do, by limiting comparisons to just the potential "maybe-in" hits. $\endgroup$ – Alex Reynolds Jun 27 '18 at 0:08
  • $\begingroup$ If you were working with a smaller range of k, you could use a bitset to store a hash table of booleans for all degenerate strings (e.g. github.com/alexpreynolds/kmer-boolean for toy example). For k=20-40, though, the space requirements for a bitset are simply too much. $\endgroup$ – Alex Reynolds Jun 27 '18 at 0:13

12 Answers 12

12
$\begingroup$

It's possible to achieve $O(nk \log k)$ worst-case running time.

Let's start simple. If you care about an easy to implement solution that will be efficient on many inputs, but not all, here is a simple, pragmatic, easy to implement solution that many suffice in practice for many situations. It does fall back to quadratic running time in the worst case, though.

Take each string and store it in a hashtable, keyed on the first half of the string. Then, iterate over the hashtable buckets. For each pair of strings in the same bucket, check whether they differ in 1 character (i.e., check whether their second half differs in 1 character).

Then, take each string and store it in a hashtable, this time keyed on the second half of the string. Again check each pair of strings in the same bucket.

Assuming the strings are well-distributed, the running time will likely be about $O(nk)$. Moreover, if there exists a pair of strings that differ by 1 character, it will be found during one of the two passes (since they differ by only 1 character, that differing character must be in either the first or second half of the string, so the second or first half of the string must be the same). However, in the worst case (e.g., if all strings start or end with the same $k/2$ characters), this degrades to $O(n^2 k)$ running time, so its worst-case running time is not an improvement on brute force.

As a performance optimization, if any bucket has too many strings in it, you can repeat the same process recursively to look for a pair that differ by one character. The recursive invocation will be on strings of length $k/2$.

If you care about worst-case running time:

With the above performance optimization I believe the worst-case running time is $O(nk \log k)$.

$\endgroup$
  • 3
    $\begingroup$ If $\Omega(n)$ strings share the same first half, which may very well happen in real life, then you haven't improved the complexity. $\endgroup$ – einpoklum Jun 25 '18 at 10:30
  • $\begingroup$ @einpoklum, sure! That's why I wrote the statement in my second sentence that it falls back to quadratic running time in the worst case, as well as the statement in my last sentence describing how to achieve $O(nk \log k)$ worst-case complexity if you care about the worst case. But I guess maybe I didn't express that very clearly -- so I've edited my answer accordingly. Is it better now? $\endgroup$ – D.W. Jun 25 '18 at 14:18
15
$\begingroup$

My solution is similar to j_random_hacker's but uses only a single hash set.

I would create a hash set of strings. For each string in the input, add to the set $k$ strings. In each of these strings replace one of the letters with a special character, not found in any of the strings. While you add them, check that they are not already in the set. If they are then you have two strings that only differ by (at most) one character.

An example with strings 'abc', 'adc'

For abc we add '*bc', 'a*c' and 'ab*'

For adc we add '*dc', 'a*c' and 'ad*'

When we add 'a*c' the second time we notice it is already in the set, so we know that there are two strings that only differ by one letter.

The total running time of this algorithm is $O(n*k^2)$. This is because we create $k$ new strings for all $n$ strings in the input. For each of those strings we need to calculate the hash, which typically takes $O(k)$ time.

Storing all the strings takes $O(n*k^2)$ space.

Further improvements

We can improve the algorithm further by not storing the modified strings directly but instead storing an object with a reference to the original string and the index of the character that is masked. This way we do not need to create all of the strings and we only need $O(n*k)$ space to store all of the objects.

You will need to implement a custom hash function for the objects. We can take the Java implementation as an example, see the java documentation. The java hashCode multiplies the unicode value of each character with $31^{k-i}$ (with $k$ the string length and $i$ the one-based index of the character. Note that each altered string only differs by one character from the original. We can easily compute the contribution of that character to the hash code. We can subtract that and add our masking character instead. This takes $O(1)$ to compute. This allows us to bring the total running time down to $O(n*k)$

$\endgroup$
  • 4
    $\begingroup$ @JollyJoker Yeah, space is something of a concern with this method. You could reduce space by not storing the modified strings, but instead storing an object with a reference to the string and the masked index. That should leave you with O(nk) space. $\endgroup$ – Simon Prins Jun 25 '18 at 12:29
  • $\begingroup$ To compute the $k$ hashes for each string in $O(k)$ time, I think you will need a special homemade hash function (e.g., compute the hash of the original string in $O(k)$ time, then XOR it with each of the deleted characters in $O(1)$ time each (though this is probably a pretty bad hash function in other ways)). BTW, this is pretty similar to my solution, but with a single hashtable instead of $k$ separate ones, and replacing a character with "*" instead of deleting it. $\endgroup$ – j_random_hacker Jun 25 '18 at 13:29
  • $\begingroup$ @SimonPrins With custom equals and hashCode methods that could work. Just creating the a*b-style string in those methods should make it bulletproof; I suspect some of the other answers here will have hash collision problems. $\endgroup$ – JollyJoker Jun 25 '18 at 14:04
  • 1
    $\begingroup$ @D.W. I modified my post to reflect the fact that calculating the hashes takes $O(k)$ time and added a solution to bring the total running time back down to $O(n*k)$. $\endgroup$ – Simon Prins Jun 25 '18 at 15:01
  • 1
    $\begingroup$ @SimonPrins Worst case could maybe be nk^2 due to String equality checking in hashset.contains when hashes collide. Of course, the worst case is when every string has the same exact hash, which would require a pretty much handcrafted set of strings, especially to get the same hash for *bc, a*c, ab*. I wonder if it could be shown impossible? $\endgroup$ – JollyJoker Jun 26 '18 at 7:15
7
$\begingroup$

I would make $k$ hashtables $H_1, \dots, H_k$, each of which has a $(k-1)$-length string as the key and a list of numbers (string IDs) as the value. The hashtable $H_i$ will contain all strings processed so far but with the character at position $i$ deleted. For example, if $k=6$, then $H_3[ABDEF]$ will contain a list of all strings seen so far that have the pattern $AB\cdot DEF$, where $\cdot$ means "any character". Then to process the $j$-th input string $s_j$:

  1. For each $i$ in the range 1 to $k$:
    • Form string $s_j'$ by deleting the $i$-th character from $s_j$.
    • Look up $H_i[s_j']$. Every string ID here identifies an original string that is either equal to $s$, or differs at position $i$ only. Output these as matches for string $s_j$. (If you wish to exclude exact duplicates, make the value type of the hashtables a (string ID, deleted character) pair, so that you can test for those that have had the same character deleted as we just deleted from $s_j$.)
    • Insert $j$ into $H_i$ for future queries to use.

If we store each hash key explicitly, then we must use $O(nk^2)$ space and thus have time complexity at least that. But as described by Simon Prins, it's possible to represent a series of modifications to a string (in his case described as changing single characters to *, in mine as deletions) implicitly in such a way that all $k$ hash keys for a particular string need just $O(k)$ space, leading to $O(nk)$ space overall, and opening the possibility of $O(nk)$ time too. To achieve this time complexity, we need a way to compute the hashes for all $k$ variations of a length-$k$ string in $O(k)$ time: for example, this can be done using polynomial hashes, as suggested by D.W. (and this is likely much better than simply XORing the deleted character with the hash for the original string).

Simon Prins's implicit representation trick also means that the "deletion" of each character is not actually performed, so we can use the usual array-based representation of a string without a performance penalty (rather than linked lists as I had originally suggested).

$\endgroup$
  • 2
    $\begingroup$ Nice solution. An example of a suitable bespoke hash function would be a polynomial hash. $\endgroup$ – D.W. Jun 25 '18 at 14:33
  • $\begingroup$ Thanks @D.W. Could you perhaps clarify a bit what you mean by "polynomial hash"? Googling the term didn't get me anything that seemed definitive. (Please feel free to edit my post directly if you want.) $\endgroup$ – j_random_hacker Jun 26 '18 at 11:45
  • 1
    $\begingroup$ Just read the string as a base $q$ number modulo $p$, where $p$ is some prime less than your hashmap size, and $q$ is a primitive root of $p$, and $q$ is more than the alphabet size. It's called "polynomial hash" because it is like evaluating the polynomial whose coefficients are given by the string at $q$. I'll leave it as an exercise to figure out how to compute all the desired hashes in $O(k)$ time. Note that this approach is not immune to an adversary, unless you randomly choose both $p,q$ satisfying the desired conditions. $\endgroup$ – user21820 Jun 26 '18 at 14:45
  • 1
    $\begingroup$ I think this solution can be further refined by observing that only one of the k hash tables needs to exist at any one time, thus reducing the memory requirement. $\endgroup$ – Michael Kay Jun 26 '18 at 14:48
  • 1
    $\begingroup$ @MichaelKay: That won't work if you want to compute the $k$ hashes of the possible alterations of a string in $O(k)$ time. You still need to store them somewhere. So if you only check one position at a time, you will take $k$ times as long as if you check all positions together using $k$ times as many hashtable entries. $\endgroup$ – user21820 Jun 26 '18 at 15:23
2
$\begingroup$

Here is a more robust hashtable approach than the polynomial-hash method. First generate $k$ random positive integers $r_{1..k}$ that are coprime to the hashtable size $M$. Namely, $0 \le r_i < M$. Then hash each string $x_{1..k}$ to $(\sum_{i=1}^k x_i r_i ) \bmod M$. There is almost nothing an adversary can do to cause very uneven collisions, since you generate $r_{1..k}$ on run-time and so as $k$ increases the maximum probability of collision of any given pair of distinct strings goes quickly to $1/M$. It is also obvious how to compute in $O(k)$ time all the possible hashes for each string with one character changed.

If you really want to guarantee uniform hashing, you can generate one random natural number $r(i,c)$ less than $M$ for each pair $(i,c)$ for $i$ from $1$ to $k$ and for each character $c$, and then hash each string $x_{1..k}$ to $(\sum_{i=1}^k r(i,x_i) ) \bmod M$. Then the probability of collision of any given pair of distinct strings is exactly $1/M$. This approach is better if your character set is relatively small compared to $n$.

$\endgroup$
2
$\begingroup$

A lot of the algorithms posted here use quite a bit of space on hash tables. Here's an $O(1)$ auxiliary storage $O((n \lg n) \cdot k^2)$ runtime simple algorithm.

The trick is to use $C_k(a, b)$, which is a comparator between two values $a$ and $b$ that returns true if $a < b$ (lexicographically) while ignoring the $k$th character. Then algorithm is as follows.

First, simply sort the strings regularly and do a linear scan to remove any duplicates.

Then, for each $k$:

  1. Sort the strings with $C_k$ as comparator.

  2. Strings that differ only in $k$ are now adjacent and can be detected in a linear scan.

$\endgroup$
1
$\begingroup$

Two strings of length k, differing in one character, share a prefix of length l and a suffix of length m such that k=l+m+1.

The answer by Simon Prins encodes this by storing all prefix/suffix combinations explicitly, i.e. abc becomes *bc, a*c and ab*. That's k=3, l=0,1,2 and m=2,1,0.

As valarMorghulis points out, you can organize words in a prefix tree. There's also the very similar suffix tree. It's fairly easy to augment the tree with the number of leaf nodes below each prefix or suffix; this can be updated in O(k) when inserting a new word.

The reason you want these sibling counts is so you know, given a new word, whether you want to enumerate all strings with the same prefix or whether to enumerate all strings with the same suffix. E.g. for "abc" as input, the possible prefixes are "", "a" and "ab", while the corresponding suffixes are "bc", "c" and "". As it obvious, for short suffixes it's better to enumerate siblings in the prefix tree and vice versa.

As @einpoklum points out, it's certainly possible that all strings share the same k/2 prefix. That's not a problem for this approach; the prefix tree will be linear up to depth k/2 with each node up to k/2 depth being the ancestor of 100.000 leaf nodes. As a result, the suffix tree will be used up to (k/2-1) depth, which is good because the strings have to differ in their suffixes given that they share prefixes.

[edit] As an optimization, once you've determined the shortest unique prefix of a string, you know that if there's one different character, it must be the last character of the prefix, and you'd have found the near-duplicate when checking a prefix that was one shorter. So if "abcde" has a shortest unique prefix "abc", that means there are other strings that start with "ab?" but not with "abc". I.e. if they'd differ in only one character, that would be that third character. You don't need to check for "abc?e" anymore.

By the same logic, if you would find that "cde" is a unique shortest suffix, then you know you need to check only the length-2 "ab" prefix and not length 1 or 3 prefixes.

Note that this method works only for exactly one character differences and does not generalize to 2 character differences, it relies one one character being the separation between identical prefixes and identical suffixes.

$\endgroup$
  • $\begingroup$ Are you suggesting that for each string $s$ and each $1 \le i \le k$, we find the node $P[s_1, \dots, s_{i-1}]$ corresponding to the length-$(i-1)$ prefix in the prefix trie, and the node $S[s_{i+1}, \dots, s_k]$ corresponding to the length-$(k-i-1)$ suffix in the suffix trie (each takes amortised $O(1)$ time), and compare the number of descendants of each, choosing whichever has fewer descendants, and then "probing" for the rest of the string in that trie? $\endgroup$ – j_random_hacker Jun 25 '18 at 14:05
  • 1
    $\begingroup$ What is the running time of your approach? It looks to me like in the worst case it might be quadratic: consider what happens if every string starts and ends with the same $k/4$ characters. $\endgroup$ – D.W. Jun 25 '18 at 14:24
  • $\begingroup$ The optimization idea is clever and interesting. Did you have in mind a particular way to do the check for mtaches? If the "abcde" has the shortest unique prefix "abc", that means we should check for some other string of the form "ab?de". Did you have in mind a particular way to do that, that will be efficient? What's the resulting running time? $\endgroup$ – D.W. Jun 25 '18 at 20:38
  • $\begingroup$ @D.W.: The idea is that to find strings in the form "ab?de", you check the prefix tree how many leaf nodes exist below "ab" and in the suffix tree how many nodes exist under "de", then choose the smallest of the two to enumerate. When all strings begin and end with the same k/4 characters; that means the first k/4 nodes in both trees have one child each. And yes, every time you need those trees, those have to be traversed which is an O(n*k) step. $\endgroup$ – MSalters Jun 26 '18 at 7:08
  • $\begingroup$ To check for a string of the form "ab?de" in the prefix trie, it suffices to get to the node for "ab", then for each of its children $v$, check whether the path "de" exists below $v$. That is, don't bother enumerating any other nodes in these subtries. This takes $O(ah)$ time, where $a$ is the alphabet size and $h$ is the height of the initial node in the trie. $h$ is $O(k)$, so if the alphabet size is $O(n)$ then it is indeed $O(nk)$ time overall, but smaller alphabets are common. The number of children (not descendants) is important, as well as the height. $\endgroup$ – j_random_hacker Jun 27 '18 at 11:53
1
$\begingroup$

Storing strings in buckets is a good way (there are already different answers outlining this).

An alternative solution could be to store strings in a sorted list. The trick is to sort by a locality-sensitive hashing algorithm. This is a hash algorithm which yields similar results when the input is similar[1].

Each time you want to investigate a string, you could calculate its hash and lookup the position of that hash in your sorted list (taking $O(log(n))$ for arrays or $O(n)$ for linked lists). If you find that the neighbours (considering all close neighbours, not only those with an index of +/- 1) of that position are similar (off by one character) you found your match. If there are no similar strings, you can insert the new string at the position you found (which takes $O(1)$ for linked lists and $O(n)$ for arrays).

One possible locality-sensitive hashing algorithm could be Nilsimsa (with open source implementation available for example in python).

[1]: Note that often hash algorithms, like SHA1, are designed for the opposite: producing greatly differing hashes for similar, but not equal inputs.

Disclaimer: To be honest, I would personally implement one of the nested/tree-organized bucket-solutions for a production application. However, the sorted list idea struck me as an interesting alternative. Note that this algorithm highly depends on the choosen hash algorithm. Nilsimsa is one algorithm I found - there are many more though (for example TLSH, Ssdeep and Sdhash). I haven't verified that Nilsimsa works with my outlined algorithm.

$\endgroup$
  • 1
    $\begingroup$ Interesting idea, but I think we would need to have some bounds on how far apart two hash values can be when their inputs differ by just 1 character -- then scan everything within that range of hash values, instead of just neighbours. (It's impossible to have a hash function that produces adjacent hash values for all possible pairs of strings that differ by 1 character. Consider the length-2 strings in a binary alphabet: 00, 01, 10 and 11. If h(00) is adjacent to both h(10) and h(01) then it must be between them, in which case h(11) can't be adjacent to them both, and vice versa.) $\endgroup$ – j_random_hacker Jun 25 '18 at 14:18
  • $\begingroup$ Looking at neighbors isn't sufficient. Consider the list abcd, acef, agcd. There exists a matching pair, but your procedure will not find it, as abcd is not a neighbor of agcd. $\endgroup$ – D.W. Jun 25 '18 at 14:27
  • $\begingroup$ You both are right! With neighbours I didn't mean only "direct neighbours" but thought of "a neighbourhood" of close positions. I didn't specify how many neighbours need to be looked at since that depends on the hash algorithm. But you're right, I should probably note this down in my answer. thanks :) $\endgroup$ – tessi Jun 25 '18 at 15:05
  • 1
    $\begingroup$ "LSH... similar items map to the same “buckets” with high probability" - since it's probability algorithm, result isn't guaranteed. So it depends on TS whether he needs 100% solution or 99.9% is enough. $\endgroup$ – Bulat Jun 26 '18 at 23:48
1
$\begingroup$

One could achieve the solution in $O(nk+ n^2)$ time and $O(nk)$ space using enhanced suffix arrays (Suffix array along with the LCP array) that allows constant time LCP (Longest Common Prefix) query (i.e. Given two indices of a string, what is the length of the longest prefix of the suffixes starting at those indices). Here, we could take advantage of the fact that all strings are of equal length. Specifically,

  1. Build the enhanced suffix array of all the $n$ strings concatenated together. Let $X = x_1.x_2.x_3 .... x_n$ where $x_i, \forall 1 \le i \le n $ is a string in the collection. Build the suffix array and LCP array for $X$.

  2. Now each $x_i$ starts at position $(i-1)k$ in the zero-based indexing. For each string $x_i$, take LCP with each of the string $x_j$ such that $j<i$. If LCP goes beyond the end of $x_j$ then $x_i = x_j$. Otherwise, there is a mismatch (say $x_i[p] \ne x_j[p]$); in this case take another LCP starting at the corresponding positions following the mismatch. If the second LCP goes beyond the end of $x_j$ then $x_i$ and $x_j$ differ by only one character; otherwise there are more than one mismatches.

    for (i=2; i<= n; ++i){
        i_pos = (i-1)k;
        for (j=1; j < i; ++j){
            j_pos = (j-1)k;
            lcp_len = LCP (i_pos, j_pos);
            if (lcp_len < k) { // mismatch
                if (lcp_len == k-1) { // mismatch at the last position
                // Output the pair (i, j)
                }
                else {
                  second_lcp_len = LCP (i_pos+lcp_len+1, j_pos+lcp_len+1);
                  if (lcp_len+second_lcp_len>=k-1) { // second lcp goes beyond
                    // Output the pair(i, j)
                  }
                }
            }
        }
    }
    

You could use SDSL library to build the suffix array in compressed form and answer the LCP queries.

Analysis: Building the enhanced suffix array is linear in the length of $X$ i.e. $O(nk)$. Each LCP query takes constant time. Thus, querying time is $O(n^2)$.

Generalisation: This approach can also be generalised to more than one mismatches. In general, running time is $O(nk + qn^2)$ where $q$ is the number of allowed mismatches.

If you wish to remove a string from the collection, instead of checking every $j<i$, you could keep a list of only 'valid' $j$.

$\endgroup$
  • $\begingroup$ Can i say that $O(kn^2)$ algo is trivial - just compare each string pair and count number of matches? And k in this formula practically can be omitted, since with SSE you can count matching bytes in 2 CPU cycles per 16 symbols (i.e. 6 cycles for k=40). $\endgroup$ – Bulat Jun 28 '18 at 2:17
  • $\begingroup$ Apologies but I could not understand your query. The above approach is $O(nk + n^2)$ and not $O(kn^2)$. Also, it is virtually alphabet-size independent. It could be used in conjunction with the hash-table approach -- Once two strings are found to have the same hashes, they could be tested if they contain a single mismatch in $O(1)$ time. $\endgroup$ – Ritu Kundu Jun 29 '18 at 7:13
  • $\begingroup$ My point is that k=20..40 for the question author and comparing such small strings require only a few CPU cycles, so practical difference between brute force and your approach probably doesn't exist. $\endgroup$ – Bulat Jun 29 '18 at 7:29
1
$\begingroup$

One improvement to all the solutions proposed. They all require $O(nk)$ memory in the worst case. You can reduce it by computing hashes of strings with * instead each character, i.e. *bcde, a*cde... and processing at each pass only variants with hash value in certain integer range. F.e. with even hash values in the first pass, and odd hash values in the second one.

You can also use this approach to split the work among multiple CPU/GPU cores.

$\endgroup$
  • $\begingroup$ Clever suggestion! In this case, the original question says $n=100,000$ and $k\approx 40$, so $O(nk)$ memory doesn't seem likely to be an issue (that might be something like 4MB). Still a good idea worth knowing if one needs to scale this up, though! $\endgroup$ – D.W. Jun 27 '18 at 0:24
0
$\begingroup$

This is a short version of @SimonPrins' answer not involving hashes.

Assuming none of your strings contain an asterisk:

  1. Create a list of size $nk$ where each of your strings occurs in $k$ variations, each having one letter replaced by an asterisk (runtime $\mathcal{O}(nk^2)$)
  2. Sort that list (runtime $\mathcal{O}(nk^2\log nk)$)
  3. Check for duplicates by comparing subsequent entries of the sorted list (runtime $\mathcal{O}(nk^2)$)

An alternative solution with implicit usage of hashes in Python (can't resist the beauty):

def has_almost_repeats(strings,k):
    variations = [s[:i-1]+'*'+s[i+1:] for s in strings for i in range(k)]
    return len(set(variations))==k*len(strings)
$\endgroup$
  • $\begingroup$ Thanks. Please also mention the $k$ copies of exact duplicates, and I'll +1. (Hmm, just noticed I made the same claim about $O(nk)$ time in my own answer... Better fix that...) $\endgroup$ – j_random_hacker Jun 26 '18 at 9:06
  • $\begingroup$ @j_random_hacker I don't know what exactly the OP wants reported, so I left step 3 vague but I think it is trivial with some extra work to report either (a) a binary any duplicate/no duplicates result or (b) a list of pairs of strings that differ in at most one position, without duplicates. If we take the OP literally ("...to see if any two strings..."), then (a) seems to be desired. Also, if (b) were desired then of course simply creating a list of pairs may take $\mathcal{O}(n^2)$ if all strings are equal $\endgroup$ – Bananach Jun 26 '18 at 9:15
0
$\begingroup$

Here is my take on 2+ mismatches finder. Note that in this post I consider each string as circular, f.e. substring of length 2 at index k-1 consists of symbol str[k-1] followed by str[0]. And substring of length 2 at index -1 is the same!

If we have M mismatches between two strings of length k, they have matching substring with length at least $mlen(k,M) = \lceil{k/M}\rceil-1$ since, in the worst case, mismatched symbols split (circular) string into M equal-sized segments. F.e. with k=20 and M=4 the "worst" match may have the pattern abcd*efgh*ijkl*mnop*.

Now, the algorithm for searching all mismatches up to M symbols among strings of k symbols:

  • for each i from 0 to k-1
    • split all strings into groups by str[i..i+L-1], where L = mlen(k,M). F.e. if L=4 and you have alphabet of only 4 symbols (from DNA), this will make 256 groups.
    • Groups smaller than ~100 strings can be checked with brute-force algorithm
    • For larger groups, we should perform secondary division:
      • Remove from every string in the group L symbols we already matched
      • for each j from i-L+1 to k-L-1
        • split all strings into groups by str[i..i+L1-1], where L1 = mlen(k-L,M). F.e. if k=20, M=4, alphabet of 4 symbols, so L=4 and L1=3, this will make 64 groups.
        • the rest is left as exercise for the reader :D

Why we don't start j from 0? Because we already made these groups with the same value of i, so job with j<=i-L will be exactly equivalent to job with i and j values swapped.

Further optimizations:

  • At every position, also consider strings str[i..i+L-2] & str[i+L]. This only doubles amount of jobs created, but allows to increase L by 1 (if my math is correct). So, f.e. instead of 256 groups, you will split data into 1024 groups.
  • If some $L[i]$ becomes too small, we can always use the * trick: for each i in in 0..k-1, remove i'th symbol from each string and create job searching for M-1 mismatches in those strings of length k-1.
$\endgroup$
0
$\begingroup$

I work everyday on inventing and optimizing algos, so if you need every last bit of performance, that is the plan:

  • Check with * in each position independently, i.e. instead of single job processing n*k string variants - start k independent jobs each checking n strings. You can spread these k jobs among multiple CPU/GPU cores. This is especially important if you are going to check 2+ char diffs. Smaller job size will also improve cache locality, which by itself can make program 10x faster.
  • If you are going to use hash tables, use your own implementation employing linear probing and ~50% load factor. It's fast and pretty easy to implement. Or use an existing implementation with open addressing. STL hash tables are slow due to use of separate chaining.
  • You may try to prefilter data using 3-state Bloom filter (distinguishing 0/1/1+ occurrences) as proposed by @AlexReynolds.
  • For each i from 0 to k-1 run the following job:
    • Generate 8-byte structs containing 4-5 byte hash of each string (with * at i-th position) and string index, and then either sort them or build hash table from these records.

For sorting, you may try the following combo:

  • first pass is MSD radix sort in 64-256 ways employing TLB trick
  • second pass is MSD radix sort in 256-1024 ways w/o TLB trick (64K ways total)
  • third pass is insertion sort to fix remaining inconsistencies
$\endgroup$

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.