Succession within the prokaryotic communities during the VAHINE mesocosms experiment in the New Caledonia lagoon

N2 fixation fuels ∼ 50 % of new primary production in the oligotrophic South Pacific Ocean. The VAHINE experiment has been designed to track the fate of diazotrophderived nitrogen (DDN) and carbon within a coastal lagoon ecosystem in a comprehensive way. For this, large-volume (∼ 50 m) mesocosms were deployed in the New Caledonian lagoon and were intentionally fertilized with dissolved inorganic phosphorus (DIP) to stimulate N2 fixation. This study examined the temporal dynamics of the prokaryotic community together with the evolution of biogeochemical parameters for 23 consecutive days in one of these mesocosms (M1) and in the Nouméa lagoon using MiSeq 16S rRNA gene sequencing and flow cytometry. Combining these methods allowed for inference of absolute cell numbers from 16S data. We observed clear successions within M1, some of which were not mirrored in the lagoon. The dominating classes in M1 were Alphaand Gammaproteobacteria, Cyanobacteria, eukaryotic microalgae, Marine Group II Euryarchaeota, Flavobacteriia, and Acidimicrobia. Enclosure led to significant changes in the M1 microbial community, probably initiated by the early decay of Synechococcus and diatoms. However, we did not detect a pronounced bottle effect with a copiotroph-dominated community. The fertilization with ∼ 0.8 μM DIP on day 4 did not have directly observable effects on the overall community within M1, as the data samples obtained from before and 4 days after fertilization clustered together, but likely influenced the development of individual populations later on, like Defluviicoccusrelated bacteria and UCYN-C-type diazotrophic cyanobacteria (Cyanothece). Growth of UCYN-C led to among the highest N2-fixation rates ever measured in this region and enhanced growth of nearly all abundant heterotrophic groups in M1. We further show that different Rhodobacteraceae were the most efficient heterotrophs in the investigated system and we observed niche partitioning within the SAR86 clade. Whereas the location inor outside the mesocosm had a significant effect on community composition, the temporal effect was significantly stronger and similar in both locations, suggesting that overarching abiotic factors were more influential than the enclosure. While temporal community changes were evident, prokaryotic diversity (Shannon index) only declined slightly from∼ 6.5 to 5.7 or 6.05 in the lagoon and M1, respectively, throughout the experiment, highlighting the importance of multiple and varying sources of organic matter maintaining competition.

some of which were not mirrored in the lagoon. The dominating classes in M1 were alpha-and gammaproteobacteria, cyanobacteria (mainly Synechococcus), eukaryotic microalgae, on days 10 and 14 Marine Group II euryarchaea, on days 12-23 also Flavobacteriia. Enclosure led to significant changes in the M1 microbial community, probably initiated by the early decay of Synechococcus and diatoms. However, we did 15 not detect a pronounced bottle effect with a copiotroph-dominated community. The fertilization with ∼ 0.8 µM DIP on day 4 did not have directly observable effects on the overall community within M1, as the data samples obtained from before and four days after fertilization clustered together, but likely influenced the development of individual populations later on, like Defluviicoccus-related bacteria and UCYN-C type dia- 20 zotrophic cyanobacteria. Growth of UCYN-C led to among the highest N 2 fixation rates ever measured in this region and enhanced growth of nearly all abundant heterotrophic groups in M1. We further show that different Rhodobacteraceae were the most efficient heterotrophs in the investigated system and we observed niche partitioning within the SAR86 clade. Whereas the location in-or outside the mesocosm had a significant effect

