The Performance of Several Multiple-Sequence Alignment Programs in Relation to Secondary-Structure Features for an rRNA Sequence

Robert E. Hickson1,4,*, Chris Simon{dagger} and Soren W. Perrey2,{ddagger}

*Department of Genetics and Molecular Biology, University of Hawaii at Manoa;
{dagger}Department of Ecology and Evolutionary Biology, University of Connecticut; and
{ddagger}Department of Mathematics, Massey University, Palmerston North, New Zealand

Abstract

The performances of five global multiple-sequence alignment programs (CLUSTAL W, Divide and Conquer, Malign, PileUp, and TreeAlign) were evaluated using part of the animal mitochondrial small subunit (12S) rRNA molecule. Conserved sequence motifs derived from an alignment based on secondary structural information were used to score how well each program aligned a data set of five vertebrate and five invertebrate taxa over a range of parameter values. All of the programs could align the motifs with reasonable accuracy for at least one set of parameter conditions, although if the whole sequence was considered, similarity to the structural alignment was only 25%–34%. Use of small gap costs generally gave more accurate results, although Malign and TreeAlign generated longer alignments when gap costs were low. The programs differed in the consistency of the alignments when gap cost was varied; CLUSTAL W, Divide and Conquer, and TreeAlign were the most accurate and robust, while PileUp performed poorly as gap cost values increased, and the accuracy of Malign fluctuated. Default settings for the programs did not give the best results, and attempting to select similar parameter values in different programs did not always result in more similar alignments. Poor alignment of even well-conserved motifs can occur if these are near sites with insertions or deletions. Since there is no a priori way to determine gap costs and because such costs can vary over the gene, alignment of rRNA sequences, particularly the less well conserved regions, should be treated carefully and aided by secondary structure and conserved motifs. Some motifs are single bases and so are often invisible to alignment programs. Our tests involved the most conserved regions of the 12S rRNA gene, and alignment of less well conserved regions will be more problematical. None of the alignments we examined produced a fully resolved phylogeny for the data set, indicating that this portion of 12S rRNA is insufficient for resolution of distant evolutionary relationships.

Introduction

Alignment of nucleic acid or protein sequences is a critical early step in comparative molecular and evolutionary studies. With increasing sequence divergence, it becomes harder to construct alignments that reliably reflect historical relationships. In these situations, multiple-sequence alignment programs are often favored because of their more objective (mathematical) criteria. However, different alignments can be generated by different programs and by different parameter values (e.g., see Tyson 1992Citation ; McClure, Vasi, and Fitch 1994Citation ; Wheeler 1995Citation ; Morrison and Ellis 1997Citation ), making the critical issue the quality of the alignment (Lipman, Altschul, and Kececioglu 1989Citation )—how is one alignment preferred over another?

Knowledge of sequence secondary structure is a valuable aid in this regard, and the accuracy or biological reliability of alignment programs is commonly evaluated by observing how well a family of protein sequences is aligned by a given program when compared with known conserved features (e.g., Barton and Sternberg 1987Citation ; Lipman, Altschul, and Kececioglu 1989Citation ; McClure, Vasi, and Fitch 1994Citation ; Thompson, Higgins, and Gibson 1994Citation ). Use of such information also permits the objective comparison of different alignment programs. McClure, Vasi, and Fitch (1994)Citation evaluated a large number of alignment programs using sets of different protein sequences. They demonstrated that alignment programs differ in their accuracy and that different parameter conditions can produce quite different alignments. Determining the biological validity of alignments for non–protein-encoding sequences, such as ribosomal RNA (rRNA), has received less attention (but see Titus and Frost 1996Citation ; Morrison and Ellis 1997Citation ). This, in part, is because of the greater difficulty in aligning (four-character-state) DNA sequences compared with amino acid sequences (20 character states). However, secondary-structure models for rRNA sequences help identify well-conserved sequence motifs which greatly assist alignment (Gutell 1994Citation ; Hickson et al. 1996Citation ; de Rijk et al. 1999Citation ; van de Peer et al. 1999Citation ).

Morrison and Ellis (1997)Citation found that for an 18S rRNA data set, the choice of alignment parameters can have a greater effect on the alignment than the choice of an alignment program. These results highlight the problem of determining what parameter values are the most appropriate to use for a data set. In this paper, we investigate this issue further by comparing several multiple-sequence alignment programs over a range of parameter values using a mitochondrial small subunit (12S) rRNA data set for which the secondary structure is well characterized. In particular, we examine to what extent the programs are sensitive to small changes in parameter values, since without prior knowledge of what these should be, it is preferable to select a program which is comparatively robust to changes in the settings.

Materials and Methods

