Multiple Sequence Alignment

Contributed by R. Gutman

Basic Terms:

DNA Sequence:

A sequence of nucleotide bases in the DNA molecule.

Motif:

Conserved pattern of amino acids that is found in one or more proteins, and

has biological significance.

PAM (Point Accepted Mutation) / BLOSUM (BLOcks SUbstitution Matrix) Matrices

Substitution matrices used for sequence alignment of proteins. Both matrices estimate the rate at which each possible residue in a sequence changes to each other residue over time.

Phylogenetic Tree:

A binary tree showing the represnting the relationship among a group of related nucleic acid or protien sequences from various/same species or other entities that are believed to have a common ancestor. In a phylogenetic tree, each node with descendants represents the most recent common ancestor of the descendants, and edge lengths correspond to the numbers of changes between adjacent node in a tree.

Protien Sequence:

A sequence of amino acids which make up given segment of a protein molecule.

Introduction:

What is Sequence Alignment?

A sequence alignment is a way of arranging the sequences of DNA, RNA, or protein to identify regions of similarity that may be a consequence of functional, structural, or evolutionary relationships between the sequences. Those similar sequences that can emerge from different organism, often perform similar or even identical functions and at other times, mutating or rearranging to perform an altered function through the forces of natural selection.

Global vs Local alignment:

Global alignment is a form of global optimization that "forces" the alignment to span the entire length of all query sequences. By contrast, local alignments identify regions of similarity within long sequences that are often widely divergent overall. Local alignments are often preferable, but can be more difficult to calculate because of the additional challenge of identifying the regions of similarity.

Why do we need Multiple sequence alignment (msa)?

Through simultaneous alignment of DNA or protein sequences, one can reveal whether there is a functional and evolutionary relationship between family of sequences. In addition, msa is often used to asess sequence conservation of amino acids in protien sequences and nucleotides in DNA sequences. It also helps to find acceptable substitution for amino acids in a protein signatures.

Main Description:

The difficulty in aligning a group of sequences varies considerably, being much greater when the degree of sequence similarity decreases. If the sequences are similar it is pretty straight forward to align them, even without using computer program. However if the sequences variation is large it is almost impossible to find an optimal alignment due to the large number of substitution, deletion and insertion. Each of these actions may create a different possible alignment. Manual multiple sequence alignment is a tedious job. Automatic msa methods had been an area of extensive research in computational biology. When one deals with automatic method - two problems need to be dealt with. The first is defining a scoring function that will assign to each alignment a score so that they will be comparable, and thus decide which one is better. The second is constructing a method that will enable to find the alignment that will retrieve the alignment with the optimal score. In the following sections we will deal with those two issues in regard to global multiple sequence alignment. In addition, we will discuss several widely used methods for local multiple sequence alignment.

Scoring Function

The msa scoring system should take into account two important features: 1) The fact that some positions are more conserved than others. 2) The fact that the sequences are not independent, and related on the phylogenetic tree. An ideal way to score msa would be to specify a complete probabilistic model of molecular sequence evolution. This model will take into account the phylogenetic tree of all the sequences and will calculate from it the probability of a certain msa. This probability will be the product of all probabilities of evolutionary events necessary to produce such an alignment through ancestral intermediate sequences times the probability of the root ancestral sequence. The probability for evolutionary event will have to depend on the evolutionary times along each branch of the tree as well as position specific structural and functional constraints imposed by natural selection. Constructing such a model requires a non-existent large amount of data. Thus relaxing notions need to be taken. Almost all alignment method assumes that the individual columns in the alignment are independent. The most wildly used method for scoring msa is the "sum of pairs" (SP). This method calculates for each column the number of matches / substitution between each pair of sequences, using either PAM or BLOSUM matrix, and applying either affine or constant penalties for gaps. The SP methodology seems pretty intuitive and reasonable; however, there exists some cases where this method is unreasonable. Examining the following example: suppose we have K sequences each having the amino acid N in column A, then BLOSUM62 score for that is 6, and the SP score will be 6 * K(K-1)/2 (number of pairs). If there is one amino acid C in this column, then the K - 1 matched N-N pairs will be replaced by K - 1 mismatched pairs. Giving a score of (6-(-3))*(K -1) = 9 * (K-1) less. Thus, the score for having one C in the column divided by no Cs in the columns is equal to 3/K. Thus, creating a counter-intuitive situation: the relative differene between the correct alignment and the incorrect alignment decreases with the increase of sequences in the alignment. The relative difference should increase with more evidence we have for conserved N. There exists other methods such as the tree-based method (Altschul 1989), trace method (Kececioglu 1993), and the information content. However, these last methods haven't been implemented in by msa programs.