Introduction
The South West Pacific ocean is recognized as an area with one of the highest dini-5 trogen (N 2 ) fixation rates in the global ocean (Garcia et al., 2007;Luo et al., 2012). In this region, plankton taxa capable of N 2 fixation (N 2 -fixing or diazotrophic organisms) are very diverse (Moisander et al., 2010) and fuel up to 60 % of the primary production (Bonnet et al., 2015a), yet their interactions with surrounding planktonic communities are rarely studied. Within this vast oceanic region, most of the studies regarding the coupling between phytoplankton and heterotrophic bacteria have targeted the Nouméa lagoon, New Caledonia. Phytoplankton and bacterial production shows seasonal patterns, with maxima in December and January and an annual heterotrophic bacterial production representing 21 to 34 % of primary production (Torréton et al., 2010). High N 2 fixation rates are a recurrent feature in the Nouméa lagoon (Bie- 15 gala and Raimbault, 2008;Garcia et al., 2007), and diazotroph-derived N (DDN) was shown to be significantly channeled through the heterotrophic bacterial compartment in this environment (Berthelot et al., 2015;Bonnet et al., 2015a). Yet, no information is available on the potential effects of this N 2 fixation on bacterial successions.
The production and quality of varying sources of organic matter is expected to vary 20 as blooms of diazotrophic organisms are transient events. Such variation is expected to influence not only biogeochemical fluxes, but also the succession of different heterotrophic bacteria, according to their metabolic capabilities. Indeed, the importance of quality and quantity of dissolved organic matter (DOM) in structuring bacterioplankton communities has been well established (Alonso-Sáez and Gasol, 2007;Beier 1988;Murray et al., 2007;Nagata, 2000;Riemann et al., 2000). The VAHINE experiment (Bonnet et al., 2015b) provided a unique opportunity to study the fate of DDN in the marine planktonic food web and the interactions between phytoplankton and heterotrophic bacteria by simultaneously addressing biogeochemical parameters, stocks and fluxes, as well as biodiversity based on 16S rRNA gene 10 sequencing and flow cytometry cell counts.
After fertilization with dissolved inorganic phosphorus (DIP), ∼ 0.8 µmol L −1 on the evening of day 4, to alleviate any potential limitation often observed in the region (Moutin et al., 2007), two periods of about 10 days each (P1 and P2, see description in Sect. 2.1) were clearly identified in terms of diazotroph succession. During the 15 first period, diatom-diazotroph associations (DDAs) were dominating  the diazotrophic community with N 2 fixation rates decreasing from 17.9 ± 2.5 to 10.1 ± 1.3 nmol N L −1 d −1 (Bonnet et al., 2015a). Unicellular N 2 -fixing cyanobacteria of the UCYN-C type (Cyanothece-like) dominated the diazotroph community in the mesocosms during the second period, but did not appear in Nouméa lagoon 20 waters . N 2 -fixation rates increased during days 15-23 and reached > 60 nmol N L −1 d −1 , which are among the highest rates measured in marine waters (Bonnet et al., 2015a;Luo et al., 2012). These two periods corresponded to characteristic successions in phytoplankton taxa, chlorophyll stocks, primary and heterotrophic bacterial production, discussed in detail in companion studies in this issue and sampling strategy are described in (Bonnet et al., 2015b). The three triplicate mesocosms were fertilized with ∼ 0.8 µM KH 2 PO 4 between day 4 and day 5 of the experiment to promote a diazotroph bloom. Samples were collected every morning at 07:00 a.m. over a period of 23 days from 3 selected depths (1, 6, 12 m) in each mesocosm (M1, M2 and M3) and in surrounding waters (hereafter referred to as Nouméa lagoon waters). Further details regarding the sampling are given separately (Pfreundt et al., 2015;Van Wambeke et al., 2015). Subsampling and analysis procedures for stocks and fluxes are not given in detail here as they were mainly used for statistical analyses (see Sect. 2.3) but they are described in full in companion papers from this issue: Inorganic nutrients and chloro- 15 phyll a (chl) concentrations, DIP turnover time, N 2 fixation rates and primary production (PP) were analyzed according to procedures detailed in (Berthelot et al., 2015). Synechococcus and Prochlorococcus cell numbers were counted by flow cytometry (Leblanc et al., 2015). Heterotrophic bacterioplankton abundances were also determined by flow cytometry and corresponded to the sum of cells with high and low nu-20 cleic acid content, discarding autotrophic cell counts (Van Wambeke et al., 2015). Heterotrophic bacterial production (BP) was determined by the 3 H leucine technique, and alkaline phosphatase activity (APA) was determined on the total fraction (i.e., unfiltered samples) using MUF-P substrate (Van Wambeke et al., 2015). Three periods of the experiment were defined. P0 was defined as days 2-4, before DIP fertilization of the 25 mesocosms. P1 (days 5-14) was characterized by DIP availability and correspondingly lower APA inside the mesocosms (Van Wambeke et al., 2015). P2 (days [15][16][17][18][19][20][21][22][23]  the mesocosms and also in the Nouméa lagoon, albeit to a lower extent. Further, APA increased steeply corresponding to DIP becoming limiting in the mesocosms during P2 (Berthelot et al., 2015;Van Wambeke et al., 2015).

16S rRNA gene amplicon sequencing and OTU clustering
DNA was extracted from seawater samples collected as described in Sect. 2.1 ev-5 ery two days at 1 and 12 m depth in M1 and in Nouméa lagoon waters. A volume of 10 L was pre-filtered through a 1 mm mesh to keep out large eukaryotes and filtered on 0.45 µm pore size, 47 mm diameter polyethersulfone filters (Pall Supor), immersed in RNA resuspension buffer (10 mM NaAc pH 5.2, 200 mM D(+)-sucrose, 100 mM NaCl, 5 mM EDTA) and snap frozen in liquid N 2 . Tubes with filters were vortexed, then agi-10 tated in a Precellys bead beater (Peqlab, Erlangen, Germany) 2×(2×15 s) at 6500 rpm after adding 0.25 mL glass beads (0.10-0.25 mm, Retsch, Frimley, UK) and 1 mL PGTX (39.6 g phenol, 6.9 mL glycerol, 0.1 g 8-hydroxyquinoline, 0.58 g EDTA, 0.8 g NaAc, 9.5 g guanidine thiocyanate, 4.6 g guanidine hydrochloride, H 2 O to 100 mL; Pinto et al., 2009). RNA and DNA was extracted simultaneously by adding 0.7 mL chloroform, vig-15 orous shaking, incubation at 24 • C for 10 min and subsequent phase separation by centrifugation. RNA and DNA was retained in the aqueous phase, precipitated together and stored at −80 • C for further use for metatranscriptomics analysis (Pfreundt et al., 2015) and 16S tag sequencing. For the latter, ∼ 100 ng of RNase treated total DNA was sent to a commercial provider (LGC Genomics, Berlin, Germany) for amplicon 20 sequencing with primer pair S-D-Bact-0341-b-S-17/S-D- per sample. Merged reads were then quality filtered, reads shorter than 350 nt discarded, and converted to fasta. All fasta files were concatenated for OTU clustering, with the following steps: dereplication at 100 % identity and keeping information on the weight of each unique sequence, sorting of unique sequences by decreasing weight and discarding of singletons (weight = 1), clustering into OTUs with the -cluster_otus 5 command using a maximum dissimilarity of 2 %. The last step includes chimera filtering (Edgar, 2013). Finally, the merged reads from the different samples were mapped back onto the generated OTUs using vsearch v1. To calculate pseudo-absolute cell numbers from these relative OTU abundances, for each sample the cumulated relative abundances of Synechococcus and Prochlorococcus were divided by the respective sum of Synechococcus and Prochlorococcus cytometric cell counts, producing a sample-specific ratio (%16S/cell). We used these two groups, because flow cytometry is well-established for them and we could expect 20 them to be reasonably precise. This ratio was then used to calculate cell numbers for all investigated phylogenetic groups by dividing their relative abundance in each sample by the sample-specific ratio.
The Shannon index as a measure of diversity was calculated separately for M1 and Nouméa lagoon waters using the diversity function in the R package vegan (Oksanen Introduction

Multivariate statistics of 16S community data, samples and environmental data
A canonical correspondence analysis (CCA) (XLSTATS Excel software) was used to test for significant effects of biogeochemical parameters on variations in the 16Scommunity structure. These parameters were PP, BP, Chl, APA, DIP, N 2 fixation rates, 5 particulate organic carbon concentrations (POC) and temperature (T ). CCA generates ordination axes which are linear combinations of the environmental variables that best explain microbial diversity composition data (ter Braak, 1986). The variables were the log-transformed biogeochemical stocks and fluxes, the objects were the 34 dominating OTUs in terms of abundance (expressed in relative OTU abundance), and the sites 10 were sampling time and origin of the samples. Non Metric Multidimensional Scaling (NMDS) was used to find the best display of differences in 16S community structures. The full OTU table containing relative abundances of each OTU per sample, with samples as rows and species as columns, was used for unconstrained ordination with the R package vegan (Oksanen et al., 2015). 15 Additionally, the same procedure was run with the subset tables containing only "heterotroph" OTUs (i.e. all that were not classified as Cyanobacteria or chloroplasts) or only "autotroph" OTUs (i.e. all classified as chloroplasts or Cyanobacteria without Melainabacteria). First, the abundances were standardized using Wisconsin double standardization, which divides abundances by the species maximum and then stan-20 dardizes samples to equal totals. This converts abundances per species into relative values and therefore reduces the weight of dominant species and increases the weight of rare species. Then, Bray-Curtis distance matrices were computed and an NMDS ordination calculated using metaMDS. MetaMDS did not use any further transformations or standardizations. Assuming that samples of day 2-8, day 10-16, and  Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | alpha = 0.05. Permutational Multivariate Analysis of Variance (PERMANOVA) was conducted using the adonis function in the R package vegan on the Bray-Curtis distance matrices described above. We tested the dependence of sample distances on the factors location (M1 or lagoon waters), depth (1 or 12 m) and time (separation of samples into three phases: day 2-8, [10][11][12][13][14][15][16][18][19][20][21].

Biogeochemical context
Here, we present data on the prokaryotic community in M1 and the Nouméa lagoon from depths of 1 and 12 m, which was close to the bottom of the mesocosm. As a foundation for understanding and discussing these data, selected biogeochemical param-10 eters are shown in Figs. 1 and 2c and d for the corresponding sampling locations. On the evening of day 4, M1 was fertilized with ∼ 0.8 µM DIP, leading to lower APA than in the lagoon until day 19, when APA in M1 increased steeply (Figs. 1c and 2c) and DIP turnover times dropped below initial levels (Van Wambeke et al., 2015). Of the three mesocosms, M1 was the one with the lowest increase of chl a during P2, still reaching 15 higher levels than in Nouméa lagoon waters, whereas PP in M1 did not exceed the values measured in the lagoon ( Fig. 1a and b). BP exhibited two prominent peaks in M1 with a 4-fold and 2-fold increase, on day 4 and on day 21, and was otherwise close to values measured for Nouméa lagoon waters, only slightly exceeding them on most days after day 11 (Fig. 1d).

An overview
We identified a total of 3600 OTUs belonging to bacteria, archaea, and chloroplasts using 16S rRNA gene sequencing data from all samples combined ( M1, 20 samples from the Nouméa lagoon, Supplement). In the following, the term "16S community" will be used for the entirety of all 16S reads belonging to these OTUs. On this basis, the Shannon index as a measure for overall diversity was calculated. With values between 5.6 and 7 at both investigated depths over the 23 days experiment, it was relatively high. Values above 5 have been reported for the edge of the South 5 Pacific gyre (Yin et al., 2013), but based on clone libraries, which might underestimate diversity. Over the course of the experiment we observed a slight decrease in diversity both in M1 and in the lagoon (Fig. S1 in the Supplement). The dominating classes in M1 were alpha-and gammaproteobacteria, cyanobacteria, chloroplasts, Marine Group II (MGII) euryarchaea (Thermoplasmata, mainly days 10 and 14), and on days 12-23 also Flavobacteriia (Fig. 3a). Of these, alpha-and gammaproteobacteria exhibited clearly different dynamics between M1 and the lagoon with completely different peaking times, whereas cyanobacteria, chloroplasts and MGII euryarchaea showed similar trends (Figs. 3a and 4a, correlations in Fig. S2a in the Supplement). 15 Amongst the classes contributing less than 5 % to the total 16S community, dynamics were generally different between M1 and the lagoon (Fig. S2b). Inside M1, Deferribacteres and Deltaproteobacteria abundance decreased from > 1.2 % to < 0.7 % of the 16S community after the DIP fertilization and increased again when DIP concentrations fell below 0.1 nM on day 12 (Fig. 3c). Cytophagia and Flavobacteriia also 20 increased at the same time, while Actinobacteria and Phycisphaerae increased on day 20 ( Fig. 3b). Except for the increase in Flavobacteriia and cyanobacteria (predominantly Synechococcus), none of these late-phase increases were observed in the lagoon, although several of the less abundant classes (Cytophagia, Betaproteobacteria, Actinobacteria) more than doubled from day 16 to 18 in the lagoon, but also decreased 25 just as quickly after that (Fig. 4b

Heterotrophic bacteria
As alpha-and gammaproteobacterial 16S abundances showed close to zero correlation between M1 and the lagoon (Fig. S2a), we examined the dynamics within these classes at family level. We observed clear successions within M1 that were not mirrored in the lagoon (Figs. 3c and d, 4c and d). Already on day 4 (two days after the meso-5 cosms had been closed around the water column), the relative share of Rhodobacteraceae 16S tags increased from 3.5 to 12.5 % in M1, exhibiting the most extreme dynamics of all families (Fig. 3c). One particular OTU (297), classified as Thalassobium sp., contributed 79 % of Rhodobacteraceae sequences on that day, but only 2-20 % in other samples. This peak corresponded to increased transcript accumulation from 10 Rhodobacteraceae seen in the metatranscriptome analysis (Pfreundt et al., 2015) and a 4-fold increase in bacterial production on that day (Fig. 1d). This situation occurred the morning before the DIP fertilization and was accompanied by a relative increase in SAR86 and Rhodospirillaceae by 35 and 45 % to > 13 and 9.7 % of the 16S community, respectively. In addition, a 60 % drop in Alteromonadaceae down to 2 %, and 15 a severe drop of Synechococcus from 20 to < 8 % of the 16S community were also observed ( Fig. 3c, d, and f). While Rhodobacteraceae dropped to initial values similar to those in the lagoon on day 6 and contributed around 3 % until day 16, relative abundance of SAR11, dominated by the Surface 1 group, increased from 8 to 25 % of the 16S community on day 8 and outcompeted all other heterotrophic bacteria in M1 20 ( Fig. 3c and d). SAR11 was replaced by SAR116, Rhodospirillaceae (mainly AEGEAN-169), and SAR86 two days later, on day 10, which constituted the three almost equally dominant groups until day 18 with 7 to 12 % of the 16S community each. On day 21, Rhodobacteraceae other than Thalassobium became the dominant heterotrophs in M1 and the lagoon (Figs. 3c and 4c).
Within alphaproteobacteria, SAR11 abundance dynamics were completely uncorrelated between M1 and the lagoon, whereas other groups were slightly positively correlated (Pearson correlation ∼ 0.5, Fig. S2c). Within gammaproteobacteria, Alteromon- adaceae and SAR86 together constituted between 75 % and almost 90 % (Fig. 3d) and whereas the former showed very similar dynamics between M1 and the lagoon waters, SAR86 did not (Fig. S2d). A second abundance peak of SAR86 occurred only in M1 on days 16 and 18 and coincided with similar peaks of the less abundant heterotrophic classes Deferribacteres, Betaproteobacteria, and Cytophagia. While SAR86 never de-5 creased below 7.5 % in M1, these bacteria constituted below 6 % of the 16S community in the lagoon throughout P2 and dropped to 2.5 % on day 21. Within Flavobacteriia, the dominant families exhibited highly correlated temporal dynamics in M1 and the lagoon and also between the different families (Fig. S2e). The Flavobacteriaceae family dominated, increasing more than three-fold between day 10 and 12, in the period of DIP consumption and coinciding with higher chl concentrations from day 12 onward (Figs. 1b and 2c). The maximum relative abundance of almost 15 % was reached in the lagoon on day 18 (Fig. 4e). In M1, Flavobacteriaceae constituted at most 9 %. Two other families, the NS7 and NS9 marine groups, had very low abundances until day 18, then increased afterwards.

Cyanobacteria
Synechococcus sp. (OTU-9) was the dominant OTU throughout the experiment in M1 and the lagoon, being almost an order of magnitude more abundant than the second ranked OTU at some time points. The steep drop of Synechococcus relative abundance on day 4 in M1, two days after the mesocosm closure, was equally seen in absolute 20 counts by flow cytometry (from ∼ 78 000 to ∼ 44 000 cells mL −1 , Figs. 2a and 3f). At the same time in M1, Prochlorococcus gained in relative and absolute abundance reaching a maximum of 8 % of the 16S community on day 6, corresponding to ∼ 16 000 and ∼ 18 000 cells mL −1 on day 6 and 8, respectively ( 4, owing mainly to Trichodesmium with a maximum of 1.5 % of the 16S community on day 2 (Fig. 3f). After the DIP fertilization, Trichodesmium 16S tags dropped quickly by a factor of ten, leaving Calothrix as the most abundant detected diazotroph until day 12. Important to mention, the SILVA database used for 16S classification lists all Richelia intracellularis 16S sequences as descendants of the Calothrix node. This is discussed 5 in Sect. 4. Candidatus Atelocyanobacterium thalassa (UCYN-A) constituted the third most abundant diazotroph until day 14. UCYN-A may be seen as a heterotrophic diazotroph, as it misses photosystem II, cannot fix CO 2 into biomass (Tripp et al., 2010), and lives in symbiosis with unicellular algae (Thompson et al., 2012). Between day 14 and 16, there was a steep increase of Cyanothece (UCYN-C), which had been increas-10 ing steadily after day 10, but at generally low levels. Nevertheless, UCYN-C was the most abundant diazotroph on day 16 (Fig. 3f). Heterotrophic diazotrophs like Bradyrhizobium and Mesorhizobium were, with a maximum of 0.01 % of the community for the former, not present in numbers anywhere close to those for N 2 -fixing cyanobacteria.
In the lagoon samples, the three dominant groups were the same, but with different 15 dynamics. Notably, Cyanothece OTUs did not increase in abundance in the lagoon. After day 4, the temporal evolution of relative (16S community) and absolute (flow cytometry) Synechococcus abundances in M1 diverged substantially (compare Figs. 3f and 2a). As Synechococcus and Prochlorococcus together constituted over 90 % of the cyanobacterial 16S community and thus up to 30 % of all 16S reads, changes 20 in absolute cell counts for these groups can highly influence the temporal evolution of relative 16S data. To elude possible ambiguous interpretation of these data, we used the combined Synechococcus and Prochlorococcus cytometric counts to calculate a %16S/cell ratio for each sample (see Sect. 2.2). This ratio was then used to calculate absolute cell numbers from relative 16S abundance in M1 for all other groups 25 investigated here (Fig. S3 in the Supplement). We assumed equal 16S gene copy numbers.
Interestingly, this approach revealed a steep increase (> 2-fold) of total nonphotosynthetic cells (sum of all abundant classes except cyanobacteria and chloro-BGD 12,2015  plasts per sample) from day 12 to day 16 in M1 that was not seen in flow cytometry counts of total heterotrophic bacteria ( Fig. S3) but was accompanied by a moderate increase (50 %) in BP (Fig. 1d), supporting the former. Spearman rank correlation between BP and 16S-derived numbers for non-photosynthetic bacteria was 0.61 (0.73 when omitting day 21 with its prominent BP peak). N 2 -fixation rates were also substan-5 tially higher throughout this period (days 12 to 16) than before and stayed high until the end of the experiment on day 23 (Fig. 1d). Preceding and overlapping this reaction of heterotrophs was a constant increase in Synechococcus from its minimum of ∼ 20 000 cells mL −1 on day 9 to maximally 103 000-104 000 cells mL −1 on day 19-20. 10 Using all OTUs as a basis, we did unconstrained ordination of all samples on Bray-Curtis dissimilarities calculated from standardized relative abundances. Non-metric multidimensional scaling (NMDS) showed that inside M1 as well as in lagoon waters, the microbial communities defined three distinct clusters (Fig. 5a). Days 2-8, days 10-16, and days 18-21 (the last day 16S rRNA gene sequencing was done for) were 15 grouped together and the clusters separated along the first MDS axis both in lagoon waters and in M1, whereas the locations lagoon waters and M1 (color coded in blue and orange, respectively) were also clearly separated along the second MDS axis. These clusters-in-time will be called early, mid-time, and late from here on and were time-delayed from the periods P0, P1, and P2 defined by biogeochemical parameters, 20 production rates, and N 2 fixation (see Sect. 2.1). Thus, communities from M1 and the lagoon were clearly different from each other but followed a similar trajectory with time.

Temporal dynamics of the 16S community
When we separated the autotrophic community (547 OTUs) from the heterotrophic community (2981 OTUs, all OTUs other than cyanobacteria and chloroplasts) for ordination, it became evident that the heterotrophs showed a pattern very similar to the full 25 community (Fig. 5b). The autotrophic community also clearly separated early, mid-time, and late phases from each other in M1 and the lagoon, but not as far, and separation by location (M1 vs. lagoon waters) only really became evident in the late phase (Fig. 5c) Analysis of variance (PERMANOVA) showed that indeed the community structures of the full OTU set, as well as both separate OTU sets, were significantly different between the three phases (predictor time) and were also different according to their origin, M1 or lagoon waters (predictor location, p values < 0.001, for autotrophs p = 0.003). Using the interaction of the predictor location with time additionally showed that time still had 5 a significant effect after taking out the location effect but not the other way around. This shows that the effect of time was stronger and included the effect of location. Depth had no significant effect.

10
The 34 dominant OTUs (in terms of summed relative abundances of all samples) were subjected to constrained ordination within the frame of eight environmental variables (BP, PP, APA, Chl, POC, N 2 fixation, DIP, and temperature (T )), taking all samples from M1 and the lagoon from 1 and 12 m depth into account (Fig. 6). The sum of all unconstrained eigenvalues was 0.304, the sum of canonical eigenvalues 0.130. Thus, 15 the amount of total variation in the species data that could be explained by the employed environmental variables was 43 % (permutation test, p < 0.001 after 999 permutations). The canonical axes F1 and F2 accounted for 50.5 and 22.3 % of this total explainable variance and 21.6 % and 9.5 % of the total variance, respectively. Both axes were significant and displayed strong species-environment correlations (p < 0.001 for 20 both axes). Samples (location and time) were mainly separated along F1 with day 6 to 14 towards negative values and day 16 to 21 towards positive values. Interestingly, samples from M1 were separated further apart from each other than those from the lagoon. The main parameters separating samples and OTUs along the F1 axis were DIP, which was clearly higher in M1 samples following the DIP fertilization (days 6-10), 25 and BP, Chl, and PP in opposition to DIP. Highly associated with BP and Chl were the late phase samples from M1 and lagoon, together with Rhodobacteraceae OTUs 1 and 37 and the flavobacterial NS4 marine group OTU 44, suggesting these directly profited 20193 Introduction

Synechococcus as the most numerous primary producer
In the world's oceans, the unicellular cyanobacteria Prochlorococcus and Synechococcus constitute the most abundant marine picophytoplankton, amounting together up to 10 6 cells mL −1 , and are responsible for up to 50 % of the total CO 2 fixation in some 5 regions (Liu et al., 1997;Li, 1994;Veldhuis et al., 1997). Previous studies (Zwirglmaier et al., 2008) showed that Synechococcus abundance is low in the oceanic oligotrophic gyres, reaching only 4 × 10 3 cells mL −1 , but up to 10 5 cells mL −1 in nutrient-rich upwelling regions (see also the data compilation in Flombaum et al., 2013). Therefore, the cell counts reported here for Synechococcus in the New Caledonian lagoon were 10 with up to 10 5 cells mL −1 among the highest reported thus far for marine waters. In contrast, the maximum values for Prochlorococcus were relatively low. Several studies indicated that Prochlorococcus is adapted to ultra-oligotrophic conditions and is outcompeted by Synechococcus in waters of higher nutrient content, such as upwelling and coastal regions (Partensky et al., 1999a, b), which is consistent with the findings 15 in this study (Leblanc et al., 2015). The Synechococcus group can be further delineated into several major clades with different ecophysiologies (Scanlan et al., 2009). Due to the insufficient phylogenetic resolution, this was not possible based on 16S rRNA gene sequencing. However, in the accompanying metatranscriptome analysis (Pfreundt et al., 2015), we found that the Synechococcus population was dominated by 20 representatives of clade II and to a lower extent, clade IX, within the picophytoplankton subcluster 5.1A. These may correspond to the Synechococcus cells with low and high orange fluorescence, detected during cell sorting by flow cytometry (Van Wambeke et al., 2015), because the phycourobilin : phycoerythrobilin ratio differs between cultured representatives of these clades (Six et al., 2007)  suggesting clade II ecotypes are competitive in warmer, oligotrophic conditions (Sohm et al., 2015), which is consistent with the results of our metatranscriptomic analysis (Pfreundt et al., 2015). Measurements of BP and PP revealed that the organic matter produced by non-diazotrophic phytoplankton was the major carbon source fueling heterotrophic prokaryotes (Van Wambeke et al., 2015). Our community analysis identifies 5 among these primary producers Synechococcus as the numerically most abundant contributor.

Effects of enclosure and DIP fertilization on the microbial community
Apart from investigating the effect of changing diazotroph communities on general community structure, we judged to which extent the enclosure and the DIP fertilization 10 in the deployed large mesocosms changed community composition and development compared to the surrounding lagoon waters. Roughly, the 16S community composition partitioned the samples into three periods, an early period of DIP availability (including samples before the spike) until day 8, a mid-time period of DIP consumption, and a late period of Chl accumulation/P depletion starting with day 18. Remarkably, although we 15 enriched the mesocosms with ∼ 0.8 µM DIP and noticed a significant difference between the communities in M1 and in lagoon waters (Fig. 5a, separation along the second ordination axis, p < 0.001), the evolution of the both communities with time was mostly parallel, shown by their parallel separation along the first ordination axis. This indicates an overarching influence of abiotic factors acting on both locations equally, 20 like temperature, which increased during the course of the study. Despite the parallel evolution of the full community (rare and abundant OTUs had similar weights in NMDS), community composition and dynamics of single abundant groups in M1 were different to those in the lagoon from the beginning. These differences first became evident with the drop in Synechococcus cell counts and 16S reads on days 2 to 4 ( Figs. 2a and 3a). 25 In this process, cyanophages might have played a role as suggested by the detected high gene expression from Myoviridae in M1 compared to the lagoon (Pfreundt et al., 2015). Based on the known host association of these  Ma et al., 2014;Sabehi et al., 2012), these were likely targeting Synechococcus. This event and its putative consequences on the availability of organic material coincided with diatoms sinking out of the water column in the mesocosms, roughly halving in numbers (Leblanc et al., 2015), and overlapped with DIP fertilization on day 4. Thus, it cannot be conclusively determined which had a higher influence on the community 5 differences we observed in the early phase of the experiment, decaying phytoplankton or increased DIP concentrations. Van Wambeke et al. (Van Wambeke et al., 2015) showed that BP was not P-limited on day 4, suggesting that phytoplankton decay had the greater effect. Indeed, POC concentrations in M1, but not in the lagoon, increased from ∼ 10 µmol L −1 on day 2 to 15 µmol L −1 on day 4 (Berthelot et al., 2015). 10 We correlated the relative abundances in M1 and the Nouméa lagoon for all investigated phylogenetic groups (Fig. S2). This approach revealed alpha-and gammaproteobacteria as the least correlated classes between M1 and the lagoon. Investigating this aspect at higher taxonomic resolution showed that SAR11 and SAR86 temporal dynamics differed the most, suggesting that these groups were the most sensitive to 15 changing nutrient sources. SAR11 was present at comparatively low numbers compared to other reports (Morris et al., 2002). We used filters with 0.45 µm pore size, so it cannot be excluded that SAR11 (average size 0.2 µm × 0.4 µm) and the even smaller Candidatus Actinomarina were lost at the beginning of filtering and retained when more cells were on the filter. Thus, we have possibly missed a substantial fraction of these 20 groups in this study. Nevertheless, we can expect the loss rate to be very similar in all our samples, due to the following reasons: 16S and metatranscriptomics analyses were done from DNA and RNA isolated together from the same filter. Transcriptome abundances were tightly correlated (Pearson correlation 0.88 in M1, and 0.96 in the lagoon) between SAR86 and SAR11 (Pfreundt et al., 2015)   . Oppositely, SAR86 genomes encode a disproportionally high number of Ton-B dependent outer membrane receptors (TBDRs) that allow transport of > 600 Da molecules, including carbohydrates and lipids, combined with lipase genes and an expanded sugar utilization metabolism, whereas SAR11 genomes are not known to contain any TBDR genes and have strongly reduced sugar utilization capabilities (Dupont 5 et al., 2012). In our accompanying metatranscriptome study, transcripts for TBDRs and phospholipase (days 14-16) from SAR86 were highly abundant, suggesting membrane lipids were used as an energy source (Pfreundt et al., 2015). Cultivation experiments of SAR11 representatives have shown their ability to utilize a wide variety of small organic molecules, including peptides, to generate energy (Tripp, 2013). The SAR11 subgroup Ia (surface 1) representative Pelagibacter ubique was able to supply one of its obligate required molecules (glycine) by glycolate (Carini et al., 2013), which is one of the main excretion by-products of phytoplankton and is assumed to be abundant in phytoplankton release. Finally, both SAR11 and SAR86 genomes contain proteorhodopsin genes, potentially providing an additional energy source, but the corresponding tran-15 script levels were much higher in SAR11 than in SAR86 in this study (Pfreundt et al., 2015). It appears striking that for both, SAR11 and SAR86, 16S abundance was not well correlated with the respective total transcript accumulation (Pfreundt et al., 2015), but that the two clades were correlated among each other both in transcript abundance 20 and abundance of dominant OTUs. Dupont et al. (2012) hypothesized that, due to their different resource specializations, SAR11 and SAR86 compete only slightly for dissolved organic carbon (DOC), which would allow for their co-existence, albeit not explaining their correlation.
Interestingly, we saw niche partitioning among dominant SAR86 OTUs in constrained 25 ordination (Fig. 6) 86.3 % (OTU 26 and 1229) and 94.9 % (OTU 7 and 98). Thus, those SAR86 OTUs sharing a niche were not more similar on 16S level than in between niches. An aspect of special interest was to recognize community changes which could be attributed to the DIP fertilization and to possibly explain why. CCA suggested that the taxon profiting the most from DIP fertilization was Defluviicoccus (Rhodospirillaceae).
Interestingly, Defluviicoccus-related bacteria were initially described from enhanced biological phosphorus removal systems (Maszenan et al., 2005;McIlroy and Seviour, 2009;Mielczarek et al., 2013). In such systems, Defluviicoccus are known as glycogenaccumulating bacteria, competing with polyphosphate accumulating organisms in enhanced biological phosphorus removal from waste water (Burow et al., 2007;McIlroy 10 and Seviour, 2009;Wong et al., 2004;Wong and Liu, 2007). It is therefore interesting that in constrained ordination, the dominant Defluviicoccus OTU (40) was associated with DIP for M1 samples on days 6, 8 and 10, immediately following the phosphorus fertilization (Fig. 6). While not showing this apparent reaction to DIP, this OTU had similar abundances also in the Nouméa lagoon. There are only few studies reporting 15 Defluviicoccus-related bacteria in marine samples, and if so, these were associated with marine invertebrates (Enomoto et al., 2012;Fan et al., 2013). Based on our data it is likely that there was a minor population of Defluviicoccus-related bacteria that could benefit from the addition of fresh phosphorus in M1. Note however, that this OTU was only 83 % identical to Defluviicoccus vanus 16S, but 99 % identical to environmental se-20 quences that, according to SILVA taxonomy for SSU release 119 and 123, are leaves of the same node as D. vanus and were thus classified as Defluviicoccus (name of the node). Therefore, additional work towards the more precise characterization of these Defluviicoccus-related bacteria is required. 25 Some aspects of the taxonomic succession after day 13, when DIP was mainly depleted, can be related to the strategies used for adaptation to P limitation. Among autotrophs, the main trend observed in M1 and the lagoon waters was the opposition 20199 Introduction of Prochlorococcus and Synechococcus with Synechococcus reaching higher abundances in M1 at the end of the experiment compared to the beginning, when DIP was comparatively low. Both expressed alkaline phosphatase (AP) mainly two days after their respective 16S abundance peak (Pfreundt et al., 2015), suggesting that the increase in cell numbers led to a higher phosphorous demand that was subsequently 5 adapted to. Incorporation of sulfolipids instead of phospholipids into membrane lipids is an ecological strategy developed by cyanobacteria for an advantage in P deficient environments (Van Mooy et al., 2006. The ability to synthetize sulphoquinovosyldiacylgycerol in place of phosphatidylglycerol in P-limited cultures can save 5 to 43 % of total cellular phosphorus demand and is highly variable according to the strain tested 10 (Van Mooy et al., 2009). SqdB, a gene encoding a sulfolipid synthesis protein, was mainly transcribed by Synechococcus and less by Prochlorococcus (Pfreundt et al., 2015) and corresponded to their changes in 16S abundance. In the North Atlantic Ocean, Synechococcus are also good competitors with Prochlorococcus for DIP and ATP uptake (Michelou et al., 2011). A variety of heterotrophs expressed AP mainly 15 after day 10, possibly to the benefit of other populations. Indeed, metagenomic analyses showed that 59 % of marine bacterial APs are extracellular, periplasmic or located at the outer membrane (Luo et al., 2009). Possibly, DIP release through extracellular APA may provide P to neighboring cells that cannot synthesize this enzyme, and such effects could possibly be emphasized in macro-aggregate gels seen during P2 20 .

Possible functional roles of bacterial subcommunities in the degradation of organic matter
The concerted increase of Cytophagia and Flavobacteriia towards the end of the experiment is a recurrent feature seen in older studies based on plateable strains (Fukami 25 et al., 1981b;Painting et al., 1989), and in molecular studies allowing access to the uncultivable community (Teeling et al., 2012 1993;Fontanez et al., 2015;Teeling et al., 2012). Members of the Flavobacteriia are known to degrade organic matter, some even specializing on algal-derived organic matter (Bauer et al., 2006;Bowman, 2006). In the lagoon, members of the Flavobacteriaceae reached even higher relative abundances than in M1. The occurrence of different Flavobacteriaceae phylotypes seems to depend on the nature of the dominant phyto-5 plankton species (Pinhassi et al., 2004). Information retrieved from analysis of cultured representatives show a high capacity to degrade high molecular weight molecules, which is in line with the high ectoenzymatic activities generally associated at the end of blooms (Chróst, 1992;Riemann et al., 2000). The quality of organic matter issued from specific phytoplankton detritus influences bacterial community succession, and 10 while the potential for polymer degradation appears to be widespread among different phyla it is restricted to some very specialized species (Landa et al., 2014;Murray et al., 2007). Our data indicate that both in M1 and the lagoon, there were organic components available favoring their growth, most pronounced during P2.
Our data further indicate that Rhodobacteraceae, when abundant, assimilate much 15 more organic matter than other heterotrophs seen in this study, as dynamics in BP were closely linked to this group. This was not only shown by CCA (Fig. 5) but is also evident from directly comparing BP (Van Wambeke et al., 2015) to Rhodobacteraceae 16S reads. Both striking BP peaks on day 4 and day 21 in M1 correspond to concomitant peaks in Rhodobacteraceae 16S, while no other heterotrophic groups showed these 20 peaks. Rhodobacteraceae are considered ecological generalists, with large gene inventories, but they cannot be easily condensed in an ecologically differentiated cluster (Newton et al., 2010). Roseobacter, a common marine member of the Rhodobacteraceae, have a variety of metabolic pathways that facilitate interaction and attachment with phytoplankton detritus or give them advantages among other bacteria. They can 25 produce quorum sensing molecules, vitamins and antimicrobial compounds (Buchan et al., 2014). CCA showed two dominating Rhodobacteriaceae OTUs linked with chl concentrations (Fig. 5), which argues for their ability to access phytoplankton detritus, while one OTU was closely linked with N 2 fixation. Notably, none of these three OTUs BGD 12,2015 Succession within the prokaryotic communities during the VAHINE mesocosms U. Pfreundt et al. represented the Thalassobius sp. OTU that was responsible for the Rhodobacteraceae peak on day 4.

The diazotroph community and its impact on the heterotroph community structure
Diazotroph community structure (Figs. 3f and 4f) probably changed in M1 due to DIP 5 fertilization and reduced turbulence, with Cyanothece-like cells profiting and forming large aggregates during P2 (Bonnet et al., 2015a). Our 16S based diazotroph abundances compare reasonably well with those found by quantitative PCR of the nifH gene , if we assume that those 16S reads classified as Calothrix are indeed mostly Richelia (Het-1). The SILVA database positions Richelia below the 10 Calothrix node, and the metatranscriptome data published in this issue show a substantial Richelia read fraction and only little Calothrix (Pfreundt et al., 2015), supporting this assumption. Further, UCYN-C nifH gene abundance was at maximum 100 copies mL −1 in M1 , which matches well the maximal 500 cells mL −1 that we calculated from relative 16S values for Cyanothece using flow cytometry counts 15 of Synechococcus and Prochlorococcus (Fig. S3f). Note that there are also nondiazotrophic Cyanothece representatives. This is remarkable given the different methods used to produce these numbers. The higher Trichodesmium abundance inside M1 after its deployment might be attributed to entrapment of filaments that were not present in the lagoon anymore when 20 the sampling started, likely due to high tide currents. The reduction in the number of 16S tags belonging to Trichodesmium following the DIP fertilization on the evening of day 4 was unexpected. To the contrary, the combined effects of reduced turbulence and P limitation was thought to trigger a Trichodesmium bloom. However, a Trichodesmium bloom was observed at the end of the VAHINE experiment in the lagoon (Spungin et al., 25 2015), demonstrating that the small numbers of existing filaments were competent to build up blooms rapidly. While we did see cyanophage gene expression in M1 in the accompanying metatranscriptome study (Pfreundt et al., 2015), the known host ranges 20202 Introduction

Tables Figures
Back Close

Full Screen / Esc
Printer-friendly Version

Interactive Discussion
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | of those phages suggest that they rather would act on Synechococcus, excluding viral lysis as the reason for the Trichodesmium decline. Thus, there were likely other factors involved that deserve further investigation. Trichodesmium has comparatively low growth rates , so enhanced competition after DIP fertilization might be such a factor.

5
Average N 2 fixation rates of 27.3 ± 1.0 nmol N L −1 d −1 were measured inside the mesocosms during P2 (days [15][16][17][18][19][20][21][22][23], and with 69.7 nmol N L −1 d −1 very high maximum rates were reached in M1 between days 18 and 21. A short-term experiment performed by Bonnet et al. (2015a) indicated that DDN was efficiently transferred to non-diazotrophic microbial communities during P2: after 24 h, 16 ± 6 % of the fixed 10 N 2 was released in the dissolved nitrogen pool and 18 ± 4 % was transferred to nondiazotrophic picoplankton. We show that in terms of relative abundance, SAR86, Deferribacteres, Deltaproteobacteria, Cytophagia, and Flavobacteriaceae increased during this period (Fig. 3). In terms of pseudo-absolute abundances, calculated from 16S abundances employing the cytometry cell counts of Synechococcus and Prochloro- 15 coccus, it became clear that indeed almost all heterotrophic groups analyzed here (> 0.5 % relative abundance) increased during this period (Fig. S3). Another publication in this issue showed that BP was always stimulated after short-term ammonium plus nitrate enrichments on day 4 and 20 (Van Wambeke et al., 2015). This finding speaks for a positive effect of released DDN on heterotrophic growth during P2, but 20 effects through phytoplankton growth and decay probably also played a role. Indeed, nitrogen coming from N 2 fixation alone was not sufficient to sustain BP (Van Wambeke et al., 2015). Constrained ordination showed that generally, temporal variation of dominant Rhodobacteraceae, Halieaceae (OM60(NOR5), in SILVA release 119 within the Alteromonadaceae), and Flavobacteriia could be partially be explained with changes 25 in N 2 fixation, whereas variation in other alphaproteobacteria (mainly SAR11, SAR116, and Rhodospirillaceae), SAR86, euryarchaea, and Alteromonas abundances were not explainable by N 2 fixation (Fig. 6), indicating the involvement of additional factors. However, it appears very clear that the shift within the diazotroph community from the dom-

Tables Figures
Back Close

Full Screen / Esc
Printer-friendly Version

Interactive Discussion
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | inant diazotrophs Richelia and Trichodesmium to UCYN-C after day 14, together with very high N 2 fixation rates (Bonnet et al., 2015a), led to an increase in abundance of almost all heterotrophic groups investigated here, consistent with N limitation of BP (Van Wambeke et al., 2015).

5
We show that despite their large size (> 50 m 3 ), the VAHINE mesocosms led to significant changes in the microbial community structure, probably initiated by decay of Synechococcus and diatoms (Leblanc et al., 2015), albeit not exhibiting a typical bottle effect with a copiotroph-dominated community. Interestingly, the accompanying metatranscriptome study for M1 showed gene expression from cyanophages, which might 10 have contributed to Synechococcus decay. Fertilization with 0.8 µM DIP on day 4 did not have directly observable effects on the overall community, as samples from before and four days after fertilization clustered together, but likely influenced development of individual populations, like Defluviicoccus-related bacteria or UCYN-C diazotrophs. The latter increased after day 10 only in the mesocosms, and were probably limited by 15 temperature before that . The shift of the dominant diazotroph from Richelia and Trichodesmium to UCYN-C after day 14 , together with extremely high N 2 fixation rates (Bonnet et al., 2015a), led to an increase in abundance of almost all heterotrophic groups investigated here, consistent with Nlimited BP (Van Wambeke et al., 2015). 20 Linked with the observation of tight correlation between SAR11 and SAR86 transcription during a diel cycle (Aylward et al., 2015) and over three weeks (Pfreundt et al., 2015), we show that the temporal dynamics of individual OTUs of these clades over three weeks can be explained by very similar combinations of environmental variables. With the results of this work we present for the first time in this ecosystem an in-depth 25 analysis of the prokaryotic diversity and community structure and of their changes over time with a focus on diazotrophic organisms. We correlate the prokaryotic community structure with effects of DIP fertilization, of enclosure vs. open waters in the lagoon, and with the fate of DDN and DOC. The determined microbial diversity was very high, and we provide evidence for previously unknown niche partitioning, e.g., among dominant SAR86 OTUs. These results suggest a plethora of still unknown metabolic and regulatory interactions to occur within these marine microbial communities that are worth to 5 be explored further.

Data availability
All raw sequencing data can be downloaded from NCBI's BioProject database under the accession number BioProject PRJNA304389.

Tables Figures
Back Close

Full Screen / Esc
Printer-friendly Version