Algorithms for normalized multiple sequence alignments

Sequence alignment lies at the foundation of bioinformatics. Several procedures rely on alignment methods for a range of distinct purposes, such as detection of sequence homology, secondary structure prediction, phylogenetic analysis, identification of conserved motifs or genome assembly. On the other hand, alignment techniques have also been reshaped and found applications in other fields, such as natural language processing, pattern recognition, or social sciences [AT00, AG97, BL02, MV93].

Given its range of applications in bioinformatics, extensive efforts have been made to improve existing or developing novel methods for sequence alignment. The simpler ones compare a pair of sequences in polynomial time on their lengths, usually trying to find editing operations (insertions, deletions, and substitutions of symbols) that transform one sequence into another while maximizing or minimizing some objective function called edit distance [HAR09]. This concept can naturally be generalized to align multiple sequences [WLXZ15], adding another new layer of algorithmic complexity, though. In this case, most multiple sequence alignment (MSA) formulations lead to NP-hard problems [Eli06]. Nevertheless, a variety of methods suitable for aligning multiple sequences have been developed, as they can outperform pairwise alignment methods on tasks such as phylogenetic inference [OR06], secondary structure prediction [CB99] or identification of conserved regions [SWD⁺11].

In order to overcome the cost of exact solutions, a number of MSA heuristics have been developed in recent years, most of them using the so-called progressive or iterative methods [HTHI95, SH14, TPP99, WOHN06]. Experimental data suggest that the robustness and accuracy of heuristics can still be improved, however [WLXZ15].

Most approaches for pairwise sequence alignment define edit distances as absolute values, lacking some normalization that would result in edit distances relative to the lengths of the sequences. However, some applications may require sequence lengths to be taken into account. For instance, a difference of one symbol between sequences of length 5 is more significant than between sequences of length 1000. In addition, experiments suggest that normalized edit distances can provide better results than both unnormalized or post-normalized edit distances [MV93]. While normalized edit distances have been developed for pairwise sequence alignment [AE99, MV93], none have been proposed for MSA to the best of our knowledge.

In this work, we propose exact and approximation algorithms for normalized MSA (NMSA). This is a first step towards the development of methods that take into account the lengths of sequences for computing edit distances when multiple sequences are compared. The remainder of this paper is organized as follows. Section 2 introduces concepts related to sequence alignment and presents normalized scores for NMSA, followed by the complexity analysis of NMSA using those scores in Section 3. Next, Sections 4 and 5 describe exact and approximation algorithms, respectively. Section 6 closes the paper with the conclusion and prospects for future work.

An alphabet $\Sigma$ is a finite non-empty set of symbols. A finite sequence $s$ with $n$ symbols in $\Sigma$ is seen as $s(1)\cdots s(n)$. We say that the length of $s$, denoted by $\left|s\right|$, is $n$. The sequence $s(p)\cdots s(q)$, with $1\leq p\leq q\leq n$, is denoted by $s(p\!:\!q)$. If $p>q$, $s(p\!:\!q)$ is the empty sequence which length is equal to zero and it is denoted by $\varepsilon$. We denote the sequence resulting from the concatenation of sequences $s$ and $t$ by $s\cdot t=st$. A sequence of $n$ symbols $\mathtt{a}$ is denoted by $\mathtt{a}^{n}$. A $k$-tuple $S$ over $\Sigma^{*}$ is called a $k$-sequence and we write $s_{1},\ldots,s_{k}$ to refer to $S$, where $s_{i}$ is the $i$-th sequence in $S$. We denote $S=\emptyset$ if every sequence of $S$ is a empty sequence.

Let $\Sigma_{\texttt{-}}:=\Sigma\cup\{\texttt{-}\}$, where $\texttt{-}\not\in\Sigma$ and the symbol - is called a space. Let $S=s_{1},\ldots,s_{k}$ be a $k$-sequence. An alignment of $S$ is a $k$-tuple $A=[s^{\prime}_{1},\ldots,s^{\prime}_{k}]$ over $\Sigma_{\texttt{-}}^{*}$, where

