Ocean acidification impacts bacteria–phytoplankton coupling at low-nutrient conditions

Abstract

The oceans absorb about a quarter of the annually produced anthropogenic atmospheric carbon dioxide (CO${}_{\mathrm{2}})$, resulting in a decrease in surface water pH, a process termed ocean acidification (OA). Surprisingly little is known about how OA affects the physiology of heterotrophic bacteria or the coupling of heterotrophic bacteria to phytoplankton when nutrients are limited. Previous experiments were, for the most part, undertaken during productive phases or following nutrient additions designed to stimulate algal blooms. Therefore, we performed an in situ large-volume mesocosm ($\sim$ 55 m${}^{\mathrm{3}})$ experiment in the Baltic Sea by simulating different fugacities of CO${}_{\mathrm{2}}$ ($f$CO${}_{\mathrm{2}})$ extending from present to future conditions. The study was conducted in July–August after the nominal spring bloom, in order to maintain low-nutrient conditions throughout the experiment. This resulted in phytoplankton communities dominated by small-sized functional groups (picophytoplankton). There was no consistent$f$CO${}_{\mathrm{2}}$-induced effect on bacterial protein production (BPP), cell-specific BPP (csBPP) or biovolumes (BVs) of either free-living (FL) or particle-associated (PA) heterotrophic bacteria, when considered as individual components (univariate analyses). Permutational Multivariate Analysis of Variance (PERMANOVA) revealed a significant effect of the $f$CO${}_{\mathrm{2}}$ treatment on entire assemblages of dissolved and particulate nutrients, metabolic parameters and the bacteria–phytoplankton community. However, distance-based linear modelling only identified $f$CO${}_{\mathrm{2}}$ as a factor explaining the variability observed amongst the microbial community composition, but not for explaining variability within the metabolic parameters. This suggests that$f$CO${}_{\mathrm{2}}$ impacts on microbial metabolic parameters occurred indirectly through varying physicochemical parameters and microbial species composition. Cluster analyses examining the co-occurrence of different functional groups of bacteria and phytoplankton further revealed a separation of the four $f$CO${}_{\mathrm{2}}$-treated mesocosms from both control mesocosms, indicating that complex trophic interactions might be altered in a future acidified ocean. Possible consequences for nutrient cycling and carbon export are still largely unknown, in particular in a nutrient-limited ocean.

Introduction

Since the industrial revolution the oceans have absorbed ca. one half of anthropogenic carbon dioxide (CO${}_{\mathrm{2}})$. This has resulted in a shift in carbonate equilibria and pH (e.g. Caldeira and Wickett, 2003) with potential consequences for organismal physiology (Taylor et al., 2012). In principal, autotrophs should be fertilised by an enhanced CO${}_{\mathrm{2}}$ availability, increasing the production of particulate (POM) and dissolved organic matter (DOM; Hein and Sand-Jensen, 1997; Riebesell et al., 2007). However, most CO${}_{\mathrm{2}}$ enrichment experiments studying natural plankton assemblages under variable nutrient conditions do not reveal a consistent response of primary production to elevated CO${}_{\mathrm{2}}$ (e.g. Engel et al., 2005; Riebesell et al., 2007; Bach et al., 2016). Both amount and stoichiometric composition of algal DOM and POM can be affected by changes in $f$CO${}_{\mathrm{2}}$. For example, Riebesell et al. (2007) and Maat et al. (2014) reported an increased stoichiometric drawdown of inorganic carbon (C) to nitrogen (N) at higher levels of $f$CO${}_{\mathrm{2}}$, most likely as a result from C overconsumption (Toggweiler, 1993).

Heterotrophic bacteria, in oligotrophic systems, are largely dependent on phytoplankton-derived organic carbon (e.g. Azam, 1998), and as such respond to alterations in both the quantity and quality of phytoplankton-derived DOM and POM (e.g. Allgaier et al., 2008; Grossart et al., 2006a). Availability and competition for nutrients, however, can substantially impact $f$CO${}_{\mathrm{2}}$-induced changes in activity and biomass of phytoplankton and subsequently of heterotrophic bacteria. In nutrient-depleted or nutrient-limited systems, bacteria are restricted in their utilization of phytoplankton-derived organic carbon (Hoikkala et al., 2009; Lignell et al., 2008). Consequently, $f$CO${}_{\mathrm{2}}$-dependent increases in inorganic C availability for autotrophs may not stimulate heterotrophic activity, causing a decoupling of heterotrophic and autotrophic processes (Thingstad et al., 2008). The accumulation of bioavailable dissolved organic carbon (DOC) and particulate organic carbon (POC), as a consequence of this decoupling in nutrient-limited oceanic surface waters, may have profound consequences for nutrient cycling and the nature of the oceanic carbon pump (e.g. Thingstad et al., 1997). Given that various studies have reported on limitation of bacterial growth by inorganic nutrients in several parts of the Baltic Sea (e.g. Hoikkala et al., 2009; Kuparinen and Heinänen, 1993), we sought to evaluate the effects of enhanced $f$CO${}_{\mathrm{2}}$ on activity and biomass of free-living (FL) as well as particle-associated (PA) bacteria during a period characterized by low nutrients and low productivity.

Methods

Experimental setup, CO${}_{\mathrm{2}}$ manipulation and sampling

Nine floating, pelagic KOSMOS (Kiel Off-Shore Mesocosms for future Ocean Simulations) mesocosms (cylindrical, 2 m diameter, 17 m long with conical sediment trap extending to 19 m depth) were moored (day $-$10 $=$ $t-10$; 10 days before CO${}_{\mathrm{2}}$ manipulation) at 59${}^{\circ }$51.5${}^{\prime }$ N, 23${}^{\circ }$15.5${}^{\prime }$ E, in the Baltic Sea at Tvärminne Storfjärden on the southwestern coast of Finland on 12 June 2012. Mesocosm bags were rinsed for a period of 5 days, covered on the top and bottom with a 3 mm net to exclude larger organisms. Thereby, the containing water was fully exchanged with the surrounding water masses. Sediment traps were attached to the bottom of each mesocosm at 17 m depth 5 days prior the start of the experiment ($t-\mathrm{5}$). In addition, submerged mesocosm bags were drawn 1.5 m above the water surface, enclosing and separating $\sim$ 55 m${}^{\mathrm{3}}$ of water from the surrounding Baltic Sea, and meshes were removed. Mesocosms were covered by a photosynthetic active radiation (PAR) transparent roof to prevent nutrient addition from birds and freshwater input from rain. Additionally, existing haloclines were removed in each mesocosm as described in Paul et al. (2015), thereby creating a fully homogeneous water body.

The experiment was conducted between 17 June ($t-\mathrm{5}$) and 4 August ($t43$) 2012. To minimize environmental stress on enclosed organisms, CO${}_{\mathrm{2}}$ addition was performed stepwise over 3 days commencing on day $t\mathrm{0}$. CO${}_{\mathrm{2}}$ addition was repeated on $t15$ in the upper mixed 7 m of water to compensate for outgassing. Different $f$CO${}_{\mathrm{2}}$ treatments were achieved by equally distributing filtered (50 $\mathrm{\mu }$m), CO${}_{\mathrm{2}}$-saturated seawater into the treated mesocosms with a water distributor as described by Paul et al. (2015). Control mesocosms were also manipulated with the water distributor and 50 $\mathrm{\mu }$M pre-filtered water without CO${}_{\mathrm{2}}$. CO${}_{\mathrm{2}}$ amendments resulted in ca. 0.04–0.35 % increases in the total water volume across mesocosms (Paul et al., 2015). Integrated water samples (0–17 m) were collected from each mesocosm and the surrounding seawater using depth-integrated water samplers (IWS, HYDRO-BIOS, Kiel). Samples for activity measurements were directly subsampled from the IWS on the sampling boat without headspace to maintain in situ $f$CO${}_{\mathrm{2}}$ levels during incubation.

