8. Moving Beyond OTU Methods

8.1.4.1 Assumptions of OTU Methods in Numerical Taxonomy and Microbiome

OTU methods used in numerical taxonomy and microbiome have different aims. Adopting OTU methods into microbiome studies is somewhat inappropriate and thus is controversial.

Although at the beginning the concept of OTU is convenient for bioinformatic analysis of microbiome sequencing data, the use of the OTU methods in microbiome literature deviates from the considerations of the original authors. As reviewed in Chap. 7, the original authors of numerical taxonomy strictly separated taxonomic classifications from phylogenetic classifications and treated the latter as “phylogenetic speculation” and the former as “taxonomic procedure.” The original authors considered OTUs as operational not functional; the functionality remained to taxa; cluster analysis of similarity is merely on datasets at hand and requires that the characters are from larger and randomly chosen samples; taxonomic similarity does not represent the phylogenetic similarity. In contrast, by assigning sequences to a particular OTU, actually it assumes that the similarities among sequences are closely related to the phylogeny, and therefore the similarities of sequences can most likely be used to derive phylogenetic divergence or ecologically rank the similarity of taxa. However, the assumption that sequence similarity of 16S hypervariable regions is a good proxy for phylogenetic and therefore ecological similarity is problematic (Prosser et al. 2007; Tikhonov et al. 2015; Preheim et al. 2013) because sequence similarity does not imply ecological similarity. It was demonstrated that sequence similarity very poorly predicted ecological similarity (Tikhonov et al. 2015). The application of OTU methods in microbiome studies actually forces OTUs to play the functional role of taxa.

8.1.4.2 Challenges of Defining Species Using OTU Methods

OTU methods have been used in microbiome studies to classify taxa in classical Linnaean taxonomy at different taxonomic levels with similar DNA sequences.

Because 16S rRNA genes evolve at different rates in different organisms, determining taxonomic rank by sequence identity thresholds is approximate at the best (Woese 1987). Hierarchical clustering methods generally use a predefined clustering threshold to the hierarchical tree and then group sequences within the threshold into one OTU. For example, because a natural entity “species” cannot be identified as a group of strains that is genetically well separated from its phylogenetic neighbors, sequences are often clustered into OTUs as proxies for species, although sequence divergence is not evenly distributed in the 16S rRNA region. Thus, defining a species by a polyphasic approach is only a pragmatic approach.

Previously, the threshold value of 70% was used in the phylogenetic definition of a species (“approximately 70% or greater DNA-DNA relatedness and with 5°C or less”) (Wayne et al. 1987). This approach of applying DNA similarity to study bacterial species was well received and well proven and has been acknowledged since the 1970s (Stackebrandt and Goebel 1994; Johnson 1985; Steigerwalt et al. 1976). In 1994 Stackebrandt and Goebel (1994) proposed to use the “similarities of 97% and higher” sequence identity as the canonical clustering threshold when few 16S rRNA sequences were available. In practice the 3% dissimilarity is often chosen as the cutoff value to define bacteria species (Sogin et al. 2006; Sun et al. 2009; Schloss 2010; Huse et al. 2007) although different cutoff values of similarity have been studied such as 99.6–95.6% similarity (equating 1–13 bp) (Hathaway et al. 2017).

However, the cutoff value of similarities of 97% was based on the observation that species having 70% or greater DNA similarity usually have more than 97% sequence identity with these 3% or 45-nucleotide differences being concentrated mainly in certain hypervariable regions (Stackebrandt and Goebel 1994) given the primary structure of the 16S rRNA is highly conserved. Thus, the criterion of 97% similarities or 3% dissimilarities (Schloss and Handelsman 2005) is observed and experienced. It is not really able to define a natural entity “species.” Especially, it was shown that 97% threshold is not optimal for an approximation to species because it is too low. To define a species, at least 99% of OTUs should be used (Edgar 2018). In summary, this practical approach of using OTUs to define a species is some kinds of arbitrary, subjective, and therefore inaccurate.

8.1.4.3 Criteria of Defining the Levels of Taxonomy

Like the 3% dissimilarities of sequences that are used to define a species (Borneman and Triplett 1997; Stackebrandt and Goebel 1994; Hugenholtz et al. 1998; Konstantinidis et al. 2006; Goris et al. 2007), accordingly in practice the 95% sequence similarity (95% homology) or the dissimilarities of 5% by species are used to define a genus (Wolfgang Ludwig et al. 1998; Everett et al. 1999; Sait et al. 2002; Schloss and Handelsman 2005); 10% dissimilarities (90% homology) are used to define a family (Schloss et al. 2009); 20% dissimilarities (80% homology) are used to define a phylum (Schloss and Handelsman 2005); an interkingdom identity range of 70–85% is used to define a kingdom (Woese 1987; Borneman and Triplett 1997); and 45% dissimilarities (55% homology) are used to define a domain (Stackebrandt and Goebel 1994).