1. ($a$)

   each sequence $s^{\prime}_{h}$ is obtained by inserting spaces in $s_{h}$;

2. ($b$)

   $\left|s^{\prime}_{h}\right|=\left|s^{\prime}_{i}\right|$ for each pair $h,i$, with $1\leq h,i\leq k$;

3. ($c$)

   there is no $j$ in $\{1,\ldots,k\}$ such that $s^{\prime}_{1}(j)=\ldots=s^{\prime}_{k}(j)=\texttt{-}$.

Notice that alignments, which are $k$-tuples over $\Sigma_{\texttt{-}}^{*}$, are written enclosed by square brackets “$[~]$”. The sequence $A(j)=[s^{\prime}_{1}(j),\ldots,s^{\prime}_{k}(j)]\in\Sigma_{\texttt{-}}^{k}$ is the column $j$ of the alignment $A=[s^{\prime}_{1},\ldots,s^{\prime}_{k}]$ and $A[j_{1}\!:\!j_{2}]$ is the alignment defined by the columns $j_{1},j_{1}+1,\ldots,j_{2}$ of $A$. We say that the pair $[s^{\prime}_{h}(j),s^{\prime}_{i}(j)]$ aligns in $A$ or, simply, that $s^{\prime}_{h}(j)$ and $s^{\prime}_{i}(j)$ are aligned in $A$, and $\left|A\right|=\left|s^{\prime}_{i}\right|$ is the length of the alignment $A$. It is easy to check that $\max_{i}\{\left|s_{i}\right|\}\leq\left|A\right|\leq\sum_{i}\left|s_{i}\right|$. We denote by $\mathcal{A}_{S}$ the set of all alignments of $S$.

An alignment can be used to represent editing operations of insertions, deletions and substitutions of symbols in sequences, where the symbol - represents insertions or deletions. An alignment can also be represented in the matrix format by placing one sequence above another. Thus, the alignments

of 3-sequence $\mathtt{a}\mathtt{a}\mathtt{a},\mathtt{a}\mathtt{b},\mathtt{c}\mathtt{a}\mathtt{c}$ can be represented respectively as

Let $I=\{i_{1},\ldots,i_{m}\}\subseteq\{1,\ldots,k\}$ be a set of indices such that $i_{1}<\cdots<i_{m}$ and let $A=[s^{\prime}_{1},\ldots,s^{\prime}_{k}]$ be an alignment of a $k$-sequence $S=s_{1},\ldots,s_{k}$. We write $S_{I}$ to denote the $m$-sequence $s_{i_{1}},\ldots,s_{i_{m}}$. The alignment of $S_{I}$ induced by $A$ is the alignment $A_{I}$ obtained from the alignment $A$, considering only the corresponding sequences in $S_{I}$ and, from the resulting structure, removing columns where all symbols are -. In the previous example, being $A=[\mathtt{a}\mathtt{a}\mathtt{a}\texttt{-},\mathtt{a}\mathtt{b}\texttt{-}\texttt{-},\texttt{-}\mathtt{c}\mathtt{a}\mathtt{c}]$ an alignment of $\mathtt{a}\mathtt{a}\mathtt{a},\mathtt{a}\mathtt{b},\mathtt{c}\mathtt{a}\mathtt{c}$, we have

is an alignment of $\mathtt{a}\mathtt{a}\mathtt{a},\mathtt{a}\mathtt{b}$ induced by $A$.