Unfortunately, three mesocosms failed during the experiment as a consequence of welding faults, resulting in unquantifiable water exchanges with the surrounding waters. Therefore, with reference to the six remaining mesocosms, CO${}_{\mathrm{2}}$ concentrations defining each treatment are reported as the mean $f$CO${}_{\mathrm{2}}$ determined over the initial 43 days ($t\mathrm{1}-t43$) as described in Paul et al. (2015). The control mesocosms (two replicates) had 365 and 368 $\mathrm{\mu }$atm $f$CO${}_{\mathrm{2}}$. The four treatment mesocosms each had 497, 821, 1007 and 1231 $\mathrm{\mu }$atm $f$CO${}_{\mathrm{2}}$. Detailed descriptions on the study site, mesocosm deployment, system performance of the mesocosm facility throughout the experiment, CO${}_{\mathrm{2}}$ addition, carbonate chemistry, cleaning of the mesocosm bags, as well as sampling frequencies of single parameters are given in Paul et al. (2015).

Physical and chemical parameters

Temperature and salinity were measured using a CTC60M memory probe (Sea and Sun Technology GmbH, Trappenkamp, Germany) and are calculated as the mean integrated over the total depth. Photosynthetic active radiation (PAR) was measured with a PAR sensor (LI-COR LI-192) on the roof of Tvärminne Zoological Station.

Samples for dissolved inorganic carbon (DIC) and total pH were gently pressure-filtered (Sarstedt Filtropur PES, 0.2 $\mathrm{\mu }$m pore size) using a membrane pump (Stepdos). As described in Dickson et al. (2007), total pH was determined on a Cary 100 (Varian) spectrophotometer in a temperature-controlled 10 cm cuvette using $m$-cresol indicator dye. DIC concentrations were determined by infrared absorption using a LI-COR LI-7000 on an AIRICA system (MARIANDA e.K., Kiel). Total pH and DIC were used to calculate carbonate chemistry speciation using the stoichiometric equilibrium constants for carbonic acid of Mehrbach et al. (1973) as refitted by Lueker et al. (2000).

Samples for dissolved organic carbon (DOC), total dissolved nitrogen (TDN), dissolved silica (DSi) and dissolved inorganic phosphate (DIP) were filtered through pre-combusted (450 ${}^{\circ }$C, 6 h) GF/F filters (Whatman, nominal pore size of 0.7 $\mathrm{\mu }$m). Concentrations of DOC and TDN were determined using a high-temperature catalytic combustion technique with a Shimadzu TOC-TN V analyser following Badr et al. (2003). DSi concentrations were determined using standard colorimetric techniques (Grasshoff et al., 1983) at the micromolar level with a nutrient autoanalyser (Seal Analytical, QuAAtro). DIP concentrations were determined with a colorimetric method using a 2 m liquid waveguide capillary cell (Patey et al., 2008) with a miniaturized detector (Ocean Optics Ltd).

Total particulate carbon (TPC), particulate organic nitrogen (PON) and total particulate phosphorus (TPP) samples were collected onto pre-combusted (450 ${}^{\circ }$C, 6 h) GF/F filters (Whatman, nominal pore size of 0.7 $\mathrm{\mu }$m) using gentle vacuum filtration and stored in glass Petri dishes at $-$20 ${}^{\circ }$C. Biogenic silica (BSi) samples were collected on cellulose acetate filters (0.65 $\mathrm{\mu }$m, Whatman) using gentle vacuum filtration (< 200 mbar) and stored in glass Petri dishes at $-$20 ${}^{\circ }$C. Filters for TPC$/$PON analyses were dried at 60 ${}^{\circ }$C, packed into tin capsules and measured on an elemental analyser (EuroEA) according to Sharp (1974), coupled by either a Conflo II to a Finnigan Delta${}^{\text{Plus}}$ isotope ratio mass spectrometer or a Conflo III to a Thermo Finnigan Delta${}^{\text{Plus}}$ XP isotope ratio mass spectrometer. Filters for TPP were treated with oxidizing decomposition reagent (MERCK, catalogue no. 112936) to oxidize organic phosphorus to orthophosphate. Particulate silica was leached from filtered material. Concentrations of dissolved inorganic phosphate as well as dissolved silica were determined spectrophotometrically according to Hansen and Koroleff (1999).

Samples for chlorophyll $a$ (Chl $a)$ were filtered on GF/F filters (Whatman, nominal pore size of 0.7 $\mathrm{\mu }$m) and stored at $-$20 ${}^{\circ }$C. Chl $a$ was extracted in acetone (90 %) and samples were homogenized. After centrifugation (10 min, 800 $\times$ $g$, 4 ${}^{\circ }$C) the supernatant was analysed on a fluorometer (TURNER 10-AU) to determine concentrations of Chl $a$ (Welschmeyer, 1994).

Further details on the determination of physical parameters, concentration of Chl $a$, as well as dissolved and particulate nutrients can be obtained from Paul et al. (2015).

Microbial standing stock

Abundance of free-living (FL) heterotrophic prokaryotes (HP) and photoautotrophic prokaryotic (Synechococcus spp.) as well as eukaryotic cells (< 20 $\mathrm{\mu }$m) were determined by flow cytometry (Crawfurd et al., 2016). Briefly, phytoplankton were discriminated based on their chlorophyll red and/or phycoerythrin orange autofluorescence (Marie et al., 1999). In combination with their side scatter signal and size fractionation, the phytoplankton community could be divided into six clusters, varying in size from 1 to 8.8 $\mathrm{\mu }$m average cell diameter (Crawfurd et al., 2016). Three groups of picoeukaryotic phytoplankton (Pico I-III), 1 picoprokaryotic photoautotroph (Synechococcus spp., SYN) and 2 nanoeukaryotic phytoplankton groups (Nano I–II) were detected. Biovolume (BV) estimations were based on cell abundance and average cell diameters by assuming a spherical cell shape. The BV sum of Synechococcus and Pico I-III is expressed as BV${}_{\text{Pico}}$. The BV sum of Nano I and II will be referred to as BV${}_{\text{Nano}}$.

Abundances of FL prokaryotes were determined from 0.5 % glutaraldehyde fixed samples after staining with the nucleic acid-specific dye SYBR green I (Crawfurd et al., 2016). Unicellular cyanobacteria (Synechococcus spp.) contributed at maximum 10 % of total counts. Two additional groups were identified based on their low (LDNA) and high (HDNA) fluorescence. This identification was based on gating of SYBR green I fluorescence against the side scatter signal (Crawfurd et al., 2016). Particle-associated (PA) prokaryotes were enumerated by epifluorescence-microscopy on a Leica Leitz DMRB fluorescence microscope with UV- and blue light excitation filters (Leica Microsystems, Wetzlar, Germany). Fresh samples were gently mixed to prevent particle settling and a 15 mL subsample was filtered on a 0.1 % Irgalan Black coloured 5.0 $\mathrm{\mu }$m polycarbonate filter (Whatman, Maidstone, UK; Hobbie et al., 1977). Filters were fixed with glutaraldehyde (Carl Roth, Karlsruhe, Germany; final conc. 2 %) and stained for 15 min with 4-6-diamidino-2-phenylindole (DAPI, final conc. 1 $\mathrm{\mu }$g mL${}^{-\mathrm{1}}$; Porter and Feig, 1980) directly on the filtration device and rinsed twice with sterile filtered habitat water before air-drying and embedding in Citifluor AF1 (Citifluor Ltd, London, UK) on a microscopic slide (Rieck et al., 2015). Counts were made based on 15 random unique squares as observed at a magnification of 1000$\times$. The total number of heterotrophic PA prokaryotes was enumerated by subtracting Chl $a$ autofluorescent cells from DAPI-stained cells (Rieck et al., 2015).