The Model and Motifs
Based on comparison of over 180 sequences from six vertebrate and nine invertebrate classes (five phyla overall), we refined the general small-subunit secondary-structure model of Gutell (1994)Citation for domain III of animal 12S rRNA (Hickson et al. 1996Citation ), a widely and frequently sequenced region of animal rRNA. This model identified well-conserved sequence motifs around which the alignment can be based (Hickson et al. 1996Citation ).

For the present evaluation, we selected 10 domain III sequences (5 vertebrate and 5 invertebrate) and 10 motifs (fig. 1 and table 1 ) to compare a range of alignment programs and alignment parameters. The sequences were those of cow (Bos taurus; GenBank accession number J01394), bird (Megalapteryx didinus; X67634), lizard (Oligosoma nigriplantare polychroma), Xenopus laevis (M10217), fish (Crossostoma lacustre; M91245), sea urchin (Paracentrotus lividus; J04815), Drosophila yakuba (X03240), honeybee (Apis mellifera; L06178), cicada (Magicicada tredecim), and earthworm (Aporrectodea rosea; L02392). Secondary-structure diagrams for some of these taxa can be found in Hickson et al. (1996)Citation .



View larger version (49K):
[in this window]
[in a new window]
 
Fig. 1.—Alignment of domain III of animal 12S rRNA derived from the alignment of Hickson et al. (1996). Note, however, that the alignment of the honeybee sequence from stem 47 onward has been altered from the original figure in Hickson et al. (1996). The approximate locations of stems are indicated above the alignment, and sets of motifs are indicated below (r = a purine; y = a pyrimidine; uppercase characters indicate highly conserved nucleotides). The 10 motifs used for evaluating the alignment programs are numbered and highlighted in black

 

View this table:
[in this window]
[in a new window]
 
Table 1 Motifs in Domain III of Animal 12S rRNA Used for the Evaluation of Alignment Programs

 
The average degree of similarity between the vertebrate sequences (based on the secondary-structure alignment) was 75%, while the invertebrate sequences had 64% similarity with each other and, on average, 62% nucleotide identity with the vertebrate sequences. The motifs used for scoring were distributed throughout domain III and vary in size (from 3 to 14 bases) and degree of conservation (fig. 1 ). Seventy-one percent of the motif nucleotides were in stem (primarily base-paired) regions; the remainder were in loops (single-stranded regions) of the rRNA molecule. We devised a simple alignment scoring scheme in the manner of McClure, Vasi, and Fitch (1994)Citation to evaluate the alignment programs and parameters. This simply counts for each motif how many of the sequences are aligned as in the secondary-structure model, with the overall score expressed as proportion correct relative to the model.

Alignment Programs
We did not compare all available multiple-sequence alignment programs but selected five which are widely used or incorporate different criteria for aligning. The five programs, which all attempt to identify a global optimal alignment, were:

  1. PileUp, a progressive pairwise alignment program using a simplified version of the Feng and Doolittle (1987)Citation algorithm. Sequence similarity is used to cluster sequences, and a dendrogram is constructed using the unweighted pair grouping method with arithmetic means (UPGMA). This dendrogram directs the ordering of subsequent pairwise alignments based on the Needleman and Wunsch (1970) method. PileUp is part of the Wisconsin Sequence Analysis package (version 8.1; Genetics Computer Group 1995).
  2. CLUSTAL W (Thompson, Higgins, and Gibson 1994Citation ), which is also based on the Feng and Doolittle (1987)Citation algorithm, but with improved choice of alignment parameters using dynamic assignment of penalties. CLUSTAL W varies gap penalties in relation to sites along the sequence and the relative evolutionary distance between sequences (and subsets of sequences), so the specified penalties are used as initial values only. As in PileUp, a "guide tree" for the order of alignment is initially constructed based on pairwise differences. However, CLUSTAL W uses neighbor joining rather than UPGMA to build the dendrogram. In CLUSTAL W, the gap and gap-length costs can be specified separately for the pairwise and multiple alignments. Transitions can either be weighted or unweighted (treated as mismatches). The "delay divergent sequences" option (see table 2 ) delays the inclusion of sequences with less than the specified sequence identity until the more similar sequences are aligned. For our analysis, we used the "slow/accurate" option in CLUSTAL W to calculate the distances between sequences, although the "fast approximate" method gave similar scores for some of the parameter values.
  3. Divide and Conquer (D&C; Tönges et al. 1996Citation ), which makes use of the MSA algorithm (Carrillo and Lipman 1988Citation ; Lipman, Altschul, and Kececioglu 1989Citation ) as a subroutine. MSA uses the (weighted) sum-of-pair score with a bounded score criterion that significantly reduces the search space for the optimal alignment. D&C recursively cuts the sequences under consideration into sets of smaller subsequences. If these subsequences are sufficiently small (specified by "stopsize" in table 2 ), they are aligned by the MSA algorithm. Eventually, the resulting subalignments are concatenated to obtain an alignment of the full-length sequence. The only heuristic step in the D&C procedure is the calculation of the cutting points, which is done according to a measure closely related to the sum-of-pair score (Tönges et al. 1996Citation ). Dividing the sequences up into shorter lengths for alignment makes it possible to align many more sequences than would be possible using MSA alone. While MSA always calculates an optimal alignment, D&C approximates such an optimum for much larger sets (up to 20) of sequences. MSA and D&C are the only algorithms which essentially align the sequences simultaneously.
  4. Malign (version 2.3; Wheeler and Gladstein 1994Citation ), which uses parsimony as an optimality criterion for aligning sequences. It examines the topologies of many alignment orders and chooses the order giving the shortest cladogram. A "minimum" cladogram can be found in a variety of ways, specified by the "score" option (0 = nontopological; 1 = branch and bound; 2 = random cladograms; 3 = rapid heuristic; 4 = exhaustive search). We used the branch-and- bound option because exhaustive searches took inordinate amounts of time to run. A large range of other options can also be specified in the Malign parameter file (see table 2 ).
  5. TreeAlign (Hein 1990, 1994Citation ), a fast approximation algorithm which infers an alignment and the phylogeny. The objective of TreeAlign is to minimize the weighted sum of insertion/deletions (indels) and substitutions. It infers ancestral sequences from the clustering of pairwise comparisons by use of string and graph comparisons (Hein 1990Citation ). These ancestral sequences are aligned, and the minimum path is selected. TreeAlign defines distance functions for transitions (a score of 2) and transversions (a score of 5).