This practice of using an arbitrary fixed global clustering threshold to group amplicons into molecular OTUs has been criticized (Mahé et al. 2015b) because global clustering thresholds have rarely been justified and are not applicable to all taxa and marker lengths (Caron et al. 2009; Nebel et al. 2011; Dunthorn et al. 2012; Brown et al. 2015). Actually it is impossible to use an accurate distance-based threshold to define taxonomic levels (Schloss and Westcott 2011). Therefore, the criteria of defining the levels of taxonomy are arbitrary, subjective, and inaccurate.

8.1.4.4 Clustering Algorithms and the Accuracy of OTU Clustering Results

As reviewed in Chap. 7, the original OTU methods are based on cluster analysis. However, when the clustering-based OTU methods were adopted in microbiome studies, several limitations have been detected that are associated with the clustering algorithms.

The complete-linkage was originally recommended in numerical ecology because it is often desirable in ecology when one wishes to delineate clusters with clear discontinuities (Legendre and Legendre 1998), while the single-linkage clustering is implemented (Wei et al. 2012) in BlastClust (BLASTLab 2004) and GeneRage (Enright and Ouzounis 2000). Other researchers demonstrated that the average-linkage algorithm produces more robust OTUs than other hierarchical and heuristic clustering algorithms (i.e., CD-HIT, UCLUST, ESPRIT, and BlastClust) (Schloss and Westcott 2011; Sun et al. 2012; Quince et al. 2009). The average-linkage algorithm was further evaluated as having higher clustering quality compared to the other common clustering algorithms and even compared to greedy heuristic clustering algorithms, including distance-based greedy clustering (DGC) and abundance-based greedy clustering (AGC) (Westcott and Schloss 2015; Ye 2011). In such analysis, an OTU is generally considered presenting a bacterial species (Kuczynski et al. 2012). Therefore, using the average-linkage algorithm for clustering 16S rRNA gene sequences into OTUs has been suggested (Huse et al. 2010; Schloss and Westcott 2011). The suggestion of using the average-linkage clustering in microbiome data is in line with the suggestion of the original authors of numerical taxonomy (see Chap. 7 for details).

The controversial and challenges of clustering approaches are due to using different criteria to assess clustering quality or test accuracy of clustering, such as simulation, information theoretic-based normalized distance measurement (Chen et al. 2013b; Baldi et al. 2000; van Rijsbergen 1979), stability of OTU assignments (He et al. 2015), and Matthew’s correlation coefficient (MCC) (Schloss and Westcott 2011; Westcott and Schloss 2015), or combining the MCC and the F-score (Preheim et al. 2013; Baldi et al. 2000).

Different criteria often reach totally different conclusions. For example, many hierarchical and greedy clustering methods are unable to stably assign or produce OTUs. Thus, some researchers, such as He et al. (2015), thought that clustering all sequences into one OTU is completely useless for downstream analysis. By using MCC, Westcott and Schloss (2015) further demonstrated that the stability of OTUs actually suffered a lack of quality because the stability concept focused on the precision of the assignments and did not reflect the quality of the OTU assignments. Actually, this point was also raised in the numerical taxonomy.

In summary, the commonly used hierarchical and heuristic-based OTU methods suffered several fundamental problems, such as (1) the assumption of sequence similarity implying ecological similarity, (2) using arbitrary criteria to define the level of taxonomy (e.g., the dissimilarity of 3% OTUs is used to define species), (3) seed sensitivity, (4) input-order dependency, and (5) read preprocessing.

8.1.5.1 Defining Taxonomy

Initially when the OTU concept was proposed to be used in microbiome sequencing, the purpose is for the taxonomic grouping, i.e., identifying bacterial population boundaries (Callahan et al. 2017). However, as we reviewed above and others such as Murat Eren et al. (2016) and Benjamin Callahan et al. (2017), the original purpose of using OTUs has not been achieved by the clustering OTU methods. For example, it is problematic to define OTUs using 97% sequence similarity threshold in the closed-reference approach because two sequences might be 97% similar to the same reference sequence, maybe less similarity (e.g., only 94%) shared by each other (de Queiroz and Good 1997; Westcott and Schloss 2015). Thus, the use of 97% OTUs to construct “species” has not been largely approved. In practice a 16S rRNA gene sequence similarity threshold range of 98.7–99% was recommended for testing the genomic uniqueness of a novel isolate(s) (Erko Stackebrandt and Ebers 2006). Similarly, it has not been largely approved to use the 95%, 90% of sequence similarity to define a genus and a family, respectively, and so on.