BV of FL and PA prokaryotes was calculated separately. For FL prokaryotes we estimated BVs on the basis of an average cell volume of 0.06 $\mathrm{\mu }$m${}^{\mathrm{3}}$ (Hagström et al., 1979). BV of PA prokaryotes was calculated from measurements of 1600 cells across three different mesocosms (346, 868, 1333 $\mathrm{\mu }$atm) and three time points ($t\mathrm{0}$, $t20$, $t39$) throughout the experiment (Massana et al., 1997). A resulting average BV of 0.16 $\mathrm{\mu }$m${}^{\mathrm{3}}$ per cell was used to calculate BV of PA prokaryotes derived from cell abundances. We subsequently adopted the term “heterotrophic bacteria”, since bacteria account for the majority of non-photosynthetic prokaryotes in surface waters (Kirchman et al., 2007).

Metabolic parameters

Rates of bacterial protein production (BPP) were determined by incorporation of ${}^{14}$[C]-leucine (${}^{14}$C-Leu; Simon and Azam, 1989) according to Grossart et al. (2006a). Triplicates and a formalin-killed control were incubated with ${}^{14}$C-Leu (213 mCi mmol${}^{-\mathrm{1}}$; Hartmann Analytic GmbH, Germany) at a final concentration of 165 nM, which ensured saturation of the uptake systems of both FL and PA bacteria. Incubation was performed in the dark at in situ temperature (between 7.8 and 15.8 ${}^{\circ }$C) for 1.5 h. After fixation with 2 % formalin, samples were filtered onto 5.0 $\mathrm{\mu }$m (PA bacteria) nitrocellulose filters (Sartorius, Germany) and extracted with ice-cold 5 % trichloroacetic acid (TCA) for 5 min. Thereafter, filters were rinsed twice with ice-cold 5 % TCA, once with ethanol (50 % $v/v$) and dissolved in ethyl acetate for measurement by liquid scintillation counting (Wallac 1414, Perkin Elmer). Afterwards, the collected filtrate was filtered on 0.2 $\mathrm{\mu }$m (FL bacteria) nitrocellulose filters (Sartorius, Germany) and processed in the same way as the 5.0 $\mathrm{\mu }$m filters. Standard deviation of triplicate measurements was usually < 15 %. The amount of incorporated ${}^{14}$C-Leu was converted into BPP by using an intracellular isotope dilution factor of 2. A conversion factor of 0.86 was used to convert the produced protein into carbon (Simon and Azam, 1989). Cell-specific BPP rates (csBPP) were calculated by dividing BPP rates by abundances of FL prokaryotes and PA HP.

Community respiration (CR) rates were calculated from oxygen consumption during an incubation period of 48 h at in situ temperature in the dark by assuming a respiratory quotient of 1 (Berggren et al., 2012). Thereby, oxygen concentrations were measured in triplicate in 120 mL O${}_{\mathrm{2}}$ bottles without headspace, using a fiber optical dipping probe (PreSens, Fibox 3), which was calibrated against anoxic and air saturated water.

Primary production (PP) was measured using radio-labelled NaH${}^{14}$CO${}_{\mathrm{3}}$ (Steeman-Nielsen, 1952) in 0–10 m depth-integrated samples. After incubation of duplicate samples with 10 $\mathrm{\mu }$L of ${}^{14}$C bicarbonate solution (DHI Lab, 20 $\mathrm{\mu }$Ci mL${}^{-\mathrm{1}})$ in 8 mL vials at 2, 4, 6, 8 and 10 m for 24 h, samples were acidified with 1 M HCl to remove remaining inorganic ${}^{14}$C. Radioactivity was determined by using a scintillation counter (Wallac 1414, Perkin Elmer). PP was calculated knowing the dark-control corrected ${}^{14}$C incorporation and the fraction of the ${}^{14}$C addition to the total inorganic carbon pool according to Gargas (1975). Further details on the measurement of CR and PP are given by Spilling et al. (2016a).

Statistical analyses

Permutational multivariate analysis of variance – PERMANOVA (Anderson, 2001) was used to determine associations between physical/chemical variables and biotic variables. PERMANOVA (perm $=$ 9999) was performed to test for significant differences in variance over time and between $f$CO${}_{\mathrm{2}}$-treated mesocosms (Anderson et al., 2008). Environmental data were normalized according to Clarke and Gorley (2006). Biotic abundance data were log($x+\mathrm{1}$) transformed (Clarke and Gorley, 2006). PERMANOVA partitions the total sum of squares based on the experimental design and calculates a distance-based pseudo-$F$ statistic for each term in the model. Distance-based linear modelling (DistLM) was implemented to relate physical/chemical predictor variables and the multivariate assemblage of biotic variables (Supplement Table S1; Legendre and Anderson, 1999; Anderson et al., 2008). The DistLM routine was based on the AIC model selection criterion (see Anderson et al., 2008) using a stepwise selection procedure. In case of equally AIC-ranked models (difference < 1), a model with fewer parameters was preferred. A Principal Component Analysis (PCA) was performed on normalized chemical data to identify chemical gradients and patterns between the differently $f$CO${}_{\mathrm{2}}$-treated mesocosms over time (Clarke and Gorley, 2006). Distance-based redundancy analysis (dbRDA) was used for visual interpretation of the DistLM in multidimensional space (Anderson et al., 2008). Multivariate analyses of physicochemical, metabolic and community data were performed on a reduced data set comprising 10 time points ($t\mathrm{5}-t29$, every 3 days; $t31$), containing all measured activity variables (BPP, areal PP and CR). Missing values of nutrient data or abundance data (based on every-other-day measurements) were estimated as means of the preceding and following measurement day. In general, no activity data were interpolated or extrapolated.

Cluster analyses were performed based on Spearman's rank correlation coefficient calculated for each mesocosm between all possible combinations of LDNA, HDNA, pico- and nanophytoplankton abundances as well as total Chl $a$. Thereafter, $p$ values were corrected for multiple testing according to Benjamini and Hochberg (1995). The R package pvclust was used to assess the uncertainty in hierarchical cluster analysis (Suzuki and Shimodeira, 2015). For each cluster, AU (approximately unbiased) $p$ values (between 0 and 1) were calculated via multiscale bootstrap resampling (Suzuki and Shimodaira, 2015).