A $k$-vector $\vec{\jmath}=[j_{1},\ldots,j_{k}]$ is a $k$-tuple, where $j_{i}\in\mathbb{N}=\{0,1,2,\ldots\}$. We say that $j_{i}$ is the $i$-th element of $\vec{\jmath}$. The $k$-vector $\vec{0}$ is such that all its elements are zero. If $\vec{\jmath}$ and $\vec{h}$ are $k$-vectors, we write $\vec{\jmath}\leq\vec{h}$ if $j_{i}\leq h_{i}$ for each $i$; and $\vec{\jmath}<\vec{h}$ if $\vec{\jmath}\leq\vec{h}$ and $\vec{\jmath}\not=\vec{h}$. A sequence of $k$-vectors $\vec{\jmath}_{1},\vec{\jmath}_{2},\ldots$ is in topological order if $\vec{\jmath}_{i}<\vec{\jmath}_{h}$ implies $i<h$. Given two $k$-vectors, say $\vec{a}=[a_{1},a_{2}\ldots a_{k}]$ $\vec{b}=[b_{1},b_{2}\ldots b_{k}]$, we say that $\vec{a}$ precedes $\vec{b}$ when there exists $l\in\mathbb{N}$ such that $a_{i}=b_{i}$ for each $i>l$ and $a_{l}<b_{l}$. A sequence of $k$-vectors $\vec{\jmath}_{1},\vec{\jmath}_{2},\ldots$ is in lexicographical order if $\vec{\jmath}_{i}$ precedes $\vec{\jmath}_{h}$ for each $i<h$. Clearly, a lexicographical order is a special case of topological order.

Consider $S=s_{1},\ldots,s_{k}$ a $k$-sequence with $n_{i}=\left|s_{i}\right|$ for each $i$ and we call $\vec{n}=[n_{1},\ldots,n_{k}]$ the length of $S$. Let $V_{\vec{n}}=\{\vec{\jmath}:\vec{\jmath}\leq\vec{n}\}$. Therefore, $\left|V_{\vec{n}}\right|=\prod_{i}(\left|s_{i}\right|+1)$ which implies that if $n_{i}=n$ for all $i$, then $\left|V_{\vec{n}}\right|=(n+1)^{k}$. Define $S(\vec{\jmath})=s_{1}(j_{1}),\ldots,s_{k}(j_{k})$ a column $\vec{\jmath}$ in $S$ and we say that $S(1\!:\!\vec{\jmath})=s_{1}(1\!:\!j_{1}),\ldots,s_{k}(1\!:\!j_{k})$ is the prefix of $S$ ending in $\vec{\jmath}$. Thus, $S=S(1\!:\!\vec{n})$. Besides that, if $A$ is an alignment and $\vec{j}=[j,j,\ldots,j]$, then $A[\vec{j}]=A(j)$.

Denote by $\mathcal{B}^{{k}}$ the set of $k$-vectors $[b_{1},\ldots,b_{k}]$, where $b_{i}\in\{0,1\}$ for each $i$. Now, for $\vec{b}\leq\vec{j}$, where $\vec{j}=[j_{1},\ldots,j_{k}]$, define

such that

Therefore, given an alignment $A$ of $S(1\!:\!\vec{n})$, there exists $\vec{b}\in\mathcal{B}^{{k}}$, with $\vec{b}\leq\vec{n}$, such that $A(\left|A\right|)=\vec{b}\cdot S(\vec{n})$. In this case, notice that $b_{i}=1$ if and only if $s_{i}(n_{i})$ is in the $i$-th row of the last column of $A$ and we say that $\vec{b}$ defines the column $\left|A\right|$ of alignment $A$. For $\vec{b}\leq\vec{\jmath}$, we also define the operation

Notice that $\left|\mathcal{B}_{k}\right|=2^{k}$. Figure 1 shows an example of alignment using vectors to define columns and operations.