8.1.5.2 Reducing the Impact of Sequencing Error

The second purpose of using OTU concept and methods is to differentiate true diversity from sequencing errors (Callahan et al. 2017). This purpose is different from and orthogonal to the first one. When the OTU concept was initially proposed and used in microbiome study, OTUs were not expected to reduce the impact of amplicon sequencing error. The task of correcting errors has been repurposed to OTU (Rosen et al. 2012). In other words, this task is an expectation or was weighted as having the function of reducing sequence errors because the use of 97% sequence similarity threshold has successfully reduced the impact of erroneous OTUs on diversity estimations (Eren et al. 2016). After the second role was repurposed to OTUs, clustering reads into OTUs by sequence similarity became a standard approach to filter the noise (Tikhonov et al. 2015).

However, the second purpose has not be achieved either. The use of OTUs for differentiating true diversity from sequencing errors is only appropriate for higher levels of taxonomy (e.g., phylum). When OTU methods were used for probing finer-scale diversity in lower levels of taxonomy, they have intrinsically high false-positive and false-negative rates because using arbitrary cutoff values not only overestimate diversity when there exist errors larger than the OTU-defining cutoff but also cannot resolve real diversity at a scale finer (Rosen et al. 2012). It was also shown that clustering reads into OTUs greatly underestimates ecological richness (Tikhonov et al. 2015). In other words, the use of OTUs has not reduced the controversy on diversity estimations and measures of community composition.

Taken together, both purposes of using OTU-based methods have problems to be achieved and the two major challenges remained: one is to identify bacterial population boundaries and another is to differentiate true diversity from sequencing errors (Preheim et al. 2013).

The third purpose of using OTUs to play the functionality of species in microbiome research is more complicated. To better understand whether this purpose is achieved, we write Sect. 8.1.6 to provide a theoretical background to discuss definitions of species and species-level analysis.

8.1.6.1 Eukaryote Species Concepts

The species concept has been discussed among microbiologists and eukaryote taxonomists (Claridge et al. 1997a). The species concepts for prokaryotes and eukaryotes are different (May 1986, 1988); briefly review of species concept for eukaryotes will definitely help to delineate the species concept for prokaryotes.

8.1.6.2 Prokaryote or Bacterial Species Concepts

The ability to isolate organisms in pure cultures started in 1872 cultivating pure colonies of chromogenic bacteria by Joseph Schroeter (Logan 2009). Cultivating microorganisms facilitates bacteria classification based on the phenotypic description of these organisms. However, early bacterial species was often defined based on monothetic groups (a unique set of features). This is a phenotypic definition of species, which lies on a database with an accurate morphologic and phenotypic (e.g., biochemical) description of type strains or typical strains and comparing the isolate to be identified to the database using standard methods to determine these characteristics for the isolate (Clarridge 2004). The limitation is that these sets of phenotypic properties were subjectively selected (Goodfellow et al. 1997).

Among the more than 20 species concepts, Mayr’s BSC is the most widely known and most controversial concept (Rosselló-Mora and Amann 2001; Sokal 1973). BSC attempts to unify genetics, systematics, and evolutionary biology (Claridge et al. 1997a). However, it lacks of practicability and hence generally it was agreed that Mayr’s BSC should be abandoned (Rosselló-Mora and Amann 2001; Sokal 1973). Therefore, BSC cannot presently be applied to prokaryotes (Stackebrandt and Goebel 1994).

The development of numerical taxonomy in the late 1950s in parallel to application of computer for multivariate analyses makes bacteria classification toward an objective approach (Sneath 1989). The numerical taxonomy began to establish when modern biochemical analytical techniques advanced to study the distributions of specific chemical constituents (i.e., amino acids, proteins, sugars, and lipids in bacteria) (Rosselló-Mora and Amann 2001; Logan 2009).

In the late 1970s, with the advent of cataloging ribosomal ribonucleic acids (rRNAs) (Stackebrandt et al. 1985), the rRNA and especially 16S rRNA sequences were shown to be a very useful molecular marker for phylogenetic analyses (Ludwig and Schleifer 1994).

