Sediment Stocks of Carbon, Nitrogen, and Phosphorus in Danish Eelgrass Meadows

Seagrass ecosystems provide an array of ecosystem services ranging from habitat provision to erosion control. From a climate change and eutrophication mitigation perspective, the ecosystem services include burial and storage of carbon and nutrients in the sediments. Eelgrass (Zostera marina) is the most abundant seagrass species along the Danish coasts, and while its function as a carbon and nutrient sink has been documented in some areas, the spatial variability of these functions, and the drivers behind them, are not well understood. Here we present the first nationwide study on eelgrass sediment stock of carbon (Cstock), nitrogen (Nstock), and phosphorus (Pstock). Stocks were measured in the top 10 cm of eelgrass meadows spanning semi-enclosed estuaries (inner and outer fjords) to open coasts. Further, we assessed environmental factors (level of exposure, sediment properties, level of eutrophication) from each area to evaluate their relative importance as drivers of the spatial pattern in the respective stocks. We found large spatial variability in sediment stocks, representing 155-4413 g C m-2, 24-448 g N m-2, and 7-34 g P m-2. Cstock and Nstock were significantly higher in inner fjords compared to outer fjords and open coasts. Cstock, Nstock, and Pstock showed a significantly positive relationship with the silt-clay content in the sediments. Moreover, Cstock was also significantly higher in more eutrophied areas with high concentrations of nutrients and chlorophyll a (chl a) in the water column. Conversely, silt-clay content was not related to nutrients or chl a, suggesting a spatial dependence of the importance of these factors in driving stock sizes and implying that local differences in sediment properties and eutrophication level should be included when evaluating the storage capacity of carbon, nitrogen, and phosphorus in Danish eelgrass meadows. These insights provide guidance to managers in selecting priority areas for carbon and nutrient storage for climate- and eutrophication mitigation initiatives. (Less)