Multiple Sequence Alignment Algorithms

Keeping in mind the scoring issues, we turn now to the algorithms for constructing multiple sequence alignments.

Dynamic Programming

The dynamic programming algorithm described for pair-wise sequence alignment between two protein / DNA sequences can be extended to an alignment of K sequences. However, this program is very limited because the amount of steps and memory used by the algorithm grows exponentially with the number of sequences analyzed. Thus if each sequence was n residues long the amount of time that is required is proportional to O(n^K) (K is the number of sequences). When a dynamic programming method of sequence alignment is used for pair-wise alignment, it builds a scoring matrix where each position provides the best alignment up to that point in the sequence comparison. This scoring matrix size is proportional to the length of the residues in one sequence times the number of residues in the other sequence. For a three sequence alignment, we need to fill a cube matrix, with scoring positions on three surfaces of the cube representing pair-wise alignment between each of the three sequences. A program that implements this algorithm is called MSA and since its, its memory requirements and the number of step required were substantially reduced. The program can be found at MSA.

The Progressive Method

This is probably the most used method in the field. It was first suggested by Feng & Doolittle 1987. The main idea of this method is that it does not compare all of the sequences together, but it compares only a pair of sequences at a time using the following method:

1. Calculate N(N-1)/2 distances between all pairs of sequences using the standard dynamic programming algorithm on each pair.
2. Construct a guide tree from the distance matrix using a tree clustering algorithm.
3. Starting from the first node added to the tree; align it to the second node. Then keep this alignment fixed and add another sequence / alignment. The procedure continues iteratively until all the sequences are aligned. The iterative step always chooses the next node added to the tree. A sequence is aligned to a group by aligning to the sequence in the with the highest alignment score to the new sequence. When aligning two groups of sequences, the alignment is performed using the pair (one sequence from each group) that had the highest alignment score. After an alignment between two sequences is completed, gap symbols are replaced with a neutral character X). This "gap property" allows pair-wise sequence alignment to be used to guide the alignment of sequence - alignment, and alignment - alignment. Otherwise, any given pair-wise sequence would not necessarily be consistent with the preexisting alignment. This strategy was used by several other novel algorithms. These algorithms differ by three main points:
4. The order by which the sequences are aligned.
5. Building up a single sequence alignment, or creating a group of alignments using a tree structure.
6. The use of already known alignments for alignment and scoring new sequences.

Progressive alignment is heuristic method, it does not promise optimal results. The main advantage of the algorithm is that it is a lot faster and more efficient than the generalized dynamic programming algorithm. In addition the algorithm returns easonable results. The main disadvantages of the algorithm are "the gap parameters" problem, since there is no cost for aligning a X (even with a gap), there is an encouraging effect to have gaps in subsequent pair-wise alignments. Another problem is the fact that an alignment that is created is declared "frozen". This way there is no chance to throw out any mistaken alignments. Mistakes made in an intermediate alignment may spread to further steps of the iteration.

Profile Alignment

The problem with the naïve progressive approach is that it performs pair-wise sequence alignment to determine all of the alignments. Once you already aligned sequences, it is advantageous to use position-specific information from the aligned sequences when trying to align a new sequence to the group. For example the degree of sequence conservation should be taken into account, or highly conserved position should be highly penalized for mismatches. Gap penalties in certain position should be reduced in those areas that are abundant with gaps, and increase when no gaps usually occur.

CLUSTAL