Clustering-based OTU methods have been combined with the 16S rRNA sequencing. The 16S rRNA sequencing method has been widely used for the prokaryotic classification (Olsen et al. 1994), and the datasets generated from 16S rRNA sequencing studies have also been increasingly used to propose new bacterial species (Stackebrandt and Goebel 1994).

The determination of relationships between distantly related bacteria has been evidenced a remarkable breakthrough in the late 1970s by cataloging 16S ribosomal ribonucleic acids (rRNAs) (Stackebrandt et al. 1985) and DNA-RNA hybridization (De Ley and De Smedt 1975) and in the mid-1980s by the full sequence analysis of rRNA (Rosselló-Mora and Amann 2001). In the 1980s, determining phylogenetic relationships of bacteria and all life forms by comparing a stable part of the genetic code (commonly the 16S rRNA gene) became a new standard for identifying bacteria (Woese et al. 1985; Woese 1987).

Due to above two important properties, bacterial phylogeny (the genealogical trees among the prokaryotes) can be established by comparing the 16S rRNA gene sequences (Ludwig et al. 1998). PSC includes monophyletic species concept (MSC) and diagnostic species concept (DSC) (Hull 1997). Both are defined as phylogenetic (or genealogical) concepts with a minimal time dimension (Rosselló-Mora and Amann 2001).

When DNA was discovered, prokaryote classification was based solely on phenotypic characteristics. With the development of numerical taxonomy (Sneath and Sokal 1973), the individuals are treated as operational taxonomic units that are polythetic (they can be defined only in terms of statistically covarying characteristics), resulting in a more objective circumscription of prokaryotic units. The discovery of genetic information gave a new dimension to the species concept for microorganisms, which enable at least a first rough insight into phylogenetic relationships. Thus, the species concept for prokaryotes evolved into a mostly phenetic or polythetic. This means that species are defined by a combination of independent, covarying characters, each of which may occur also outside the given class, thus not being exclusive of the class (Van Regenmortel 1997).

Compared to DNA-RNA hybridization approach, the 16S rRNA gene sequence approach is much simpler and thus has replaced DNA-RNA hybridization and become the new gold standard to define a species (Fournier et al. 2003; Harmsen 2004). However, as we will describe in Sect. 8.1.6.3, using the 16S rRNA method to define a species still remains challenging. Before we move forward to definition of species in 16S rRNA method, we quote some microbiological definitions of species.

8.1.6.3 16S rRNA Method and Definition of Species

Using the 16S rRNA similarities to define bacteria species remains challenging. Much early in the 1870s, Ferdinand Cohn investigated whether the similarity exits between bacteria and whether bacteria, like animals and plants, can be arranged in distinct taxa (Cohn 1972; Schlegel and Köhler 1999). Cohn considered the form genera as natural entities but considered species as largely artificial. Cohn emphasized that character of the classification system is artificial. The original authors of OTU method Sokal and Sneath (1963) also criticized their then current taxonomy made little increase in understanding the nature and evolution of the higher categories (Sokal and Sneath 1963) (p. 5).

When clustering-based OTU methods emerged into the 16S rRNA sequencing analysis, an important concept, the phylogenetic species concept (PSC), was developed. PSC compares the 16S rRNA similarities with DNA similarity (DNA-DNA similarities), assuming that nucleic acid (DNA) reassociation can determine whether taxa are phylogenetically homogeneous. PSC defines a species that would generally include strains with “approximately 70% or greater DNA-DNA relatedness and with 5°C or less ΔT_m” (Wayne et al. 1987). Here, T_m is the melting temperature or the thermal denaturation midpoint, and ΔT_m is the difference between the homoduplex DNA T_m and the heteroduplex DNA T_m. Thus, ΔT_m is a reflection of the thermal stability of the DNA duplexes.

However, as reviewed in Sect. 8.1.4.2, the rationale for using DNA reassociation as the gold standard to delineate species originates from the observations that a high degree of correlation existed between DNA similarity and chemotaxonomic, genomic, serological, and numerical phenetic similarity (Stackebrandt and Goebel 1994). The 16S rRNA sequencing method does not have sufficient power to correctly define bacterial species, which was reviewed in several articles including Ash et al. (1991); Amann et al. (1992); Fox et al. (1992); Martinez-Murcia et al. (1992); Stackebrandt and Goebel (1994); and Johnson et al. (2019).