I: Free-living (0.2–5.0 $\mathrm{\mu }$m) and II: particle-associated ($>5.0$ $\mathrm{\mu }$m), size-fractionated rates of (a) bacterial protein production (BPP) [ $\mathrm{\mu }$g C L${}^{-\mathrm{1}}$ d${}^{-\mathrm{1}}$] as well as (b) cell-specific bacterial protein production (csBPP) [fg C cell${}^{-\mathrm{1}}$ d${}^{-\mathrm{1}}$] and (c) heterotrophic prokaryotic biovolume (BV${}_{\text{HP}}$) [$\times$10${}^{\mathrm{5}}$ $\mathrm{\mu }$m${}^{\mathrm{3}}$ mL${}^{-\mathrm{1}}$] during the course of the experiment. I: (d) Ratio of high vs. low nucleic acid-stained prokaryotic heterotrophs (HDNA : LDNA). Colours and symbols indicate average $f$CO${}_{\mathrm{2}}$ [$\mathrm{\mu }$atm] between $t\mathrm{1}-t43$.

[Figure omitted. See PDF]

Results of two-factor permutational multivariate analysis of variance (PERMANOVA)${}^{\text{a}}$ on a resemblance matrix (Euclidian distance) of normalized chemical variables (Phosphate, DOC, TDN, DSi, TPC, PON, POP, BSi). Degrees of freedom (d$f$); sum of squares (SS); mean square (MS).

Source of variation	d$f$	SS	MS	Pseudo-$F$	$p$ (perm)	Unique perms
Time	9	309.93	34.436	11.118	0.0001	9920
$f$CO${}_{\mathrm{2}}^{\text{b}}$	4	31.974	7.9936	2.5808	0.0246	9936
Time $x$ $f$CO${}_{\mathrm{2}}$	36	80.177	2.2271	0.71906	0.8794	9904
Residuals	10	30.973	3.0973
Total	59	472

${}^{\text{a}}$ Permutation was performed with unrestricted permutation of raw data. ${}^{\text{b}}$ Pairwise test could only be performed for control mesocosms ($n=$ 2) with each $f$CO${}_{\mathrm{2}}$ treatment ($n=$ 1), due to missing replication for each $f$CO${}_{\mathrm{2}}$ treatment. Pairwise comparison was only significant between controls and the highest $f$CO${}_{\mathrm{2}}$ treatment ($p_{\text{perm}}=$ 0.029).

PERMANOVA, distLM and dbRDA were carried out using Primer 6.0 and PERMANOVA $+$ for PRIMER software (Clarke and Gorley, 2006; Anderson et al., 2008). All other analyses and the visualization of results were performed with R 3.2.5 (R Core Team, 2016) using packages Hmisc (Harrell et al., 2016), vegan (Oksanen et al., 2016), pvclust (Suzuki and Shimodeira, 2015), gplots (Warnes et al., 2016) and ggplot2 (Wickham, 2009).

Eigenvectors and explained variation of the first four axes of a PCA on normalized variables of dissolved and particulate nutrients. Ordination of the PCA is visualized in Fig. 6.

	Variable	PC1	PC2	PC3	PC4
	DOC	$-$0.4	$-$0.23	0.04	0.68
	TDN	0.39	0.21	0.21	0.47
	Phosphate	$-$0.1	0.48	$-$0.74	0.35
	DSi	0.3	0.52	$-$0.03	$-$0.24
	TPC	0.48	$-$0.06	0.03	0.13
	PON	0.46	$-$0.05	$-$0.05	0.16
	POP	0.36	$-$0.39	$-$0.04	0.21
	BSi	0.17	$-$0.51	$-$0.63	$-$0.22
% variation		49.2	19.7	11.4	7.2
cum. % variation		49.2	68.9	80.4	87.6

Results

Bacterial dynamics

Heterotrophic bacterial BV was predominantly comprised of FL bacteria. PA bacteria contributed at maximum 2 $\pm$ 0.7–10 $\pm$ 0.7 % (mean 4.8 $\pm$ 0.6 %) of total bacterial BV. PA bacteria, however, accounted for a substantial fraction of overall BPP (27 $\pm$ 1–59 $\pm$ 7 %, mean 39 $\pm$ 4 %). There was no significant effect of $f$CO${}_{\mathrm{2}}$ on BPP, csBPP or BV of neither FL nor PA heterotrophic bacteria ($p_{\text{perm}}$ > 0.05); however, a significant temporal effect was observed ($p_{\text{perm}}$ < 0.05). Both bacterial size-fractions had distinct dynamics in abundance, BPP and csBPP during the course of the experiment (Fig. 1). BPP and bacterial abundances were closely related to Chl $a$ and BV of nano- and picophytoplankton, trending along with Chl $a$ until $t$10 and then continuing to increase with BVs of nano- and photoautotrophs and Chl $a$. The period between $t16$ and $t26$, following a sharp decrease in Chl $a$ at $t16$, revealed the highest BPP rates across the experiment with lower rates at higher $f$CO${}_{\mathrm{2}}$ for PA as well as FL bacteria. CsBPP rates were lower at elevated $f$CO${}_{\mathrm{2}}$ for only the FL bacteria during this period. Additionally, BVs of FL and PA bacteria revealed contrasting dynamics (Figs. 1, S1 in the Supplement). PA bacterial BVs declined with the decay of Chl $a$, whereas FL bacterial BVs increased strongly associated with an increase in BV of picophotoautotrophs during this period. The ratio of HDNA : LDNA prokaryotes, which both make up FL bacteria, also showed differences between the experimental treatments. Between $t14-t25$ the ratio of HDNA : LDNA was lower at higher $f$CO${}_{\mathrm{2}}$.

Phytoplankton dynamics

Chl $a$ concentration exhibited distinct maxima at two time periods ($t\mathrm{5}$ and $t16$). The second maximum was associated with an increase in the BV of nanophotoautotrophs (BV${}_{\text{Nano}}$, Fig. 2). This increase was reduced in mesocosms containing higher levels of $f$CO${}_{\mathrm{2}}$ between $t13-t17$. The differences in BV${}_{\text{Nano}}$ between the treatments were reflected in lower concentrations of Chl $a$ in the three highest $f$CO${}_{\mathrm{2}}$-treated mesocosm at $t16$. Chl $a$ and BV${}_{\text{Nano}}$ concentrations declined after $t16$. In contrast, BV of picophotoautotrophs (BV${}_{\text{Pico}})$ increased after $t11$, associated with an increase in BV of Synechococcus spp., which accounted for 31 $\pm$ 2 to 59 $\pm$ 2 % of BV${}_{\text{Pico}}$ across the period of this study (Fig. S2). All four groups of picoautotrophs distinguished by flow cytometry, exhibited time-dependent positive or negative relationships with $f$CO${}_{\mathrm{2}}$ (Figs. 3, S2, S3). The Pico I ($\sim$ 1 $\mathrm{\mu }$m) and Pico II taxa infrequently exhibited strong fertilization effects in response to the $f$CO${}_{\mathrm{2}}$ treatment. In contrast, Synechococcus spp. and Pico III were infrequently negatively affected by the $f$CO${}_{\mathrm{2}}$ treatment.

Relation between functional heterotrophic and autotrophic groups