The CLUSTAL algorithm is widely used profile-based progressive alignment algorithm. The CLUSTAL method relies on the fact that the similarity between sequences is probable to evolutionary. CLUSTALW (Higgins D. et al, 1994) succeeded the CLUSTEL version by adding the 'W' standing for weighted. This represents the program's ability to provide weights to the sequences and to the program parameters. CLUSTALW works using the same progressive alignment algorithm described before, except for its carefully tuned use of the profile alignment method. The initial alignments (first step of the progressive algorithm) used to produce the guide tree may be obtained by either fast k-tuple or pattern-finding approach that is similar to FASTA, or by the slower dynamic programming method.

CLUSTALW is an ad-hoc method in regards to its alignment construction and scoring stage. In addition it has added heuristic methods of profile construction and sequence alignment that helps improve its accuracy:

1. Sequences are weighted in order to compensate for defects that may arise from use of sum-of pairs scoring function.
2. Different substitution matrices are chosen on the basis of the similarity expected of the alignment; closely related sequences and distant related matrices are aligned using different matrices (e.g. BLOSUM 80, BLOSUM 50).
3. It gives different penalties for position-specific gaps. In general, hydrophobic residues (highly conserved) give higher gap penalties from hydrophilic or flexible residues.
4. Gap-open penalties decrease if the position is spanned by consecutive stretch of five or more gaps.
5. Gaps open and gap-extend penalties are increased if there is no gaps in the alignment, but there is gaps in nearby position. This rule tries to force all the gaps to occur in the same columns.
6. In the progressive alignment phase if the alignment score is low the guide tree may be adjusted on the fly to defer the low-scoring alignment, until more profile information is accumulated.

A related algorithm that also uses the progressive method is called T-COFFEE. This method begins by making an optimal global sequence alignment of all possible sequence pairs using CLUSTAL program, and a set of possible local alignments. The best alignments found are then used by the program to devise a weighting system for creating the guide tree of the progressive method.

An example for multiple sequence alignment result by CLUSTALW

Hidden Markov Model (HMM)

The Hidden Markov Models (HMM) is a probabilistic statistical model that considers all possible combinations of matches, mismatches and gaps when generating msa. HMMs can produce a single highest-scoring output but can also generate a family of possible alignments that can then be evaluated for biological significance. Because HMMs are probabilistic, they do not produce the same solution every time they are run on the same dataset; thus they cannot be guaranteed to converge to an optimal alignment. A HMM can be visualized as a finite state machine. Finite state machines move through a series of states and produce some kind of output, either when the machine has reached a particular state or when it is moving from state to state. The HMM generates a protein sequence by emitting amino acids as it progresses through a series of states. Each state has a table of amino acid emission probabilities, and transition probabilities for moving from state to state. Transition probabilities define a distribution over the possible next states. Any sequence can be represented by a path through the model. This path follows the Markov assumption, that is, the choice of the next state is only dependent on the choice of the current state. However, the state sequence is not known; it is hidden. Finally, the alignment probability, given the model, is the product of the emission and transition probabilities along the path.

Taken from Bioinformatics Sequence and Genome Analysis, David W. Mount

A simple HMM example for proteins can be seen in the following example. Each column in the model represents the possibility of matches, deletion, or insertion. Sequences can be generated from the model starting with the BEG state. Each diamond shape represents an insertion of amino acids between two states, and each circle represents a deletion in the alignment. For example one possible way to create the N K Y L T sequence is BEG->M1->I1->M3->M4->END. For each of the possible ways to create an alignment the HMM assign a probability based on the emission and transition probabilities associated with it. Thus, the main problem in building an HMM is to define those probabilities. Usually, the "training" (finding the emission and transition probabilities) of the HMM requires 20 or more sequences. Then the trained model can be used to produce the most probable alignment of any new sequence. There are several ways to train such a model (e.g. Baum-Welch, Gibbs sampler, etc.) that will be described later in this online book. This model advantages are that it often provides an msa as good as, if not better than, other alignments. It is also well grounded in theory, and does not require and sequence ordering or insertion / deletion penalty guesses. The main disadvantage of the model is its initial requirement of 20 or more sequences for accommodating the evolutionary history of the model.

Local Multiple Sequence Alignment