INTRODUCTION
The coastal ocean is a highly dynamic component of the Earth's system and features many complex interactions between oceanic, terrestrial, and atmospheric processes. The central role of coastal ecosystems in global and local biogeochemical cycles is in part due to disproportionately high biological activity in conjunction with extensive inputs of carbon and nutrients from rivers, terrestrial runoff, and upwelling (Gattuso et al., 1998;Duarte et al., 2005). Seagrasses are pervasive Frontiers in Marine Science | www.frontiersin.org 1 December 2018 | Volume 5 | Article 474 benthic macrophytes in coastal and estuarine waters spanning all continents except Antarctica (Short et al., 2007). In addition to providing important habitat for a wide range of species, seagrasses also play important roles in biogeochemical cycles including carbon and nutrient cycling (Duarte et al., 2005;McGlathery et al., 2007). Seagrasses, along with other marine vegetated habitats such as mangroves, saltmarshes, and macroalgae, have been recognized as a potential solution to combat anthropogenic perturbations such as excessive emissions of greenhouse gasses and nutrients to the environment (Nordlund et al., 2016;Gattuso et al., 2018;Himes-Cornell et al., 2018). The ability of marine vegetated habitats to sequester and store atmospheric CO 2 was first introduced by Smith (1981), and the term Blue Carbon was defined in the Blue Carbon report by Nellemann et al. (2009). This research field has propelled over the last decade and Blue Carbon has been proposed as a viable tool in both climate mitigation-and adaptation strategies (McLeod et al., 2011;Duarte et al., 2013;Greiner et al., 2013;Luisetti et al., 2013). Similarly, the potential of these habitats to act as nutrient filters has been suggested to alleviate effects of local eutrophication (McGlathery et al., 2007;Piehler and Smyth, 2011). Due to a combination of high productivity and the ability to trap sestonic particles, seagrasses produce and bury large amounts of organic matter (OM) of which a portion is stored in the often anoxic underlying sediment and effectively sequestered from the ocean-atmosphere carbon pool (Fourqurean et al., 2012a). The sediment stock has been shown to be the predominant contributor to the total pool of organic carbon in seagrass meadows (Kennedy et al., 2010;Fourqurean et al., 2012a;Macreadie et al., 2014). In addition, a fraction of the seagrassderived detritus and dissolved organic matter (DOM) is exported from the meadow and may either be buried in depositional basins and the deep sea, or enter the microbial loop as refractory DOM (Barrón et al., 2014;Duarte and Krause-Jensen, 2017).
Nutrient cycling and retention in seagrass meadows is regulated both directly through the uptake and assimilation in leaves, roots, and rhizomes and indirectly through trapping of sestonic particles containing variable amounts of nutrients (Touchette and Burkholder, 2000;Eyre and Ferguson, 2002;Romero et al., 2006). In addition, elevated rates of denitrification carried out by microbes associated with the rhizosphere have been shown to contribute to additional export of nitrogen at the community and ecosystem scale (Eyre et al., 2016;Reynolds et al., 2016), although nitrogen fixation in seagrass meadows may also be significant Welsh, 2000). However, the size and variability of the nutrient sink capacity of seagrass meadows, the mechanisms governing the variability, and the temporal aspects of the retention are not well constrained Nielsen et al., 2004;McGlathery et al., 2012).
Several studies have investigated the drivers behind Blue Carbon storage, and it has been shown that they consist of a complex set of physical, chemical, and biological factors that vary considerably between and within systems and are highly dependent on which functional scale is assessed (e.g., spatial, temporal, species, community, ecosystem) (see e.g., Ricart et al., 2015;Röhr et al., 2016;Samper-Villarreal et al., 2016;Belshe et al., 2017a,b;Dahl, 2017;Oreska et al., 2017a;Mazarrasa et al., 2018;Röhr et al., 2018). Sediment properties such as grain size, porosity, degree of sorting, and density are largely controlled by the hydrodynamic regime where, for instance, grain size typically increases with increasing wave-and current exposure, as the higher hydrodynamic energy and short residence times do not allow finer grains to settle (Folk and Ward, 1957;Fonseca and Bell, 1998;van Keulen and Borowitzka, 2003;Mazarrasa et al., 2017). Contrarily, low-energy environments typically observed in sheltered locations, tend to have a larger fraction of fine-grained material such as silt and clay (Jankowska et al., 2016;Röhr et al., 2016). Thus, information on sediment characteristics can serve as a useful proxy for the level of exposure in unvegetated areas, but the relationship may be confounded by water depth and in vegetated sediments by the seagrass canopy and rhizosphere's effect on sedimentation and resuspension (Fonseca and Bell, 1998;Madsen et al., 2001;Yang et al., 2008;Linders et al., 2018). Moreover, sediment grain size is a strong predictor of sediment OM pools and seagrass sediments exhibit increases in both the amount of carbon stored Röhr et al., 2016;Samper−Villarreal et al., 2016;Miyajima et al., 2017) and carbon burial rates (Mazarrasa et al., 2017) with increasing silt-clay content. This correlation may be due to the larger surface area for adsorbing organic molecules at high silt-clay content (Keil et al., 1994;Bergamaschi et al., 1997), and low oxygen availability reducing the remineralization of OM (Dauwe et al., 2001;Serrano et al., 2016), as well as the relatively higher sedimentation of sestonic particles in these low-energy environments (Mazarrasa et al., 2017(Mazarrasa et al., , 2018. However, plant traits such as shoot density, above-and belowground biomass, and net primary productivity also explain a large part of the variability in carbon stock and sink capacity in local areas Gillis et al., 2017;Oreska et al., 2017a). The relative importance of biological and environmental factors seem to largely depend on which spatial scale is assessed, with biological factors being more important on smaller scales (meadow, patch) and environmental factors more important on larger scales (regional, global) (Dahl, 2017).
In the present study, we made a concerted effort to assess sediment stocks of carbon (C stock ), nitrogen (N stock ), and phosphorus (P stock ) in eelgrass (Zostera marina L.) meadows across Denmark and evaluated environmental properties as explanatory variables for the variability in sediment stocks. Eelgrass is a prevalent seagrass species which, although it has been in drastic decline, creates vast meadows in areas along the Danish coasts (Boström et al., 2014). Despite receiving extensive attention and research efforts in Denmark for over a century, the role of eelgrass in carbon and nutrient storage is not well known. To our knowledge, this study is the first to assess species-specific seagrass sediment storage of carbon, nitrogen, and phosphorus on a nationwide scale ranging from semi-enclosed estuaries to the open coast. Given the large spatial range in sampling areas, we proposed the following research questions: (i) what is the spatial variability in sediment C stock , N stock , and P stock in Danish eelgrass meadows? (ii) To what extent is stock size related to sediment characteristics and/or level of eutrophication? Based on these questions we formulated the following hypotheses: (1) C stock , N stock , and P stock are higher in sheltered areas with less Frontiers in Marine Science | www.frontiersin.org 2 December 2018 | Volume 5 | Article 474 exposure and higher silt-clay content and (2) stocks are higher in areas with higher nutrient and chl a concentrations in the water column.

Assessing C, N, and P Stocks
To assess the sediment stocks of each element (C stock , N stock , and P stock ), particulate organic carbon [POC % dry weight (% dw)], total nitrogen (TN % dw), and total phosphorus (TP % dw) concentrations in the top 10 cm sediment were compiled from six studies conducted across Denmark in the period 1999-2014 (Supplementary Table S1). In all studies, sediment cores (1-4 cores per site) were collected in eelgrass meadows using similar techniques and subsequently analyzed for elemental composition and sediment properties. Sediment carbon density (g cm −3 ) was calculated by multiplying dry bulk density [DBD (g cm −3 ) with POC (% dw)] according to Howard et al. (2014). The sediment density of nitrogen and phosphorus was calculated in the same way but using TN and TP, respectively. In studies where OM [defined as loss on ignition, LOI (%)] was measured instead of POC, a conversion was performed using the empirical linear relationship between POC and OM found in studies measuring both parameters (Supplementary Figure S1). For sites where a full depth profile was obtained (0-10 cm), the accumulated amount of C, N, and P was assessed by multiplying the average (n = 1-4 cores) sediment density of carbon, nitrogen, and phosphorus by the length of each core section (cm) and integrating across the 10 cm depth to express as mass per area (g m −2 ). The total number of sites, where all information needed to assess sediment stock was available, equaled n = 50 for C stock , followed by n = 47 for N stock , and n = 36 for P stock . Maps of C stock , N stock , and P stock were created using the MATLAB R R2017b mapping toolbox (Mathworks, Inc.) and projected using the WGS1984 ellipsoid UTM zone 32N.