View this table:
[in this window]
[in a new window]
 
Table 2 Parameters Manipulated in the Alignment Programs

 
The major parameters which are modifiable in each of these programs are listed in table 2 . The two key parameters for all of the programs are the costs assigned to the opening of a gap (gap cost) and the cost for extending a gap (gap-length cost) (total gap cost = gap cost + [gap-length cost x gap length]). These parameters were varied to evaluate the performance of the programs and find the conditions which gave alignments most similar to the structurally based alignment. We did not undertake a comprehensive evaluation of all combinations of gap and gap-length costs (see Wheeler 1995Citation ; Morrison and Ellis 1997Citation ), but evaluated a range of parameter values in a linear fashion (see table 2 ) to provide an indication of the programs' sensitivity to changes. In this paper, we use the term "accuracy" to refer to how similar automatic alignments are to the alignment based on secondary-structure information.

Results

All of the programs aligned the 10 motifs in a manner similar to the structural alignment for at least one set of parameter conditions (table 1 ), with small gap and gap-length costs generally producing more accurate motif alignments. However, the programs differed in their robustness to changes in parameters, particularly with respect to gap cost (fig. 2a ). Over the range of gap costs tested, CLUSTAL W, D&C, and TreeAlign showed only minor changes in alignment scores, whereas the accuracy of alignments from PileUp declined as gap cost values increased, and those for Malign fluctuated (fig. 2a ). Larger gap costs did not produce improved alignment scores. Gap-length cost had a smaller influence on motif alignment scores than did gap cost (data not shown). If the whole alignment from each program was considered, then overall similarity to the structural alignment (fig. 1 ) was low, ranging between 25% and 34% among the different programs for the alignment giving the highest motif scores.



View larger version (18K):
[in this window]
[in a new window]
 
Fig. 2.—The effect of changing gap cost on the quality of alignments. a, Influence on accuracy of motif alignments. Gap-length costs for each program were 0.1 for PileUp, 1 for Malign, and 2 for the rest of the programs. b, Influence on alignment length. Gap-length costs were the same as for a. Note that the length of the alignment based on secondary structure (fig. 1 ) was 317 nt

 
The length of the alignment for PileUp, CLUSTAL W, or D&C was less affected by changes in gap cost, but many more gaps were introduced by Malign and TreeAlign when gap cost was <5 (fig. 2b ). For Malign, the parameter "extragaps" must be explicitly set as 1; otherwise, many gaps of length 1 can be incorporated into the alignment.

With D&C, we found that accuracy improved (alignment scores increased up to 8%) if we included a transition/transversion matrix with values set to 5/2 (as in TreeAlign). However, if the same transition/transversion matrix was used in Malign, the alignment became substantially worse—the best score declined from 0.89 to 0.76, and the alignment length increased substantially to 673. Different relative costs for transitions and transversions could not be assigned in all the programs, so we did not evaluate the effect of a large range of different transition and transversion costs. However, it is clear that attempting to select similar parameter values in different programs does not necessarily result in more similar alignments.