A cluster analysis of pairwise Spearman correlations between functional bacterial and phytoplankton groups revealed a separation based on $f$CO${}_{\mathrm{2}}$ treatment. Specifically the four CO${}_{\mathrm{2}}$-amended mesocosms were readily distinguishable from the control treatments. Multiple bootstrap resampling (Suzuki and Shimodaira, 2015) supported this, but only significantly for the three highest $f$CO${}_{\mathrm{2}}$-treated mesocosms. The two highest $f$CO${}_{\mathrm{2}}$ treatments revealed a positive correlation between LDNA bacteria and Pico I, which could not be observed in any other experimental treatment. In all CO${}_{\mathrm{2}}$-treated mesocosms we observed positive correlations between Synechococcus spp. and Pico III as well as Synechococcus spp. and Pico I, which were not present in both control mesocosms. In contrast, positive correlations between LDNA and HDNA were not detected in any $f$CO${}_{\mathrm{2}}$ treatment. Additionally, positive correlations between Pico I and Nano II as well as HDNA and Cyanobacteria were only present in both controls and the lowest $f$CO${}_{\mathrm{2}}$ treatment (Fig. 4).

After $t10$, the ratio between heterotrophic prokaryotic BV and Chl $a$ varied between the $f$CO${}_{\mathrm{2}}$ treatments, but did not show a consistent pattern. After $t17$, however, the control mesocosms revealed a higher ratio compared to all $f$CO${}_{\mathrm{2}}$-treated mesocosms (Fig. 5).

(a) Concentration of chlorophyll $a$ [$\mathrm{\mu }$g L${}^{-\mathrm{1}}$], (b) biovolume of nanophytoplankton (Nano I and Nano II) [$\times$10${}^{\mathrm{5}}$ $\mathrm{\mu }$m${}^{\mathrm{3}}$ mL${}^{-\mathrm{1}}$] and (c) biovolume of picophytoplankton (Synechococcus spp., Pico I-III) [$\times$10${}^{\mathrm{5}}$ $\mathrm{\mu }$m${}^{\mathrm{3}}$ mL${}^{-\mathrm{1}}$] during the course of the experiment. Colours and symbols indicate average $f$CO${}_{\mathrm{2}}$ [$\mathrm{\mu }$atm] between $t\mathrm{1}$ and $t43$.

[Figure omitted. See PDF]

(a) Biovolume of Synechococcus spp. [$\times$10${}^{\mathrm{5}}$ $\mathrm{\mu }$m${}^{\mathrm{3}}$ mL${}^{-\mathrm{1}}$] and (b–d) biovolume of picoeukaryote groups I–III (Pico I-III) [$\times$10${}^{\mathrm{5}}$ $\mathrm{\mu }$m${}^{\mathrm{3}}$ mL${}^{-\mathrm{1}}$] during the course of the experiment. Colours and symbols indicate average $f$CO${}_{\mathrm{2}}$ [$\mathrm{\mu }$atm] between $t\mathrm{1}$ and $t43$.

[Figure omitted. See PDF]

Heatmap and cluster analysis based on significant Spearman's rank correlation coefficients calculated for each mesocosm between all possible combinations of abundances between different functional heterotrophic prokaryotic and phytoplanktonic groups (high and low nucleic acid-stained prokaryotic heterotrophs (HDNA : LDNA), Synechococcus spp. (SYN), picoeukaryotes I-III (Pico I-III), nanophytoplankton I-II (Nano I-II)) and Chl $a$ based on daily measurements between $t\mathrm{1}$ and $t39$. Colours indicate the Spearman'r rank coefficient ($r_{\text{s}})$ between two variables. $P$ values of correlations were corrected for multiple testing according to Benjamini and Hochberg (1995). Uncertainty in hierarchical clustering was assessed with multiscale bootstrap resampling using approximately unbiased (AU) $p$ values (between 0 and 1; Suzuki and Shimodeira, 2015). Clusters of the three highest $f$CO${}_{\mathrm{2}}$ treatments are significantly different at the 0.05 level. Numbers indicate the $f$CO${}_{\mathrm{2}}$ treatment with average $f$CO${}_{\mathrm{2}}$ [$\mathrm{\mu }$atm] between $t\mathrm{1}$ and $t43$.

[Figure omitted. See PDF]

Results of two-factor permutational multivariate analysis of variance (PERMANOVA)${}^{\text{a}}$ on a resemblance matrix (Euclidian distance) based on normalized metabolic variables (bacterial protein production (BPP), areal primary production (PP) and community respiration (CR)). Degrees of freedom (d$f$); sum of squares (SS); mean square (MS).

Source of variation	d$f$	SS	MS	Pseudo-$F$	$p$ (perm)	Unique perms
Time	9	92.128	10.236	6.73	0.001	9931
$f$CO${}_{\mathrm{2}}^{\text{b}}$	4	16.044	4.011	2.637	0.023	9944
Time $x$ $f$CO${}_{\mathrm{2}}$	36	42.721	1.1867	0.78018	0.792	9904
Residuals	10	15.21	1.521
Total	59	182.46

${}^{\text{a}}$ Permutation was performed with unrestricted permutation of raw data. ${}^{\text{b}}$ Pairwise test could only be performed for control mesocosms ($n=$ 2) with each $f$CO${}_{\mathrm{2}}$ treatment ($n=$ 1), due to missing replication for each $f$CO${}_{\mathrm{2}}$ treatment. Pairwise comparisons were significant between control and all $f$CO${}_{\mathrm{2}}$ treatments ($p_{\text{perm}}$ < 0.04).

Multivariate physicochemical characterization

Integrated water temperature and PAR ranged between 8.0–15.9 ${}^{\circ }$C and 11.2–66.8 mol m${}^{-\mathrm{2}}$ day${}^{-\mathrm{1}}$ during the experiment respectively. Integrated water temperature reached the maximum at $t15$ and dropped again to 8.2 ${}^{\circ }$C at $t31$.

PERMANOVA results (Table 1) on a multivariate assemblage of dissolved (DOC, TDN, Phosphate, BSi) and particulate (TPC, PON, POP, BSi) nutrients showed significant temporal (Time-$F_{\mathrm{9},10}=$ 11.1, $p=$ 0.0001) and spatial variations along the $f$CO${}_{\mathrm{2}}$ gradient ($f$CO${}_{\mathrm{2}}$-$F_{\mathrm{4},10}=$ 2.6, $p=$ 0.02). PCA ordination of the same chemical dataset strongly reflects the temporal pattern, separating the initial time points before $t11$ from other time points of the experiment along the first PCA axis (Fig. 6). Thereby, eigenvectors of TPC and PON were highest on PCA axis 1 (Table 2). PCA axis two was mainly characterized by high eigenvectors of dissolved phosphate as well as dissolved and particulate silica. The first two PCA axes explained 69 % of variation, and cumulatively 80 % including axis three (Table 2).

Ratio of heterotrophic prokaryotic biovolume and Chl $a$ (BV${}_{\text{HP}}$ : Chl $a)$ during the course of the experiment. PA BV${}_{\text{HP}}$ was interpolated using splines with software R (R Core Team, 2016) for time-points, where no data were available. Colours and symbols indicate average $f$CO${}_{\mathrm{2}}$ [$\mathrm{\mu }$atm] between $t\mathrm{1}$ and $t43$.

[Figure omitted. See PDF]

Multivariate characterization of metabolic parameters

PERMANOVA on the resemblance matrix of normalized metabolic variables (BPP, areal PP, CR) revealed significant temporal (Time-$F_{\mathrm{9},10}=$ 6.7, $p=$ 0.0002) and spatial variations along the $f$CO${}_{\mathrm{2}}$-gradient ($f$CO${}_{\mathrm{2}}$-$F_{\mathrm{4},10}=$ 2.64, $p$ < 0.03) (Table 3). DistLM identified significant effects of temperature ($p$ < 0.03), phosphate ($p$ < 0.02), DOC ($p$ < 0.05) and BSi ($p$ < 0.02) on the multivariate assemblage of metabolic variables (Table 4). The stepwise procedure selects PAR, temperature, DOC and phosphate as determining factors (AIC $=$ 59.6, $R^{\mathrm{2}}=$ 0.26, number of variables $=$ 4). The dbRDA ordination separates the temporal development. Thereby, 92 % of the variability in the fitted model and 24 % of the total variation is explained by the first two dbRDA axes (Fig. 6).

