<p><strong>Upd:</strong> Space complexity improved to just $O(2^LML)$ which makes the solution much more practical.</p>

<p>The problem can be solved at least in $O(M^2NL + 2^L MNL)$ time and $O(2^L ML)$ additional memory as follows.</p>

<p>Let's denote by $S^j$ the prefix $s_1 s_2 \ldots s_j$ of $S$. First, for each binary string $s$ of length $L$ precompute mismatch profile $mp&#39;(s, i, j, k) = mp(s, S_i^j, k)$ of that string with all prefixes of length $j$ of all strings $S_1, S_2, \ldots, S_M$. That will require $O(2^L MNL)$ memory to store and can be done in $O(2^L MNL)$ in straightforward way.</p>

<p>Then iteratively build mismatch profiles $mp_j(i,i&#39;,k)$ for each pair of $j$-prefixes of $S_1, S_2, \ldots, S_M$ for $j=L,L+1,\ldots,N$, so $mp_N(i,i&#39;,k)$ will be your desired answer. Let's denote by $end(S)$ the suffix of $S$ of the length $L$ (i.e. $S[(length(S)-L+1) \ldots length(S)]$). Notice that $$
\begin{align*}
mp_j(i,i&#39;,k) = &amp;mp(S_i^{j-1},S_{i&#39;}^{j-1},k) + mp(end(S_i^j), S_{i&#39;}^{j-1}, k) + \\&amp;mp(S_i^{j-1}, end(S_{i&#39;}^j),  k) + mp(end(S_i^j), end(S_{i&#39;}^j),k)
\end{align*}
$$ which is indeed
$$
\begin{align*}
&amp;mp_{j-1}(i,i&#39;,k) + mp&#39;(end(S_i^j), i&#39;, j-1, k) + \\ &amp;mp&#39;(end(S_{i&#39;}^j), i, j-1, k) + mp(end(S_i^j), end(S_{i&#39;}^j),k)
\end{align*}

$$ </p>

<p>The first 3 terms of that sum are already computed and the last one can be computed in amortized $O(1)$ time (it is $O(L)$ but you need to compute it only once for all values of $k$ with fixed $j,i,i'$). Thus you can compute $mp_j$ using $mp_{j-1}$ and $mp&#39;$ in $O(M^2L)$ time, which leads to $O(M^2NL)$ term in the total complexity.</p>

<p>Notice that you need values of $mp&#39;$ only for a fixed prefix length on each iteration. Thus you can compute $mp&#39;$ not all of a piece at the beginning but rather iteratively for $j=L,L+1, \ldots, N$ in between of $mp_j$ computations. That can be easily done in $O(2^LML)$ time per iteration (so total complexity does not change) but requires only $O(2^LML)$ additional memory.</p>

<p>The main issue of this solution is exponential in $L$ time.</p>

<p>P.S. If $L$ is too small ($2^L &lt; N$) then the naive method can be improved to $O(M^2 4^L)$ time by calculating the number of occurrences of each string $s$ of length $L$ in $S_i$. This allows you to process all different substrings of $S$ of length $L$ in $O(2^L)$ time instead of $O(N)$ (you will need to precompute pairwise differences between all such substrings first).</p>

CC BY-SA 3.0
2012-02-13T09:50:38.133Z
Actually it can be done in $O(2^LML)$ additional memory!
Dmytro Korduban
CC BY-SA 3.0
2012-02-12T11:30:04.737Z
2012-02-13T10:09:05.733Z
2012-02-13T10:09:05.733Z
Dmytro Korduban
<p>the details would be slightly tricky to work out but (as a rough sketch) there appears to be a recursive simplification that if you have computed all the mismatches of substrings of size n-1 then you can use that to compute mismatches of strings of size n. </p>

<p>your problem also seems to have connections to boolean convolutions (ie using a sliding window and the mod 2 sum) &amp; also matrix rigidity (originally proposed by Valiant)</p>

<p>if you are actually working with DNA bases as you indicate with your comment about ACGT, then why are you reformulating it in binary? it would probably be better to keep the ACGT formulation. the xor problem is not exactly equivalent and it would give wrong answers for subsequences off by one digit if you're really interested in a 4-symbol ACGT matching system.</p>

<p>you do not state this directly but from my knowledge of bioinformatics it would seem you might be interested in finding long subsequences such that there are low mismatches. if you add that criteria it significantly changes the optimization problem because one does not have to consider all mismatches, only those that are "better" than the current best record.</p>

<p>in that case there may be some similarity to this problem and the well studied problem of finding the <a href="http://en.wikipedia.org/wiki/Longest_common_substring_problem" rel="nofollow">longest matching substring</a> in a string, with your case interested in the "longest <em>near</em> matching substring". similar mechanisms for that problem like the suffix tree may prove useful</p>

CC BY-SA 3.0
2012-02-11T19:30:03.557Z
2012-02-11T19:35:53.457Z
2012-02-11T19:35:53.457Z
vzn
<p>We have M=10000 binary sequences of length N=1000.</p>

<p>given length L=15, for each pair of sequences, $S_1$ and $S_2$, we define the mismatch profile, mp($S_1$,$S_2$,$L$)[$m$], for m=0,1,...,L as following: </p>

<p>mp($S_1$,$S_2$,$L$)[$m$] = number of all substrings of length $L$ in $S_1$ and $S_2$ that exactly differ in m positions. </p>

<p>for example:</p>

<p>$\begin{array}{ccc}
L=5    \\
N=7    \\
S1=&#39;1010101&#39; \\
S2=&#39;1011000&#39; 
 
\end{array}$</p>

<p>set of all L-mers in $S_1$ is { 10101, 01010, 10101} and 
set of all L-mers in $S_2$ is { 10110, 01100, 11000}</p>

<p>All the pairs of L-mers (one from S1 and the other from S2) are (with the number of mismatches for each pair in the parenthesis): </p>

<p>$\begin{array}{ccc}
10101 &amp; 10110 &amp; (2) \\
10101 &amp; 01100 &amp; (3) \\
10101 &amp; 11000 &amp; (3) \\

01010 &amp; 10110 &amp; (3) \\
01010 &amp; 01100 &amp; (2) \\
01010 &amp; 11000 &amp; (2) \\

10101 &amp; 10110 &amp; (2) \\
10101 &amp; 01100 &amp; (3) \\
10101 &amp; 11000 &amp; (3) 
\end{array}$</p>

<p>hence the mismatch profile for $S_1$ and $S_2$ for this example is: 
mp($S_1$,$S_2$,$L$)[0]= 0<br>
mp($S_1$,$S_2$,$L$)[1]= 0
mp($S_1$,$S_2$,$L$)[2]= 4
mp($S_1$,$S_2$,$L$)[3]= 5
mp($S_1$,$S_2$,$L$)[4]= 0
mp($S_1$,$S_2$,$L$)[5]= 0</p>

<p>what is a fast algorithm to calculate the mismatch profile for every pair of sequences. The trivial algorithm would be $O(N^2 M^2)$ . can we make this better ? </p>

CC BY-SA 3.0
2012-02-10T08:57:25.650Z
I think you mean "substring" not "subsequence".  11111 is a subsequence of 0010010001100000100, but not a substring.
Jeffε
CC BY-SA 3.0
2012-02-10T08:59:31.233Z
Also, this looks like a reasonable algorithms homework problem.  Can you give some motivation/background?
Jeffε
CC BY-SA 3.0
2012-02-10T12:12:53.383Z
yes I meant substring.  no it is not a hw problem. We use this in gapped-kmer kernel calculation of an SVM classification which has a lot of applications. The simplification I made here is that the sequences are not binary, they are usually in base four defined over an alphabet of {A,C,G,T}
mghandi
CC BY-SA 3.0
2012-02-10T14:34:39.103Z
I assume you have at least $O(M^2 L)$ memory (at least for the answer), right? Is $L$ typically small? And on somewhat unrelated topic: was my answer on your previous question helpful?
Dmytro Korduban
CC BY-SA 3.0
2012-02-10T16:19:01.227Z
there appears to be a recursive simplication that if you have computed all the mismatches of substrings of size n-1 then you can use that to compute mismatches of strings of size n.
vzn
CC BY-SA 3.0
2012-02-10T17:12:01.193Z
@maldini.ua: yes L is topically small. Your answer to the previous question was great. But as you mentioned for small L, it does not have big effect. In fact I used xor to calculate the number of mismatches and it is pretty fast. Thanks for your answer anyways. It was very interesting. I am wondering if we can apply similar idea to get rid of N^2 and make it NlogN here
mghandi
<ds.algorithms><ds.data-structures><bioinformatics>
an efficient algorithm for mismatch profile
mghandi
# an efficient algorithm for mismatch profile

(Asked by: username_0 on 2012-02-10)

We have M=10000 binary sequences of length N=1000.

given length L=15, for each pair of sequences, $S\_1$ and $S\_2$, we define the mismatch profile, mp($S\_1$,$S\_2$,$L$)[$m$], for m=0,1,...,L as following:

mp($S\_1$,$S\_2$,$L$)[$m$] = number of all substrings of length $L$ in $S\_1$ and $S\_2$ that exactly differ in m positions.

for example:

$\begin{array}{ccc}
L=5 \\
N=7 \\
S1='1010101' \\
S2='1011000'
\end{array}$

set of all L-mers in $S\_1$ is { 10101, 01010, 10101} and
set of all L-mers in $S\_2$ is { 10110, 01100, 11000}

All the pairs of L-mers (one from S1 and the other from S2) are (with the number of mismatches for each pair in the parenthesis):

$\begin{array}{ccc}
10101 & 10110 & (2) \\
10101 & 01100 & (3) \\
10101 & 11000 & (3) \\
01010 & 10110 & (3) \\
01010 & 01100 & (2) \\
01010 & 11000 & (2) \\
10101 & 10110 & (2) \\
10101 & 01100 & (3) \\
10101 & 11000 & (3)
\end{array}$

hence the mismatch profile for $S\_1$ and $S\_2$ for this example is:
mp($S\_1$,$S\_2$,$L$)[0]= 0  
mp($S\_1$,$S\_2$,$L$)[1]= 0
mp($S\_1$,$S\_2$,$L$)[2]= 4
mp($S\_1$,$S\_2$,$L$)[3]= 5
mp($S\_1$,$S\_2$,$L$)[4]= 0
mp($S\_1$,$S\_2$,$L$)[5]= 0

what is a fast algorithm to calculate the mismatch profile for every pair of sequences. The trivial algorithm would be $O(N^2 M^2)$ . can we make this better ?

**Tags:** ds.algorithms, ds.data-structures, bioinformatics

> **Comments:**
> **username_1**: I think you mean "substring" not "subsequence".  11111 is a subsequence of 0010010001100000100, but not a substring.  
> **username_1**: Also, this looks like a reasonable algorithms homework problem.  Can you give some motivation/background?  
> **username_0**: yes I meant substring.  no it is not a hw problem. We use this in gapped-kmer kernel calculation of an SVM classification which has a lot of applications. The simplification I made here is that the sequences are not binary, they are usually in base four defined over an alphabet of {A,C,G,T}  
> **username_2**: I assume you have at least $O(M^2 L)$ memory (at least for the answer), right? Is $L$ typically small? And on somewhat unrelated topic: was my answer on your previous question helpful?  
> **username_3**: there appears to be a recursive simplication that if you have computed all the mismatches of substrings of size n-1 then you can use that to compute mismatches of strings of size n.  
> **username_0**: @maldini.ua: yes L is topically small. Your answer to the previous question was great. But as you mentioned for small L, it does not have big effect. In fact I used xor to calculate the number of mismatches and it is pretty fast. Thanks for your answer anyways. It was very interesting. I am wondering if we can apply similar idea to get rid of N^2 and make it NlogN here  

---

## 2 answers

### Answer by username_2 on 2012-02-12. Score: 4 | Accepted answer

**Upd:** Space complexity improved to just $O(2^LML)$ which makes the solution much more practical.

The problem can be solved at least in $O(M^2NL + 2^L MNL)$ time and $O(2^L ML)$ additional memory as follows.

Let's denote by $S^j$ the prefix $s\_1 s\_2 \ldots s\_j$ of $S$. First, for each binary string $s$ of length $L$ precompute mismatch profile $mp'(s, i, j, k) = mp(s, S\_i^j, k)$ of that string with all prefixes of length $j$ of all strings $S\_1, S\_2, \ldots, S\_M$. That will require $O(2^L MNL)$ memory to store and can be done in $O(2^L MNL)$ in straightforward way.

Then iteratively build mismatch profiles $mp\_j(i,i',k)$ for each pair of $j$-prefixes of $S\_1, S\_2, \ldots, S\_M$ for $j=L,L+1,\ldots,N$, so $mp\_N(i,i',k)$ will be your desired answer. Let's denote by $end(S)$ the suffix of $S$ of the length $L$ (i.e. $S[(length(S)-L+1) \ldots length(S)]$). Notice that $$
\begin{align\*}
mp\_j(i,i',k) = ∓(S\_i^{j-1},S\_{i'}^{j-1},k) + mp(end(S\_i^j), S\_{i'}^{j-1}, k) + \\∓(S\_i^{j-1}, end(S\_{i'}^j), k) + mp(end(S\_i^j), end(S\_{i'}^j),k)
\end{align\*}
$$ which is indeed
$$
\begin{align\*}
∓\_{j-1}(i,i',k) + mp'(end(S\_i^j), i', j-1, k) + \\ ∓'(end(S\_{i'}^j), i, j-1, k) + mp(end(S\_i^j), end(S\_{i'}^j),k)
\end{align\*}
$$

The first 3 terms of that sum are already computed and the last one can be computed in amortized $O(1)$ time (it is $O(L)$ but you need to compute it only once for all values of $k$ with fixed $j,i,i'$). Thus you can compute $mp\_j$ using $mp\_{j-1}$ and $mp'$ in $O(M^2L)$ time, which leads to $O(M^2NL)$ term in the total complexity.

Notice that you need values of $mp'$ only for a fixed prefix length on each iteration. Thus you can compute $mp'$ not all of a piece at the beginning but rather iteratively for $j=L,L+1, \ldots, N$ in between of $mp\_j$ computations. That can be easily done in $O(2^LML)$ time per iteration (so total complexity does not change) but requires only $O(2^LML)$ additional memory.

The main issue of this solution is exponential in $L$ time.

P.S. If $L$ is too small ($2^L < N$) then the naive method can be improved to $O(M^2 4^L)$ time by calculating the number of occurrences of each string $s$ of length $L$ in $S\_i$. This allows you to process all different substrings of $S$ of length $L$ in $O(2^L)$ time instead of $O(N)$ (you will need to precompute pairwise differences between all such substrings first).

> **Comments:**
> **username_2**: Actually it can be done in $O(2^LML)$ additional memory!  

---

### Answer by username_3 on 2012-02-11. Score: 0

the details would be slightly tricky to work out but (as a rough sketch) there appears to be a recursive simplification that if you have computed all the mismatches of substrings of size n-1 then you can use that to compute mismatches of strings of size n.

your problem also seems to have connections to boolean convolutions (ie using a sliding window and the mod 2 sum) & also matrix rigidity (originally proposed by Valiant)

if you are actually working with DNA bases as you indicate with your comment about ACGT, then why are you reformulating it in binary? it would probably be better to keep the ACGT formulation. the xor problem is not exactly equivalent and it would give wrong answers for subsequences off by one digit if you're really interested in a 4-symbol ACGT matching system.

you do not state this directly but from my knowledge of bioinformatics it would seem you might be interested in finding long subsequences such that there are low mismatches. if you add that criteria it significantly changes the optimization problem because one does not have to consider all mismatches, only those that are "better" than the current best record.

in that case there may be some similarity to this problem and the well studied problem of finding the [longest matching substring](http://en.wikipedia.org/wiki/Longest_common_substring_problem) in a string, with your case interested in the "longest *near* matching substring". similar mechanisms for that problem like the suffix tree may prove useful


---