Some of the motifs were well aligned by all of the programs over the range of gap costs tested, but this was not necessarily a simple reflection of their length or degree of conservation. For example, the long and well-conserved motif 6 and the short motif 4 (fig. 1 and table 1 ) were generally aligned quite accurately, while another long motif (number 5) was sometimes misaligned by PileUp, Malign, and TreeAlign (fig. 3a ). This variability seems to be the result of having to include a large alignment gap before it (see fig. 1 ), but there is not a simple relationship between gap cost value and alignment accuracy (fig. 3a ).



View larger version (16K):
[in this window]
[in a new window]
 
Fig. 3.—Quality of alignment for two different motifs situated near large indels, derived from the different programs for a range of gap cost values. The "score" of the Y axis refers to how many of the 10 taxa had the motif aligned as it was in the secondary-structure model. a, The well-conserved motif 5. b, The less strongly conserved motif 8. Gap-length cost values are as in figure 2

 
Motif 8, being both in a region where several gaps are required and not absolutely conserved (see fig. 1 ), presented the greatest challenge for the programs, yielding a range of differing alignments both within and between programs (fig. 3b ). For gap costs of >6, a single large gap tended to be introduced by the programs. Motif 8 was more accurately aligned when gap cost values were moderate (about six—see fig. 3b ), but another similarly conserved motif nearby ("rYgrr" adjacent to stem 47; see fig. 1 ) was then poorly aligned by many of the programs (data not shown). As illustrated in figure 3 for Malign, changing the gap cost can improve the alignment of one motif while degrading the alignment of another.

The number of sequences and degree of similarity between them can also affect alignment accuracy, with the degree of similarity having a greater influence. This is illustrated in figure 4 , where alignments of three sets of 10 taxa of increasing dissimilarity (10 vertebrates, 5 vertebrates and 5 invertebrates, and 10 invertebrates). For simplicity, we only show results for the CLUSTAL W and Malign programs to investigate whether lack of alignment robustness in the latter program can be a consequence of using sequences that are too dissimilar. If just vertebrate sequences were aligned, then only motifs 8 and 10 presented difficulties, while as more invertebrate sequences were included, other motifs became misaligned. Note that even with just vertebrate sequences, Malign still shows more variability in alignment accuracy than does CLUSTAL W (fig. 4 ). Increasing the number of sequences did not necessarily improve alignments, especially if more dissimilar sequences were added. However, D&C, which simultaneously aligns, performed better than the other programs as more sequences were added (data not shown).



View larger version (16K):
[in this window]
[in a new window]
 
Fig. 4.—The effect of including more divergent sequences in an alignment using the programs CLUSTAL W (Cw) and Malign (Ma). The 10-vertebrate data set (10 Verts) included the five sequences used in the previous analyses, as well as the fin whale (Balaenoptera physalus; GenBank accession number X61145), the opossum (Didelphis virginiana; Z29573), the kiwi (Apteryx australis; X67626), the tuatara (Sphenodon punctatus; L28076), and the alligator (Alligator mississippiensis; Y13113). The data set including five vertebrate and five invertebrate (5+5) sequences was the same one used in the previous analyses. The 10- inveterbrate (10 Inverts) data set included the three insects and the earthworm from the original data set, as well as the damselfly (Ischnura cervula), silverfish (Ctenolepisma longicaudata; L02381), spider (Tetragnatha mandibulata; U00118), brineshrimp (Artemia franciscana; X69067), onychophoran (Euperipatoides leuckartii; L02380), and chiton (Ischnochiton australis; L02388) sequences

 
All of the programs evaluated had trouble aligning the honeybee sequence, which deviates from some of the motifs (notably, numbers 3 and 7; fig. 1 ). Relative scores improved if this sequence was excluded, but differences in accuracy and consistency between the programs remained (data not shown).

The degrees of treelike structure in the alignments are shown in figure 5 using split decomposition, as implemented in SplitsTree (Huson 1998Citation ). None of the alignments, including the one based on secondary structure, generated fully resolved phylogenies (fig. 5 ), irrespective of whether just the sites with gaps were omitted or more extensive removal of alignment-ambiguous sites occurred prior to phylogenetic reconstruction. We attribute this poor phylogenetic performance to the short length of the sequence examined and to site saturation becoming a problem when domain III of the 12S rRNA is used for comparing different animal orders and phyla (unpublished data).



View larger version (12K):
[in this window]
[in a new window]
 