Summary of a DistLM procedure for modelling the relationship between physicochemical variables and a resemblance matrix based on a multivariate assemblage comprising normalized data of bacterial protein production (BPP), areal primary production (PP) and community respiration (CR). Redundant physicochemical variables were removed prior to analysis. Therefore, PON and pH were excluded from the subsequent analysis due to high correlations ($r_{\text{s}}$ > 0.9) to TPC and $f$CO${}_{\mathrm{2}}$ respectively.

Variable	SS (trace)	Pseudo-$F$	$p$	Prop.
$f$CO${}_{\mathrm{2}}$	5.0551	1.6527	0.1759	0.03
Temp${}^{*}$	10.209	3.4376	0.0229	0.055
PAR${}^{*}$	6.2466	2.056	0.1067	0.034
DOC${}^{*}$	8.6228	2.8769	0.0474	0.047
TDN	4.7628	1.5545	0.1984	0.026
Phosphate${}^{*}$	12.319	4.1994	0.0111	0.068
DSi	0.26167	0.083	0.9648	0.001
TPC	7.7827	2.5842	0.0613	0.004
POP	5.0171	1.6399	0.1818	0.027
BSi	11.688	3.9696	0.0111	0.064

${}^{*}$ Variables selected in stepwise procedure based on AIC.

Multivariate characterization of the bacterioplankton and phytoplankton community

PERMANOVA on the resemblance matrix of a multivariate assemblage comprising variables of bacterial and phytoplankton communities (abundances of Pico I-III, Nano I-II, FL bacteria (HDNA, LDNA), PA bacteria, Synechococcus spp. and Chl $a)$ revealed significant temporal (Time-$F_{\mathrm{9},10}=$ 56.8, $p=$ 0.0001) and spatial variations along the $f$CO${}_{\mathrm{2}}$-gradient ($f$CO${}_{\mathrm{2}}$-$F_{\mathrm{4},10}=$ 14.9, $p=$ 0.0001) (Table 5). DistLM identified significant effects of $f$CO${}_{\mathrm{2}}$ ($p$ < 0.02), temperature ($p$ < 0.001), phosphate ($p$ < 0.003), TPC ($p$ < 0.001), BSi ($p$ < 0.001) and POP ($p$ < 0.001) on the multivariate assemblage of bacterial and phytoplankton communities (Table 6). The stepwise procedure selects $f$CO${}_{\mathrm{2}}$, temperature, TPC and phosphate as determining factors (AIC$=$ 67.2; $R^{\mathrm{2}}=$ 0.44; number of variables $=$ 4). The dbRDA reveals a separation along the gradient of $f$CO${}_{\mathrm{2}}$ on the second dbRDA axis. The first dbRDA axis represents the overall temporal development. Thereby, the first two dbRDA axes capture 74 % of the variability in the fitted model and 32 % of the total variation.

Discussion

Although OA and its ecological consequences have received growing recognition during the last decade (Riebesell and Gattuso, 2015), surprisingly little is known about the ecological effects on heterotrophic bacterial biomass or the production and coupling of bacterio- and phytoplankton at nutrient-limited conditions. Previous experiments were, for the most part, conducted during productive phases of the year (e.g. phytoplankton blooms), under eutrophic conditions (e.g. coastal areas) or with nutrient additions (e.g. Grossart et al., 2006a; Allgaier et al., 2008; Brussaard et al., 2013; Bach et al., 2016). However, large parts of the oceans are nutrient limited or experience extended nutrient-limited periods during the year (Moore et al., 2013). Thus, we conducted our experiment in July–August, when nutrients and phytoplankton production were relatively low in the northeastern Baltic Sea (Hoikkala et al., 2009; Lignell et al., 2008), and exposed a natural plankton community to different levels of CO${}_{\mathrm{2}}$.

Results of two-factor permutational multivariate analysis of variance (PERMANOVA)${}^{\text{a}}$ on a resemblance matrix (Bray–Curtis similarity coefficient) based on log($X+\mathrm{1}$) transformed abundances of Pico I-III, Nano I-II, FL bacteria (HDNA, LDNA), PA bacteria, SYN and Chl $a$. Degrees of freedom (d$f$); sum of squares (SS); mean square (MS).

Source of variation	d$f$	SS	MS	Pseudo-$F$	$p$ (perm)	Unique perms
Time	9	201.83	22.426	56.754	0.0001	9923
$f$CO${}_{\mathrm{2}}^{\text{b}}$	4	23.631	5.9077	14.951	0.0001	9940
Time $x$ $f$CO${}_{\mathrm{2}}$	36	19.859	0.55164	1.396	0.151	9915
Residuals	10	3.9515	0.39515
Total	59	271.01

${}^{\text{a}}$ Permutation was performed with unrestricted permutation of raw data. ${}^{\text{b}}$ Pairwise test could only be performed for control mesocosms ($n=$ 2) with each $f$CO${}_{\mathrm{2}}$ treatment ($n=$ 1), due to missing replication for each $f$CO${}_{\mathrm{2}}$ treatment. Pairwise comparisons were significant between control and all $f$CO${}_{\mathrm{2}}$ treatments ($p_{\text{perm}}$ < 0.01).

Summary of a DistLM procedure for modelling the relationship between physicochemical variables and a multivariate assemblage, comprising variables of the bacterial and phytoplankton community. The resemblance matrix (Bray–Curtis similarity coefficient) was based on log($X+\mathrm{1}$) transformed abundances of Pico I-III, Nano I-II, FL bacteria (HDNA, LDNA), PA bacteria, Synechococcus spp. and Chl $a$. Redundant physicochemical variables were removed prior analysis. Therefore, PON and pH were excluded from the subsequent analysis due to high correlations ($r_{\text{s}}$ > 0.9) to TPC and $f$CO${}_{\mathrm{2}}$ respectively.

Variable	SS (trace)	Pseudo-$F$	$p$	Prop.
$f$CO${}_{\mathrm{2}}^{*}$	20.469	4.7386	0.0119	0.075
Temp${}^{*}$	51.838	13.718	0.0001	0.191
PAR	10.791	2.4051	0.0813	0.039
DOC	11.14	2.4864	0.0769	0.041
TDN	9.4456	2.0945	0.1078	0.034
Phosphate${}^{*}$	25.649	6.063	0.0029	0.095
DSi	9.5766	2.1246	0.103	0.035
TPC${}^{*}$	36.038	8.8955	0.0002	0.133
POP	52.171	13.827	0.0001	0.193
BSi	36.439	9.01	0.0005	0.134

${}^{*}$ Variables selected in stepwise procedure based on AIC.

Phytoplankton–bacterioplankton coupling at low nutrient conditions

