RESEARCH Open Access Disentangling the mechanisms shaping the surface ocean microbiota Ramiro Logares 1,2* , Ina M. Deutschmann 1 , Pedro C. Junger 3 , Caterina R. Giner 1,4 , Anders K. Krabberød 2 , Thomas S. B. Schmidt 5 , Laura Rubinat-Ripoll 6 , Mireia Mestre 1,7,8 , Guillem Salazar 1,9 , Clara Ruiz-González 1 , Marta Sebastián 1,10 , Colomban de Vargas 6 , Silvia G. Acinas 1 , Carlos M. Duarte 11 , Josep M. Gasol 1,12 and Ramon Massana 1 Abstract Background: The ocean microbiota modulates global biogeochemical cycles and changes in its configuration may have large-scale consequences. Yet, the underlying ecological mechanisms structuring it are unclear. Here, we investigate how fundamental ecological mechanisms (selection,dispersal and ecological drift) shape the smallest members of the tropical and subtropical surface-ocean microbiota: prokaryotes and minute eukaryotes (picoeukaryotes). Furthermore, we investigate the agents exerting abiotic selection on this assemblage as well as the spatial patterns emerging from the action of ecological mechanisms. To explore this, we analysed the composition of surface-ocean prokaryotic and picoeukaryotic communities using DNA-sequence data (16Sand 18S-rRNA genes) collected during the circumglobal expeditions Malaspina-2010 and TARA-Oceans. Results: We found that the two main components of the tropical and subtropical surface-ocean microbiota, prokaryotes and picoeukaryotes, appear to be structured by different ecological mechanisms. Picoeukaryotic communities were predominantly structured by dispersal-limitation, while prokaryotic counterparts appeared to be shaped by the combined action of dispersal-limitation, selection and drift. Temperature-driven selection appeared as a major factor, out of a few selected factors, influencing species co-occurrence networks in prokaryotes but not in picoeukaryotes, indicating that association patterns may contribute to understand ocean microbiota structure and response to selection. Other measured abiotic variables seemed to have limited selective effects on community structure in the tropical and subtropical ocean. Picoeukaryotes displayed a higher spatial differentiation between communities and a higher distance decay when compared to prokaryotes, consistent with a scenario of higher dispersal limitation in the former after considering environmental heterogeneity. Lastly, random dynamics or drift seemed to have a more important role in structuring prokaryotic communities than picoeukaryotic counterparts. (Continued on next page) © The Author(s). 2020 Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data. * Correspondence:
[email protected] 1 Institute of Marine Sciences (ICM), CSIC, 08003 Barcelona, Catalonia, Spain 2 Department of Biosciences, Section for Genetics and Evolutionary Biology, University of Oslo, 0316 Oslo, Norway Full list of author information is available at the end of the article Logares et al. Microbiome (2020) 8:55 https://doi.org/10.1186/s40168-020-00827-8
(Continued from previous page) Conclusions: The differential action of ecological mechanisms seems to cause contrasting biogeography, in the tropical and subtropical ocean, among the smallest surface plankton, prokaryotes and picoeukaryotes. This suggests that the idiosyncrasy of the main constituents of the ocean microbiota should be considered in order to understand its current and future configuration, which is especially relevant in a context of global change, where the reaction of surface ocean plankton to temperature increase is still unclear. Keywords: Ocean, Plankton, Microbiota, Picoeukaryotes, Prokaryotes, Community structure, Ecological processes, Selection, Dispersal, Drift Background The surface ocean microbiota is a pivotal underpinning of global biogeochemical cycles [1,2]. The smallest ocean microbes, the picoplankton, have a key role in the global carbon cycle, being responsible for an important fraction of the total atmospheric carbon and nitrogen fixation in the ocean [3–5], which supports ≈46% of the global primary productivity [6]. Oceanic picoplankton plays a fundamental role in processing organic matter by recycling nutrients and carbon to support additional production as well as by channelling organic carbon to upper trophic levels through food webs [5,7,8]. The ocean picoplankton includes prokaryotes (both bacteria and archaea) and tiny unicellular eukaryotes (hereafter picoeukaryotes), which feature fundamental differences in terms of cellular structure, feeding habits, metabolic diversity, growth rates and behaviour [9]. Even though marine picoeukaryotes and prokaryotes are usually investigated separately, they are intimately connected through biogeochemical and food web networks [10–12]. The underlying ecological mechanisms determining the biogeography of prokaryotes and picoeukaryotes in the global ocean are unclear [13,14]. In particular, we do not know whether these crucial components of the ocean microbiota are structured by the action of the same or different ecological processes. Comprehending such processes is fundamental, as their differential action can produce changes in the ocean microbiota composition that could impact global ecosystem function [15– 17]. A recent ecological synthesis explains the structure of communities and the emergence of biogeography as a consequence of the action of four main processes: selection,dispersal,ecological drift and speciation [18]. Selection involves deterministic reproductive differences among individuals from different or the same species as a response to biotic or abiotic conditions. Selection can act in two opposite directions; it can constrain (homogeneous selection) or promote (heterogeneous selection) the divergence of communities [19]. Dispersal is the movement of organisms across space, and rates can be high (homogenising dispersal), moderate, or low (dispersal limitation)[19]. Dispersal limitation occurs when species are absent from suitable habitats because potential colonizers are too far away [20], and the significance of dispersal limitation increases as geographic scale increases [21]. Ecological drift (hereafter drift) in a local community refers to random changes in species’relative abundances derived from stochastic birth, death, offspring production, immigration and emigration [18]. The action of drift in a metacommunity, that is, local communities that are connected via dispersal of multiple species [22], may lead to neutral dynamics [21], where random dispersal is the main mechanism of community assembly. Finally, speciation is the evolution of new species [18], and it will not be considered hereafter as it is expected to have a small impact in the turnover of communities that are connected via dispersal [23], being also difficult to measure this ecological process in the wild. The action of the previous ecological processes is typically manifested as different taxonomic or phylogenetic patterns of community turnover, that is, β-diversity. At the moment, there are several estimators of β-diversity which capture different aspects of community turnover [24]. Most of these indices consider taxonomic or phylogenetic aspects of communities, but not speciesassociation patterns, which can also manifest the action of ecological processes. For example, selection exerted by an environmental variable can drive species cooccurrences generating groups of highly associated species or modules in association networks that correspond with specific environmental conditions [25]. Different members of these modules may be more abundant in specific regions of the ocean, contributing to increase βdiversity estimates between these regions when based on standard compositional or phylogenetic β-diversity metrics. Yet, β-diversity estimates based on associationaware metrics may point to higher similarity between these regions, as taxa belong to the same modules. Furthermore, modules may display correlations with environmental heterogeneity. Thus, association aware metrics of β-diversity may allow unveiling community patterns and their relationships with environmental variables (i.e. selection), which would be missed by standard approaches [26]. So far, most studies investigating the structure of the ocean microbiota have not considered species associations in their analyses of β-diversity. Logares et al. Microbiome (2020) 8:55 Page 2 of 17
The differential action of selection, dispersal and drift may generate different microbial assemblages that could feature diverse metabolisms and ecologies [16,17]. Moderate or high selection together with moderate dispersal rates may couple environmental heterogeneity with combinations of species, leading to a spatial pattern known as species sorting [27]. In contrast, high or low levels of dispersal may decouple environmental heterogeneity (i.e. selection) from the composition of species assemblages. High dispersal rates may maintain populations in habitats to which they are maladapted [16,22]. Inversely, low dispersal rates may promote microbial assemblages that become more different as the geographic distance between them increases (distance decay). If environmental heterogeneity and geographic distance covary, then distance decay could reflect both selection and dispersal limitation [28]. Drift is expected to cause important random effects in local community composition in cases where selection is weak and populations are small [15,29]. Here, we investigate the mechanisms that shape the smallest members of the surface-ocean microbiota by using DNA-sequence data collected in two of the largest circumglobal oceanographic expeditions to date, Malaspina 2010 [30] and TARA Oceans [31]. Specifically, we ask: What is the relative importance of selection, dispersal and drift in structuring the sunlit ocean microbiota? Do these processes act similarly on main components of this microbiota (prokaryotes and picoeukaryotes)? What are the main agents that exert abiotic selection? Do species association networks reflect the action of selection in the upper ocean microbiota? What are the main spatial-structure patterns that emerge due to the action of selection, dispersal and drift? Results Quantifying the mechanisms that structure the surface ocean picoplankton We analysed 16S and 18S rRNA-genes from prokaryotes and picoeukaryotes in 120 globally distributed tropical and subtropical stations sampled during the Malaspina 2010 expedition [30](Fig.1a; Figure S1, Additional file 1). TARA Oceans data were not included in these analyses as the type of generated DNA fragments could not be used for phylogenetic reconstructions (see details in ‘Methods’section). Operational taxonomic units were delineated at 99% similarity (OTUs -99% ) and as unique sequence variants (OTUs -ASVs , the maximum resolution for the 18S and 16S rRNAgene). Analyses using both, OTUs -99% and OTUs -ASVs indicated that dispersal limitation was the dominant factor structuring picoeukaryotic communities, explaining ≈76– 67% of community turnover, while this process had a lower importance in prokaryotes (≈35–25%; Fig. 1b). Note that percentage refers to the percentage of pairs of communities that appear to be driven by dispersal limitation. In contrast, homogenising dispersal had a very limited role in the structuring of the tropical and subtropical upper-ocean microbiota (< 3% for both picoeukaryotes and prokaryotes). Drift had a limited role in the structuring of picoeukaryotic communities as indicated by both OTUs -99% and OTUs -ASVs , representing ≈21–6% of community turnover (Fig. 1b). In contrast, drift appeared as a relevant factor structuring prokaryotic communities, explaining ≈44–31% of the community turnover according to OTUs -99% and OTUs -ASVs (Fig. 1b). The role of selection was higher in prokaryotes compared to picoeukaryotes according to both OTUs -99% and OTUs -ASVs , explaining ≈34–27% of the turnover of prokaryotic communities, and ≈17–11% of that in picoeukaryotes (Fig. 1b). Heterogeneous selection had a relatively higher importance in structuring picoeukaryotes as compared to prokaryotes (≈16–7% vs. ≈9–4%, respectively). Instead, homogeneous selection appeared more important in structuring prokaryotic (≈24–23%) than picoeukaryotic (≈ 1–4%) communities (Fig. 1b). Our quantifications indicated different roles of ecological processes in structuring communities of marine prokaryotes and picoeukaryotes populating the tropical and subtropical surface-ocean (Fig. 1b). We then aimed at confirming these results using other more traditional approaches. In these analyses, considering Malaspina data, we used OTUs -99% , given that these likely correspond to well-defined lineages, while OTUs -ASVs may reflect, in some cases, intraspecific variation [32]. We found moderate correlations between picoeukaryotic and prokaryotic β-diversity (Bray-Curtis: ρ= 0.58, gUniFrac: ρ= 0.61, p= 0.01, Mantel tests; Figure S2, Additional file 2). Given that rare species tend to occupy less sites than more abundant ones [33], communities featuring different proportions of abundant or rare species may display different spatial turnover. We found that picoeukaryotes had proportionally more regionally rare (i.e. mean abundances across all samples < 0.001%) species than prokaryotes (71% vs. 48% respectively) (Table S1, Additional file 3). This is consistent with the observation that picoeukaryotes had more restricted species distributions (i.e. occurring in < 20% of the stations) than prokaryotes (95% vs. 88% of the species respectively) (Figure S3, Additional file 4, Table S2, Additional file 5). Selection acting on the microbiota We investigated the agents exerting abiotic selection on the tropical and subtropical surface-ocean microbiota by analysing β-diversity together with the environmental variables included in the Meta-119 Malaspina dataset (temperature (°C), conductivity (S m −1 ), fluorescence, salinity and dissolved oxygen (mL L −1 )). We used different indices that capture distinct facets of β-diversity (Bray-Curtis, TINA w ,PINA w , gUniFrac; see ‘Methods’ Logares et al. Microbiome (2020) 8:55 Page 3 of 17
section). Water temperature was the most important driver of selection on prokaryotes (Fig. 2), ranging between 15.7 and 29.3 °C, with a mean of 24.5 °C and a standard deviation of 3.2 °C across the whole Meta-119 Malaspina dataset (Fig. 1a). Furthermore, water temperature appeared to affect prokaryotic association networks, given that TINA w [26] explained ≈50% of community variance (ADONIS R 2 ) (Fig. 2), while other used β-diversity indices that do not consider species associations explained considerably lower proportions (Fig. 2). In contrast, temperature had limited effects on picoeukaryotic community turnover (Fig. 2). Analyses using both the Malaspina and TARA Oceans datasets indicated stronger positive correlations between TINA w and watertemperature differences in prokaryotes (Mantel r=0.8–0.5, p< 0.01) than in picoeukaryotes [Mantel r=0.3,p<0.05] (Fig. 3). In particular, TARA Oceans samples displayed a higher correlation with water temperature than Malaspina samples (Fig. 3). Overall, TINA w results indicate that locations with similar temperatures include Fig. 1 Ecological mechanisms shaping the tropical and subtropical surface-ocean picoplankton. aPosition of the 120 stations included in this work that were sampled as part of the Malaspina-2010 expedition (green dots) in the tropical and subtropical ocean. A snapshot of the global sea surface temperature, a main environmental driver affecting microbial distributions, is shown as a general representation of the temperature gradients in the surface ocean (as inferred using the ‘optimum interpolation sea surface temperature’dataset from the NOAA corresponding to the 17 of March of 2018). Note that temperatures measured in situ were used in all analyses, not the ones displayed here. bPercentage of the community turnover associated to different ecological processes in prokaryotes and picoeukaryotes in the tropical and subtropical upper ocean as calculated using OTUs -99% and OTUs -ASVs . Note that percentage refers to the percentage of pairs of communities that appear to be driven by a given process Logares et al. Microbiome (2020) 8:55 Page 4 of 17
prokaryotic species that tend to co-occur, with this pattern disappearing as the temperature difference between stations increases. The previous pattern was either weak or non-existent in microbial eukaryotes (Fig. 3). We expanded the exploration of the role of abiotic selection on microbiota structuring by analysing a larger number of environmental variables (total 17) that were available for only 57 globally distributed Malaspina stations (see details in Supplementary Methods, Additional file 6; Figure S4, Additional file 7). Results supported the importance of temperature-driven selection for prokaryotic community structuring (Figure S5, Additional file 8)andindicatedthat fluorescence (a proxy for Chlorophyll aconcentration) explained 31% of PINA w -based prokaryotic community variance (ADONIS R 2 ), being non-significant for picoeukaryotes (Figure S5, Additional file 8). The remaining tested abiotic variables explained a minor fraction of community variance, suggesting that abiotic selection, at the whole oceanmicrobiota level, operates via few agents, mainly temperature, although we cannot rule out that other unmeasured abiotic variables may also be exerting selection. The different correlations between temperature and βdiversity as measured by TINA w in prokaryotes and picoeukaryotes suggest that they may feature different species association networks. We found that prokaryotes sampled in both Malaspina and TARA Oceans were more associated between themselves than protists (Figure S6, Additional file 9; Table S3, Additional file 10; Table S4, Additional file 11; Table S5, Additional file 12). Furthermore, the prokaryotic networks were more modular (in terms of cliques) than the picoeukaryotic counterparts (Table S3, Additional file 10), which may reflect to certain extent, temperature-driven selection [25]. Given that selection exerted by variables that lack phylogenetic signal, typically biotic variables, could inflate estimates of dispersal limitation, we have checked whether the high dispersal limitation we estimated for picoeukaryotes could reflect zooplankton grazing. For that, we have analysed globally distributed surface TARA Oceans stations for which we could estimate both the community composition of picoeukaryotes (here defined as the 0.8–5μm size-fraction; 36 or 38 stations) as well as that of microzooplankton (20–180 μm size-fraction; 36 stations) or mesozooplankton (180–2,000 μm sizefraction; 38 stations) based on 18S-rRNA genes [34]. Analyses considering abiotic (total 6, see Supplementary Methods, Additional file 6) and biotic (estimated zooplankton abundance) variables indicated that microand mesozooplankton had a minor influence on picoeukaryotic community structure (≈5% of the variance explained, ADONIS R 2 ). In addition, the correlation between picoeukaryotic and zooplankton β-diversity was either weak (microzooplankton, ρ= 0.34) or absent (mesozooplankton) [p< 0.01, Mantel tests]. Thus, zooplankton grazing does not appear to influence βdiversity in picoeukaryotes. Fig. 2 Main variables influencing the structure of the surface-ocean microbiota as captured by different β-diversity metrics. Percentage of variance in picoeukaryotic and prokaryotic community composition (ADONIS R 2 ) explained by water temperature and Longhurst Provinces when using different β-diversity metrics. Figure based on the Malaspina Meta-119 dataset (see ‘Methods’section). TINA w TINA weighted, gUniFrac generalized Unifrac, PINAw PINA weighted, N.S. non-significant. Note that TINA w , which considers species association networks, captures a significantly higher proportion of community variance associated to temperature than Bray-Curtis, a compositional index, in prokaryotes Logares et al. Microbiome (2020) 8:55 Page 5 of 17
Selection acting on single species The previous analyses investigated how selection may operate on the entire assemblage of species, without considering the different responses to selection that are expected in individual species. We therefore evaluated the potential action of selection on single species by determining their individual correlations with multiple abiotic environmental variables using the maximal information coefficient (MIC). In the Malaspina dataset (Fig. 1a), temperature was the variable with the highest number of associated prokaryotic species (1.7%), representing ≈17% of the 16S rRNA genesequence abundance, while picoeukaryotic species displayed limited associations with temperature (≈0.3% of the species representing ≈5% of the 18S rRNA gene-sequence Fig. 3 Temperature-driven selection seems to affect species association networks in prokaryotes but not in pico-/nano-eukaryotes. Differences in community composition (as 1-[TINA-weighted] = TINA w dissimilarities) vs. temperature differences (as Euclidean distances based on dimensionless z-scores) for both small unicellular eukaryotes and prokaryotes sampled during the Malaspina and TARA Oceans expeditions. Note that, in contrast to other indices, TINA w considers species-association patterns (i.e. co-occurrences and co-exclusions ) when estimating β-diversity [26]. NB: While only picoeukaryotes were included in Malaspina (cell sizes < 3 μm), TARA Oceans data included picoand nano-eukaryotes (cell sizes < 5 μm). Picoand nanoeukaryotes from both expeditions (left panels) displayed low or no correlations between TINA w distances and temperature differences (Mantel test results included in the panels). On the contrary, prokaryotes (right panels) displayed high to moderate correlations between TINA w distances and temperature differences. These differences in the correlations are likely due to the wider temperature ranges covered by TARA Oceans compared to Malaspina (see Discussion).The regression line is shown in red (Malaspina microbial eukaryotes N.S., Malaspina Prokaryotes R 2 = 0.3, TARA Oceans microbial eukaryotes R 2 = 0.1, TARA Oceans Prokaryotes R 2 = 0.7; p< 0.05). The maps at the bottom indicate the surface stations from the expeditions Malaspina (119 stations for both prokaryotes and picoeukaryotes) and TARA Oceans (63 stations for prokaryotes and 40 stations for small unicellular eukaryotes) that were used to calculate TINA w Logares et al. Microbiome (2020) 8:55 Page 6 of 17
abundance) (Figure S7, Additional file 13). Picoeukaryotic and prokaryotic species were also associated with oxygen, conductivity and salinity (Figure S7, Additional file 13), which covary with temperature. The remaining variables displayed limited associations with individual prokaryotic or picoeukaryotic species (Figure S7, Additional file 13), thus agreeing with our previous results suggesting that abiotic selection on the tropical and subtropical surface-ocean microbiota operates via few variables, with a dominant role for temperature among prokaryotes. Overall, prokaryotes featured proportionally more individual-species associations with environmental parameters than picoeukaryotes (Figure S7, Additional file 13), suggesting that environmental heterogeneity in the tropical and subtropical surface-ocean has a stronger effect on prokaryotic assemblages than on picoeukaryotic counterparts. Analyses of TARA Oceans data supported this by indicating that prokaryotic species were associated predominantly with temperature and oxygen in the upper global ocean, while unicellular eukaryotes had weak associations to multiple variables (Table S6, Additional file 14). Dispersal Abiotic environmental conditions in adjacent stations over the trajectory of the Malaspina cruise, typically separated by 250–500 km, in the tropical and subtropical ocean (Fig. 1a) are generally comparable [35]. Therefore, compositional differences between pairs of neighbouring communities could manifest the differential capability of distinct microbial assemblages to disperse. Following these premises, we analysed the change in picoeukaryotic and prokaryotic community composition along the trajectory of the Malaspina cruise by comparing each community to the one sampled immediately before in a sequential manner (i.e. sequential βdiversity) (Fig. 4a–c). Both picoeukaryotic and prokaryotic communities displayed variable amounts of sequential βdiversity (Fig. 4a, b), although picoeukaryotes featured, on average, a higher sequential β-diversity than prokaryotes (Fig. 4c). This agrees with the overall mean β-diversity, which was significantly higher for picoeukaryotes than for prokaryotes (Figure S8, Additional file 15). Tests by subsampling the number of picoeukaryotic OTUs -99% to the Fig. 4 Picoeukaryotic communities display a higher spatial differentiation than prokaryotic counterparts in the tropical-subtropical surface-ocean. a–cSequential change in community composition across space (sequential β-diversity). Communities were sampled along the Malaspina expedition (a, b black arrows), and the composition of each community was compared against its immediate predecessor. In panels a,b, the size of each bubble represents the Bray-Curtis dissimilarity between a given community and the community sampled previously. Blue squares in panels a, b represent the stations where β-diversity displayed abrupt changes (Bray-Curtis values > 0.8 for picoeukaryotes and > 0.7 for prokaryotes). Abrupt changes coincided in a total of 11 out of 14 stations for both picoeukaryotes and prokaryotes, while one station displayed marked changes only for picoeukaryotes and two only for prokaryotes. Panel csummarizes the sequential Bray-Curtis values for prokaryotes and picoeukaryotes (Means were significantly different between domains [Wilcoxon text, p< 0.05]). Panel dindicates the differences in distance-decay between prokaryotes and picoeukaryotes in the tropical and subtropical surface-ocean. Mantel correlograms between geographic distance and βdiversity featuring distance classes of 1000 km for both picoeukaryotes and prokaryotes are shown. Coloured squares indicate statistically significant correlations (p< 0.05). Note that β-diversity in picoeukaryotes displayed positive correlations with increasing distances up to ≈3000 km, while prokaryotes had positive correlations with distances up to ≈2000 km. Correlations tended to be smaller in prokaryotes than in picoeukaryotes, indicating smaller distance decay in the former compared to the latter Logares et al. Microbiome (2020) 8:55 Page 7 of 17
same number of prokaryotic ones (7025) indicated that different numbers of OTUs -99% in these groups did not affect mean Bray-Curtis estimates of β-diversity displayed in Figure S8, Additional file 15 [36]. When geographic distance covaries with environmental heterogeneity, spatial community variance may be the manifestation of both selection and/or dispersal limitation. β-diversity in picoeukaryotes and prokaryotes displayed positive correlations with geographic distance (i.e. distance decay) predominantly within 1000 km (Fig. 4d). Yet, correlations were weaker in prokaryotes than in picoeukaryotes, pointing to stronger dispersal limitation or selection in the latter. Variance partitioning analyses considering both environmental [temperature (°C), conductivity (S m −1 ), fluorescence, salinity and dissolved oxygen (mL L −1 )] and geographic variables (ocean basin and subdivisions, as well as Longhurst biogeographic provinces [37], Figure S1, Additional file 1) indicated that in prokaryotes, geographic variables explained most of the variance (24%), while environmental variables explained 10%, and 13% was explained by both variables; 53% of the variance remained unexplained. In contrast, picoeukaryotes displayed non-significant results in the same analyses. Still, after controlling for the effects of the most important environmental variables, Longhurst provinces (but not ocean basins nor subdivisions) accounted for ≈20–25% of community variance in both picoeukaryotes and prokaryotes (ADONIS R 2 ) (Fig. 2). All in all, the previous analyses seem coherent with our quantifications of ecological processes (Fig. 1b), in the sense that they indicate that both selection and dispersal limitation (represented by geographic variables such as distance or ocean provinces), do seem to have a role in the structuring of the surface ocean picoplankton. Selection and dispersal limitation may operate more strongly in geographic areas that constitute ecological boundaries, leading to abrupt changes in microbiota composition. We identified 14 communities where sequential β-diversity displayed abrupt changes, with 11 of them coinciding for both picoeukaryotes and prokaryotes (Fig. 4a, b). The Local Contributions to Beta Diversity (LCBD) index [38] (Figure S9, Additional file 16) indicated that ≈22% of both picoeukaryotic and prokaryotic communities (26 stations each, totaling 36 different stations) contributed the most to the β-diversity, with 16 communities coinciding for both prokaryotes and picoeukaryotes (Figure S9, Additional file 16; Table S7, Additional file 17). In addition, eight of the 36 stations featuring a significant LCBD were also identified as zones of abrupt community change in sequential βdiversity analyses (Table S7, Additional file 17). These zones point to selection or dispersal operating simultaneously and strongly upon both prokaryotic and picoeukaryotic communities in the surface ocean. Discussion Applying an innovative ecological framework [23] allowed us to quantify the mechanisms that shape the tropical and subtropical upper-ocean microbiota. Yet, this approach has limitations (summarised by Zhou and Ning [19]) that need to be considered in the context of our results. First, our results represent the overall action of ecological processes at the whole microbiota level, and not their operation on every taxonomic group or lineage (for example, different taxonomic classes may be structured by different processes). In addition, our results reflect the action of ecological mechanisms at the global ocean level, and we expect that other spatial scales (ocean basin for example) may lead to other results. Furthermore, our results provide a snapshot of the importance of ecological processes at the global-ocean scale, and future studies should investigate how the relative importance of these mechanisms change over time [39]. Second, the measured ecological mechanisms are associated with the evolutionary diversification that is reflected by the variation in the chosen molecular markers. OTUs -99% and OTUs -ASVs based on the 16S and 18S rRNA genes likely reflect defined species (or gene flow units [40]) or in some cases population variation [32], and therefore, the measured ecological mechanisms in the tropical and subtropical ocean apply to those evolutionary levels. Hence, our results do not reflect the mechanisms shaping intra-population variation or those shaping taxonomic ranks above the species level. Furthermore, our results indicate that delineating OTUs based on sequence clustering (OTUs -99% ) or sequence variants (OTUs -ASVs ) can affect measurements of ecological mechanisms, although in our study, main trends were maintained. It could be hypothesized that OTUs -99% and OTUs -ASVs may represent different taxonomic units in prokaryotes or picoeukaryotes, especially if one group was evolving faster than the other. Yet, both prokaryotes and picoeukaryotes show a wide range of evolutionary rates [41,42], including lineages evolving slow or fast, therefore potential differences in unit definitions associated to different evolutionary rates will likely compensate when analysing complex assemblages of species. Third, failure to detect selection could inflate estimates of dispersal limitation. We consider that our estimates indicating substantial dispersal limitation in picoeukaryotes were not inflated, as picoeukaryotes displayed more restricted spatial distributions than prokaryotes and important biotic variables, such as potential zooplankton grazing, did not seem to affect the structure of picoeukaryotic assemblages. Furthermore, another study also suggests that dispersal limitation influences protist distributions in the global ocean [34]. Altogether, the used framework [23] can be considered as a guide that can provide important insights on the ecological Logares et al. Microbiome (2020) 8:55 Page 8 of 17
mechanisms structuring the global ocean microbiota, while more data (e.g. single nucleotide variants in genes or genomes) and experiments are necessary to understand such mechanisms in further detail. Our results indicated that the differential action of ecological processes may promote different biogeographic patterns in prokaryotic and picoeukaryotic assemblages in the upper global-ocean. This is consistent with other works using similar approaches to ours indicating that protistan and bacterial assemblages are shaped by different ecological processes [39,43–45]. In particular, selection, which is known to have an important role in structuring prokaryotic communities [27,28], explained a higher proportion of community turnover in surface-ocean prokaryotes (≈34–27% of the turnover) than in picoeukaryotes (≈17–11%). This modest role of selection in structuring the tropical and subtropical sunlit-ocean microbiota is consistent with the moderate environmental gradients characterizing this habitat. In other habitats featuring a higher selective pressure, the role of selection in structuring microbiotas was, as expected, higher [43]. The quantifications of the importance of selection are also associated to the global scale of our survey. Thus, for example, at smaller geographic scales, where dispersal limitation is expected to have a lower impact than at global scales [20], the relative importance of selection could increase. Congruently, in surface waters of the East China Sea, it was found that selection was ~ 40% more important than dispersal limitation in structuring bacterial communities [44], while in our global study, selection and dispersal limitation had a similar importance in structuring prokaryotes. Furthermore, the previous study [44] found that selection was considerably more important than dispersal limitation in structuring communities of microbial eukaryotes. In contrast, our global assessment yields dispersal limitation to be ≈5 times more important than selection in structuring picoeukaryotic communities. We found that heterogeneous selection was more important in structuring picoeukaryotic than prokaryotic communities, while homogeneous selection was more important in structuring prokaryotic than picoeukaryotic communities. This suggests that prokaryotes and picoeukaryotes respond differently to the same environmental heterogeneity, which in the tropical and subtropical surface-ocean would be preventing community divergence in prokaryotes while promoting it in picoeukaryotes. Different adaptations in prokaryotes and picoeukaryotes [9] may determine such contrasting responses to the same environmental heterogeneity. For example, a given environmental heterogeneity could select for a few species featuring wide environmental tolerance or several species that are adapted to narrow environmental conditions. Several studies have indicated that water temperature is one of the main abiotic variables affecting the structure and diversity of the ocean microbiota [46–52]. Furthermore, temperature is known to structure microbial assemblages in seasonal time-series, pointing also to the importance of this variable at local scales over yearly cycles [53–55]. In our study, the higher correlation between TARA Oceans communities with temperature as compared to Malaspina (Fig. 3) is coherent with the importance of this variable, as TARA Oceans sampled a wider temperature range (range ≈ 0–30 °C, mean ≈21 °C, SD ≈7°C)thanMalaspina (range ≈15–30 °C, mean ≈24 °C, SD ≈3°C).Furthermore,and consistent with our results, recent global-scale studies reported strong correlations between ocean-microbiota composition (predominantly prokaryotic) and temperature, and weak correlations with nutrients [56,57]. In sum, the previous agrees with our results indicating that temperature is one of the most important agents exerting abiotic selection on the surface-ocean microbiota, although we cannot rule out the selective action of other unmeasured abiotic factors. Our analyses also unveiled an additional layer of information by indicating that temperature-driven selection affects prokaryotic taxa co-occurrences, a pattern not observed in picoeukaryotes. Such β-diversity related to species associations is typically not captured by classic compositional indices like Bray-Curtis, possibly due to variations in the relative abundance of the co-occurring species [58]. In contrast to prokaryotes, less is known about the effects of temperature on the community structure of ocean picoeukaryotes, which according to our results are modest. Yet, specific picoeukaryotic lineages, such as MAST-4, do seem to be affected by temperature [59], pointing to taxonomic-group specific responses to selection. One of the possible reasons why picoeukaryotes do not show co-occurrence patterns comparable to those observed in prokaryotes is dispersal limitation, which precludes picoeukaryotic species with similar niches to share the same geographic zone. Overall, our work indicates that species association patterns are informative on the β-diversity of marine prokaryotes, therefore taxa association networks should be contemplated in future analyses of the ocean microbiota. To what extent dispersal limitation affects the distribution of ocean microbes is a matter of debate. The impact of dispersal limitation is expected to increase with increasing body size [60]; therefore, larger protists are expected to be more limited by dispersal than smaller prokaryotes. Ocean protists seem to follow the previous tenet, as it has been observed that dispersal limitation appears to increase with increasing cell size [34]. Furthermore, in surface open-ocean waters, prokaryotes typically display abundances of 10 6 cells/mL, while picoeukaryotes normally have abundances of 10 3 cells/ mL [61]. Due to random dispersal alone, the more Logares et al. Microbiome (2020) 8:55 Page 9 of 17
53. Giner CR, Balague V, Krabberod AK, Ferrera I, Rene A, Garces E, Gasol JM, Logares R, Massana R. Quantifying long-term recurrence in planktonic microbial eukaryotes. Mol Ecol. 2019;28(5):923–35. 54. Lambert S, Tragin M, Lozano J-C, Ghiglione J-F, Vaulot D, Bouget F-Y, Galand PE. Rhythmicity of coastal marine picoeukaryotes, bacteria and archaea despite irregular environmental perturbations. ISME J. 2019;13(2): 388–401. 55. Bunse C, Pinhassi J. Marine bacterioplankton seasonal succession dynamics. Trends Microbiol. 2017;25(6):494–505. 56. Sunagawa S, Coelho LP, Chaffron S, Kultima JR, Labadie K, Salazar G, Djahanschiri B, Zeller G, Mende DR, Alberti A, et al. Structure and function of the global ocean microbiome. Science. 2015;348(6237):1261359. 57. Salazar G, Paoli L, Alberti A, Huerta-Cepas J, Ruscheweyh HJ, Cuenca M, Field CM, Coelho LP, Cruaud C, Engelen S, et al. Gene expression changes and community turnover differentially shape the global ocean metatranscriptome. Cell. 2019;179(5):1068–83 e1021. 58. Chase JM. Community assembly: when should history matter? Oecologia. 2003;136(4):489–98. 59. Rodriguez-Martinez R, Rocap G, Salazar G, Massana R. Biogeography of the uncultured marine picoeukaryote MAST-4: temperature-driven distribution patterns. ISME J. 2013;7(8):1531–43. 60. De Bie T, De Meester L, Brendonck L, Martens K, Goddeeris B, Ercken D, Hampel H, Denys L, Vanhecke L, Van der Gucht K, et al. Body size and dispersal mode as key traits determining metacommunity structure of aquatic organisms. Ecology Letters. 2012;15(7):740–7. 61. Kirchman DL. Microbial Ecology of the Oceans. Hoboken, New Jersey: John Wiley & Sons; 2008. 62. Foissner W. Biogeography and dispersal of micro-organisms: a review emphasizing protists. Acta Protozoologica. 2006;45:111–36. 63. Casteleyn G, Leliaert F, Backeljau T, Debeer AE, Kotaki Y, Rhodes L, Lundholm N, Sabbe K, Vyverman W. Limits to gene flow in a cosmopolitan marine planktonic diatom. Proc Natl Acad Sci U S A. 2010;107(29):12952–7. 64. Cermeno P, Falkowski PG. Controls on diatom biogeography in the ocean. Science. 2009;325(5947):1539–41. 65. Whittaker KA, Rynearson TA. Evidence for environmental and ecological selection in a microbe with no geographic limits to gene flow. Proc Natl Acad Sci U S A. 2017;114(10):2651–6. 66. Bass D, Richards TA, Matthai L, Marsh V, Cavalier-Smith T. DNA evidence for global dispersal and probable endemicity of protozoa. BMC Evol Biol. 2007; 7(1):162. 67. Lewis J, Harris ASD, Jones KJ, Edmonds RL. Long-term survival of marine planktonic diatoms and dinoflagellates in stored sediment samples. J Plankton Res. 1999;21(2):343–54. 68. Billard C, Inouye I. What is new in coccolithophore biology? In: Thierstein HR, Young JR, editors. Coccolithophores: From Molecular Processes to Global Impact. Berlin, Heidelberg: Springer Berlin Heidelberg; 2004. p. 1–29. 69. Milici M, Tomasch J, Wos-Oxley ML, Decelle J, Jauregui R, Wang H, Deng ZL, Plumeier I, Giebel HA, Badewien TH, et al. Bacterioplankton biogeography of the Atlantic Ocean: a case study of the distance-decay relationship. Front Microbiol. 2016;7:590. 70. Sintes E, De Corte D, Ouillon N, Herndl GJ. Macroecological patterns of archaeal ammonia oxidizers in the Atlantic Ocean. Mol Ecol. 2015;24(19):4931–42. 71. Louca S, Parfrey LW, Doebeli M. Decoupling function and taxonomy in the global ocean microbiome. Science. 2016;353(6305):1272–7. 72. Jones SE, Lennon JT. Dormancy contributes to the maintenance of microbial diversity. Proc Natl Acad Sci U S A. 2010;107(13):5881–6. 73. Locey KJ. Synthesizing traditional biogeography with microbial ecology: the importance of dormancy. J Biogeography. 2010;37(10):1835–41. 74. Louca S, Polz MF, Mazel F, Albright MBN, Huber JA, O’Connor MI, Ackermann M, Hahn AS, Srivastava DS, Crowe SA, et al. Function and functional redundancy in microbial systems. Nat Ecol Evol. 2018;2:936–43. 75. Östman Ö, Drakare S, Kritzberg ES, Langenheder S, Logue JB, Lindström ES. Regional invariance among microbial communities. Ecology letters. 2010; 13(1):118–27. 76. Salazar G, Cornejo-Castillo FM, Benitez-Barrios V, Fraile-Nuez E, AlvarezSalgado XA, Duarte CM, Gasol JM, Acinas SG. Global diversity and biogeography of deep-sea pelagic prokaryotes. ISME J. 2016;10(3):596–608. 77. Zinger L, Boetius A, Ramette A. Bacterial taxa–area and distance–decay relationships in marine environments. Mol Ecol. 2014;23(4):954–64. 78. Díez B, Massana R, Estrada M, Pedrós-Alió C. Distribution of eukaryotic picoplankton assemblages across hydrographic fronts in the Southern Ocean, studied by denaturing gradient gel electrophoresis. Limnol Oceanography. 2004;49(4):1022–34. 79. Flaviani F, Schroeder D, Lebret K, Balestreri C, Schroeder J, Moore K, Paszkiewicz K, Pfaff M, Rybicki E. Distinct oceanic microbiomes (from viruses to protists) found either side of the Antarctic Polar Front. Front Microbiol. 2018;9. 80. Grasshoff K, Ehrhardt M, Kremling K. Methods of seawater analysis. Weinheim: Verlag Chemie; 1983. 81. Estrada M, Delgado M, Blasco D, Latasa M, Cabello AM, Benitez-Barrios V, Fraile-Nuez E, Mozetic P, Vidal M. Phytoplankton across tropical and subtropical regions of the Atlantic, Indian and Pacific Oceans. PLoS One. 2016;11(3):e0151699. 82. Boyer TP, Antonov JI, Baranova OK, Coleman C, Garcia HE, Grodsky A, Johnson DR, Locarnini RA, Mishonov AV, O'Brien TD, et al. In: Levitus S, Mishonov A, editors. World Ocean Database 2013. In: NOAA Atlas NESDIS 72. Silver Spring, MD: NOAA; 2013. 83. Massana R, Murray AE, Preston CM, DeLong EF. Vertical distribution and phylogenetic characterization of marine planktonic Archaea in the Santa Barbara Channel. Appl Environ Microbiol. 1997;63(1):50–6. 84. Stoeck T, Bass D, Nebel M, Christen R, Jones MD, Breiner HW, Richards TA. Multiple marker parallel tag environmental DNA sequencing reveals a highly complex eukaryotic community in marine anoxic water. Mol Ecol. 2010;19(Suppl 1):21–31. 85. Parada AE, Needham DM, Fuhrman JA. Every base matters: assessing small subunit rRNA primers for marine microbiomes with mock communities, time series and global field samples. Environ Microbiol. 2016;18(5):1403–14. 86. Logares R. Workflow for Analysing MiSeq Amplicons based on Uparse v1.5. In.:https://doi.org/10.5281/zenodo.259579; 2017. 87. Nikolenko SI, Korobeynikov AI, Alekseyev MA. BayesHammer: Bayesian clustering for error correction in single-cell sequencing. BMC Genomics. 2013; 14 Suppl 1:S7. 88. Schirmer M, Ijaz UZ, D'Amore R, Hall N, Sloan WT, Quince C. Insight into biases and sequencing errors for amplicon sequencing with the Illumina MiSeq platform. Nucleic Acids Res. 2015;43(6):e37. 89. Zhang J, Kobert K, Flouri T, Stamatakis A. PEAR: a fast and accurate Illumina Paired-End reAd mergeR. Bioinformatics. 2014;30(5):614–20. 90. Edgar RC. Search and clustering orders of magnitude faster than BLAST. Bioinformatics. 2010;26(19):2460–1. 91. Edgar RC. UPARSE: highly accurate OTU sequences from microbial amplicon reads. Nat Methods. 2013;10(10):996–8. 92. Callahan BJ, McMurdie PJ, Rosen MJ, Han AW, Johnson AJ, Holmes SP. DADA2: High-resolution sample inference from Illumina amplicon data. Nat Methods. 2016;13(7):581–3. 93. Wang Q, Garrity GM, Tiedje JM, Cole JR. Naive Bayesian classifier for rapid assignment of rRNA sequences into the new bacterial taxonomy. Appl Environ Microbiol. 2007;73(16):5261–7. 94. Quast C, Pruesse E, Yilmaz P, Gerken J, Schweer T, Yarza P, Peplies J, Glockner FO. The SILVA ribosomal RNA gene database project: improved data processing and web-based tools. Nucleic Acids Res. 2013;41(Database issue):D590–6. 95. Altschul SF, Gish W, Miller W, Myers EW, Lipman DJ. Basic local alignment search tool. Journal of molecular biology. 1990;215(3):403–10. 96. Guillou L, Bachar D, Audic S, Bass D, Berney C, Bittner L, Boutte C, Burgaud G, de Vargas C, Decelle J, et al. The protist ribosomal reference database (PR2): a catalog of unicellular eukaryote small sub-unit rRNA sequences with curated taxonomy. Nucleic Acids Res. 2013;41(Database issue):D597–604. 97. Logares R, Sunagawa S, Salazar G, Cornejo-Castillo FM, Ferrera I, Sarmento H, Hingamp P, Ogata H, de Vargas C, Lima-Mendez G, et al. Metagenomic 16S rDNA Illumina tags are a powerful alternative to amplicon sequencing to explore diversity and structure of microbial communities. Environ Microbiol. 2014;16(9):2659–71. 98. Oksanen J, Kindt R, Legendre P, O'Hara B, Simpson GL, Solymos P, Stevens MHH, Wagner H. vegan: Community ecology package. R package version 1. 15-0. In.; 2008. 99. Logares R, Audic S, Bass D, Bittner L, Boutte C, Christen R, Claverie JM, Decelle J, Dolan JR, Dunthorn M, et al. Patterns of rare and abundant marine microbial eukaryotes. Curr Biol. 2014;24(8):813–21. 100. Schloss PD, Westcott SL, Ryabin T, Hall JR, Hartmann M, Hollister EB, Lesniewski RA, Oakley BB, Parks DH, Robinson CJ, et al. Introducing mothur: open-source, platform-independent, community-supported software for describing and comparing microbial communities. Appl Environ Microbiol. 2009;75(23):7537–41. Logares et al. Microbiome (2020) 8:55 Page 16 of 17
101. Capella-Gutierrez S, Silla-Martinez JM, Gabaldon T. trimAl: a tool for automated alignment trimming in large-scale phylogenetic analyses. Bioinformatics. 2009;25(15):1972–3. 102. Price MN, Dehal PS, Arkin AP. FastTree: computing large minimum evolution trees with profiles instead of a distance matrix. Mol Biol Evol. 2009;26(7): 1641–50. 103. R-Development-Core-Team. R: A language and environment for statistical computing. Vienna, Austria: R Foundation for Statistical Computing; 2008. 104. Paradis E, Claude J, Strimmer K. APE: Analyses of Phylogenetics and Evolution in R language. Bioinformatics. 2004;20(2):289–90. 105. Wickham H. ggplot2: Elegant Graphics for Data Analysis: Springer-Verlag; 2009. 106. Chen J, Bittinger K, Charlson ES, Hoffmann C, Lewis J, Wu GD, Collman RG, Bushman FD, Li H. Associating microbiome composition with environmental covariates using generalized UniFrac distances. Bioinformatics. 2012;28(16):2106–13. 107. Kembel SW, Cowan PD, Helmus MR, Cornwell WK, Morlon H, Ackerly DD, Blomberg SP, Webb CO. Picante: R tools for integrating phylogenies and ecology. Bioinformatics. 2010;26(11):1463–4. 108. Dray S, Blanchet G, Borcard D, Clappe S, Guenard G, Jombart T, Larocque G, Legendre P, Madi N, Wagner HH. adespatial: Multivariate multiscale spatial analysis. In.; 2017. 109. Cavender-Bares J, Kozak KH, Fine PV, Kembel SW. The merging of community ecology and phylogenetic biology. Ecology letters. 2009;12(7): 693–715. 110. Losos JB. Phylogenetic niche conservatism, phylogenetic signal and the relationship between phylogenetic relatedness and ecological similarity among species. Ecol Lett. 2008;11(10):995–1003. 111. Stegen JC, Lin X, Konopka AE, Fredrickson JK. Stochastic and deterministic assembly processes in subsurface microbial communities. ISME J. 2012;6(9): 1653–64. 112. Andersson AF, Riemann L, Bertilsson S. Pyrosequencing reveals contrasting seasonal dynamics of taxa within Baltic Sea bacterioplankton communities. ISME J. 2010;4(2):171–81. 113. Chase JM, Kraft NJB, Smith KG, Vellend M, Inouye BD. Using null models to disentangle variation in community dissimilarity from variation in α-diversity. Ecosphere. 2011;2(2):1–11. 114. Reshef DN, Reshef YA, Finucane HK, Grossman SR, McVean G, Turnbaugh PJ, Lander ES, Mitzenmacher M, Sabeti PC. Detecting novel associations in large data sets. Science. 2011;334(6062):1518–24. 115. Friedman J, Alm EJ. Inferring correlation networks from genomic survey data. PLoS Comput Biol. 2012;8(9):e1002687. 116. Watts SC, Ritchie SC, Inouye M, Holt KE. FastSpar: rapid and scalable correlation estimation for compositional data. Bioinformatics. 2018;35(6):1064-6. 117. Shannon P, Markiel A, Ozier O, Baliga NS, Wang JT, Ramage D, Amin N, Schwikowski B, Ideker T. Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res. 2003;13(11): 2498–504. 118. Csardi G, Nepusz T. The igraph software package for complex network research. InterJournal. 2006; Complex Systems:1695. Publisher’sNote Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. Logares et al. Microbiome (2020) 8:55 Page 17 of 17