Fig. 5.—Phylogenies derived from the alignments indicated in table 1 . Split decomposition, performed using the program SplitsTree2.4 (Huson 1997Citation ), was used to generate the trees using the LogDet correction (Lockhart et al. 1994Citation ) to compensate for differences in base composition. Areas of uncertain alignment were excluded. Conflicting phylogenetic signals are indicated by parallelograms between branches. Results from parsimony and neighbor-joining reconstruction methods were similar. a, Alignment based on secondary structure (228 characters used). b, Alignment derived from PileUp (200 characters). c, Alignment derived from CLUSTAL W (228 characters). d, Alignment derived from Divide and Conquer (207 characters). e, Alignment derived from Malign (225 characters). f, Alignment derived from TreeAlign (211 characters)

 
Although the phylogenies are incompletely resolved, some of the alignments do produce alternative clusterings of taxa. For example, the clustering of the three insect sequences is less well supported using the PileUp alignment (fig. 5b ), and the Malign alignment places the earthworm within the insects (fig. 5e ). Alignments based on secondary structure (fig. 5a ) and the CLUSTAL W output (fig. 5c ) show more ambiguity in the relationships of the five vertebrates relative to the urchin sequence.

Discussion

Automated sequence alignment, as with other steps in phylogenetic reconstruction, requires selection of parameter values whose most appropriate settings can be difficult to determine in advance. Use of secondary-structure information is valuable for helping determine what values generate realistic alignments, and it also provides a way of objectively evaluating the performances of different programs (McClure, Vasi, and Fitch 1994Citation ).

All of the programs we evaluated produced similar alignment scores for at least one set of parameter values (table 1 ). However, the programs differed in their robustness to changes in the cost assigned to gaps (fig. 2 ) and in how well they aligned specific motifs (fig. 3 ), as well as in phylogenetic resolution (fig. 5 ). CLUSTAL W, D&C, and TreeAlign produced more consistent alignments than Malign or PileUp over the range of gap costs tested. This is an important observation, since without prior knowledge of what the parameter values should be, programs that are less sensitive to relatively small changes in settings are preferable. Default settings of the programs did not give the best results. More structurally accurate alignments were usually obtained when lower gap (fig. 2a ) and gap-length costs were used, although Malign and TreeAlign generated longer alignments when gap costs were low (fig. 2b ). We concentrated on the alignment of a few conserved motifs, largely ignoring the rest of the sequence alignment, but the alignment of these latter regions showed a similar pattern to that for the motifs, with outputs from programs such as CLUSTAL W showing little difference as gap costs changed, while alignments from Malign exhibited variation (data not shown).

Morrison and Ellis (1997)Citation found that the choice of parameter values in CLUSTAL W had a larger influence on alignment than the choice of an alignment program (which also included PileUp, Malign, and TreeAlign in their study) for a data set of 43 parasitic protozoan 18S rRNA sequences. This differs from our results and appears to be due to the different data sets used. Morrison and Ellis (1997)Citation examined a longer sequence length and a larger number of sequences than we did, and, although more closely related than the taxa in our data set, their 18S sequences show a bigger difference in the number and size of indels, which may make alignments more sensitive to parameter values. Morrison and Ellis (1997)Citation tested a larger range of gap and gap-length costs for CLUSTAL W than we did here, but we did not find that increasing the gap costs improved alignments. The robustness of programs other than CLUSTAL W to changes in parameter values were not evaluated by Morrison and Ellis (1997)Citation , but considering our findings, there may be a larger effect of program type on alignment than they concluded.

Both our results and those of Morrison and Ellis (1997)Citation demonstrate that for rRNA data, as with alignment of protein-encoding sequences, alignment accuracy can be influenced by the choice of alignment program and the parameter values. As can be expected, differences between alignments increase as more insertion/deletion events have to be incorporated (see fig. 4 ). Use of more robust programs such as CLUSTAL W and D&C would appear to be preferable in these situations, but it is desirable to compare alignments from any program with information about secondary structure to ensure that the alignments are biologically reasonable.

Determining why some programs are more robust than others is difficult. Differences in alignment scores between some of the programs were sometimes due to misalignment of well-conserved motifs, and the major factor influencing this is the presence of large indels near them (see fig 3b ). PileUp uses the same basic algorithm as CLUSTAL W (Feng and Doolittle 1987Citation ), but the latter has more sophisticated features, such as neighbor joining to guide alignments and dynamic penalty assignments, which will contribute to CLUSTAL W's greater precision and more consistent performance. Nonetheless, PileUp still aligned most of the motifs accurately when the gap cost was <=5, and it provides a good starting point from which to improve the alignment. The variability in alignments produced by Malign are likely to be due in part to its use of a parsimonious phylogenetic tree in the alignment process. With our data set, phylogenetic resolution is poor (see fig. 5 ), which will, therefore, influence the robustness of alignments. Other factors must also be involved, however, since TreeAlign, which also employs a tree-based alignment procedure, was not as variable under the conditions we tried (see fig. 2 ).