Heterotrophic bacteria are important recyclers of autochthonous DOM in aquatic systems and play an important role in nutrient remineralization in natural plankton assemblages (Kirchman, 1994). BV and production of heterotrophic bacteria are highly dependent on quantity and quality of phytoplankton-derived organic carbon and are usually tightly related to phytoplankton development (e.g. Grossart et al., 2006b; Allgaier et al., 2008). During this study, low nitrogen availability limited overall autotrophic production (Paul et al., 2015; Nausch et al., 2016). This resulted in a post spring bloom phytoplankton community, dominated by picophytoplankton (Paul et al., 2015). This is consistent with previous reports of picophytoplankton accounting for a large fraction of total phytoplankton biomass in oligotrophic, nutrient-poor systems (e.g. Agawin et al., 2000). Chl $a$ dynamics indicated two minor blooms of larger phytoplankton during the first half of the experiment, although picophytoplankton still accounted for mostly > 50 % of the total Chl $a$ during this period (Paul et al., 2015; Spilling et al., 2016b). The phytoplankton development was also reflected in the PCA ordination of dissolved and particulate nutrients, clearly separating the preceding period before $t11$, including the first peak of Chl $a$, from the other observations during the experiment (Fig. 6). The separation was primarily driven by concentrations of particulate matter (Table 2), which decreased until $t11$ and subsequently sank out of the water column (Paul et al., 2015).

Bacterial BV and BPP paralleled phytoplankton development during this period. With the decay of the initial phytoplankton bloom, a second bloom event resulted, comprised primarily of nanophytoplankton and picophytoplankton (Crawfurd et al., 2016). A decrease in nanophytoplankton BV and Chl $a$ concentrations after $t16/t17$ benefitted both FL heterotrophic bacteria and picophotoautotrophs. The increased availability of DOM, resulting from cell lysis and remineralization of POM, was associated with increases in the BV of both groups and bacterial production levels (Figs. 1, S1). We attributed these increases to the cells of Picoplankton, which, due to their high volume-to-surface ratio as well as a small boundary layer surrounding these cells, are generally favoured compared to larger cells in terms of resource acquisition at low nutrient conditions (Raven, 1998; Moore et al., 2013). If cell size is the major factor determining the access to dissolved nitrogen and phosphorous, bacteria should be able to compete equally or better with picophytoplankton at low concentrations (Suttle et al., 1990; Drakare et al., 2003). However, when phytoplankton is restricted in growth due to the lack of mineral nutrients, a strong commensalistic relationship between phytoplanktonic DOM production and bacterioplanktonic DOM utilization may evolve (Azam et al., 1983; Bratbak and Thingstad, 1985; Joint et al., 2002). Although heterotrophic microbes may indirectly limit primary production by depriving phytoplankton of nutrients, they would not be able to outcompete autotrophs completely since this would remove their source of substrates for carbon and energy (Bratbak and Thingstad, 1985; Joint et al., 2002). Such a relationship might explain the paralleled increase in FL bacterial and picophytoplankton BV.

PA bacteria are typically impacted to a lesser extent by nutrient limitation due to consistently higher nutrient availability at particle surfaces (e.g. Grossart and Simon, 1993). This was reflected in this study by the maintenance of high csBPP rates associated with PA heterotrophic bacteria throughout the experiment. Overall, PA bacteria contributed only a minor fraction (at maximum 10 $\pm$ 0.7 %) to the overall bacterial BV, which is typical for oligotrophic or mesotrophic ecosystems (Lapoussière et al., 2011). Nevertheless, their substantial contribution to overall BPP emphasizes their importance, especially during such low productive periods (e.g. Grossart, 2010). PA heterotrophic bacteria are essential for the remineralization of nutrients from autotrophic biomass, which would otherwise sink down from surface waters (Grossart, 2010). Leakage of hydrolysis products and the attachment and detachment of bacteria to and from particles stimulate production amongst free-living bacteria (Smith et al., 1992; Grossart, 2010) and picophytoplankton.

Effects of $f$CO${}_{\mathrm{2}}$/pH on phytoplankton–bacterioplankton coupling at low nutrient conditions

(a) First and second axis of a principal component analysis (PCA) calculated on normalized variables of dissolved and particulate nutrients ($n=$ 60). The set of variables and the eigenvectors of the first four axes can be derived from Table 2. (b) Ordination of a distance-based redundancy analysis (dbRDA) for visual interpretation of distance-based linear modelling (DistLM) between physical/chemical predictor variables and metabolic variables as well as (c–d) abundances of functional bacterial and phytoplankton groups. A table comprising the set of variables used for DistLM can be derived from the Supplement (Table S1).

[Figure omitted. See PDF]

The response of heterotrophic bacteria to changes in $f$CO${}_{\mathrm{2}}$ has been previously shown to be related to phytoplankton rather than being a direct effect of pH or CO${}_{\mathrm{2}}$ (e.g. Allgaier et al., 2008; Grossart et al., 2006a). Here, neither BPP nor BV of neither FL nor PA bacteria suggested a direct effect of CO${}_{\mathrm{2}}$ on these variables (PERMANOVA). Differences in FL bacterial BV, BPP and the ratio of HDNA : LDNA occurred along the gradient of $f$CO${}_{\mathrm{2}}$, but were limited to short time periods. Furthermore, these changes were not consistent with $f$CO${}_{\mathrm{2}}$, resulting in both increases and decreases of a particular variable at specific times (Fig. 1). Periods where $f$CO${}_{\mathrm{2}}$-related effects were apparent comprised periods with high organic matter turnover (e.g. breakdown of Chl $a$ maximum). However, Paul et al. (2015) could not reveal any effect of $f$CO${}_{\mathrm{2}}$ on the export of carbon, neither across the study period nor at individual time points. Thus, it is reasonable to assume that these small $f$CO${}_{\mathrm{2}}$-related differences in bacterial variables are a consequence of other altered components of the aquatic food web, and do not necessarily manifest as changes in carbon export.

Given the inability to relate individual aspects of microbial metabolism or community composition to $f$CO${}_{\mathrm{2}}$ concentrations, we sought to determine whether an impact was evident using a multivariate approach. Chemical, metabolic and community matrices exhibited large variations in relation to a strong temporal effect throughout the whole sampling period ($p$ < 0.01; Tables 1, 3, 5). In addition, an effect of the $f$CO${}_{\mathrm{2}}$ treatment was also evident in all three multivariate assemblages, albeit explaining far less of the observed variability in chemical and metabolic variables ($p$ < 0.03; Tables 1, 3, 5). However, when relating physicochemical to metabolic variables (DistLM, Table 4), neither $f$CO${}_{\mathrm{2}}$ nor pH were suitable to explaining the observed variability. In contrast, $f$CO${}_{\mathrm{2}}$ contributed to explaining the variability amongst the bacterioplankton–phytoplankton community dynamics (DistLM, Table 6). Taken together, this suggests that effects of $f$CO${}_{\mathrm{2}}$ treatments manifest indirectly, through either altering physicochemical parameters or more likely the composition of the microbial community with possible but so far hidden consequences for microbial metabolism.

$f$CO${}_{\mathrm{2}}$/pH effects on phytoplankton indirectly alter phytoplankton–bacterioplankton coupling at low-nutrient conditions

