Warming and ocean acidification may decrease estuarine dissolved organic carbon export to the ocean

Relative to their surface area, estuaries make a disproportionately large contribution of dissolved organic carbon (DOC) to the global carbon cycle, but it is unknown how this will change under a future climate. As such, the response of DOC fluxes from microbially dominated unvegetated sediments to individual and combined future climate stressors of temperature change (from Δ−3 to Δ+5 ^∘C compared to ambient mean temperatures) and ocean acidification (OA, ∼ 2× current CO₂ partial pressure, pCO₂) was investigated ex situ. Warming alone increased sediment heterotrophy, resulting in a proportional increase in sediment DOC uptake; sediments became net sinks of DOC (3.5 to 8.8 mmol C m⁻² d⁻¹) at warmer temperatures (Δ+3 and Δ+5 ^∘C, respectively). This temperature response changed under OA conditions, with sediments becoming more autotrophic and a greater sink of DOC (up to 4× greater than under current pCO₂ conditions). This response was attributed to the stimulation of heterotrophic bacteria with the autochthonous production of labile organic matter by microphytobenthos. Extrapolating these results to the global area of unvegetated subtidal estuarine sediments, we find that the future climate of warming (Δ+3 ^∘C) and OA may decrease estuarine export of DOC by ∼ 80 % (∼ 150 Tg C yr⁻¹) and have a disproportionately large impact on the global DOC budget.

1 Introduction

The aquatic dissolved organic carbon (DOC) pool is one of the largest pools of organic carbon on earth (Hedges, 1987), roughly equivalent in size to the atmospheric CO₂ reservoir (Siegenthaler and Sarmiento, 1993). The role of DOC in long-term carbon storage in the ocean has been a focus of research for decades (Siegenthaler and Sarmiento, 1993; Hansell et al., 2009; Bauer and Bianchi, 2011; Wagner et al., 2020), with DOC reaching the ocean interior being effectively stored for millennia (Hansell et al., 2009). Although phytoplankton in the surface ocean are the main source of DOC globally, with an estimated production of around 50 Pg C yr⁻¹, only 0.3 % of the DOC they produce reaches the ocean interior (Hansell et al., 2009). Most of the DOC produced by phytoplankton is rapidly remineralized in the water column by heterotrophic bacteria (Azam, 1998). Only more recently has the coastal zone been considered a major source of DOC export to the open ocean and deep sea (Duarte et al., 2005; Maher and Eyre, 2010; Krause-Jensen and Duarte, 2016). The shallow coastal zone accounts for 1 to 10 % of global net primary production (NPP) (Duarte and Cebrián, 1996), with up to 33 % of the associated DOC exported offshore and stored in the ocean interior (Krause-Jensen and Duarte, 2016).

Although shallow estuaries and fringing wetlands make up only ∼ 22 % of the world's coastal area (Costanza et al., 1997) and 8.5 % of the total marine area (Costanza et al., 1997), they are quantitatively significant in terms of DOC processing and offshore transport (Smith and Hollibaugh, 1993). In 1998, Bauer and Druffel used radioisotopic carbon (¹⁴C) to identify the source and age of DOC and particulate organic carbon inputs into the open ocean interior. They found that ocean margins accounted for greater organic carbon inputs into the ocean interior than the surface ocean by more than an order of magnitude. Assuming one-third of the DOC produced in the coastal zone (100–1900 Tg C yr⁻¹; Duarte, 2017) is subducted and reaches the ocean interior (Krause-Jensen and Duarte, 2016), 30 to 630 Tg C yr⁻¹, or up to 3.5× more DOC, could reach the ocean interior from coastal areas than from the open ocean (180 Tg C yr⁻¹, Hansell et al., 2009). This is despite coastal areas having a DOC production rate only 0.2 to 3.9 % of that of the open ocean (Duarte, 2017). As such, small changes to the coastal production and export of DOC may have a disproportionate influence on the global DOC budget.

Euphotic estuarine sediments occupy the coastal boundary between terrestrial and marine ecosystems. Microalgal communities (microphytobenthos, or MPBs) are ubiquitous in these sediments, occupying ∼ 40 to 48 % of the coastal surface area (Gattuso et al., 2020) and generating up to 50 % of total estuarine primary productivity (Heip et al., 1995; MacIntyre et al., 1996; Underwood and Kromkamp, 1999). MPBs exude some of the carbon they fix as extracellular substances, including carbohydrates (Oakes et al., 2010), and can therefore be a source of relatively labile DOC in net-autotrophic sediments (Cook et al., 2004; Oakes and Eyre, 2014; Maher and Eyre, 2010). The dominant sink of DOC in estuarine sediments, however, is uptake by heterotrophic bacteria (Azam, 1998). These heterotrophic bacteria not only consume autochthonous DOC from upstream (Boto et al., 1989), but their biomass is influenced by the lability of sediment organic matter (OM) (Hardison et al., 2013), which can be altered by MPB production (Hardison et al., 2013; Cook et al., 2007). Estuarine sediments are therefore a potentially important sink for DOC.

Unvegetated estuarine sediments can affect the quantity and quality of DOC input to the ocean by (1) acting as a source of autochthonous DOC, through MPB production (Duarte, 2017; Krause-Jensen and Duarte, 2016; Maher and Eyre, 2010), or (2) modifying allochthonous and terrigenous DOC inputs (Fichot and Benner, 2014). Through efficient mineralization of DOC (Opsahl and Benner, 1997), estuaries can act as a sink for DOC and a source of CO₂ to the ocean (Frankignoulle et al., 1998; Fichot and Benner, 2014; Sandberg et al., 2004). Given the disproportionate contribution of estuaries to the export of DOC to offshore marine ecosystems, relative to their surface area, it is important to understand how this balance of DOC sources and sinks within estuaries may change with future shifts in climate, particularly expected increases in temperature and ocean acidification (OA) associated with elevated atmospheric CO₂ concentrations.

Climate projection models assuming a high-emission scenario suggest that atmospheric CO₂ concentrations could more than double by the end of the century, increasing the partial pressure of CO₂ (pCO₂) in surface waters to 1000 µatm and decreasing pH by 0.3 units, together termed OA (RCP8.5; IPCC, 2019). There is also expected to be an increase in mean surface ocean temperature by 2–4 ^∘C (RCP8.5; IPCC, 2019) and increased frequency of unseasonably warm days (Morak et al., 2013; Fischer and Knutti, 2015).