McClure, Vasi, and Fitch (1994)Citation found that global alignment programs performed better than those employing local optimization. We evaluated only global alignment programs, but greater accuracy for some of the individual motifs could have been obtained using local alignment programs. D&C aligns small fragments of the sequence data, so it may be able to align the more length-variable regions with greater precision than the other programs tested. However, this did not appear to be the case. For gap costs of <=5, D&C was more accurate than the other programs in aligning the "rYgrr" motif near stem 47, but at the cost of misaligning motif 8.

Accuracy of alignment can also be influenced by the number of sequences used and their degree of similarity (see Gatesy, DeSalle, and Wheeler 1993Citation ; McClure, Vasi, and Fitch 1994Citation ). However, as illustrated in figure 4 , the latter is the more important factor, particularly with respect to indels. Even well-conserved motifs (such as motif 5 in our data set) may become misaligned if they are near regions which require insertions or deletions of several nucleotides in some sequences. We deliberately tested a data set which included some quite divergent sequences (such as the honeybee) to reflect situations in which alignment programs are more frequently relied on. PileUp, CLUSTAL W, and Malign were more sensitive to input order than were D&C and TreeAlign, so it is also advisable to try different orders of sequence input when aligning divergent sequences.

In our comparisons, the alignment of motifs between the vertebrate sequences was generally accurate, reflecting their high degree of similarity (75% over the whole alignment, 94% between motifs). The three insect sequences (77% similarity between their motifs) often had some of their motifs aligned correctly to each other but not to the other sequences. McClure, Vasi, and Fitch (1994)Citation drew attention to this problem, which can often be corrected by visual inspection.

In any phylogenetic study, alignment of the less well-conserved and the variable-length regions is often problematical, even when sequence divergence is moderate (e.g., between vertebrate orders). For regions of ambiguous alignment, a critical consideration is how these sites should be treated in phylogenetic analyses. Inclusion of ambiguously aligned regions can support erroneous patterns of branching, while removal of such regions can reduce resolution (Gatesy, DeSalle, and Wheeler 1993Citation ; Wheeler, Gatesy, and DeSalle 1995Citation ). Wheeler, Gatesy, and DeSalle (1995)Citation advocated a strategy called "elision." Elision involves aligning sequences using a range of parameter values to identify alignment-ambiguous sites and combining these multiple multiple alignments into a single alignment with alignment-ambiguous sites downweighted. Such a technique has its advantages but, as our analyses indicate, some programs (e.g., PileUp and Malign) can misalign well-conserved motifs when even small changes are made to parameter values. Knowledge of secondary structure is helpful to delimit the range of parameter values which preserve the alignment of well-conserved motifs. Alignment ambiguity can be reduced when secondary-structure features (such as motifs and paired regions) are considered, although this may still not help fully resolve a phylogeny (fig. 5 ).

An added difficulty for the determination of the most suitable alignment parameter values is that different costs may be appropriate for different regions, especially where large insertions or deletions occur. In this regard, CLUSTAL W, with its dynamic penalty assignments, offers a potential advantage over the other programs, although it still results in inaccuracies (e.g., for motif 8; fig. 3b ), and the performances of D&C and TreeAlign (which require greater computational power) in aligning motifs were as good or better (fig. 2 ).

An important area for improving alignment programs is the ability to constrain alignment of specific regions, as is possible for MSA (Lipman, Altschul, and Kececioglu 1989Citation ). Bell, Coggins, and Milner-White (1993)Citation used a "mix-and- match" strategy to first identify strongly conserved regions and then use these as anchors to align adjacent less well conserved regions. We attempted this in our data set by aligning only the region between motifs 7 and 8, but we did not obtain alignments that more closely matched that based on secondary structure. Where ambiguity remains, it is preferable to adopt a more conservative analysis and exclude these regions, although the presence or size of indels may be able to provide additional phylogenetic information.

van de Peer et al. (1993)Citation used secondary-structure information and a matrix method for calculating the variability of each site in a sequence as a framework to produce more meaningful alignments of rRNA genes. They developed an alignment program (called DCSE) which uses secondary-structure information (de Rijk and de Wachter 1994). We did not evaluate DCSE in the present study, but because it requires prior knowledge of the secondary structure, its main advantage lies in refining rather than initiating an alignment.