In summary, clustering sequence similarity of 16S rRNA can never be solely used to define bacterial species. On the one hand, it was recognized that the DNA hybridization (DNA reassociation) remains the optimal method for measuring the degree of relatedness between highly related organisms and hence cannot be replaced by 16S rRNA method for defining species (Stackebrandt and Goebel 1994; Rosselló-Mora and Amann 2001). On the other hand, it has been accepted that a prokaryotic species classification should be defined by integrating both phenotypic and genomic parameters, as well as should analyze and compare as many phenotypic and genomic parameters as possible (Wayne et al. 1987). This approach is known as “polyphasic taxonomy” (Vandamme et al. 1996). Based on this polyphasic approach, in 2001 a more pragmatic prokaryote species concept (called a phylo-phenetic species concept) has been described as “a monophyletic and genomically coherent cluster of individual organisms that show a high degree of overall similarity in many independent characteristics, and is diagnosable by a discriminative phenotypic property” (Rosselló-Mora and Amann 2001). Obviously, this is far beyond the clustering-based OTU method can achieve.

8.1.6.4 16S rRNA Method and Physiological Characteristics

Whether the 16s rRNA method can inference the physiological characteristics is debatable. On the one hand, the 16S rRNA approach has the practical advantage for species identification and has been used frequently in bacterial phylogenetics, species delineation, and microbiome studies; for example, it is recommended to include the ribosomal sequence for describing new prokaryotic species (Ludwig 1999; Chakravorty et al. 2007). Some researchers considered that comparative analysis of 16S rRNA is a very good method to find a first phylogenetic affiliation in both potentially novel and poorly classified organisms (Goodfellow et al. 1997). Thus, although the 16S rRNA sequence analysis is not sufficient for numerically delineating borders of the prokaryotic species, it can indicate the ancestry pattern of the taxon studied, as well as confirm the monophyletic nature in grouping the members (Rosselló-Mora and Amann 2001).

On the other hand, other researchers have questioned the validity of 16S rRNA as a marker for phylogenetic inferences and thought that it is very difficult to infer the physiological characteristics of an organism by clustering the rRNA sequence method (Rosselló-Mora and Amann 2001; Cohan 1994; Gupta 1998; Gribaldo et al. 1999). Thus, using arbitrary threshold of 97% similarity (or 3% dissimilarity) to define OTUs that could be interpreted as a proxy for bacterial species is computationally convenient, but at the expense of accurate ecological inference (Eren et al. 2016) because 3% OTUs are often phylogenetically mixed and inconsistent (Koeppel and Wu 2013; Nguyen et al. 2016; Eren et al. 2014).

Recently, Hassler et al. (2022) found that phylogenies of the 16S rRNA gene and its hypervariable regions lack concordance with core genome phylogenies, and they concluded that the 16S rRNA gene has a poor phylogenetic performance and has far-reaching consequences.

By comparing core gene phylogenies to phylogenies constructed using core gene concatenations to estimate the strength of signal for the 16S rRNA gene, its hypervariable regions, and all core genes at the intra- and inter-genus levels, Hassler et al. (2022) showed that at both intra- and inter-genus taxonomic levels, the 16S rRNA gene was recombinant and subject to horizontal gene transfer, suffering from intragenomic heterogeneity, accumulating recombination and an unreliable phylogenetic signal. The poor phylogenetic performance using the 16S rRNA gene not only results in incorrect species/strain delineation and phylogenetic inference but also has the potential to confound community diversity metrics if incorporating phylogenetic information in the analysis such as using Faith’s phylogenetic diversity and UniFrac (thus, their use to measure the diversity is not recommended) (Hassler et al. 2022). Additionally, taxonomic abundance and hence measures of microbiome diversity may also be inflated and confounded by the multiple copies within a genome and wide range in 16S rRNA gene among genomes (De la Cuesta-Zuluaga and Escobar 2016; Louca et al. 2018).

8.2.2.1 Distribution-Based Clustering (DBC)

In 2013, Preheim et al. (2013) used DBC to refine the OTU for single-nucleotide resolution but without relying on sequence data alone to serve the dual purpose for the use of OTUs: identifying taxonomic groups and eliminating sequencing errors. Alternatively, the DBC method compares the distribution of sequences across samples. In other words, it uses ecological information (distribution of abundance across multiple biological samples) to supplement sequence information.

