4 | Discussion

In this study, our main objectives were to assess the level of correspondence between the taxonomic and functional profiles derived from amplicon-based 16S rRNA metabarcoding data and shotgun metagenomics, evaluate the strengths and weaknesses of both approaches, and determine whether functional profiles from HSP methods can be used as a substitute to metagenomics for monitoring functional changes associated with fish farming in marine environments.
The overall taxonomic richness recovered by 16S rRNA metabarcoding was up to 20-fold higher (depending on the taxonomic level) than metagenomics. This is due, in part, to the differences in the ability of identifying sequences of different lengths (~ 450 bp [metabarcoding] versus ~ 150 bp [metagenomics]) and by differences in the taxonomic assignment methods used. However, it is also well recognized that metagenomics requires much more sequencing effort than metabarcoding to reach equivalent 16S rRNA coverage as it captures all DNA material (Cottier et al., 2018; Singer, Shekarriz, McCarthy, Fahner, & Hajibabaei, 2020). This is especially problematic in highly diverse communities such as the marine sediment samples assessed in this study⁠. Small differences in taxonomic composition were anticipated because 16S rRNA primer sets are never truly universal (Pollock et al., 2018). While both datasets were relatively similar in terms of dominant Phyla and Families, several taxa such as BD2-2 (Bacteroidetes), Sulfurovaceae, Thermoanaerobaculaceae and Rubritaleaceae were more predominant in the metabarcoding data. These results clearly illustrate the advantage of using targeted amplicon 16S profiling over metagenomics when it comes to providing a comprehensive overview of communities in complex environments, although biases due to preferential PCR amplification and primer specificity are inevitably introduced.
Despite the differences in taxonomic richness, functional richness between the HSP methods and metagenomics was relatively similar. This counter intuitive result occurs because most ECs/KOs are typically redundant across microbial communities (Louca et al., 2018; Starke et al., 2020), many of which performing functions that are essential to cellular activities. Additionally, several non-ubiquitous functions, and especially those occurring at shallow phylogenetic depth, are difficult to accurately predict (Martiny et al., 2015) and may therefore be omitted by HSP methods (Bowman & Ducklow, 2015)⁠. These false negatives were prominent for Paprica, which had the lowest specificity. The opposite scenario where functional richness is artificially increased due to false positive predictions is also a possibility. For example, phylogenetic plasticity and genomic variability can result in loss of functions within taxa that can’t be detected by HSP methods. However, because of the substantial differences in terms of recovered taxonomic diversity between metagenomics and 16S metabarcoding methods and limits imposed by sequencing depth, it is also possible that some functions predicted by the HSP methods were not detected by metagenomics. As such, the true sensitivity of the HSP methods, which was lowest for Tax4Fun2 and highest for paprica, was likely underestimated.
In general, we observed little correlation and/or no significant difference with the null expectation datasets based on the abundance of functions shared between HSP methods and metagenomics. Other studies have reported weak correlations of HSP derived pathways abundance with metagenomics when compared to null expectation datasets, with decreasing performance for more complex and/or less characterized environments (Douglas et al., 2019; Sun, Jones & Fodor, 2020). These low correlations could be due to preferential amplification of certain DNA sequences, primers biases, and varying gene copy numbers of 16S rRNA per taxa, although HSP methods typically try to correct this bias (Bowman & Ducklow, 2015; Wemheuer et al., 2018; Douglas et al., 2019). The increased detection sensitivity of 16S rRNA metabarcoding can also create a bias in the number of contributing taxa to certain functions, which can negatively affect correlations with functions derived from metagenomics. Considering that amplicon-based 16S rRNA metabarcoding and metagenomics uncovered a substantially different bacterial diversity, the weak correlation in functional abundance between the two methods is expected.
An alternative and possibly more appropriate approach to comparing functional profiles of HSP and metagenomics approaches is by contrasting their correlation with metadata (Sun et al., 2020) or by evaluating the correspondence between the ordination of their functional communities. Using procrustes analyses, there was a very strong and significant correlation between HSP methods and metagenomics, especially when using presence/absence data. The ASV communities showed no significant relationship with the functional profiles, suggesting that the taxonomic and functional communities were influenced differently by the biological and/or environmental conditions. We also tested the correspondence of the bacterial taxonomic and functional profiles with macrofaunal communities and physico-chemical data. While both profiles correlated relatively well to physico-chemical data (r ≥ 0.5), with the EC-based data performing best, only functional profiles correlated strongly and significantly with the macrofaunal communities. A higher association of functional versus taxonomic beta-diversity with macrofaunal data was also reported by Laroche et al. (2018), which suggest that interactions between these communities are especially driven by microbial metabolic capabilities rather than specific phylogenetic association.
The sensitivity of the different datasets in detecting the effect of fish farm activities was evaluated by comparing changes in community composition between near-field (0 m from pen) and far-field samples (>= 1,200 m from pen). In general, we found higher sensitivity for the HSP methods, and especially for Paprica and Tax4Fun2, compared to metagenomics and ASV communities. These results indicate that despite the lower accuracy and increased detection sensitivity of HSP methods, they may be more accurate in assessing how microbial communities respond to environmental changes than metagenomics. This is likely enhanced when complex microbial communities are present, such as in marine sediments. The results improved when transforming functional abundance data, including those of metagenomics, to presence/absence data, as it reduced within group variability. In addition, we observed higher stochasticity of microbial taxonomic shifts in response to a contamination gradient compared to functional community changes. This observation has also been reported by Hornick & Buschmann (2018), Laroche et al. (2018) and Ren et al. (2016) and is likely due to several taxa sharing the same metabolic capabilities and their succession in the ecosystem may have less to do with environmental changes than with biological properties (e.g. growth cycle and bacterial interactions) and geo-topographic factors (e.g. depth and geographic distance). As such, functional profiles may be slightly more robust and sensitive in detecting environmental alterations caused by fish farm activities, although further research is needed to properly test this assumption.
Benthic environments under cage fin-fish aquaculture are usually enriched in organic waste and nutrients such as phosphorus and nitrogen compounds (from faeces and uneaten fish feed for example), which can lead to eutrophic conditions, microbial anaerobic activities and the production of ammonia and hydrogen sulfide gasses (Brooks & Mahnken, 2003; Buschmann et al., 2006; Valdemarsen, Kristensen & Holmer, 2009; Wang et al., 2012). In the present study, we were particularly interested in comparing the response of classes of pathways associated to the nitrogen and sulfur cycles between the pen and reference sites, and between the metagenomics and HSP-based data. Overall, results from both approaches were very similar, with an increase near the pens of pathways associated to nitrate reduction and glycosaminoglycan degradation, and a decrease of pathways affiliated to allantoin degradation and sulfur oxidation. Additionally, the Paprica analysis showed a decrease in pathway abundance associated with sulfur reduction, sulfite and sulfide reduction, and an increase of dimethylsulfide degradation near the pens. While we expected pathways of nitrate-reduction to be in higher abundance near the fish farms, due to enriched nutrients and possibly anoxic conditions, it was somewhat surprising that pathways associated with sulfur compounds were less abundant in both the metagenomics and HSP datasets. However, sulfite respiration was found to be associated with the pens when using taxonomic information of the 16S rRNA data and literature-based functional association (Faprotax methodology). It is likely that certain pathways associated to the sulfur cycle, such as sulfite oxidation and reduction, were indeed more prominent near the pens but were not fully picked-up by metagenomics and HSP-based functional profiling, possibly due to the incompleteness of reference databases. For example, pathways associated to sulfite respiration were absent from both the Humann2 and Picrust2 datasets. Glycosaminoglycan degradation is responsible for the degradation of long linear polysaccharides made of repeating disaccharide units, also referred to as mucopolysaccharides (Ernst et al., 1995). It is probable that high quantities of mucopolysaccharides originate from mucus produced and excreted by the caged salmon (see Reverter et al., 2018; Jacobsen et al., 2019) and this is being catabolized by a specialized group of bacteria. Allantoin represents a product of uric acid, an important metabolic intermediate compound produced by both animals and bacteria. Under limited nutrient conditions, allantoin can be degraded into ammonia by some bacteria, to serve as a secondary source of nitrogen (Switzer et al., 2020). Correspondingly, our results suggest that pathways associated to allantoin degradation are less abundant near the pens, where strong nutrient enrichment occurs. Overall, these results show congruence between metagenomics and HSP methods for the classes of pathways of particular interest, and highlight both the potential and caveats of the current functional profiling methods in providing further understanding of the metabolic and environmental changes occurring in benthic ecosystems. Further research involving more samples and taking into account regional and temporal variability is needed to confidently identify potential functional indicators of fish farm ecological impacts that could be eventually integrated into benthic health indexes.
Collectively our results suggest that the lower specificity of HSP methods may be offset by the ability of amplicon-based metabarcoding to provide a much more exhaustive assessment of the taxonomic community, and hence of functions that have low genomic variability. This allows HSP methods to provide functional profiles that are relatively similar to those of metagenomics, and which respond similarly to environmental changes. Although the accuracy and sensitivity of HSP methods are still strongly affected by the incompleteness of reference databases, our results demonstrate that they provide a useful functional profiling alternative to metagenomics, and a valuable tool in detecting and evaluating the effects of salmon farming on benthic ecosystems.