Alignment of rRNA sequences is often more problematical than that of protein-encoding sequences, because for rRNA molecules some secondary-structure elements cannot be as simply predicted from the sequence data. While some nucleotides in rRNA sequences are highly conserved (see Gutell 1994;Citation de Rijk et al. 1999Citation ; van de Peer et al. 1999Citation ), other parts of the molecule can show little sequence similarity even when secondary-structure features are preserved. Base-paired regions must be identified by direct examination of the data, and even then the nucleotides involved may not be homologous, since paired and unpaired regions can increase or decrease in length and shift position. Utilization of secondary-structure information is, however, vital to help identify which regions are homologous, to provide a more reliable basis for aligning and comparing rRNA sequences, and to give the researcher a better understanding of the sequences being analyzed.

Use of alignment programs has been advocated because of their objectivity and reproducibility (Thorne and Kishino 1992Citation ). Reproducible alignments, though, may still not be accurate, and alignments should always be compared with secondary-structure information if it is available. Sequence alignment should not be a one-step process and, as with other steps in phylogenetic reconstruction, the assumptions of the alignment program must be considered and the biological accuracy evaluated. The identification of conserved motifs in protein and rRNA molecules offers an easy way to assess how an alignment program and the parameter values reflect biological accuracy, and visual inspection of the data can also provide the investigator with added insight into how taxa vary. Secondary-structure models are available for both large- and small-subunit rRNAs (Gutell 1994Citation ; de Rijk et al. 1999Citation ; van de Peer et al. 1999Citation ), so outputs generated by multiple-alignment programs can, and should, be compared with these prior to further phylogenetic analyses of the data.

A related but separate issue concerning alignment accuracy is that of how well a specific alignment reflects a known phylogeny. Alignments from the different programs showed some differences in support or conflict for some clusterings (fig. 5 ), but because of site saturation, we could not address the issue of alignment and phylogeny in detail with our data set. However, others have shown that even small differences in alignment can result in quite different phylogenies (Kjer 1995;Citation Morrison and Ellis 1997Citation ), so programs that produce alignments which are relatively insensitive to changes in parameter values offer advantages as a starting point for establishing an alignment.

Acknowledgements

We thank Peter Gogarten for running our data sets through TreeAlign. Jotun Hein and Des Higgins provided answers to queries about their programs. Comments by one of the anonymous reviewers helped improve the manuscript. During this study, R.E.H. was supported initially by an Alfred P. Sloan Postdoctoral Fellowship and later by a New Zealand Science and Technology Fellowship. This work was partially supported by grants to C.S. from the U.S. National Science Foundation and the University of Connecticut.

Footnotes

Mike Hendy, Reviewing Editor

1 Present address: ERMA NZ, Wellington, New Zealand. Back

2 Present address: Münster, Germany. Back

3 Keywords: rRNA sequence alignment phylogenetics secondary structure Back

4 Address for correspondence and reprints: Robert Hickson, ERMA NZ, P.O. Box 131, Wellington, New Zealand. E-mail: robert.hickson{at}ermanz.govt.nz Back