Autotrophic organisms can be fertilized by an enhanced CO${}_{\mathrm{2}}$ availability, altering growth conditions of phytoplankton and increasing the production of particulate (POM) and dissolved organic matter (DOM; Hein and Sand-Jensen, 1997; Riebesell et al., 2007). As a consequence of this increased photosynthetic fixation rate, both quantity and quality of dissolved organic matter (DOM) available for heterotrophic bacteria are impacted, with potential implications for the nature of coupling between phytoplankton and bacterioplankton at low-nutrient conditions (Azam et al., 1983; Bratbak and Thingstad, 1985). So far, CO${}_{\mathrm{2}}$ enrichment experiments examining natural plankton assemblages (e.g. Engel et al., 2005; Riebesell et al., 2007; Bach et al., 2016) have not revealed a consistent pattern of species response or primary production to elevated CO${}_{\mathrm{2}}$. During this experiment, Spilling et al. (2016a) could not detect any effect of increased CO${}_{\mathrm{2}}$ on total primary production, even though Crawfurd et al. (2016) reported effects of CO${}_{\mathrm{2}}$ on several groups of picophytoplankton. Although one larger picoeukaryote (Pico III) was negatively impacted by $f$CO${}_{\mathrm{2}}$ during our study, two small picoeukaryotes (Pico I, Pico II) benefitted from the CO${}_{\mathrm{2}}$ addition, yielding significantly higher growth rates and BVs at higher $f$CO${}_{\mathrm{2}}$ (Crawfurd et al., 2016). This is consistent with recent evidence suggesting a positive impact of enhanced $f$CO${}_{\mathrm{2}}$ on the abundance of small picoeukaryotic phytoplankton (Brussaard et al., 2013; Newbold et al., 2012; Sala et al., 2015). Both picoeukaryotic groups were identified as variables explaining the separation along the $f$CO${}_{\mathrm{2}}$ gradient on the second and third dbRDA axes in the DistLM ordination of the bacteria–phytoplankton community. Specifically, Pico I was highly negatively correlated ($r_{\text{s}}=-0.67$) to dbRDA axis two. However, dbRDA also indicated opposing effects of $f$CO${}_{\mathrm{2}}$ on Pico II ($r_{\text{s}}=$ 0.54) and HDNA prokaryotes ($r_{\text{s}}=-0.31$), being positively or negatively correlated with axis three. Indeed, sharp increases in BV${}_{\text{Pico II}}$ at high $f$CO${}_{\mathrm{2}}$ between $t14$ and $t17$ were associated with decreases in BV${}_{\text{HDNA}}$.

Although we are not able to draw solid conclusions on the interaction of these two particular groups of organisms, a cluster analysis of pairwise Spearman correlations between functional groups of bacteria and phytoplankton revealed a distinct clustering with mesocosms based on $f$CO${}_{\mathrm{2}}$ levels (Fig. 4). We also detected a change in the ratio of heterotrophic bacterial BV to Chl $a$ between the different $f$CO${}_{\mathrm{2}}$ treatments, though this change was not visible for the entire study duration and not consistent with $f$CO${}_{\mathrm{2}}$. These results strongly suggest that trophic interactions between functional groups of bacteria and phytoplankton might be changing in a future acidified ocean.

In nutrient-poor systems, variable growth rates of phytoplankton, DOM quality and quantity, and also loss of phyto- and bacterioplankton due to grazing or viral lysis may potentially contribute to this observed decoupling of phytoplankton and bacterioplankton at high $f$CO${}_{\mathrm{2}}$ (Azam et al., 1983; Bratbak and Thingstad, 1985; Sheik et al., 2014). The viral shunt or bacterivory may release phytoplankton from competition with bacteria for limited nutrients (e.g. Bratbak and Thingstad, 1985; Caron and Goldman, 1990). How increased $f$CO${}_{\mathrm{2}}$ will affect these processes (e.g. viral lysis and bacterial grazing) under nutrient-limited conditions remains uncertain so far. Bacterial grazing by mixotrophs, which would also directly benefit from increased CO${}_{\mathrm{2}}$ availability (Rose et al., 2009), may provide a mechanism for recycling of inorganic nutrients, otherwise bound in bacterial biomass, as a means for supporting phytoplankton growth (e.g. Mitra et al., 2014). However, other studies examining bacterial grazing under different nutrient conditions reported conflicting positive and negative results of increased $f$CO${}_{\mathrm{2}}$ (e.g. Brussaard et al., 2013; Rose et al., 2009). Although we are unable to draw defined conclusions on how this myriad of complex biological processes are impacted by $f$CO${}_{\mathrm{2}}$, it is very likely that there is an impact on trophic interactions which may account for the portion of unexplained variance we observed in our multivariate analyses.

Conclusions

The use of large-volume mesocosms allowed us to test for multiple $f$CO${}_{\mathrm{2}}$-related effects on dynamics of heterotrophic bacterial activity and their biovolume in a near-realistic ecosystem by including trophic interactions from microorganisms up to zooplankton. Complex interactions between various trophic levels, which can only be properly addressed at the scale of whole ecosystems, are important for understanding and predicting $f$CO${}_{\mathrm{2}}$-induced effects on aquatic food webs and biogeochemistry in a future, acidified ocean. We examined these impacts in a nutrient-depleted system, which is representative for large parts of the oceans (Moore et al., 2013). Heterotrophic bacterial productivity was, for the most part, tightly coupled to the availability of phytoplankton-derived organic matter. When accounting for temporal development and taking into account trophic interactions using multivariate statistics, changes in nutrient composition, metabolic parameters and bacteria–phytoplankton communities revealed a significant effect of the $f$CO${}_{\mathrm{2}}$ treatment. Although not consistent throughout the experiment, differences in the ratio of heterotrophic bacterial BV to Chl $a$ during the last half of the experiment suggest that a future ocean will become more autotrophic during low productive periods as a result of altered trophic interactions between functional groups of bacteria and phytoplankton. There is additional support for this conclusion from examining the atmospheric exchange of CO${}_{\mathrm{2}}$ (Spilling et al., 2016b). During the limited timescale of this study, the observed effects of $f$CO${}_{\mathrm{2}}$ did not manifest as altered carbon export (Paul et al., 2015). However, over several years, maintained changes in nutrient cycling, as a consequence of a permanent decoupling between bacteria and phytoplankton, are likely to arise and impact the nature of the carbon pump.

Data availability

Data of primary production and respiration can be obtained from Spilling et al. (2016c). Other variables from the experiment (e.g. total particulate and dissolved nutrients) can be found in Paul et al. (2016). Flow cytometry data can be obtained from Crawfurd et al. (2016). Bacterial protein production rates and abundances of particle-associated bacteria can be obtained from Hornick et al. (2016).

The Supplement related to this article is available online at doi:10.5194/bg-14-1-2017-supplement.

Acknowledgements

We thank the KOSMOS team and all of the participants in the mesocosm campaign for organization, maintenance and support during the experiment. In particular, we would like to thank Andrea Ludwig for coordinating the campaign logistics and assistance with CTD operations, the diving team and Allanah Paul for her help in data acquisition. Furthermore, we thank the Tvärminne Zoological Station for the opportunity to carry out such a big mesocosm experiment at their research station and for technical support on site. Additionally, we acknowledge the captain and crew of R/V ALKOR for their work transporting, deploying (AL394) and recovering (AL397) the mesocosms. We also thank the two reviewers for their comments to improve this article. The collaborative mesocosm campaign was funded by BMBF projects BIOACID II (FKZ 03F06550) and SOPRAN Phase II (FKZ 03F0611). Thomas Hornick was further financially supported by DFG grant GR 1540/23-1, given to Hans-Peter Grossart. Corina P. D. Brussaard was financially supported by the Darwin project, the Royal Netherlands Institute for Sea Research (NIOZ) and the EU project MESOAQUA (grant agreement number 228224). Edited by: J. Engström-Öst Reviewed by: L. Rhodes and one anonymous referee