ARK: Aggregation of Reads by K-Means for Estimation of Bacterial Community Composition

Abstract

Motivation: Estimation of bacterial community composition from high-throughput sequenced 16S rRNA gene amplicons is a key task in microbial ecology. Since the sequence data from each sample typically consist of a large number of reads and are adversely impacted by different levels of biological and technical noise, accurate analysis of such large datasets is challenging.

Results: There has been a recent surge of interest in using compressed sensing inspired and convex-optimization based methods to solve the estimation problem for bacterial community composition. These methods typically rely on summarizing the sequence data by frequencies of low-order k-mers and matching this information statistically with a taxonomically structured database. Here we show that the accuracy of the resulting community composition estimates can be substantially improved by aggregating the reads from a sample with an unsupervised machine learning approach prior to the estimation phase. The aggregation of reads is a pre-processing approach where we use a standard K-means clustering algorithm that partitions a large set of reads into subsets with reasonable computational cost to provide several vectors of first order statistics instead of only single statistical summarization in terms of k-mer frequencies. The output of the clustering is then processed further to obtain the final estimate for each sample. The resulting method is called Aggregation of Reads by K-means (ARK), and it is based on a statistical argument via mixture density formulation. ARK is found to improve the fidelity and robustness of several recently introduced methods, with only a modest increase in computational complexity.

Availability: An open source, platform-independent implementation of the method in the Julia programming language is freely available at https://github.com/dkoslicki/ARK. A Matlab implementation is available at http://www.ee.kth.se/ctsoftware.

Fig 1. A flow-chart of the ARK method.

Fig 2. Results for the random K-means clustering on the simulated data.

Mean VD error at the genus level as a function of the number of clusters. Note the improvement that ARK contributes to each method.

Fig 3. Results for the random K-means clustering on the simulated data.

Mean execution time increase (factor given in comparison to running SEK or Quikr in the absence of ARK) as a function of number of clusters. The dashed line represents a line with slope 1.

Fig 4. Comparison of the underlying algorithms with and without ARK.

Results are for the random K-means clustering on the simulated data when fixing the number of clusters to 75. Mean VD error at the genus level. Included for comparison are results for RDP’s NBC (compare to Fig 2(b) of [3]).

Fig 5. Comparison of the underlying algorithms with and without ARK.

Results are for the random K-means clustering on the simulated data when fixing the number of clusters to 75. Boxplot of the individual simulated sample execution times. Mean execution times for Quikr and ARK Quikr were 1.75 seconds and 4.71 minutes, while for SEK and ARK SEK they were 21.26 seconds and 19.21 minutes respectively. Mean execution time for RDP’s NBC was 38.19 minutes.

Fig 6. Total execution time for each method on the 28 samples of real biological data.

Fig 7. PCoA plots using the Jensen-Shannon divergence for RDP’s NBC.

Fig 8. PCoA plots using the Jensen-Shannon divergence for ARK SEK.

Fig 9. ARK Quikr PCoA plots (using the Jensen-Shannon divergence) on the real biological data.

In this case, we have labeling by body site. Note the clustering.

Fig 10. ARK Quikr PCoA plots (using the Jensen-Shannon divergence) on the real biological data.

In this case, we have labeling by variable region. Note the clustering.