TAXONDC: CALCULATING THE SIMILARITY VALUE OF THE 16S RRNA GENE SEQUENCES OF PROKARYOTES OR ITS REGIONS OF FUNGI

1. Introduction

Analysis of 16S rRNA gene sequences has an important role in systematics of prokaryotes (Kämpfer, Glaeser, 2013; Kim, Chun, 2014; Tindall et al., 2010). There are two ways to use 16S rRNA gene sequences: for phylogenetic analyses following multiple sequence alignments, and for calculating pairwise sequence similarities. The first approach is used to find out evolutionary relationships between related taxa and is important for generic or suprageneric classification. While the second approach provides a simple way for identification and delineation of novel isolates, and is a critical checkpoint at the species level (Kim, Chun, 2014; Stackebrandt, Ebers, 2006; Stackebrandt, Goebel, 1994).

The similarity value between sequences of the prokaryotic isolate and the closest species is calculated using a proper and robust global pairwise sequence alignment algorithm, not using local sequence alignment or fast searching algorithm (Tindall et al., 2010; Kim, Chun, 2014). Suchlike capability is available on a popular web-based tool EzBioCloud. However, in contrast to the first version (Chun et al., 2007), in which it was possible to compare any two sequences, the most recent versions of this tool (Kim et al., 2012; Yoon et al., 2016) are designed to identify isolates using their own database. Therefore, it is not suitable for routine comparison of two or more desirable sequences.

There are a lot of programs for the alignment of the sequences and determination of their similarity values. However, these programs have a number of shortcomings. First, the similarity values calculated with these software do not always match with the values calculated using the widely accepted taxonomy-oriented EzBioCloud. Second, the software usually makes multiple alignment of all the entered sequences and is not suitable for fully independent pairwise comparisons of sequences. For this reason, these programs might produce different similarity values depending on the sequences included in the analysis. Thus, there exists a necessity for software, which uses EzBioCloud algorithm, and allows a comparison of two or more desirable sequences between themselves rather than with pre-set database.

At present, the fungal taxonomy is based primarily on analysis of the internal transcribed spacer (ITS) region of the nuclear ribosomal repeat unit (Peay et al., 2008; Schoch et al., 2012). To facilitate ITS-based molecular identification of fungi, specific databases were created like UNITE (Kõljalg et al., 2004) and EzBioCloud. These tools also provide the possibility of determining the similarity value between sequences of the fungal isolate and the closest species using pre-set database, but do not allow the comparison of a desired set of sequences.

2. Implementation

Here, we introduce a new program TaxonDC (Taxon Distance Calculator) for calculating the 16S rRNA gene sequence similarities of prokaryotes or ITS regions of fungi. TaxonDC is a desktop software written in C++ and using Qt libraries. It could be run under Windows and Linux. This program has easy-to-use interface and could operate in two modes.

In the first mode the sequences to compare could be entered in the ‘Pair Input’ tab (Supplementary Fig. S1a) for single comparison of two sequences. Resulting alignment, similarity value, and statistics are shown in ‘Pair Output’ tab (Supplementary Fig. S1b).

The second mode is the key feature of proposed TaxonDC that allows performing the independent pairwise alignment of more than two sequences with subsequent calculation of similarity values. Our program is based on pairwise alignment only, but not on multiple sequence alignment. Therefore, pairs of sequences have their own results of comparison, which do not depend on each other. In this mode initial sequences are entered in the ‘Fasta Input’ tab (Supplementary Fig. S1c) in FASTA format. Results of calculation are shown in ‘Fasta Output’ tab (Supplementary Fig. S1d). Summary table with sequence similarity values is located at the top of the ‘Fasta Output’ tab. If the cell with the similarity value is chosen in this table, then alignment and statistics are displayed in the bottom area for the visual inspection. If the sequence name is selected in the table, all alignments with this sequence will be displayed. The results could be exported in TSV format for Microsoft Excel or OpenOffice/LibreOffice Calc.

Selecting the checkbox in the bottom of input tabs allows using pre-aligned sequences in both comparison modes.

TaxonDC uses CLUSTALW (Thompson et al., 1994) for global sequence alignment. The alignment gaps and ambiguous bases (e.g., N) are not considered during calculation of similarity (Supplementary Fig. S2). If ones are considered, similarity values will be lower in most cases, especially in low quality sequences.

The order of entered sequences has an effect on the result of the alignment, thereby affecting the calculated similarity value. This happens due to specific feature of the algorithm of the global sequence alignment. This specific feature is inherent to the EzBioCloud too, and hard to exclude completely keeping the results comparable with the previous ones. For instance, the 16S rRNA gene sequence similarity values of the Ignicoccus islandicus (CP006867) and Ignicoccus pacificus (AJ271794) pair compared in both directions are 98.35% or 98.27% (Supplementary Fig. S3–S4, Supplementary Table S2), and the ITS region similarity values of the Geosmithia flava (HF546291) and Geosmithia morbida (FN434081) pair compared in both directions are 98.06% or 98.25% (Supplementary Fig. S5–S6, Supplementary Table S1). Thus, this phenomenon may occur already at the species level. To consider this issue in TaxonDC, each pair of sequences is compared in both directions (A vs B and B vs A) to reveal any possible values difference. In a case, when such difference is found, two alignments and corresponding statistics are displayed. If the order of entered sequences gives different alignment results in the second mode, the corresponding similarity values are marked in bold in the summary table. Furthermore, different values are highlighted by color. For such situation, we recommend to choose a larger value since it corresponds to potentially a smaller number of evolution events. We believe it could help prevent discrepancies among different researchers.

To verify the TaxonDC capabilities, the sequences of type strains from the previous version of EzBioCloud (Kim et al., 2012) were used. Total of 180 sequence pairs were downloaded from the domains of Bacteria (120), Archaea (30), and the kingdom of Fungi (30) (Supplementary Table S1). Each pair was randomly selected. Further, total of 60 sequence pairs were downloaded from the domains of Bacteria (30) and Archaea (30) from the latest version of tool (Yoon et al., 2016) (Supplementary Table S2). The 16S rRNA gene sequence and ITS region similarities compared were identical in both cases calculated using the TaxonDC and the EzBioCloud (Supplementary Table S1, S2).

3. Conclusion

TaxonDC was developed to provide an improved tool for calculation of pairwise similarity values between desired sequences. To our knowledge, among publicly available programs, TaxonDC is the only which combines all previously described features: comparability with EzBioCloud, fully independent pairwise alignment and revealing of the effect of the order of sequence input on calculated values. Our program could be useful in phylogenetics, ecology and systematics of prokaryotes and eukaryotes.