For a problem $\mathbf{P}$, we call $\mathbb{I}_{\mathbf{P}}$ the set of instances of $\mathbf{P}$. If $\mathbf{P}$ is a decision problem, then $\mathbf{P}(I)\in\{\texttt{Yes},\texttt{No}\}$ is the image of an instance $I$. If $\mathbf{P}$ is an optimization (minimization) problem, there is a set $\mathrm{Sol}(I)$ for each instance $I$, a function $v$ defining a non-negative rational number for each $X\in\mathrm{Sol}(I)$, and a function $\mathrm{opt}_{v}(I)=\min_{X\in\mathrm{Sol}(I)}\{v(X)\}$. We use $\mathrm{opt}$ instead of $\mathrm{opt}_{v}$ if $v$ is obvious. Let $\mathbf{A}(I)=v(X\in\mathrm{Sol}(I))$ a solution computed by an algorithm $\mathbf{A}$ with input $I$. We say that $\mathbf{A}$ is an $\alpha$-approximation for $\mathbf{P}$ if $\mathbf{A}(I)\leq\alpha\,\mathrm{opt}(I)$ for each $I\in\mathbb{I}_{\mathbf{P}}$. We say that $\alpha$ is an approximation factor for $\mathbf{P}$.

The alignment problem is a collection of decision and optimization problems whose instances are finite subsets of $\Sigma^{*}$ and $\mathrm{Sol}(S)=\mathcal{A}_{S}$ for each instance $S$. Function $v$, used for scoring alignments, is called criterion for $\mathbf{P}$ and we call $v[A]$ the cost of the alignment $A$. The $v$-optimal alignment $A$ of $S$ is such that $v[A]=\mathrm{opt}(S)$. Thus, we state the following general optimization problems using the criterion $v$ and integer $k$:

We also need the decision version of the alignment problem with criterion $v$ where we are given a $k$-sequence $S$ and a number $d\in\mathbb{Q}_{\geq}$, and we want to decide whether there exists an alignment $A$ of $S$ such that $v[A]\leq d$.

Usually the cost of an alignment $v$ is defined from a scoring matrix. A scoring matrix $\gamma$ is a matrix whose elements are rational numbers. The matrix has rows and columns indexed by symbols in $\Sigma_{\texttt{-}}$. For $\mathtt{a},\mathtt{b}\in\Sigma_{\texttt{-}}$ and a scoring matrix $\gamma$, we denote by ${\gamma}_{{\mathtt{a}}\rightarrow{\mathtt{b}}}$ the entry of $\gamma$ in row $\mathtt{a}$ and column $\mathtt{b}$. The value ${\gamma}_{{\mathtt{a}}\rightarrow{\mathtt{b}}}$ defines the score for a substitution if $\mathtt{a},\mathtt{b}\in\Sigma$, for an insertion if $\mathtt{a}=\texttt{-}$, and for a deletion if $\mathtt{b}=\texttt{-}$. The entry ${\gamma}_{{\texttt{-}}\rightarrow{\texttt{-}}}$ is not defined.

We study now the complexity of the multiple sequence alignment problem for each new criterion defined in Section 2. We consider the decision version of the computational problems and we prove $\mathbf{NMSA}\text{-}z$ is NP-complete for each $z$ even though the following additional restrictions for the scoring matrix $\gamma$ hold: ${\gamma}_{{\mathtt{a}}\rightarrow{\mathtt{b}}}={\gamma}_{{\mathtt{b}}\rightarrow{\mathtt{a}}}$ and ${\gamma}_{{\mathtt{a}}\rightarrow{\mathtt{b}}}=0$ if and only if $\mathtt{a}=\mathtt{b}$ for each pair $\mathtt{a},\mathtt{b}\in\Sigma_{\texttt{-}}$. Elias [Eli06] shows that, even considering such restrictions, $\mathbf{MSA}$ is NP-complete. We start showing a polynomial time reduction from $\mathbf{MSA}$ to $\mathbf{NMSA}\text{-}z$.