DBC algorithm is implemented under the framework of clustering OTU methods. It first uses Jukes-Cantor-corrected genetic distances or Jensen-Shannon divergence (JSD) as appropriate to measure the genetic distances among sequences and then implements either complete algorithm (all sequences are analyzed together in the analysis) or parallel algorithm (sequences are pre-clustered with a heuristic approach) or different algorithms (e.g., nearest-neighbor single-linkage clustering). DBC first uses the chi-squared test to determine whether two sequences have similar distributions across libraries. Then, based on ecological distribution to form OTUs, so 16S rRNA sequences that even differ by only 1 base but that are found in different samples are still put into different OTUs. Conversely, the sequences drawn from the same underlying distribution across samples are grouped together into the same OTU regardless of interoperon variation or sequence variation or random sequencing errors. The underlying assumption is that bacteria in different populations are often highly correlated in their abundance across different samples (Preheim et al. 2013). Thus, 16S rRNA sequences derived from the same population regardless of their variations or random sequencing error will have the same underlying distribution across sampled environments.

In summary, DBC is a clustering OTU method or a refined OTU method. The difference from other OTU methods is that it uses ecology to model multiple samples to refine the OTU. As a refined OTU method, DBC is still within traditional framework of OTU methods.

8.2.2.2 Swarm2

The de novo amplicon clustering methods share two fundamental problems (reviewed in Sect. 6.​4). One uses arbitrary fixed global clustering thresholds; another is input-order dependency induced by centroid selection, i.e., the clustering results are strongly influenced by the input order of amplicons. To solve these two fundamental problems, in 2014, Mahé et al. (2014) proposed an amplicon clustering algorithm Swarm to fine-scale OTUs without relying on arbitrary global clustering thresholds and input-order dependency.

However, Swarm has the disadvantage: it cannot resolve at the single-base level (Hathaway et al. 2017).

Like distribution-based clustering (Sect. 8.2.2.1), Swarm and Swarm2 aim to provide single-nucleotide resolution (Mahé et al. 2014, 2015b). Other bioinformatic software that aim to single-nucleotide resolution include oligotyping (Section 8.3.1), denoising-based methods (Sect. 8.3.2) such as cluster-free filtering (Sect. 8.3.2.3), DADA2 (Sect. 8.3.2.4), UNOISE2 (Sect. 8.3.2.5), Deblur (Sect. 8.3.2.6), and SeekDeep (Sect. 8.3.2.7).

8.3.2.1 Denoising-Based Methods Versus Clustering-Based Methods

8.3.2.2 Pyrosequencing Flowgrams

Compared to denoising pyrosequencing flowgrams, Illumina denoisers have been developed more recently. Some Illumina methods can also be used with the 454 sequencing platform, such as cluster-free filtering (CFF) and SeekDeep. Below we will introduce CFF in Sect. 8.3.2.3, DADA2 in Sect. 8.3.2.4, UNOISE2 in Sect. 8.3.2.5, Deblur in Sect. 8.3.2.6, and SeekDeep in Sect. 8.3.2.7.

8.3.2.3 Cluster-Free Filtering (CFF)

In 2015, Tikhonov et al. (2015) proposed a cluster-free filtering (CFF) denoiser approach to achieve sub-OTU resolution with ecologically differing by as little as one nucleotide (nt) (99.2% similarity). The CFF denoiser was proposed in the contexts: on the one hand, to improve estimates of ecological diversity by reducing the impact of amplicon sequencing error, researchers have tried various denoising algorithms to reassign or eliminate noisy reads; however, the challenging still remains because no approach can fully address the denoiser issues. On the other hand, the alternative DBC approach is necessarily unreliable due to computational cost and severe limitation for low-count sequences in the cross-sample comparisons. Under this background, the focuses of the CFF denoiser are not to further try to improve the existing OTU clustering approaches and also not to attempt to identify rare species. Instead, the CFF denoiser combines error-model-based denoising and systematic cross-sample comparisons to resolve the fine (sub-OTU) structure of moderate-to-high-abundant taxonomy data (Tikhonov et al. 2015).

CFF employs cross-sample comparisons to achieve sub-OTU resolution. CFF was developed based on the observation that sequence similarity need not imply dynamical similarity and vice versa in time-series data. In contrast, the maximum correlation between the time traces (normalized counts versus observation day) of two sequences depends on their abundance measured by Poisson model. Thus CFF denoiser (Tikhonov et al. 2015) first defines the “dynamical similarity” of two traces as the Pearson correlation of their abundance, normalized by their maximum possible correlation. Then it models the abundance traces of these two sequences by an additive Poisson counting noise, setting an upper bound on the correlation coefficient for low-abundance sequences because Poisson sampling noise becomes non-negligible for low-abundance sequences. Finally it uses the Hamming distance metric to measure the pairwise sequence distances based on normalized counts.