Environmental Parameters
To assess the level of exposure, a wave energy layer covering the sampling sites was obtained from the EUSeaMap 1 with a resolution of 300 m × 300 m. The wave energy layer is based on mean wind data in the period (Wijkmark and Isaeus, 2010. The process of eutrophication is defined as "an increase in the rate of supply of organic matter to an ecosystem" (Nixon, 1995). However, several alternative definitions exist that concern the various proxies from which the level of eutrophication can be inferred (Cloern, 2001;Andersen et al., 2006). These include water column concentrations of macronutrients, chl a, and Secchi depth . In order to assess relationships between sediment stocks and level of eutrophication we extracted data from the Danish National Monitoring and Assessment Program for the Aquatic and Terrestrial Environments (DNAMAP). Due to the long temporal range of our sediment data (1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014), we used annual averages 1 www.EMODnet.eu of water column concentrations of TN, dissolved inorganic nitrogen (DIN), TP, dissolved inorganic phosphorus (DIP), chl a, salinity, and Secchi depth from each site for the 5 years prior to the respective sampling events (Wulff et al., 1990;Christiansen and Emelyanov, 1995) (Supplementary Data Sheet S1). Based on levels of chl a and Secchi depth at the Danish intercalibration sites of the Water Framework Directive (WFD), all sites in this study can be considered eutrophic (Henriksen, 2009). Thus, we discuss various levels of eutrophication when referring to water column TN, TP, chl a, and Secchi depth.
To obtain information on the effects of sediment fine-grained material on sediment stocks, and to explore the relationship between level of exposure and fine-grained material, we used siltclay content (% dw, < 62.5 µm) reported in each study. We used the average silt-clay content of the top 10 cm at all sites except for 10 sites, where only surface samples (0-2 cm) were available.

Statistical Analyses
The sediment cores were sampled in different waterbodies, which we qualitatively categorized into three types (inner fjord, outer fjord, and open coast) based on visual inspection of the sites' respective location (Figure 1). Out of the 60 sampling sites, 34 were categorized as inner fjord, 18 as outer fjord, and 8 as open coast. The sites differed in water depth ranging from 0.5 to 12 m, and we performed a one-way ANOVA to assess whether mean water depth differed between waterbodies. Similarly, we performed one-way ANOVAs to assess whether the silt-clay content and exposure differed between waterbodies.
To assess whether sediment stocks in Danish eelgrass meadows depend on waterbody type, we performed oneway ANOVAs. C stock and N stock were log-transformed to meet assumptions of normality and homogeneity of variances, obtained using Shapiro-Wilk's test of normality and Levene's test of homoscedasticity, respectively. Tukey's HSD post hoc tests were performed to test the pairwise difference in sediment stocks between waterbodies.
To examine environmental effects on sediment stock size, we constructed generalized linear models (GLM) for each sediment stock (C stock , N stock , and P stock ). Stocks were modeled as a function of silt-clay, TN, TP, chl a, Secchi depth, salinity, and exposure, using the Gaussian family with identity links. The model setup was chosen following evaluation of Akaike Information Criterion (AIC), likelihood-ratio test, and residual deviance. We also assessed the multicollinearity between explanatory variables by using variance inflation factors (VIF). In order to constrain the predictive power of each explanatory variable, we standardized all variables to z-scores by subtracting the mean from each observation and dividing by its standard deviation.
A caveat of using GLMs was that due to missing values at some sites, the model runs were applied to a reduced number of observations and thus did not cover the entire range of sites. Specifically, the GLMs for C stock , N stock , and P stock were only applied to 27, 29, and 24 sites, respectively, and the relative number of waterbodies (inner/outer fjord, open coast) was not always proportional to the full dataset. Consequently, we also assessed the relationships between sediment stocks and environmental variables through simple linear regression analyses of one explanatory variable at a time. Similarly, we assessed the linear correlation between C stock , N stock and P stock . Considering our aim to describe variation across all sites, we focus on these results and primarily use the GLMs to explore the predictive power of environmental factors for a subset of sites.
All statistics were performed with R statistical software (R Core Team, 2017) using the following packages: basic, boot, bbmle, car, and lmtest. A significance level of α < 0.05 was used for all parametric tests.

Eelgrass Sediment Stocks of C, N, P in Denmark
Concentrations (% dw) of sediment POC, TN, and TP exhibited significant heterogeneity ranging from 0.005-11.33, <0.001-1.13, and 0.007-0.17, respectively. Mean (± SE) C stock , N stock , and P stock inferred from the top 10 cm sediment equaled 1013 ± 116, 109 ± 12, and 17 ± 1 g m −2 , respectively. The spatial variability was large and ranged 155-4413, 24-448, and 7-34 g m −2 , for C stock , N stock , and P stock , respectively (Figure 2). Furthermore, median stocks of C and N were substantially lower than mean values (711 and 68 g m −2 , respectively), reflecting a lognormal distribution with a majority of lower stocks and a few very high stocks.
We identified significant relationships between waterbody (inner or outer fjord, open coast) and sediment stock for C stock (ANOVA, F 2,47 = 9.39, p < 0.001) and N stock (ANOVA, F 2,44 = 4.10, p = 0.02) but not for P stock (ANOVA, F 2,33 = 0.42, p = 0.52). C stock was significantly higher in inner fjord locations than in outer fjords (post hoc, p < 0.001) whereas N stock was significantly higher in inner fjords compared to open coast (post hoc, p = 0.03) (Figure 3 and Supplementary  Table S2). Neither silt-clay content nor level of exposure differed significantly between waterbodies, but the water depth of the sampling sites differed significantly between waterbodies (ANOVA, F 2,43 = 5.76, p = 0.006) and was on average significantly higher in outer fjords compared to the inner fjords (post hoc, p = 0.006), reflecting that the depth range of eelgrass is deeper in outer fjords (Riemann et al., 2016).

Sediment Properties
The silt-clay content in the sediment revealed a positive linear relationship with C stock , N stock , and P stock with the highest slope for C stock and highest R 2 value (0.30) for N stock (Figure 4 and Supplementary Figures S4, S5).
Silt-clay content was not significantly correlated with exposure (p = 0.29). Similarly, exposure was not significantly related to C stock (p = 0.37), N stock (p = 0.41), or P stock (p = 0.15).
C stock displayed a significant positive correlation with both water column nutrient concentrations and chl a, and a significant negative correlation with Secchi depth (Figure 4). In contrast, N stock and P stock displayed no significant  correlations with any water column parameter (Supplementary Figures S4, S5). Furthermore, no significant correlation was observed between any of the water column parameters and silt-clay content, except for Secchi depth which exhibited a significant negative correlation with silt-clay content, displaying the highest silt-clay content at low Secchi depths (Supplementary Table S3). Elemental ratios in the sediment were not related to water column nutrients, except that sediment N:P ratio was positively related to TP concentration in the water column (p = 0.03). The variation in sediment elemental ratios was highest in C:N, followed by C:P and N:P as indicated by their respective coefficients of variation (CV = standard deviation/mean) of 0.84, 0.58, and 0.57, respectively (Supplementary Table S4).

Model Predictions
The three GLMs, based on a subset of sites for which all variables were measured, overall supported the patterns of the   variable-by-variable linear regression analyses. Consistent for all three sediment stocks was that silt-clay content was significantly positively related to stock size ( Table 1). In contrast to the simple linear regression analyses, however, chl a was the only other variable that was significantly related to C stock (p = 0.01) and Secchi depth was significantly related to N stock (p = 0.005).
Overall, the GLM predictions performed well and yielded R 2 values of 0.78, 0.66, and 0.71 for C stock , N stock , and P stock , respectively (Figure 5).

DISCUSSION C, N, P Stocks in Danish Eelgrass Sediments
We observed large spatial variability in sediment stocks of carbon, nitrogen, and phosphorus where sheltered, inner fjord locations generally comprised larger stocks whereas outer fjords and open coasts exhibited stocks that were sometimes orders of magnitude lower. These findings demonstrate the necessity of considering spatial variability when assessing regional or national sediment stocks and provide guidance into which eelgrass areas may be of particular importance as reservoirs of these elements.

Sediment Characteristics and Exposure as Explanatory Variables for the Variability in Eelgrass Sediment Stocks
As hypothesized, we observed strong positive relationships between fine-grained material (silt-clay) in the sediment and C stock , N stock , and P stock , which is in line with previous studies on eelgrass sediment carbon content (e.g., Dahl et al., 2016;Röhr et al., 2016;Mazarrasa et al., 2017;Oreska et al., 2017a;Röhr et al., 2018). In general, the grain size composition of marine sediments is controlled by the hydrodynamic regime and is often considered a useful proxy for wave-and current exposure (Yang et al., 2008;Cabaço et al., 2010;Mazarrasa et al., 2017). In absence of in situ current or exposure measurements, we utilized the publicly available EUSeaMap exposure layer. Notably, we did not find a relationship between level of exposure and silt-clay content in the sediment. This could be due to the relatively coarse resolution of the exposure layer used (300 m × 300 m) compared to the point measurements that the sampled cores constitute. The exposure layer also does not address temporal differences in exposure between sites, being based on mean wind conditions for a 5-year period between 2002 and 2007. Moreover, the level of exposure is strongly depth-dependent where the generally shallower sampling depths in inner fjords in this study represent a higher exposure level than the deeper sampling depths in outer fjords and open coast locations. However, we did not observe a significant relationship between water depth and silt-clay content or sediment stocks. The lack of correlation between silt-clay content and exposure may also be due to the attenuating effect of the eelgrass canopy on hydrodynamic energy which allows for fine-grained particles to settle, in addition to the stabilizing effect of the rhizosphere on limiting resuspension (Fonseca et al., 1982;De Boer, 2007;Hansen and Reidenbach, 2012). Lastly, considering that our sites ranged from semi-enclosed estuaries to open coast, variable patterns in sediment discharge from runoff may explain the lacking correlation, where low-exposure sites are often located in areas of higher riverine and/or terrestrial sediment input compared to the open coast and thereby offsetting the role of eelgrass.

Eutrophication Effects on Eelgrass Sediment Stocks
The literature is equivocal regarding the effects of eutrophication on seagrass carbon storage and burial, reporting both higher storage and burial due to an increased input of sestonic particles (Gacia et al., 2002;Mazarrasa et al., 2017;Samper-Villarreal et al., 2017), and decreased storage and burial due to detrimental effects on seagrass productivity and survival (Macreadie et al., 2012;Jiang et al., 2018). Here, the level of eutrophication appeared to increase C stock (Figures 4C-F and Table 1), suggesting that Danish eelgrass meadows have a relatively high reliance on allochthonous carbon (Oreska et al., 2017b). In fact, for 10 of the sites used in this study, the eelgrass net primary production only explained 2.3% of the variation in C stock indicating a low contribution of autochthonous carbon (Röhr et al., 2016). Moreover, the positive relation between C stock and silt-clay content and turbidity (Secchi depth) observed in this study suggests that the sources of carbon in the eelgrass sediment are largely of allochthonous origin (Samper−Villarreal et al., 2016). This was also demonstrated for 10 of the sites where analysis of stable isotopes (δ 13 C and δ 15 N) and C:N ratios in sediments and in carbon sources was carried out (Supplementary Table S1). An isotope mixing model suggested that the carbon content in eelgrass sediments was largely comprised (up to 78%) of sestonic particles (e.g., phytoplankton and macroalgae), further pointing to the influence of allochthonous carbon on C stock (Röhr et al., 2016;Kindeberg et al., unpublished). Unexpectedly, N stock or P stock were uncorrelated with water column nutrient concentrations ( Table 1 and Supplementary  Figures S4, S5). This could be due to the lower sample size for these stocks compared to C stock which thus did not cover the entire range of eutrophication levels (TN: 248-1123 µmol L −1 ; TP: 21-75 µmol L −1 ; chl a: 2-6 µg L −1 ). Moreover, it may reflect the different burial efficiencies and remineralization rates for the three elements due to differing nutrient content of the sediment OM, where OM rich in nutrients (low C:N and C:P ratios) is decomposed faster and to a larger extent than nutrient-depleted OM (Kristensen and Hansen, 1995;Banta et al., 2004;McGlathery et al., 2007). It is difficult to ascertain the processes affecting relationships between sediment stocks and C:N:P ratios as elemental ratios of OM (e.g., eelgrass leaves and phytoplankton) are affected both by nutrient supply and light availability, and a confounding effect of eutrophication on increased nutrients and associated decreased light on productivity may also confound the effects of OM decomposition (Zimmerman et al., 1987;Grice et al., 1996). Overall, we observed a quite large variability in C:N:P ratios across the different sites (Supplementary Table S4), illustrating the spatial difference in sources and/or age of OM and nutrient availability as well as sediment biogeochemistry (Fourqurean et al., 1997;Kennedy et al., 2010). Nutrient dynamics in eelgrass sediments are comprised of many complex processes, some of which are directly governed by uptake of dissolved nutrients in the water column and sediment porewaters and others that are facilitated by the rhizosphere's microbial communities . These dynamics may explain the lack of correlation between water column nutrients and N stock and P stock . Pedersen and Borum (1993) estimated a total annual uptake rate of nitrogen of 25.3-50 g m −2 year −1 in an eelgrass meadow in Øresund, Denmark, of which 38-64% was attributed to the roots' uptake from sediment porewater. If uptake into eelgrass biomass is high and this biomass is subsequently exported out of the system (e.g., during storms), this could account for a significant loss term. Depending on the fraction of the plant detritus that is subsequently deposited in the sediment compared to exported elsewhere (Heck et al., 2008;Duarte and Krause-Jensen, 2017), this could be a factor which may explain the discrepancy between available nutrients in the water column and nutrient stocks in the sediment.
In addition, denitrification and anammox are both loss processes that convert inorganic nitrogen into nitrogen gas (N 2 , N 2 O, or NO) which escapes the sediment, and have been shown to be a major contributing process to nitrogen loss in seagrass habitats (Eyre et al., 2016;Reynolds et al., 2016). Rates of eelgrass-associated denitrification differ between studies (3.2-10.1 mg N m −2 d −1 ) but the net effect on the sediment nitrogen fluxes should also consider nitrogen fixation, which can be significant (1.2-6.5 mg N m −2 d −1 ) (see e.g., Flindt, 1994;McGlathery et al., 1998;Risgaard-Petersen and Ottosen, 2000;Welsh, 2000;Welsh et al., 2000;Cole and McGlathery, 2012;Russell et al., 2016). Although dependent on several factors and far from straightforward, denitrification rates typically increase and nitrogen fixation rates decrease in seagrass sediments with increasing eutrophication (Seitzinger and Nixon, 1985;Howarth et al., 1988;Welsh, 2000;Seitzinger et al., 2006;Murray et al., 2015;Asmala et al., 2017). This could be a potential explanation as to why N stock did not mirror the nutrient concentrations in the water column.
Phosphorus cycling in marine sediments is governed by a complex set of biogeochemical processes and is distributed into various pools, where phosphorus loosely adsorbed or bound to iron oxides (Fe-P complexes) are generally the largest pools (Ruttenberg, 1992). In carbonate-rich areas, phosphorus is also often tightly bound to calcium carbonate (CaCO 3 -P complexes) whereas in siliciclastic muddy sediments, phosphorus is considered relatively more available due to the effects of coupled sulfate-and iron-reduction on releasing ironbound phosphorus (Canfield et al., 1993;Jensen et al., 1995). The effects of eutrophication on increasing sediment organic content may stimulate sulfate reduction thus releasing ironbound phosphorus (Holmer et al., 2006). As opposed to its tropical seagrass counterparts, eelgrass growing in fine-grained, siliciclastic sediments is rarely limited by phosphorus, although contrasting evidence exists (Murray et al., 1992;Holmer et al., 2006). Eelgrass may even excrete some phosphorus into the water column (Brix and Lyngby, 1985), and depending on the rate and magnitude of this process, eelgrass may facilitate an efflux of phosphorus out of the sediment to the water column. Taken together, these two processes may mask the effects of water column nutrient concentrations on P stock and could, in addition to export of eelgrass biomass, explain the lack of correlation in our dataset.
Notably, neither water column nutrients nor chl a were significantly correlated with silt-clay, suggesting a spatial dependence in the explanatory power of these variables and it is therefore difficult to attribute the variability in stocks solely to either factor. Both eutrophication and too high silt-clay content can have negative effects on eelgrass subsistence and the positive effects on stock size observed here are likely to diminish when reaching a certain threshold (De Boer, 2007;Viaroli et al., 2008). Eelgrass is, as seagrasses in general, sensitive to turbidity and organic-rich sediments due to high light requirements and sensitivity to highly reducing sediments (Goodman et al., 1995;Holmer and Bondgaard, 2001;Kemp et al., 2004;Ochieng et al., 2010). Once a tipping point is exceeded and the eelgrass is lost, the system can enter an alternative stable state where sediments are eroded, and the stored carbon, nitrogen, and phosphorus can re-enter the water column (Pendleton et al., 2012;Arias-Ortiz et al., 2017;Moksnes et al., 2018). Depending on the bioavailability of the eroded nutrients, this efflux can result in net CO 2 emissions and fuel a positive feedback loop on water quality deterioration. Further investigation into this complex balance is warranted and in particular the observed decoupling between the effects of eutrophication and sediment silt-clay content on a larger spatial scale should be explored further.

Danish Eelgrass Sediment Stocks in a Broader Seagrass Context
Due to the relatively shallow sediment depth (10 cm) used in this study, it is difficult to directly compare our results with those from other studies which usually assess sediment stocks based on deeper cores [e.g., ≥100 cm proposed as the practicable standard for carbon stock assessments (Howard et al., 2014;Emmer et al., 2015)] (Fourqurean et al., 2012a,b; Frontiers in Marine Science | www.frontiersin.org 8 December 2018 | Volume 5 | Article 474 Miyajima et al., 2015). However, when comparing the mean POC concentration (1.1 ± 0.07%) found in our study, Danish eelgrass sediments fall within the lower end of the global range of 0-48.2% and below the global mean of 2.5 ± 0.1% for seagrasses (Fourqurean et al., 2012a). Nevertheless, the identified POC concentrations are well above the reported estimates for Z. marina in the southern Baltic Sea (0.14%, Jankowska et al., 2016) and similar to the recently reported global average for Z. marina of 1.4 ± 0.4% (Röhr et al., 2018). It should be noted, however, that mean POC concentration alone does not provide information of the total accumulated carbon in the sediment (C stock ), which requires information on sediment density and integration over sediment depth (Howard et al., 2014). Similar to our study, Jankowska et al. (2016) measured C stock in the top 10 cm, sampled in eelgrass meadows in a sheltered inner bay, an exposed outer bay, and a shallow open coast area in Puck Bay, Poland. In accordance with our results, the highest C stock values were found in the sheltered, inner location (228.0 ± 11.6 g C m −2 ) and lowest in the exposed, outer location (50.2 ± 2.2 g C m −2 ). However, these values are about an order of magnitude lower than the average C stock we observed in respective waterbodies in Denmark (Inner fjord: 1357 ± 173 g C m −2 , Outer fjord: 525 ± 117 g C m −2 ; Supplementary Table S2), which could partly be due to a higher fraction of larger grain sizes (mean ϕ = 1.68-2.39) in Puck Bay. Postlethwaite et al. (2018) measured C stock in three eelgrass meadows in Clayoquot Sound, BC, Canada, but integrated over deeper sediment depths (34 cm on average) compared to our study. They found C stock to range from 820 ± 26 g C m −2 in areas with sandier sediments closest to the coast and up to 2099 ± 365 g C m −2 near the river mouth. Considering that these C stock measurements were obtained from more than three times as deep cores than used in our study, their values are generally lower than presented here, as their carbon depth profiles were relatively stable with depth. This difference could partly be due to the more pristine conditions with low anthropogenic influence, high Secchi depth (4.8-8.8 m), and likely lower sediment discharge as compared to most of the sites in Denmark (Postlethwaite, 2018;Postlethwaite et al., 2018). A similar study was carried out south of British Columbia, in Puget Sound, Washington, DC, United States. Poppe and Rybczyk (2018) found eelgrass C stock to range from 1180-1900 g C m −2 , averaging at 1420 ± 110 g C m −2 , but also here the integrated core depth was three times as deep (30 cm) as in our study. A recent comprehensive assessment of C stock (integrated over 25 cm sediment depth) across the distributional range of Z. marina further established the spatial variability, which ranged from 318 ± 10 -26523 ± 667 g C m −2 , and averaged 2721 ± 989 g C m −2 (Röhr et al., 2018). In line with the results presented here, silt-clay content was found to be the most important predictor explaining 53% of the variation in C stock . Among the largest observed stocks were typically found in sheltered locations in the Kattegat-Skagerrak ocean margin in which the sites in this study are located. Based on the datasets in the abovementioned studies, Danish eelgrass meadows seem to hold larger C stock than other areas of comparable latitudes which highlights the influence of the local environmental setting.
Fewer studies have assessed nitrogen and phosphorus stocks in seagrass sediments and the current study is, to our knowledge, the largest-scale assessment available in the literature. McGlathery et al. (2012) measured the sediment nitrogen stock in the top 5 cm sediment in a 9-year-old restored eelgrass meadow in Virginia, United States and found an average N stock of 16.2 g N m −2 . Assuming a stable profile between 5-10 cm, and extrapolating across 10 cm, would put this restored eelgrass meadow (32.4 g m −2 ) on the lower end of the range of N stock presented here (24-448 g m −2 ). While not explicitly assessing sediment stocks,  measured the sediment nitrogen concentration in three eutrophic sites colonized by Z. marina and Z. japonica in northeastern China, where the average nitrogen content in the top 10 cm sediment was almost three times higher (0.35%) than what we found in Denmark. Fourqurean et al. (2012b) assessed nitrogen and phosphorus in the top 100 cm sediment in two seagrass systems (mainly Thalassia testudinum) in Western Australia and Florida, United States. The average nitrogen concentration in the sediment was similar between the two areas (0.15 ± 0.004 and 0.19 ± 0.01%, respectively), which is comparable to 0.12 ± 0.01% in our study. Phosphorus content exhibited much greater difference between their sites and was on average 87.7 ± 4.3 µg g −1 in Florida and 135.7 ± 3.6 µg g −1 in Australia, which is lower than what we found (310.1 ± 18.6 µg g −1 ) and was explained by the carbonate-rich sediments in these locations in which phosphorus content in the porewater is known to be low (Short, 1987;Fourqurean et al., 2012b).

Uncertainty, Limitations, and Future Direction
For best practices assessments of sediment stocks and longterm storage in seagrass meadows, sediment cores should ideally extend at least 100 cm to estimate buried carbon and nutrients below the depth of decomposition (Howard et al., 2014;Macreadie et al., 2014). One important caveat of assessing stocks in the top 10 cm is that this shallow sediment layer is considered biologically (e.g., bioturbation and -irrigation), physically (e.g., wave action, sediment resuspension) and biogeochemically (e.g., high remineralization rates, redox oscillations) active, and the sediment properties (including C, N, and P stocks) can vary over short time-scales (e.g., seasonal, annual) (Middelburg et al., 2004;Burdige, 2007;Johannessen and Macdonald, 2016). Consequently, although the top 10 cm in similar areas have been shown to represent approximately the past 60 years (Jankowska et al., 2016), any assertions regarding long-term carbon and nutrient storage cannot be adequately made (Johannessen and Macdonald, 2016;Nilsson, 2018). Normalizing stocks based on 10 cm cores to 25 or 100 cm requires untested assumptions of invariant change with sediment depth and that the 10 cm profile pattern is representative of deeper depths. However, deeper sediment cores (25 cm) from 10 sites in this study revealed that most depth profiles of carbon exhibited either a stable, or even increasing, pattern with sediment depth between 10-25 cm, especially in sites with large C stock (Supplementary Figure S6), suggesting that the risk of overestimation may be limited in this case. Furthermore, from these sites the vertically extrapolated C stock (from 10 to 25 cm) was generally lower than the measured Frontiers in Marine Science | www.frontiersin.org 9 December 2018 | Volume 5 | Article 474 C stock with a mean offset (± SE) of 895 ± 527 g C m −2 , corresponding to 64.7 ± 20.8% of the average C stock . Although uncertainties associated with deeper extrapolation persist, the estimated C stock values presented here should be considered conservative (Röhr et al., 2016).
A second limitation of this study is the lack of adjacent, unvegetated sites to allow for estimation of buried carbon and nutrients in the absence of eelgrass cover. However, Kennedy et al. (2010) reported almost twice as high median organic carbon content (0.34% compared to 0.19%) and total nitrogen (0.031% compared to 0.024%) in sediments within the seagrass meadow compared to adjacent bare. Similarly, Jankowska et al. (2016) observed 1.5-4.8 times higher organic carbon concentrations and 1.6-3.1 times higher nitrogen concentrations in eelgrass-vegetated sediments compared to unvegetated. Similar differences in OM content between vegetated and unvegetated reference sediments have also been documented elsewhere (see e.g., Miyajima et al., 1998;Ricart et al., 2015;Postlethwaite et al., 2018). This aspect should, however, be a focal point in future studies as the lack of unvegetated reference sediment is prevalent in the seagrass sediment stock literature. Yet, selecting a reference is by no means uncomplicated as both the OM produced by seagrass itself and the sestonic particles it captures have been documented to bury far away from the canopy (Kennedy et al., 2010;Duarte and Krause-Jensen, 2017).
A third important consideration concerning studies on seagrass sediment storage is the ability to scale measurements to basin-wide, regional, or even global level. Large scale extrapolation of sediment stocks warrants considerable precaution, as spatial (and temporal) variability is large not only between meadows as shown here, but can also be significant within meadows (Oreska et al., 2017a). Any upscaling of meadow point measurements therefore requires assumptions of invariable abiotic properties (e.g., hydrodynamic regime, residence time, sediment dynamics, geomorphology, water column chemistry) and biotic properties (e.g., meadow configuration, productivity, shoot density). This aspect was to some extent illustrated in the, albeit limited, discrepancy between the results of the GLM and the variable-by-variable regressions where the number of sites and environmental settings differed. Nevertheless, there is a need to provide stakeholders, managers, and policymakers with incentives to conserve and protect these valuable habitats, which are declining at rapid rates (Orth et al., 2006;Waycott et al., 2009). However, rather than producing first-order estimates of sediment storage on a larger scale, we emphasize the small-scale spatial variability and provide guidance into selecting high priority areas for carbon and nutrient storage on the local scale. Denmark does indeed comprise hotspots for carbon and nutrient storage, which can be targeted as high priority areas in the management of seagrass meadows for mitigation of climate change and eutrophication. By further pinning down the associated environmental factors, this study also indicates where additional hotspots may be located. All in all, this nationwide information on C, N, and P stocks highlights the combined climate change and eutrophication mitigation potentials of eelgrass meadows and thereby stimulates incentives for conservation efforts.

CONCLUSION
In this study, we have further established the large spatial variability in sediment storage of carbon, nitrogen, and phosphorus in seagrass systems. Our study is the first to assess these three elements in concert for a single species on a nationwide basis covering a variety of waterbodies. We have shown that carbon and nitrogen stocks are higher in sheltered inner fjords compared to more exposed locations. Furthermore, we have shown that environmental factors that are often considered detrimental to eelgrass, such as eutrophication and high sediment mud content (De Boer, 2007), to a certain extent have a positive effect on the sediment stocks and this intricate balance should be explored further. The heterogeneous nature of eelgrass sediment stocks and the varying influence of environmental parameters observed here highlights the importance of considering the environmental setting when assessing carbon and nutrient storage in sediments, constructing regional budgets, and evaluating ecosystem services in eelgrass meadows on a national scale.

AUTHOR CONTRIBUTIONS
TK, SØ, MH, and DK-J conceived and designed the study. MR collected the data. TK and SØ carried out data analysis. TK drafted the manuscript. All authors revised and approved the manuscript.

FUNDING
This work was built upon initial work presented (TK) at the 4th International Symposium on Research and Management of Eutrophication in Coastal Ecosystems (EUTRO2018). The study was part of the project Havets skove, funded by the Villum Foundation (Grant Number 18530).

ACKNOWLEDGMENTS
We would like to thank Prof. Jacob Carstensen for providing water column chemistry data from the DNAMAP database. DK-J also acknowledges support from the Ministry of Environment and Food of Denmark (contract 33010-NIFA-16651) and from the Danish Centre for Environment and Energy. We thank CG, JH, and TO for constructive comments that significantly improved the quality of this paper.