For an instance $(S,C)$ of $\mathbf{MSA}$ ($\mathbf{NMSA}\text{-}z$), where $S=s_{1},\ldots,s_{k}$ is a $k$-sequence and $C$ is an integer, $\mathbf{MSA}(S,C)$ ($\mathbf{NMSA}\text{-}z(S,C)$) is the decision problem version asking whether there exists an alignment $A$ of $S$ such that $v{\mathrm{SP}_{\!\gamma}}[A]\leq C$ (${\mathrm{V}^{i}_{\!\gamma}}[A]\leq C$).

Consider a fixed alphabet $\Sigma$ and scoring matrix $\gamma$ with the restrictions above that are ${\gamma}_{{\mathtt{a}}\rightarrow{\mathtt{b}}}={\gamma}_{{\mathtt{b}}\rightarrow{\mathtt{a}}}$ and ${\gamma}_{{\mathtt{a}}\rightarrow{\mathtt{b}}}=0$. We also assume each entry of $\gamma$ is integer. Let $\sigma\not\in\Sigma_{\texttt{-}}$ be a new symbol and $\Sigma^{\sigma}=\Sigma\cup\{\sigma\}$ and $G$ the maximum number in $\gamma$. We define a scoring matrix $\gamma^{\sigma}$ such that ${\gamma^{\sigma}\!\!\!}_{{\mathtt{a}}\rightarrow{\mathtt{b}}}={\gamma}_{{\mathtt{a}}\rightarrow{\mathtt{b}}},{\gamma^{\sigma}\!\!\!}_{{\mathtt{a}}\rightarrow{\sigma}}={\gamma^{\sigma}\!\!\!}_{{\sigma}\rightarrow{\mathtt{a}}}=G$ and ${\gamma^{\sigma}\!\!\!}_{{\sigma}\rightarrow{\sigma}}=0$ for each pair $\mathtt{a},\mathtt{b}\in\Sigma_{\texttt{-}}$. Notice that ${\gamma^{\sigma}}_{{\mathtt{a}}\rightarrow{\mathtt{b}}}={\gamma^{\sigma}}_{{\mathtt{b}}\rightarrow{\mathtt{a}}}$ and ${\gamma^{\sigma}}_{{\mathtt{a}}\rightarrow{\mathtt{b}}}=0$ if and only if $\mathtt{a}=\mathtt{b}$ for each pair $\mathtt{a},\mathtt{b}\in\Sigma^{\sigma}_{\texttt{-}}$ and therefore $G\geq 1$. Also, we denote $S^{L}=s_{1}\sigma^{L},\ldots,s_{k}\sigma^{L}$, where $L=Nk^{2}MG$, $M=\max_{i}\{\left|s_{i}\right|\}$ and $N=\binom{k}{2}M$.

Let $A$ be an alignment of $S^{L}$. A $\sigma$-column in $A$ is a column where every symbol is equal to $\sigma$. The tail of $A$ is the alignment $A[j+1:\left|A\right|]$ if each its column is $\sigma$-columns but $A(j)$ is not; in this case the tail length of $A$ is $\left|A\right|-j$ and the column $j$ is the tail base. We say that an alignment of $S^{L}$ is canonical if its tail length is $L$. If $A=[s_{1}^{\prime\prime},\ldots,s_{k}^{\prime\prime}]$ is the alignment of $S$, then we denote by $A^{L}$ the canonical alignment $[s_{1}^{\prime\prime}\sigma^{L},\ldots,s_{k}^{\prime\prime}\sigma^{L}]$ of $S^{L}$, i.e., $A^{L}[1:\left|A\right|]=A$ and $A^{L}[\left|A\right|+1:\left|A^{L}\right|]=[\sigma^{L},\sigma^{L},\ldots,\sigma^{L}]$. Notice that $\left|A^{L}\right|\geq L$. Then, we establish below a lower bound for $v{\mathrm{SP}_{\!\gamma^{\sigma}}}[A^{L}]$ and $v{\mathrm{SP}_{\!\gamma}}[A]$.