Three important characteristics that characterize CFF as a method of moving beyond the traditional OTU methods are as follows: (1) aims to differ 16S sequences by one nucleotide or a single base pair, which is similar to oligotyping and distribution-based clustering approaches; (2) does not rely on clustering similar sequences together, which is similar to oligotyping; and (3) does not assume that sequence similarity implies ecological similarity. Instead CFF actually considers sequence similarity is a very poor predictor of ecological similarity (Tikhonov et al. 2015). Thus CFF uses highly conservative filtering criteria, assuming that sequence similarity contains only true biological sequences, that is, there are no false positives. Regarding this point, CFF is similar to Sokal and Sneath’s numerical taxonomy (Sokal and Sneath 1963).

However, CFF also has limitations: (1) It is specifically designed to be run on large multi-sample datasets. (2) It is not a replacement for OTU clustering; it focuses on moderate-to-high-abundance sequences and hence discards low-abundance sequences. Thus, CFF is unsuitable for studying population-level alpha or beta diversity (Tikhonov et al. 2015), which probably is an important limitation of this method.

8.3.2.4 DADA2

DADA2 and q2-dada2 plugin were introduced in Sect. 4.​2.​1. In line with providing single-nucleotide resolution without using arbitrary taxonomic similarity thresholds to define molecular OTUs, Callahan et al. (2017) proposed to use exact “sequence variants” to replace OTUs in marker-gene data analysis and developed software package DADA2 (DADA: Divisive Amplicon Denoising Algorithm) for modeling and correcting Illumina-sequenced amplicon errors. The proposal of using ASVs instead of the traditional OTUs was due to the limitations of OTU methods, such as clustering sequence reads into OTUs often eliminates biological information present in the data, and OTUs are not species, and their construction is not necessitated by amplicon errors (Callahan et al. 2016a), and especially was motivated by a series of new computational methods (Eren et al. 2013, 2016; Tikhonov et al. 2015; Edgar 2016b; Amir et al. 2017).

DADA2 is a parametric model that relies on input read abundances and distances. It assumes that true reads are likely to be more abundant and less abundant reads are likely error-derived which may be only a few base differences away from a more abundant sequence (Prodan et al. 2020; Callahan et al. 2016a). DADA2 was developed based on the DADA (Rosen et al. 2012). DADA was built on AmpliconNoise’s error-modeling approach (Quince et al. 2011). The distinctive feature of DADA is its divisive hierarchical clustering algorithm, while previous methods, including AmpliconNoise and simple OTU clustering, use agglomerative approach.

The core denoising algorithm (Divisive Amplicon Denoising Algorithm) in the DADA2 is based on a Poisson model to estimate the errors in Illumina-sequenced amplicon reads. This error model estimates the probability for each amplicon read with each sequence is produced from sample sequence as a function of sequence composition and quality. The Poisson model assumes that the number of amplicon reads (the abundance) of each sequence is consistent with the error model and the p-value of the null hypothesis is calculated. The calculated p-values or quality scores (base Q-scores) are used as the division criteria for an iterative partitioning algorithm; sequencing reads will continue to be divided until all partitions are consistent with being produced from their central sequence.

Overall, DADA2 has been recognized as so far the best denoising-based method for focusing on the highest possible biological resolution (e.g., differentiating closely related strains) (Prodan et al. 2020).

8.3.2.5 UNOISE2 and UNOISE3

To improve error correction in Illumina 16S and ITS amplicon sequencing, distinguishing the sequence errors caused by PCR and sequencing from true biological variation, Edgar (2016b) proposed UNOISE2, an updated version of the UNOISE algorithm (Edgar and Flyvbjerg 2015) for denoising (error-correcting) Illumina amplicon reads.

UNOISE2 uses a one-pass clustering strategy without using quality (Q) scores. The clustering strategy is different from DADA2’s model-based approach, which uses quality scores in an iterative divisive partitioning clustering strategy. In UNOISE2, both algorithms of UCHIME2 (Edgar 2016a) and DADA2 are incorporated to reduce the number of incorrect sequences. However, both UNOISE2 and DADA2 were reviewed as still generating some amount of false positives with DADA2 likely to increase the number of false-positive chimera predictions (Edgar 2016b). In UNOISE2, the obtained predicted biological sequences are called as ZOTUs (zero-radius OTUs) (Edgar 2016b, 2018). The ZOTU is similar to Callahan et al.’s ASV. In DADA2, the 97% similarity of OTUs is replaced by ASVs as the standard unit of marker-gene analysis and reporting; in UNOISE2, similarly the [97%] OTUs are replaced by ZOTUs (Edgar 2018).