Primary producers fix dissolved inorganic carbon (DIC) during photosynthesis and release DOC directly through exudation and/or indirectly when they are grazed upon. Photosynthetically produced DOC is the main source of DOC in the ocean (Hansell et al., 2009) and fuels local microbial mineralization (Azam, 1998). Heterotrophic bacteria within estuarine sediments respire the carbon from DOC as CO₂, which can then be recaptured by photoautotrophs (Riekenberg et al., 2018), closing the microbial loop (Azam, 1998). DOC or DIC that is not captured is ultimately effluxed to the overlying water column and may be transported from estuaries to the coastal ocean. Individually, increased temperature and CO₂ have been reported to enhance primary productivity and DOC production in Arctic (Czerny et al., 2013) and temperate phytoplankton communities (Wohlers et al., 2009; Engel et al., 2011; Liu et al., 2017; Novak et al., 2018; Taucher et al., 2012), and in temperate stream sediments (Duan and Kaushal, 2013). However, one study in a temperate fjord reported no enhancement of DOC production despite CO₂-enhanced phytoplankton productivity (Schulz et al., 2017). This uncertainty of the response to individual climate stressors is exacerbated when considering how the combination of OA and warming will affect DOC processing. To date, only one study has considered this combined stressor effect on DOC (Sett et al., 2018), observing no difference in DOC production by temperate phytoplankton relative to current conditions under the combined stressors (Sett et al., 2018).

To understand the potential effect of future climate on DOC fluxes, it is essential that both individual and combined effects of OA and warming are considered. Here we focus on changes in DOC fluxes in unvegetated estuarine sediments, as these systems have the potential for significant uptake of DOC that is currently exported to the coastal ocean. In this study, benthic DOC responses in unvegetated estuarine sediments were investigated over an 8 ^∘C temperature range under both current and projected future high-pCO₂ conditions in an ex situ laboratory incubation.

We expected that warming would promote a stronger heterotrophic than autotrophic microbial response in shallow euphotic sediments (Patching and Rose, 1970; Vázquez-Domínguez et al., 2012; Yang et al., 2016) and, as such, there would likely be more DOC remineralization (Lønborg et al., 2018) than “new” DOC production by photoautotrophs (Wohlers et al., 2009; Engel et al., 2011; Novak et al., 2018). Moreover, despite the potential stimulation of primary productivity in unvegetated muddy sediments by OA (Vopel et al., 2018), or more likely high pCO₂, and potential enhancement of DOC production (Engel et al., 2013; Liu et al., 2017), this increase in labile DOC may promote bacterial productivity and DOC mineralization (Hardison et al., 2013). In addition, increased DOC availability alone may increase heterotrophic bacterial biomass production and activity (Engel et al., 2013). We therefore predicted that increases in DOC production from OA alone or in combination with warming may be counteracted by increased consumer activity, potentially diminishing the available DOC pool under future climate conditions.

2 Methods

2.1 Study site

A subtidal site (∼ 1.5 m below mean sea level) in the subtropical Clarence River estuary, Australia, was used for this study (29^∘24.21^′ S, 153^∘19.44^′ E; Fig. 1). Sediment at the site was unvegetated and characterized as euphotic cohesive sandy mud (31–36 % of grains within 250–500 µm, 61–65 % within 63–250 µm, and ∼ 2 % < 63 µm; Lewis and McConchie, 1994). Temperature ±0.3 ^∘C, pH ± 0.5 units, and salinity (± < 1 %) were measured over 24 h using a Hydrolab HL7 submerged at the site. The tidal cycle introduced a salinity range of 10–35, pH range of 7.92–8.15 units (min–max), and mean daily temperature of 23.9 ± 1.6 ^∘C (20–25 ^∘C). The surface sediments (0–2 cm) had a porosity of 0.43 and an organic matter content of ∼ 3.5 % (of dry weight), determined from mass loss after combustion (490 ^∘C) of dried sediment (60 ^∘C) (Luczak et al., 1997). The Clarence River estuary has low nutrient loading (Eyre and Pont, 2003) with dissolved inorganic nitrogen (DIN) concentrations < 2 µM (Eyre, 2000). This is consistent with concentrations determined at the time of this study (∼ 0.9–1.9 µM DIN, Supplement Sect. S4).

Figure 1 Study location (x; 29^∘24.21^′ S, 153^∘19.44^′ E) marked on a map of Yamba, NSW, embedded in an east coast map of Australia. © OpenStreetMap contributors 2020. Distributed under a Creative Commons BY-SA License.

2.2 Core collection

Sediment (∼ 20 cm depth) was collected and capped in acrylic cores (9 cm diameter × 47 cm length) allowing for ∼ 1.8 L of overlying water on 9 and 16 January 2018 (15 and 12 cores, respectively). Thalassinidean shrimp, Trypaea australiensis, burrows were avoided and therefore excluded from the collected cores as their occasional inclusion would result in considerable variability in sediment processes (Webb and Eyre, 2004) that would mask potential treatment effects. To ensure sediments were subtidal, cores were collected during low tide. Immediately after core collection, ∼ 700 L of site water was also collected to fill a laboratory incubation setup.

2.3 Incubation setup

Within 6 h of core collection, all cores were in the laboratory, submerged and uncapped in site water. The cores were placed in one of four temperature tanks – control (23 ^∘C), Δ−3 ^∘C (20 ^∘C), Δ+3 ^∘C (25 ^∘C), and Δ+5 ^∘C (28 ^∘C) – filled with ∼ 80 L of site water, with temperatures maintained and monitored via thermo-regulating aquarium pumps. Each tank had three cores (n= 3).

The ex situ study design allowed control of temperature, pCO₂, and light that would be difficult to achieve in situ. Due to limited space, this investigation was conducted over 2 weeks with two complementary incubations repeated back to back. The incubation in the first week (9–12 January 2018) had cores in the four temperature tanks subjected to high pCO₂ (∼ 1000 µatm), achieved with a CO₂-enriched airstream (initially adjusted and set when attached to a LICOR LI-7000) bubbled into tank water via airstones and air pumps. The incubation in the second week (16–20 January 2018) maintained current pCO₂ conditions (∼ 450  µatm) by circulating ambient laboratory air through the tank water via airstones and air pumps. An additional tank was included in week 1 alongside the high-pCO₂ incubation. This tank was a control tank equivalent to the control tank present week 2, allowing for comparison of the two separate incubations (see Table 1 for details). The temperature and pCO₂ manipulations were within 12 and 4 %, respectively, of their in situ ranges (see Sect. 2.1) to reduce any potential shock effect for the sediment community.

Table 1 Start conditions for incubations at current and high pCO₂ levels (^*) showing mean (±SD) of various carbonate parameters. CON^* is the overlapping control core present in the high-pCO₂ incubation week.

Water columns within cores were stirred at ∼ 60 rpm throughout the incubations via magnetic stir bars (∼ 5 cm above sediment surface) interacting with an external rotating magnet, ensuring water columns were well mixed whilst avoiding sediment disturbance (Ferguson et al., 2003, 2004). High-pressure sodium lamps (400 W; Philips Son-T Argo 400) were used to simulate mean daytime field conditions, providing ∼ 270–280 µmol quanta m⁻² s⁻¹ of photosynthetically active radiation at the water surface of the tanks. Lamps were turned on in the mornings in line with natural diel light cycling, following a similar in situ ∼ 12 h dark–12 h light cycle. Cores were pre-incubated at treatment conditions for 36–48 h, before solute-flux incubations. This pre-incubation period would be sufficient for three to six generations of the dominant microbial members of unvegetated estuarine sediments, diatoms, and cyanobacteria (Mori et al., 1996; Greene et al., 1992), allowing time for the microbial community to acclimatize to the new treatment conditions.