In the previous sections most of the methods produced global alignment using most or all parts of the sequences. Local msa methods align only the most similar region in a sequence, ignoring dissimilarities. Local msa is useful in identifying signatures of protein families and can often be used as tools for the prediction of protein function. The results of local msa are usually displayed in three main forms:

1. Regular expression - a string that describes or matches a set of strings, according to certain syntax rules. The main database that uses such an expression is the PROSITE database. The syntax rules are as follows:

   1. A single letter representing a specific amino acid in that place.
   2. [] - Any amino acid that is found in these brackets.
   3. x - Any amino acid.
   4. () - Number of repetition of any of any expression before it.

   A simple example to illustrate it is the Interleukins 4 and 13 signature that has the following regular expression: LI-x-E-LIVM(2)-x(4,5)-LIVM-TL-x(5,7)-C-x(4)-IVA-x-DNS- LIVMA. The first amino acid of the signature can be either a L or an I, the second can be any of the 20 amino acids, followed by amino acid E, followed by 2 of either L, I, V, M or A, followed by 4 to 5 any of the available amino acids and so forth.

2. Position Specific Scoring Matrix - Matrix representing the distribution of amino acids / DNA nucleotide for each position in an alignment of protein / DNA sequences. The matrix usually contains the number / probability of a specific amino acid in each position of the motif.
3. Direct alignment - An alignment of several sequences of the particular region of interest. An example of worm globin signature:

   GAKVLAQIGVAV

   GAKVLAQIGVAV

   GKVIIDQIGKAV

   GKVIIDQIGKAV

There are three main methods for retrieving local msa:

1. Profile analysis - The same method that was described earlier, identifies a portion of the global alignment that is highly conserved and produces a scoring matrix called profile.
2. Block analysis - Searches a global msa for regions that contains substituted regions without gaps (blocks), and use those regions in the sequence alignment.
3. Pattern searching - Search of patterns (motifs) or use of statistical methods to

find a localized region of similarity in a set of sequences. These methods don't necessarily depend on the production of a global msa.

Profile Analysis

A profile is a matrix that includes scores in a local region for matches, mismatches and gap penalties. Once this profile is produced it is used to search target sequences for possible match of the profile within the sequence, this is done using the score saved in the profile's matrix. The main disadvantages of this method are that profile is extracted from globally created msa. If an msa is only a representative of the variation in the family of sequences, the profile will be based on this variation. If the sequences in the msa are similar, the profile will be biased in favor of those sequences. Another problem with profiles is that some amino acids may not be represented in a particular column, due to small number of sequences that are included in the msa. Although some amino acids could not occur in a certain position it is usually useful to have several counts to such positions to increase the usefulness of the profile.

Block Analysis

Like profiles, blocks represent conserved regions in an msa. Blocks are different than profiles in that they only include matches and mismatches for each position in conserved region. Block can be created like profiles, searching regions that are highly conserved. However, alignable regions may also be found by searching each sequence in turn for similar patterns of the same length. These regions may include one or more similar characters followed by gaps and than continues with more similar characters, until the sequences are very distant. Usually block searching software do not use the PAM or BLOSUM substitution matrix, and they rely on finding exact matches that have the same spacing in the input sequence and that may be repeated in the given sequence.

Pattern Searching

This line of algorithms does not rely on global alignment of sequences. Instead, it searches for localized region with high similarity to predefined patterns. The pre-defined patterns can either be expressed explicitly or through a statistical model. This type of analysis was first performed in by Smith et al. (1990), using a program called MOTIF. The MOTIF program performed well in helping construct the BLOCKS database, but it suffered from few major deficiencies. The program is limited in the pattern sizes that can be found. In addition, the program distinguishes true motifs from random background patterns by requiring that motifs occur in a number of the input sequences and tend not to be internally repeated in any one sequence, which can be a problem when the length of the motif increases, and there are many possible combinations of patterns of a given length where only a few characters match. Lastly, the MOTIF program always provides a motif, even for random sequences, thus making it difficult to decide how significant the found motif really is. In order to compensate for this type of crude methods statistical based methods were incorporated to enhance the performance of such methods. An example of such algorithm are the EM algorithm(used by MEME server), Gibbs sampling (used by Phylogibbs server) and HMM and EM algorithm (Meta-MEME software). Results from using HMM models in global and local msa are stored at the Pfam database.