The UNOISE3 denoising method first uses the “unoise3” command to rank sequences in decreasing order of abundance and then discards those sequences with counts less than the specified minimum abundance threshold (the default value is 8). By using 100% identity threshold, each distinct sequence defines a separate OTU. Both OTU-level clustering with UPARSE and ASV-level denoising with UNOISE3 are implemented via USEARCH, and actually the author of USEARCH (Edgar 2010) recommended that UPARSE and UNOISE3 should be performed together.

8.3.2.6 Deblur

We illustrated Qiime2-Deblur in Chap. 4 (Sect. 4.​4). Amir et al. (2017) developed a novel sub-OTU (sOTU) method called Deblur to fast and accurately identify exact amplicon sequences and to integrate large datasets. Like DADA2 and UNOISE2, Deblur aims to identify real ecological differences between taxa with a single base pair of amplicons differentiation. Deblur in concept is similar as DADA2 and UNOISE2, but uses a different algorithm to obtain single-nucleotide resolution. Deblur performs each sample independently and compares sequence-to-sequence Hamming distances within a sample to an upper-bound error profile by combining with a greedy algorithm (Amir et al. 2017).

The Deblur algorithm is implemented via three steps (Amir et al. 2017): First, Deblur sorts sequences by abundance. Second, Deblur subtracts the number of predicted error-derived reads from the counts of neighboring reads based on their read-to-read Hamming distance using an upper bound on the error probability based on the most to least abundant sequence. Finally, Deblur removes any sequence whose abundance drops to 0 after a subtraction. Thus, after applying Deblur, the invalid (i.e., noise) sequences are removed, and only reads likely to have been presented to the sequencer are retained.

Deblur also has the disadvantages, such as Deblur has poor performance compared to DADA2 and open-reference OTU clustering when it is used to detect low abundant taxa in the extreme dataset at 97% identity (Nearing et al. 2018). For comparing Deblur and DADA2, the reader is also referred to Sect. 4.​4.​6.

8.3.2.7 SeekDeep

SeekDeep (Hathaway et al. 2017) is an error-correction method for de novo (i.e., reference-free) analysis of amplicons. SeekDeep was developed in the contexts of comparing to three pipelines that aim for single-base resolution to determine the local PCR amplicon haplotypes (simply as haplotypes for brevity), MED (Eren et al. 2014), DADA2 (Callahan et al. 2016a), and UNOISE in USEARCH (Edgar 2016b), as well as two OTU-based clustering pipelines which cannot resolve at the single-base level: USEARCH (aka UCLUST/UPARSE) (Edgar 2016b) and Swarm (Edgar 2013). SeekDeep aims for single-base resolution to provide improved resolution toward species and strains. The central algorithm of SeekDeep is the qluster (for quality clustering) that improves the correction of PCR and sequencing errors through base quality values and k-mer frequencies and other multiple key ways. SeekDeep can be utilized with Ion Torrent, 454, and Illumina sequencing technologies. SeekDeep consists of four components (Mahé et al. 2014): (1) extractor for de-multiplexing and read filtering, (2) qluster for rapid and accurate clustering based on quality, (3) processClusters for replicate and population comparisons, and (4) popClusteringViewer for viewing and manipulating final results. The first three main components are central to generating clustering results, and the fourth is an additional component to aid in viewing and sharing the results.

However, Hathaway et al. (2017) also acknowledged that SeekDeep creates more false haplotypes than DADA2 and UNOISE although the false haplotypes have much lower abundance.

8.3.2.8 Remarks on Denoising-Based Methods

Denoising-based methods are widely used for identifying low-abundance or rare species against a noisy background, often with the aim of improving estimates of ecological diversity (Hathaway et al. 2017), such as Chao1 richness.

Using 97% identity may merge phenotypically different strains with distinct sequences into a single cluster (Tikhonov et al. 2015; Callahan et al. 2016a, b, c). In contrast, the denoising usually results in a set of predicted biological sequences or valid OTUs, which often have single-base resolution. Thus, compared to OTU methods, in most cases denoising algorithms provide the maximum possible biological resolution and hence superior to conventional OTU methods.

However, the denoising-based methods cannot fully address all the issues of bioinformatic analysis of microbiome data although these methods were proposed for identifying rare species and improving the accuracy of ecological diversity estimates, because (1) all error models are necessarily approximate, and no denoising algorithm can deal with errors that are beyond their model capabilities, and (2) the denoising-based methods are particularly problematic when they are used for calling low-abundance species. As reviewed above, DADA2 had the best sensitivity and resolution, but also produced more number of spurious ASVs, while UNOISE2 was fast, but was not recommended for using singlet sequences due to decreased haplotype recovery at lower read depths.