2.3.1 Solute-flux incubation

Immediately after pre-incubation, cores were capped for a 20 h (10 h dark, 10 h light) solute-flux incubation to measure rates of O₂, DIC, and DOC production and consumption over a diel cycle. To adhere to natural diel cycling, cores were capped at dusk to start the incubation on a dark cycle. Samples were collected from three cores per tank at each of three time points in the diel cycle (dark start (dusk), dark end and light start (dawn), and light end (dusk)). Water was collected and syringe-filtered to determine concentrations of DIC (0.45 µm Minisart filter, 100 mL serum bottle; without headspace, poisoned with 50 µL of saturated HgCl₂, stored at room temperature) and DOC (GF/F filter, 40 mL glass vial with silicon septum; without headspace, poisoned with 20 µL of HgCl₂, injected with 200 µL of 85 % H₃PO₄, stored at room temperature). As water was removed for sampling, it was replaced with gravity-fed water maintained in a collapsible bag under the same atmospheric conditions and temperature. After all cores were sampled, dissolved oxygen (DO) concentrations, temperature, and pH were measured using a high-precision Hach HQ40d multiprobe meter with an LDO probe and pH probe, calibrated to 3-point NIST buffer scale (R²= 0.99). Probes were inserted into a resealable port fitted in each lid, ensuring no incubation water exchanged with tank water. After the dawn sampling, lamps were switched on.

DIC concentrations were determined with an AIRICA system (MARIANDA, Kiel, Germany) via infrared absorption using a LICOR LI-7000 and corrected for accuracy against certified reference material, batch no. 171 (Dickson, 2010). Measurements on four analytical replicates of 1.5 mL sample volume were used to calculate DIC concentration as the mean of the last three out of four measurements (typical overall uncertainty, < 1.5 µmol kg⁻¹). DIC and pH measurements were then used to calculate the remaining carbonate chemistry parameters (Table 1) using CO₂SYS (Pierrot et al., 2006). Total borate concentrations (Uppström, 1974) and boric acid (Dickson, 1990) and stoichiometric equilibrium constants for carbonic acid from Mehrbach et al. (1973), as refit by Dickson and Millero (1987), were used. Comparison of pH (free scale) measured with a Hach HQ40d multiprobe meter and pH calculated from measured total alkalinity and DIC (Table S1) indicated an uncertainty for potentiometric pH measurements of ±0.05 pH units. Propagating the uncertainty in pH measurements with the uncertainty of DIC measurements translates to a pCO₂ uncertainty of $\pm \sim$ 110 and ∼ 56 µatm under high-pCO₂ and current conditions, respectively. This uncertainty is well within the treatment variability measured among cores (Table 1) and is therefore considered unlikely to have contributed substantially to differences in treatment response. DOC concentrations were measured via continuous-flow wet oxidation using an Aurora 1030 W total organic carbon analyser (Oakes et al., 2011) (uncertainty of ∼ 3 %).

2.4 Data analysis

The dissolved oxygen and DIC measurements were used to estimate benthic microalgal production inside the cores. Net primary production and respiration (NPP and R, µmol O₂ m⁻² h⁻¹) were defined as the light and dark cycle oxygen evolution, respectively, where DIC and DOC light and dark fluxes (µmol C m⁻² h⁻¹) were defined using the evolution of DIC and DOC concentrations, respectively. Fluxes (NPP, R, DIC, or DOC) were calculated as

$$\begin{array}{}\text{(1)}& \mathrm{flux}={\displaystyle \frac{\left(\text{end}-\text{start}\right)\times V}{(T\times A)}},\end{array}$$

where end and start are the dissolved oxygen, inorganic carbon, or organic carbon concentrations (µmol O₂ L⁻¹ or µmol C L⁻¹) at the end and start of the light or dark cycle, V is the water column volume (L), T is hours of incubation, and A is surface area of the core.

Gross primary productivity (GPP, µmol O₂ m⁻² h⁻¹) was calculated using NPP and R, as follows:

$$\begin{array}{}\text{(2)}& \mathrm{GPP}=-R+\mathrm{NPP}.\end{array}$$

The production-to-respiration ratio (P $/$ R) was calculated using GPP and R scaled for a 12 h light–12 h dark diel cycle (Eyre et al., 2011).

$$\begin{array}{}\text{(3)}& P/R={\displaystyle \frac{\left(\mathrm{GPP}\right)\times \mathrm{12}\mathrm{h}}{\left(-R\times \mathrm{24}\mathrm{h}\right)}}\end{array}$$

Finally, net fluxes for DIC and DOC were calculated from the dark and light fluxes from Eq. (1) and presented as mmol C m⁻² d⁻¹ for a 12 h light–12 h dark diel cycle.

$$\begin{array}{}\text{(4)}& \mathrm{Net}\mathrm{flux}=\left(\left(\mathrm{dark}\mathrm{flux}\times \mathrm{12}\mathrm{h})+(\mathrm{light}\mathrm{flux}\times \mathrm{12}\mathrm{h}\right)\right)/\mathrm{1000}\end{array}$$

Temperature sensitivity coefficients (Q₁₀ values) were used to evaluate the temperature dependence of metabolic rates on temperature increases of 10 ^∘C. This was expressed simply as an exponential function:

$$\begin{array}{}\text{(5)}& Q_{\mathrm{10}}={\left({\displaystyle \frac{R_{\mathrm{2}}}{R_{\mathrm{1}}}}\right)}^{\mathrm{10}{}^{\circ }\mathrm{C}/\left(T_{\mathrm{opt}}-T_{\mathrm{1}}\right)},\end{array}$$

where R₁ and R₂ are the R, NPP, or GPP rates measured at temperatures of 20 ^∘C (T₁) and optimal temperatures (T_opt), where rates are highest, respectively.

2.4.1 Scaling rates

Rates in the overlapping control cores each week were checked to ensure comparability between incubations. If means (±SD) were significantly different (did not overlap), rates from individual treatment cores were scaled to the overall mean control rate of both weeks (n= 6). This was done by calculating the relative proportion of treatment rates (tProp., µmol N m⁻² h⁻¹) to the control rates present in its week (Eq. 6),

$$\begin{array}{}\text{(6)}& \mathrm{tProp}.={\displaystyle \frac{\mathrm{tRate}}{\mathrm{control}}},\end{array}$$

where tRate is the individual core rate, and control is the mean rate of control cores present during the incubation (n= 3). This proportional rate was then multiplied by the overall control mean rates (averaged across both weeks, n= 6) to scale individual core rates and calculate comparable treatment means (n= 3) across incubations (see Sect. 3.1 for details on scaled rates).