Conclusions:

Multiple sequence alignment has been proven to be a powerful tool for many fields of studies such as phylogenetic reconstruction, illumination of functionally important regions, and prediction of higher order structures of proteins and RNAs. Today's progress in multiple sequence alignment methods has been made by the multidisciplinary endeavors of mathematicians, computer scientists, statisticians, and biologists. This piece gave a brief review of most widley used global and local multiple sequence alignments. As the amount of data increase better models could be obtained by both upgrading the scoring formula, and creation of statistical models that will incorporate evolutionary dependencies as well position specific stractural and functional constraints. In the near future, generalized multiple alignment may play a central role in studies of protien / DNA sequence, structure, function, and evolution relationships.

References:

• Lecture Notes and Textbooks
  • Bioinformatics Sequence and Genome Analysis, Second Edition
  • Biological Sequence Analysis: Probabilistic Models of Proteins and Nucleic Acids
  • List of Multiple Sequence Alignment Tools
  • Lecture on Multiple Sequence Alignment
• Databases, Webservers and applications
  • List of Multiple Sequence Alignment Tools
  • T-COFFEE Multiple Sequence alignment tool
  • Multiple Em for Motif Elicitaion
  • PROSITE Database
  • Pfam Database
• Class video
• References
  • Altschul S.F. 1989. Gap costs for multiple sequence alignment. J. Theor. Biol. 138: 297–309. (Link)
  • Feng D.F. and Doolittle R.F. 1987. Progressive sequence alignment as a prerequisite to correct phylogenetic trees. J. Mol. Evol. 25: 351-360. (Link)
  • Grundy W.N., Bailey T.L., Elkan C.P., and Baker M.E. 1997. Meta-MEME: Motif-based hidden Markov models of protein families. Comput. Appl. Biosci. 13: 397-406.(Link)
  • Higgins DG, Thompson JD, Gibson TJ. 1996. Using CLUSTAL for multiple sequence alignments Methods Enzymol. 266: 383-402.(Link)
  • Kececioglu J. 1993. The maximum weight trace problem in multiple sequence alignment. In Proceedings of the 4th Symposium on Combinatorial Pattern Matching: Lecture notes in computer science, no. 684, pp. 106–119. Springer Verlag, New York.(Link)
  • Krogh A., Brown M., Mian I.S., Sjölander K., and Haussler D. 1994. Hidden Markov models in computational biology. Applications to protein modeling. J. Mol. Biol. 235: 1501-1531.(Link)
  • Lawrence C.E. and Reilly A.A. 1990. An expectation maximization (EM) algorithm for the identification and characterization of common sites in unaligned biopolymer sequences. Proteins Struct. Funct. Genet. 7: 41-51.(Link)
  • Lawrence C.E., Altschul S.F., Boguski M.S., Liu J.S., Neuwald A.F., and Wootton J.C. 1993. Detecting subtle sequence signals: A Gibbs sampling strategy for multiple alignment. Science 262: 208-214.(Link)
  • Lipman D.J., Altschul S.F., and Kececioglu J.D. 1989. A tool for multiple sequence alignment. Proc. Natl. Acad. Sci. 86: 4412-4415.(Link)
  • Notredame C. and Higgins D.G. 1996. SAGA: Sequence alignment by genetic algorithm. Nucleic Acids Res. 24: 1515-1524.(Link)
  • Smith H.O., Annau T.M., and Chandrasegaran S. 1990. Finding sequence motifs in groups of functionally related proteins. Proc. Natl. Acad. Sci. 87: 826-830.(Link)
  • Thompson JD, Plewniak F, Poch O. (1999). A comprehensive comparison of multiple sequence alignment programs. Nucleic Acids Res 7(13):2682-90.(Link)
  • Wang L, Jiang T. (1994) On the complexity of multiple sequence alignment. J Comput Biol 1:337-348. (Link)