literature cited

    Barton, G. J., and M. J. E. Sternberg. 1987. A strategy for the rapid multiple alignment of protein sequences—confidence levels from tertiary structure comparisons. J. Mol. Biol. 198:327–337.[ISI][Medline]

    Bell, L. H., J. R. Coggins, and E. J. Milner-White. 1993. Mix ‘n’ Match: an improved multiple sequence alignment procedure for distantly related proteins using secondary structure predictions, designed to be independent of the choice of gap penalty and scoring matrix. Protein Eng. 6:683–690.[Abstract]

    Carrillo, H., and D. Lipman. 1988. The multiple sequence alignment problem in biology. SIAM J. Appl. Math. 48:1073–1082.[ISI]

    de Rijk, P., and R. de Wachter. 1993. DCSE v2.54, an interactive tool for sequence alignment and secondary structure research. Comput. Appl. Biosci. 9:735–740.[Abstract]

    de Rijk, P., E. Robbrecht, S. de Hoog, A. Caers, Y. van de Peer, and R. de Wachter. 1999. Database on the structure of large subunit ribosomal RNA. Nucleic Acids Res. 27:174–178.[Abstract/Free Full Text]

    Feng, D. F., and R. F. Doolittle. 1987. Progressive sequence alignment as a prerequisite for correct phylogenetic trees. J. Mol. Evol. 25:351–360.[ISI][Medline]

    Gatesy, J., R. DeSalle, and W. C. Wheeler. 1993. Alignment-ambiguous nucleotide sites and the exclusion of systematic data. Mol. Phylogenet. Evol. 2:152–157.[Medline]

    Genetics Computer Group. 1995. Wisconsin sequence analysis package. Version 8.1. GCG, Madison, Wis.

    Gutell, R. 1994. Collection of small subunit (16S- and 16S-like) ribosomal RNA structures: 1994. Nucleic Acids Res. 22:3502–3507.[Abstract]

    Hein, J. 1990. Unified approach to alignment and phylogenies. Methods Enzymol. 183:626–645.[ISI][Medline]

    ———. 1994. TreeAlign. Pp. 349–364 in A. M. Griffin and H. G. Griffin, eds. Methods in molecular biology, Vol. 25. Computer analysis of sequence data, part II. Humana Press, Totowa, NJ.

    Hickson, R. E., C. Simon, A. C. Cooper, G. S. Spicer, J. Sullivan, and D. Penny. 1996. Conserved sequence motifs, alignment, and secondary structure for the third domain of animal 12S rRNA. Mol. Biol. Evol. 13:150–169.[Abstract]

    Huson, D. H. 1997. SplitsTree. Version 2.4. Distributed by the author. University of Bielefeld, Germany.

    ———. 1998. SplitsTree: a program for analyzing and visualizing evolutionary data. Bioinformatics 14:68–73.

    Kjer, K. 1995. Use of rRNA secondary structure in phylogenetic studies to identify homologous positions: an example of alignment and data presentation from frogs. Mol. Phylogenet. Evol. 4:314–330.[ISI][Medline]

    Lipman, D. J., S. F. Altschul, and J. D. Kececioglu. 1989. A tool for multiple sequence alignment. Proc. Natl. Acad. Sci. USA 86:4412–4415.

    Lockhart, P. J., M. A. Steel, M. D. Hendy, and D. Penny. 1994. Recovering evolutionary trees under a more realistic model of sequence evolution. Mol. Biol. Evol. 11:605–612.[Free Full Text]

    McClure, M. A., T. K. Vasi, and W. M. Fitch. 1994. Comparative analysis of multiple protein-sequence alignment methods. Mol. Biol. Evol. 11:571–592.[Abstract]

    Morrison, D. A., and J. T. Ellis. 1997. Effects of nucleotide sequence alignment on phylogeny estimation: a case study of 18S rDNAs of Apicomplexa. Mol. Biol. Evol. 14:428–441.[Abstract]

    Needleman, S. B., and C. D. Wunsch. 1970. A general method applicable to the search for similarities in the amino acid sequence of two proteins. J. Mol. Biol. 48:443–453.[ISI][Medline]

    Thompson, J. D., D. G. Higgins, and T. J. Gibson. 1994. CLUSTAL W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, positions-specific gap penalties and weight matrix choice. Nucleic Acids Res. 22:4673–4680.[Abstract]

    Thorne, J. L., and H. Kishino. 1992. Freeing phylogenies from artifacts of alignment. Mol. Biol. Evol. 9:1148–1162.[Abstract]

    Titus, T. A., and D. R. Frost. 1996. Molecular homology assessment and phylogeny in the lizard family Opluridae (Squamata: Iguania). Mol. Phylogenet. Evol. 6:49–62.[ISI][Medline]

    Tönges, U., S. W. Perrey, J. Stoye, and A. W. M. Dress. 1996. A general method for fast multiple sequence alignment. Gene 172:33–41.

    Tyson, H. 1992. Relationships between amino acid sequences determined through optimum alignments, clustering and specific distance patterns: application to a group of scorpion toxins. Genome 35:360–371.

    van de Peer, Y., J.-M. Neefs, P. de Rijk, and R. de Wachter. 1993. Reconstructing evolution from eukaryotic small-ribosomal-subunit RNA sequences: calibration of the molecular clock. J. Mol. Evol. 37:221–232.[ISI][Medline]

    van de Peer, Y., E. Robbrecht, S. de Hoog, A. Caers, P. de Rijk, and R. de Wachter. 1999. Database on the structure of small subunit ribosomal RNA. Nucleic Acids Res. 27:179–183.[Abstract/Free Full Text]

    Wheeler, W. C. 1995. Sequence alignment, parameter sensitivity, and the phylogenetic analysis of molecular data. Syst. Biol. 44:321–331.[ISI]

    Wheeler, W. C., J. Gatesy, and R. DeSalle. 1995. Elision: a method for accommodating multiple molecular sequence alignments with alignment-ambiguous sites. Mol. Phylogenet. Evol. 4:1–9.[Medline]

    Wheeler, W., and D. Gladstein. 1994. Malign. Version 2.3. American Museum of Natural History, New York, NY.

Accepted for publication November 1, 1999.