2.5 Statistical analysis

Homogeneity of variances (Levene's test) was tested before analysis to minimize the potential for type I error. All tests were run in MATLAB (Mathworks, 2011) with significance defined at a maximum alpha of < 0.05. Where Levene's test returned a significant result, either datasets were log-transformed or else, if negative values were present, an alpha of 0.01 was used for the subsequent analyses of variance (ANOVAs).

2.5.1 Net variability with temperature and CO₂

Net fluxes were compared among treatments to identify the individual and combined effects of temperature and pCO₂ on O₂, DIC, and DOC fluxes. To investigate the effect of increased pCO₂ alone, data from control temperature cores at both current and high pCO₂ levels (n= 2) were compared using a one-way ANOVA. A two-way ANOVA on each dataset identified whether there were interacting effects on O₂, DIC, and DOC fluxes of temperature (n= 4) and pCO₂ (n= 2). Finally, one-way ANOVAs were also run for each pCO₂ level to investigate differences in sediment responses across temperatures (n= 4). Post hoc Tukey's tests were then used to determine which temperatures had similar or different responses.

2.5.2 Diel variability with temperature for DIC and DOC fluxes

Differences between dark and light cycles were compared to further investigate changes observed in DIC and DOC net variability. Similar analyses to those described above were applied. To examine differences among temperatures (n= 4), light condition (n= 2), and whether light condition significantly interacted with temperature response, two-way ANOVAs were applied to current- and high-pCO₂ cores, separately. Following this, each light condition was further investigated to consider the individual temperature responses in the light and dark separately using one-way ANOVAs and post hoc Tukey's tests.

3 Results

3.1 Overlapping control scaling

Mean rates calculated from overlapping control cores present in each week were compared to establish whether the two sets of incubations were directly comparable, and whether changes attributed to high pCO₂ levels were truly due to that treatment, not just a temporal shift in how the sediments were behaving. The P $/$ R ratios were similar for incubations (0.84 ± 0.01 and 0.83 ± 0.04, Table S2); however, the magnitude of the R and NPP fluxes was ∼ 23 % greater for control cores in the high-pCO₂ week (Table S2; discussed in Sect. 4.0). As such, R and NPP rates of cores were scaled to mean control rates (n= 6) using the proportional rate difference calculated between the treatments and the individual controls present in the respective weeks (n= 3) (Eq. 6). Scaled rates were within ±13 % of actual rates. There were no significant differences between controls for light or dark production of DIC or DOC.

3.2 Productivity and respiration responses to OA

High-pCO₂ conditions alone (at mean ambient temperatures, 23 ^∘C) significantly increased P $/$ R by ∼ 20 % over control ratios (one-way interaction: $F_{\mathrm{3},\mathrm{4}}=$ 101.9, p= 0.0005; Fig. 2d). This was a result of significant increases in NPP (∼ 42 %) compared to control conditions (one-way interaction: $F_{\mathrm{3},\mathrm{4}}=$ 241.4, p < 0.0005; Fig. 2b), in concert with no significant change in R (one-way interaction: $F_{\mathrm{3},\mathrm{4}}=$ 4.5, p= 0.10; Fig. 2a). Insignificant increases of DIC uptake in the light reflected the significant increases in NPP under high-pCO₂ conditions at ambient temperature (one-way interaction: $F_{\mathrm{3},\mathrm{4}}=$ 5.9, p= 0.07; Fig. 3c). Like R, DIC in the dark did not change with pCO₂ (one-way interaction: $F_{\mathrm{3},\mathrm{4}}=$ 1.3, p= 0.33; Fig. 3b). GPP also significantly increased under high-pCO₂ conditions at ambient temperatures (one-way interaction: $F_{\mathrm{3},\mathrm{4}}=$ 65.3, p= 0.001; Fig. 2c), with net DIC significantly shifting from a slight efflux to a slight influx (one-way interaction: $F_{\mathrm{3},\mathrm{4}}=$ 24.3, p= 0.008; Fig. 3a).

Figure 2 Effect of temperature on rates of (a) respiration (R, µmol O₂ m⁻² h⁻¹), (b) net primary production (NPP, µmol O₂ m⁻² h⁻¹), (c) gross primary production (GPP, µmol O₂ m⁻² h⁻¹), and (d) P $/$ R under current (open boxes) and high-pCO₂ conditions (hatched boxes). Panels show mean values “×” at three temperature offsets from control (CON = 23^∘C). Middle horizontal line in each box is the exclusive median, with the start of the upper and lower quartiles represented by the top and bottom edges of the box, respectively. Letters identify significantly different means across temperatures under current pCO₂ conditions, and numerals identify significantly different means across temperatures under high-pCO₂ conditions, where letters or numerals that are the same indicate no significant difference, as determined by post hoc Tukey's test. Solid and dashed horizontal lines identify significant effects of OA and temperature, respectively, where double solid and dashed lines identify significant interaction of temperature and OA (two-way ANOVA). Levels of significance are denoted with “^*” for 0.05, “${}^{**}$” for 0.01, and “${}^{***}$” for 0.001. Data in Tables S3 and S4.

3.3 Productivity and respiration responses to temperature and OA

Temperature had a strong effect on R, NPP, GPP, and P $/$ R, whereas only light cycle NPP and, in turn, GPP and P $/$ R were affected by OA.

The response of R to temperature was similar at both current and high pCO₂ levels (no two-way interaction: $F_{\mathrm{3},\mathrm{16}}=$ 0.77, p= 0.53; Fig. 2a), and it was not affected by pCO₂ (CO₂ effect, two-way interaction: $F_{\mathrm{1},\mathrm{16}}=$ 0.99, p= 0.34; Fig. 2a). Accordingly, Q₁₀ values for R were similar for current (1.66) and high-pCO₂ (1.69) conditions (Table 2). R changed significantly across the 8 ^∘C temperature range, increasing by ∼ 11 and ∼ 29 % in higher-temperature cores (Δ+3 and Δ+5 ^∘C, respectively) and decreasing by ∼ 16 % in Δ−3 ^∘C cores (temperature effect, two-way interaction: $F_{\mathrm{3},\mathrm{16}}=$ 36.93, p < 0.0001; Fig. 2a).

Sediment NPP was significantly affected by the interaction of pCO₂ and temperature (two-way interaction: $F_{\mathrm{3},\mathrm{16}}=$ 8.92, p= 0.001; Fig. 2b). Under current pCO₂ conditions, NPP decreased significantly with increased temperature (one-way interaction: $F_{\mathrm{3},\mathrm{8}}=$ 41.94, p < 0.0001; Fig. 2b), shifting from net autotrophy in the light in low- and control temperature cores (efflux of 590 ± 121 and 613 ± 10 µmol O₂ m⁻² h⁻¹, respectively) to net heterotrophy in higher-temperature cores (influx of 163 ± 228 and 390 ± 97 µmol O₂ m⁻² h⁻¹ for Δ+3 and Δ+5 ^∘C, respectively). Warming alone therefore resulted in a reduction in NPP by 126 % at Δ+3 ^∘C and 164 % at Δ+5 ^∘C, compared to the control (Fig. 2a). OA in general significantly increased NPP rates over those measured under current pCO₂ conditions (CO₂ effect, two-way interaction: $F_{\mathrm{1},\mathrm{16}}=$ 21.92, p= 0.0003; Fig. 2b), and Q₁₀ of NPP increased from 1.13 to 1.92 when OA was present (Table 2). As such, under high-pCO₂ conditions NPP maintained net autotrophy in the light at Δ+3 ^∘C and only resulted in net heterotrophy in the highest-temperature treatment (one-way interaction: $F_{\mathrm{3},\mathrm{8}}=$ 53.01, p < 0.0001; Fig. 2b).

Table 2 Q₁₀ and T_opt values for current and high-pCO₂ climates.

GPP displayed a similar interactive stressor response to that described for NPP (two-way interaction: $F_{\mathrm{3},\mathrm{16}}=$ 9.39, p= 0.0008; Fig. 2c). Under current pCO₂ conditions, GPP had a slight but insignificant rate increase from lowered to control temperatures (∼ 12 %), where rates significantly decreased at temperatures higher than control (∼ 45 and ∼ 50 % for Δ+3 and Δ+5 ^∘C, respectively) (one-way interaction: $F_{\mathrm{3},\mathrm{8}}=$ 16.89, p= 0.001; Fig. 2c). OA significantly increased GPP at ambient and Δ+3 ^∘C temperatures (CO₂ effect, two-way interaction: $F_{\mathrm{1},\mathrm{16}}=$ 24.77, p= 0.0001; Fig. 2c), resulting in a stronger temperature sensitivity in GPP under high-pCO₂ conditions (one-way interaction: $F_{\mathrm{3},\mathrm{8}}=$ 40.90, p < 0.0001; Fig. 2c) than under current pCO₂ conditions (p= 0.001). This increased sensitivity of GPP to temperature was supported by GPP Q₁₀ value differences between current and high-pCO₂ conditions, increasing from 1.46 to 2.27 (Table 2).

The differences in P $/$ R among treatments further highlighted significant interaction of temperature and pCO₂ (two-way interaction: $F_{\mathrm{3},\mathrm{16}}=$ 5.86, p= 0.007; Fig. 2d), suggesting the effect of pCO₂ on primary productivity was strong enough to alter the overall productivity of the sediments. Under current pCO₂ conditions, P $/$ R reflected GPP with a clear separation between control and Δ−3 ^∘C sediments having a higher P $/$ R (0.84 ± 0.01 and 0.89 ± 0.07, respectively) than the significantly lower ratios calculated in increased-temperature cores (0.42 ± 0.11 and 0.33 ± 0.05 for Δ+3 and Δ+5 ^∘C, respectively) (one-way interaction: $F_{\mathrm{3},\mathrm{8}}=$ 49.41, p < 0.0001; Fig. 2d). Similarly, under high-pCO₂ conditions, the effect of GPP on P $/$ R was clear. The positive effect of high-pCO₂ conditions on GPP pushed the P $/$ R ratio of Δ−3 ^∘C and control temperature cores to ∼ 1 (1.09 ± 0.16 and 1.03 ± 0.03, respectively), suggesting the ecosystem shifted toward net autotrophy under those conditions. The positive effect of high-pCO₂ conditions was also highlighted at Δ+3 ^∘C, with P $/$ R (0.77 ± 0.13) remaining close to the current ecosystem ratio (0.84 ± 0.01) instead of significantly dropping like those calculated under current pCO₂ conditions or in Δ+5 ^∘C cores (0.25 ± 0.04, one-way interaction: $F_{\mathrm{3},\mathrm{8}}=$ 38.58, p<0.0001; Fig. 2d).

Figure 3 Effect of temperature on (a) net dissolved inorganic carbon (DIC) production (mmol C m⁻² d⁻¹), and (b) dark and (c) light DIC fluxes (µmol C m⁻² h⁻¹) under current (open boxes) and high-pCO₂ conditions (hatched boxes). Panels show mean values “×” at three temperature offsets from control conditions (CON = 23^∘C). Middle horizontal line in each box is the exclusive median, with the start of the upper and lower quartiles represented by the top and bottom edges of the box, respectively. Letters identify significantly different means across temperatures under current pCO₂ conditions, and numerals identify significantly different means across temperatures under high-pCO₂ conditions, where letters or numerals that are the same indicate no significant difference, as determined by a one-way ANOVA and post hoc Tukey's test. Solid and dashed horizontal lines identify significant effects of OA and temperature, respectively, where double solid and dashed lines identify significant interaction of temperature and OA (two-way ANOVA). Levels of significance are denoted with “^*” for 0.05, “${}^{**}$” for 0.01, and “${}^{***}$” for 0.001. Data in Table S5.

Figure 4 DIC : DO fluxes from sediment (µmol C m⁻² h⁻¹ or µmol O₂ m⁻² h⁻¹) for all temperatures in dark-and-light cycles subject to current and high-pCO₂ conditions (mean ± SD). Dashed line highlights the 1 : 1 ratio (±18 %; Hopkinson, 1985), with values falling on this line likely a result of aerobic respiration. Arrows indicate the position values would fall in if sediments were experiencing chemical oxidation or anaerobic respiration.

Figure 5 Effect of three temperature offsets from control (CON = 23^∘C) on (a) net dissolved organic carbon (DOC) fluxes (µmol C m⁻² d⁻¹) under current (open boxes) and high-pCO₂ conditions (hatched boxes). Light (grey boxes) and dark fluxes (spotted boxes) of DOC (µmol C m⁻² h⁻¹) for (b) current pCO₂ and (c) high-pCO₂ conditions. In panel (a) letters identify significantly different means across temperatures under current pCO₂ conditions, and numerals identify significantly different means across temperatures under high-pCO₂ conditions; in panels (b) and (c) letters identify significantly different means across temperatures in dark cycles, and numerals identify significantly different means across temperatures in light cycles. Letters or numerals that are the same indicate no significant difference, as determined by a one-way ANOVA and post hoc Tukey's test. Solid and dashed horizontal lines identify significant effects of pCO₂ or light and temperature, respectively, where double solid and dashed lines identify significant interaction of temperature and light (two-way ANOVA). Levels of significance are denoted with “^*” for 0.05, “${}^{**}$” for 0.01, and “${}^{***}$” for 0.001. Data in Table S6.

3.4 DIC fluxes

DIC fluxes mirrored those of dissolved oxygen (Figs. 2 and 3), with both light and dark DIC : DO ratios near 1 : 1 (Fig. 4). In the dark, DIC reflected R responses to temperature; like R, DIC responses to temperature did not differ with pCO₂ (two-way interaction: $F_{\mathrm{3},\mathrm{16}}=$ 0.92, p= 0.45; Fig. 3b), and rates increased with increasing temperature (temperature effect, two-way interaction: $F_{\mathrm{3},\mathrm{16}}=$ 12.66, p= 0.0002; Fig. 3b). In the light, there was a significant interactive effect of temperature and pCO₂ on DIC fluxes (two-way interaction: $F_{\mathrm{3},\mathrm{16}}=$ 12.01, p= 0.0002; Fig. 3c). Under current pCO₂ conditions, DIC reflected the significant NPP responses to temperature, with DIC taken up at Δ−3 ^∘C and control temperatures and effluxed at Δ+3 and Δ+5 ^∘C (one-way interaction: $F_{\mathrm{3},\mathrm{8}}=$ 21.33 p= 0.0004; Fig. 3c).

Net DIC responses were significantly affected by the interaction of pCO₂ and temperature (two-way interaction: $F_{\mathrm{3},\mathrm{16}}=$ 9.69, p= 0.001; Fig. 3a). Like differences in O₂, significant differences between pCO₂ conditions were also measured in the Δ+3 ^∘C temperature cores. At Δ+3 ^∘C, net DIC production in high-pCO₂ cores was ∼ 62 % lower than that measured at the same temperature under current pCO₂ conditions (one-way interaction: $F_{\mathrm{3},\mathrm{4}}=$ 17.1, p= 0.01; Fig. 3a). This again reflected changes in light cycle production, with light DIC effluxes at Δ+3 ^∘C under current pCO₂ conditions becoming influxes under high-pCO₂ conditions (132 ± 74 to −617 ± 88 µmol C m⁻² h⁻¹, Fig. 3a).

3.5 DOC fluxes

At current pCO₂ levels, increasing temperature resulted in a significant shift in net DOC fluxes, going from effluxes at the two lower temperatures (Δ−3 ^∘C and control) to uptakes at the two higher temperatures (one-way interaction: $F_{\mathrm{3},\mathrm{8}}=$ 6.96, p= 0.013; Fig. 5a). The relative light-and-dark cycle contributions of these net trends at current pCO₂ levels were also affected by temperature (two-way interaction: $F_{\mathrm{3},\mathrm{16}}=$ 13.18, p= 0.0001; Fig. 5b). Significant changes in DOC fluxes in the dark shifted from an efflux at Δ−3 ^∘C to an uptake at control temperature, with higher uptake rates at Δ+5 ^∘C (26 % higher than control rates; one-way dark: $F_{\mathrm{3},\mathrm{8}}=$ 8.64, p= 0.007; Fig. 5b). In contrast, the highest DOC effluxes in the light were at control temperatures, significantly decreasing with both increasing and decreasing temperatures to DOC fluxes around zero (one-way interaction: $F_{\mathrm{3},\mathrm{8}}=$ 16.76, p= 0.001; Fig. 5b).

High-pCO₂ conditions alone (at ambient mean temperatures, 23 ^∘C) had a significant effect on net DOC, shifting from a slight efflux at current pCO₂ levels (∼ 0.5 mmol C m⁻² d⁻¹) to a significant uptake at high pCO₂ levels (∼ 10.9 mmol C m⁻² d⁻¹; one-way interaction: $F_{\mathrm{3},\mathrm{4}}=$ 25.1, p= 0.007; Fig. 5a). The trend in temperature response was similar for current and high-pCO₂ conditions (two-way interaction: $F_{\mathrm{3},\mathrm{16}}=$ 0.88, p= 0.47; Fig. 5a), but there was a significant shift from small efflux at lower temperatures to considerable uptakes at all temperatures under high-pCO₂ conditions (two-way CO₂ effect: $F_{\mathrm{1},\mathrm{16}}=$ 61.46, p < 0.0001; Fig. 5a). Differences between dark and light DOC fluxes under high-pCO₂ conditions were independent of temperature (two-way interaction: $F_{\mathrm{3},\mathrm{16}}=$ 1.94, p= 0.16; Fig. 5c), with the overall magnitude of influxes in the dark being significantly greater than those in the light (two-way light condition: $F_{\mathrm{1},\mathrm{16}}=$ 15.83, p= 0.001; Fig. 5c). Loss of statistically different temperature responses for high-pCO₂ light and dark responses (temperature effect, two-way interaction: $F_{\mathrm{3},\mathrm{16}}=$ 1.05, p= 0.40; Fig. 5c) was in large part due to within-treatment variability in the high-pCO₂ cores.

4 Discussion

An important component of this study was testing the interaction and individual effects of warming and OA on DOC processing. This was necessarily achieved through a comparison of core incubations occurring in different weeks. As such, it is important to consider the limitations of this approach. Control treatments in different weeks would ideally be the same in all respects, but there were some differences. For instance, NPP and R were higher in the incubation week for current pCO₂ conditions (Table S2), likely due to small changes in environmental conditions, such as salinity differences (24 versus 17.7 for current and high-pCO₂ conditions, respectively; Table 1). However, these differences did not significantly affect DOC fluxes, nor the heterotrophy of the sediments (P $/$ R = 0.84 ± 0.01 and 0.83 ± 0.04; Table S2). Moreover, sediments in separate weeks maintained the same OM content (∼ 3.5 %) and molar C : N ratio (∼ 16), suggesting that differences in processing have very little short-term impact on the overall OM pool in the sediment due to the OM pool size being about 3 orders of magnitude higher than any diel flux (organic carbon pool ∼ 12 000 mM). Thus, because all conditions in the laboratory setup were the same for each incubation (with the exception of pCO₂ in treatment tanks, which was intentionally manipulated to be different), the difference in fluxes between controls was attributed to differences in when the sediments and overlying waters were collected. Therefore, the scaling of NPP and R (Table S3) was done for the sake of treatment comparison, resulting in scaled rates within 13 % of actual measured values, which had a negligible effect on P $/$ R (< 1 % across all treatments). The final NPP and R rates in comparisons across treatments should thus be considered relative to control rates and be interpreted as approximate values (±13 %).

Understanding current ecosystem functioning is of primary interest when trying to determine how disturbances in the environment may change metabolic rates and pathways of OM mineralization (Jørgensen, 1996; D'Avanzo et al., 1996; Malone and Conley, 1996). Based on unadjusted R rates, the nearly 1 : 1 ratio of DIC production to O₂ consumption in the dark (respiratory quotient of ∼ 1.13 ± 0.05; Fig. 4) suggests that aerobic respiration dominated the sediments (Eyre and Ferguson, 2002). Similarly, unadjusted NPP rates suggest that aerobic processes dominated benthic production in the light, with a 1 : 1 ratio of O₂ and DIC fluxes (Fig. 4; Eyre and Ferguson, 2002). Sediments in the current study were net heterotrophic, with a P $/$ R in control cores of ∼ 0.84 ± 0.01 and ∼ 0.83 ± 0.04 during current- and high-pCO₂ incubation weeks, respectively. Despite the undeniable range of P $/$ R ratios that unvegetated estuarine sediments may experience (1.2 (Oakes et al., 2012) to 0.01 (Ferguson and Eyre, 2013)), the ratios in the current study were similar to mean global estimates for unvegetated estuarine sediments (∼ 0.82, calculated from values in Duarte et al., 2005) and calculated from P and R values of 22 estuaries globally (∼ 0.87, compiled by Smith and Hollibaugh, 1993), suggesting that the metabolic function of sediments in the current study are representative of estuarine sediments globally and the impacts observed in this study should be broadly applicable.

4.1 DOC fuels benthic respiration

DOC appeared to be a significant driver of benthic respiration (Fig. 5b). At control temperatures (23 ^∘C) net DOC fluxes were near zero (0.47 ± 0.93 mmol C m⁻² d⁻¹), indicating that the diel production and uptake of DOC across the sediment–water interface were balanced (Fig. 5a). The control rates in the present study were close to benthic DOC flux rates reported for subtropical estuarine sediments in most seasons, ∼ 1.5 mmol C m⁻² d⁻¹, except summer (Maher and Eyre, 2010). Relative to our control (summer) rates, Maher and Eyre (2010) reported higher net DOC flux rates (∼ 10 mmol C m⁻² d⁻¹) as a result of DOC effluxes in both the light and dark (Maher and Eyre, 2010). We observed similar light DOC effluxes (610 µmol C m⁻² h⁻¹) to summer sediments in Maher and Eyre (2010) (∼ 647 µmol C m⁻² h⁻¹), whereas dark DOC uptake was observed in the current study (−571 µmol C m⁻² h⁻¹), and Maher and Eyre (2010) reported dark DOC effluxes (254 µmol C m⁻² h⁻¹). This release of DOC in the dark was attributed to enhanced microbial coupling in the sediments in warmer temperatures (Maher and Eyre, 2010). In the current study, and in previous reports, DOC uptake suggests that bacteria not only intercepted DOC produced from within the pore waters (potentially satisfying up to 60 % of total mean bacterial production; Boto et al., 1989) but also took up available DOC from the water column to satisfy its metabolic requirements (Boto et al., 1989; Brailsford et al., 2019), effectively acting as a DOC sink. Under conditions of reduced light availability and/or intensity, sediments are expected to have an amplified heterotrophic response in addition to a reduction in microalgal production of DOC.

4.2 OA increases DOC assimilation

Positive responses in primary production were associated with OA. The ∼ 72 % increase in NPP rates at ambient temperatures was consistent with general stimulation of primary production in finer sediments with increased DIC availability (Vopel et al., 2018; Oakes and Eyre, 2014). Sediments may become DIC-limited when algal demand is relatively high compared to porewater supply of CO₂ (Cook and Røy, 2006), and MPBs therefore may benefit from an increase in CO₂ availability. MPBs in fine sediments are restricted to dissolved substrates (i.e. nutrients and DIC) accessed via diffusion from deeper and adjacent sediments, and the overlying water column (Boudreau and Jørgensen, 2001). This makes them more likely to deplete accessible DIC than MPBs in permeable sediments. Primary producers in permeable sediments, like those in reef ecosystems, therefore do not often experience the same increase in primary production with increased CO₂ (Trnovsky et al., 2016; Cyronak and Eyre, 2016; Eyre et al., 2018; Cook and Røy, 2006; Vopel et al., 2018). As well as differences in diffusive versus advective modes of solute transfer between the sediment types (Cook and Røy, 2006), variable response may also be attributable to sandier sediments being limited by other factors such as nutrient and OM availability, given that coarser sediments are generally more oligotrophic (Admiraal, 1984; Heip et al., 1995). In comparison, nutrients were non-limiting in the less permeable sediments used in the current study, based on nutrient concentrations that increased during all incubations (see the Supplement and Table S7). MPB growth rates in sediments with low permeability are more likely limited by DIC availability. Primary productivity responses to pCO₂ would likely differ in permeable sediments where general access to CO₂ is greater.

Given that MPBs exude carbon (Maher and Eyre, 2010), we would expect increased GPP to correspond with increased DOC production and flux. However, although OA stimulated primary production (Fig. 2), we instead saw increased DOC uptake in the dark (Fig. 5). A likely explanation is that bacterial uptake of DOC was stimulated through the provision of labile carbon from MPBs (Morán et al., 2011; Hardison et al., 2013). As such, DOC appeared to fuel much of the dark cycle respiration, as DOC uptake in the dark reflected dark DIC production (respiration), except for sediments at Δ−3 ^∘C under current pCO₂ conditions. Under current pCO₂ conditions, uptake of DOC in the dark accounted for only ∼ 50 % of the DIC respired in the dark. The portion of DIC accounted for by dark DOC uptake increased from 50 to 100 % under the high-pCO₂ conditions. In part, this may have been due to increased availability of labile organic carbon (Moran and Hodson, 1990) arising from the increase in NPP under high-pCO₂ conditions across all temperatures (Fig. 2b), which would reduce the need for bacteria to synthesize ectoenzymes (Chróst, 1991, 1992), resulting in more rapid turnover of carbon to the water column.

4.3 Warming drives increased heterotrophy and DOC assimilation

Sediments in this study, like other manipulative studies in both permeable sands (Lantz et al., 2017; Trnovsky et al., 2016) and cohesive sediments (Apple et al., 2006), demonstrated increased heterotrophy with increased temperature. This shift to heterotrophy is often attributed to the imbalance in the thermal sensitivity of heterotrophic over autotrophic metabolism (Yang et al., 2016; Allen et al., 2005). More specifically, differences in activation energy dictated by differences in physiology and biochemical processes (Patching and Rose, 1970; Apple et al., 2006) result in increases in heterotrophic activity with increasing temperature than exceed increases in autotrophic activity (Yang et al., 2016). However, in this study, under current pCO₂ conditions, the increases in R and GPP from Δ−3 ^∘C to control temperatures were similar (∼ 16 and ∼ 11 %, respectively), whereas at higher temperatures GPP decreases far exceeded increases in R (7× and 3× for 26 and 28 ^∘C, respectively). Therefore, unlike previous studies, decreases in MPB productivity at higher temperatures appeared to be a greater driver towards heterotrophy than increases in respiration rates. In other words, temperature increases not only increased the rate of DOC uptake but also likely decreased the rate of DOC production.

4.3.1 Warming reduces GPP and DOC production under current pCO₂ conditions

Primary production is the main source of DOC in marine ecosystems (Wagner et al., 2020). Decreasing trends in GPP with warming under current pCO₂ conditions seen here have been described previously where photosynthetic growth and production decline at higher temperatures (Thomas et al., 2012). Photosynthetic productivity is often linked to seasonal temperature (Apple et al., 2006), which is also associated with differences in environmental factors such as light, nutrient concentrations, and dissolved organic matter (DOM) quality and availability (Geider, 1987; Herrig and Falkowski, 1989). Although the relative availability of light and nutrients does influence productivity rates (Kana et al., 1997) and would be expected to influence in situ seasonal production, the current study controlled light and initial nutrient concentrations in the water column to isolate the effect of temperature. Thus, decreasing GPP was driven by warming, suggesting that MPBs in these subtropical sediments likely had a temperature optimum around current mean summer temperatures of ∼ 23 ^∘C (GPP: 1515 ± 37 µmol O₂ m⁻² h⁻¹; Fig. 2c). Longer-term warming could allow for possible migration of more tolerant species to settle from lower latitudes (Hallett et al., 2018), shifting the composition of the benthic community. The introduction of more tolerant species could reduce the increase in heterotrophy and net DOC removal from the water column seen here. However, the species diversity of the estuarine sediments will ultimately decrease as they are pushed to temperature extremes (Thomas et al., 2012), reducing the functional redundancy of the microbial community. This decreased functional redundancy has the potential to make unvegetated estuarine sediments less resilient to environmental perturbations under future climate conditions.

4.3.2 Warming increases respiration and DOC assimilation

Unlike photosynthetic productivity, heterotrophic respiration often has a linear rate increase with temperature to the thermal optimum due to heterotrophs not being constrained by the same abiotic variables (e.g. nutrient and light availability) as primary producers (Apple et al., 2006, 2008; Geider, 1987; Yap et al., 1994). In this study, respiration rates under both current and high pCO₂ conditions increased from the lowest rates measured at Δ−3 ^∘C to maximum rates (> 50 % greater) at Δ+5 ^∘C (Fig. 2a). Consistent with overall lower respiration rates relative to other subtropical unvegetated sediments (∼ 900 to ∼ 1500 µmol O₂ m⁻² h⁻¹; Ferguson and Eyre, 2013), the temperature dependence of respiration under both current and high-pCO₂ conditions (Q₁₀= 1.66 and 1.69, respectively) was slightly lower than is typical for biological systems (Q₁₀= 2; Valiela, 1995) but similar to temperature dependence described in other estuarine systems (Q₁₀= 1.5–1.9; Morán et al., 2011), with values towards the lower end of this range possibly being a result of resource limitation (López-Urrutia and Morán, 2007).

A potential limiting resource for bacteria in estuarine sediments is DOM (Church, 2008), ultimately controlling the flow of carbon through the microbial loop (Kirchman and Rich, 1997). However, in the presence of sufficient DOM, warming has been associated with increased bacterial DOM incorporation (Kirchman and Rich, 1997). In line with this, increased uptake of DOC at higher temperatures and efflux at lower temperatures was observed. Although DOC is mainly produced by photoautotrophs, DOC can be produced in the dark through, for example, chemodegradation of detrital organic carbon and cell lysis by viruses and during grazing (Carlson, 2002). As such, the efflux of DOC in the dark at Δ−3 ^∘C suggests that heterotrophic bacterial productivity, and therefore DOC uptake, was reduced by lowered temperatures (Raymond and Bauer, 2000), resulting in a failure to intercept all DOC produced in the pore waters. This failure to intercept DOC may be compounded if nutrient supply is limited (Brailsford et al., 2019), as it is common for heterotrophic bacteria to rely on refractory DOC when labile sources are not readily available (Chróst, 1991), which can occur under conditions of nutrient-limited biological productivity (Allen, 1978).

4.3.3 Global estuarine loss of DOC from unvegetated sediments in the future

Up to 3.5× more DOC reaches the ocean interior from coastal areas than the open ocean (Duarte, 2017; Krause-Jensen and Duarte, 2016; Hansell et al., 2009). As such, small changes to the coastal export of DOC may have a disproportionately large influence on the global DOC budget. Our findings suggest a reduced export of DOC to the ocean from the coastal zone under high-pCO₂ conditions, across the full 8 ^∘C temperature range due to changes in carbon processing within unvegetated sediments. Despite the lack of seasonality in the study, the inclusion of an 8 ^∘C temperature range, including temperatures below current mean temperatures, suggests that seasonal temperature variation is unlikely to have a significant effect on the relative change in DOC in the future (Fig. 5). Although any upscaling of a single controlled experiment to a global scale is highly speculative, we believe it is valuable to demonstrate the potential for a high-pCO₂ climate to cause globally significant change in DOC export from coastal zones. Furthermore, putting our findings in a global context provides a guideline value for potential change. The following estimates should be considered in this context, and it should be expected that different hydrodynamic settings, sediment and/or sediment community composition, and sources of organic matter could affect the outcome. For example, the response to warming and pCO₂ may be different for pelagic communities and/or in deeper waters that are subject to stratification (Li et al., 2020), where access to nutrients and CO₂ may become limiting (Rost et al., 2008). We have applied our results to global coastal DOC exports (Maher and Eyre, 2010; Duarte, 2017) as an initial step in estimating responses of unvegetated sediment habitats to future high-pCO₂ climate. We do not assume that the responses of unvegetated sediments to the future climate found here are applicable to other ecosystems dominated by macrophytes, and we thus did not apply our findings to vegetated coastal habitats.

To estimate total DOC export from coastal zone in a future high-pCO₂ climate of Δ+3 ^∘C and OA, the sediment uptake rate of 19 ± 4 mmol C m⁻² d⁻¹ was scaled to the global surface area of unvegetated estuarine sediments (1.8×10¹² m²; Costanza et al., 1997). On this basis, an estimated 150 Tg C would be removed from the coastal zone by unvegetated estuarine sediments annually under OA conditions with an accompanying 3 ^∘C temperature increase. To then calculate the potential impact of this uptake, we applied our estimates to existing future global coastal DOC export estimates (Maher and Eyre, 2010; Duarte, 2017). Mean benthic DOC export from estuaries, including intertidal and vegetated habitats, has been estimated at 168 Tg C yr⁻¹ (90–247 Tg C yr⁻¹) (Maher and Eyre, 2010). Under this scenario, the switch to DOC uptake by sediments under future climate conditions (Fig. 5a) would result in ∼ 90 % reduction in total mean estuarine DOC export (Maher and Eyre, 2010), decreasing the load from ∼ 168 to ∼ 18 Tg C yr⁻¹. Other global estimates of DOC exported from coastal vegetated ecosystems range from 114 up to 1853 Tg C yr⁻¹ (Duarte, 2017), with scaled estimates suggesting unvegetated estuarine sediments may consume 8 to 132 % of this DOC under a future high-pCO₂ climate. As such, this basic upscaling suggests that, by impacting DOC fluxes in unvegetated sediments, future climate conditions have the potential to significantly impact global DOC export from coastal systems to the open ocean. This has implications for global marine productivity and carbon transfer to the ocean interior (Krause-Jensen and Duarte, 2016). However, to get a more accurate insight into global carbon cycling, the response of DOC export from estuarine vegetated habitats to future climate also needs to be studied.