Introduction

Drivers of ecosystem function are complex and include environmental and biological components and their interactions. For example, temperature, pH, and oxygen directly affect leaf litter decomposition in streams (Webster and Benfield 1986) but also constrain the microbial and macroinvertebrate communities involved in decomposition (Gessner et al. 1999). In forest ecosystems, climate, nutrient availability, and plant diversity influence productivity (Lieth 1975; Vicca et al. 2012; Duffy et al. 2017), where climate determines the shape of biodiversity-productivity relationships (Fei et al. 2018). Drivers of ecosystem function are also context and scale dependent: decomposition rates can vary across riffles within a single stream reach (Tiegs et al. 2009), and forest productivity can vary according to local resource availability mediated by topography (Hao et al. 2018) and canopy structure (Gough et al. 2019). Controls on carbon cycling in vegetated coastal ecosystems, such as seagrass meadows, are similarly complex, so holistically examining abiotic and biotic variables at targeted scales can identify drivers of variability in the ecosystem function of sediment carbon storage.

Seagrasses are marine angiosperms that grow along coastlines globally (except for Antarctica) and provide ecosystem services including but not limited to habitat provisioning for fish and invertebrates, sediment stabilization, water quality enhancement, and carbon storage (Nordlund et al. 2016). Seagrass ecosystems have gained recognition for their ability to trap and store organic carbon in underlying sediments, preserving stored carbon under anoxic conditions for decades to millennia (Fourqurean et al. 2012; Duarte et al. 2013). Seagrass meadows store an estimated 4.2–8.4 Pg organic carbon globally (Fourqurean et al. 2012), but storage is highly variable spatially and therefore difficult to predict (McLeod et al. 2011). Improved understanding of the mechanisms controlling spatial variability in seagrass carbon storage will inform global carbon budgets (Duarte et al. 2013) and better equip conservation practitioners and other end users to utilize seagrass restoration and management for nature-based climate change mitigation (Macreadie et al. 2017, 2019). While environmental characteristics such as coastal geomorphology may constrain carbon storage capacity at global scales (Kennedy et al. 2022), drivers of smaller scale variability in storage are unresolved.

Seagrass sediment organic carbon storage varies substantially geographically, even throughout estuaries and within meadows (Oreska et al. 2017; Ewers Lewis et al. 2018; Ricart et al. 2020). Previous work has emphasized hydrology, water depth, and landscape configuration as broad-scale controls on carbon accumulation and storage (Serrano et al. 2014, 2016b; Lima et al. 2019; Asplund et al. 2021), but at smaller spatial scales, seagrass shoot density, meadow extent, and canopy architecture can also influence carbon burial efficiency (Samper-Villarreal et al. 2016; Mazarrasa et al. 2018). Canopy complexity, an emergent property of leaf size and morphology, affects wave velocity above the canopy, which in turn affects organic matter deposition at the sediment surface (Hendriks et al. 2008; Barcelona et al. 2021). Belowground, seagrasses deliver dissolved organic carbon via root exudates directly to the sediment (Moriarty et al. 1986; Barrón and Duarte 2009), in some cases accumulating high concentrations of sugar in near-surface sediments (Sogin et al. 2022) and promote anoxic conditions that prevent carbon remineralization (Trevathan-Tackett et al. 2017). Canopy complexity, release of root exudates, and sediment redox potential may therefore be local drivers of carbon storage, and these features likely vary according to meadow characteristics such as plant density and species composition.

Temporal fluctuations in the amount of seagrass cover and species assemblage associated with succession, disturbance, and season can substantially alter canopy complexity and belowground processes that may affect carbon accumulation and burial (Williams 1990; Rasheed 2004; Kilminster et al. 2015; Dahl et al. 2020 2020; Zhu et al. 2022). Loss of seagrass cover often results in diminished carbon stocks (Arias-Ortiz et al. 2018; Moksnes et al. 2021), but natural seagrass recovery and restoration can help to return this ecosystem service (Greiner et al. 2013; Macreadie et al. 2015; Marbà et al. 2015; Aoki et al. 2021). Given seagrass meadow cover can be highly dynamic, site history, especially in terms of seagrass cover, may be an important determinant of sediment carbon storage. In addition to cover, studies have found meadows dominated by certain seagrass genera or species support greater carbon stocks than other seagrass assemblages (Lavery et al. 2013; Gillis et al. 2017; Alemu et al. 2022). However, seagrass assemblage and carbon stocks were assessed contemporaneously, and sampling sites spanned a large area (> 100 km2). Therefore, mismatched spatial and temporal scales may have confounded species effects: seagrass assemblage may differ spatially along unobserved environmental gradients and/or may have changed over the time in which the measured carbon stocks accumulated. Sampling carbon densities at the meadow-scale and considering site history may help to clarify the role of seagrass cover and species composition in carbon storage.

If species composition is a determinant of carbon storage, species diversity may also drive storage in seagrass meadows where multiple species co-exist. Mixed-species meadows are more prevalent in tropical and subtropical regions, but carbon storage research is biased toward temperate meadows that are often dominated by a single species (Piñeiro-Juncal et al. 2022). Biodiversity-ecosystem function theory predicts that biodiversity will promote ecosystem stability and functioning, and many positive relationships have been observed between biodiversity and ecosystem function in manipulative experiments and natural systems (Balvanera et al. 2006; Cardinale et al. 2011; O’Connor et al. 2017; Duffy et al. 2017). The role of biodiversity in seagrass ecosystem functioning has typically been studied in terms of phenotypic (i.e., structural complexity) and genetic diversity. Seagrass meadow functions of habitat, food, and predation refugia provisioning generally increase with structural complexity (sensu shoot density or leaf surface area; Duffy 2006). Genetic diversity can confer resistance to and enhance recovery following disturbance (Hughes and Stachowicz 2004, 2011; DuBois et al. 2021), increase shoot density (Hughes et al. 2009), and promote the recovery of ecosystem services following restoration (Reynolds et al. 2012). Species diversity has received less attention but may also drive seagrass ecosystem functioning: species-driven variation in plant physical structure and physiology affect canopy and belowground properties that ultimately affect carbon accumulation and storage. Species diversity, in terms of both richness and identity, may also reflect the community response to an ecological disturbance via succession within a seagrass meadow (Williams 1990; Martínez López et al. 2019; Correia and Smee 2022). In this way, species diversity may be correlated to meadow successional status and disturbance history, which are also likely drivers of carbon storage.

We set out to investigate historical (14-year averaged) and contemporary meadow-scale controls on surface sediment organic carbon storage in a subtropical, mixed-species seagrass meadow in the northeastern Gulf of Mexico. Surface carbon densities (here defined as 0–10 cm) are subject to community metabolism, bioturbation and other processes that can erode carbon stores and therefore represent potential quantities of carbon stored in the longer term; and for these same reasons, they are most appropriate to identify decadal-scale historical and contemporary drivers of carbon storage. We measured surface carbon densities at long-term seagrass monitoring stations to assess how storage varies by both historical and contemporary seagrass cover, species composition, and species diversity. Given the well-documented effects of environmental setting on carbon accumulation, we determined whether sediments with high carbon densities were spatially clustered across the ~ 25 km2 study area and considered the effect of relative exposure, elevation (as a proxy for depth), and distance to navigation channels (as a proxy for physical disturbance related to boating activity) on carbon densities. There are several potential environmental drivers of carbon storage, but we focused on drivers that directly affect hydrology (i.e., relative exposure, elevation) or cause physical disturbance (i.e., boating activity). We also assessed the relationship between sediment grain size and carbon densities as grain size commonly covaries with organic carbon content (Serrano et al. 2016a). Finally, we used path analysis to evaluate hypothesized causal relationships between predictor variables to shed light on historical and contemporary processes that may influence carbon storage at the meadow-scale.

We hypothesized biodiversity would affect seagrass carbon storage indirectly: historical species diversity may be correlated with meadow stability (inversely measured as variability in cover) and productivity, which create and sustain conditions favorable for carbon accrual (i.e., increased seagrass canopy complexity that promotes particle trapping and stabilized sediments or otherwise modified rhizosphere conditions that promote long-term storage); and contemporary biomass and species composition may influence the vulnerability of surface carbon densities to remineralization via species-specific metabolism, rhizosphere oxygenation, and root exudation rates. We expected environmental drivers would be directly related to carbon densities but would also affect carbon storage indirectly via historical and contemporary meadow features.

Methods

Sampling location and collection

We conducted this study in nearshore seagrass beds proximal to Cedar Key, Florida (USA) that are part of a small river delta system fed by the Suwannee River in the northeastern Gulf of Mexico. Cedar Key experiences a mean tidal range of 0.86 m and seagrass beds are predominantly subtidal, though seagrasses in some areas may briefly experience emersion during extreme low tides. Seagrass beds in Cedar Key host up to five seagrass species, growing singly or intermixed, and include Thalassia testudinum, Syringodium filiforme, Halodule wrightii, Halophila engelmanii, and Ruppia martima. The species present in our study area exhibit a spectrum of life history strategies, ranging from smaller pioneering species (e.g., H. wrightii) to larger persistent species (e.g., T. testudinum; Kilminster et al. 2015). Our study sites were co-located with 25 fixed stations in Cedar Key (Fig. 1) where the Florida Department of Environmental Protection (DEP) monitors annual total seagrass cover and species-specific cover as part of the Big Bend Seagrasses Aquatic Preserve seagrass monitoring program. In August 2020, we collected samples at each of the 25 sites to quantify surface (0–10 cm) carbon densities and contemporary meadow features.

Fig. 1
figure 1

The main map shows the 25 sampling sites, and the inset map shows the approximate study area in Cedar Key, Florida. Sites were co-located at Florida Department of Environmental Protection Big Bend Seagrasses Aquatic Preservative fixed monitoring stations and are colored according to surface carbon density quantiles

Full size image

We targeted the sediment depth that likely reflects recent accumulation coinciding with seagrass monitoring initiated in 2006; globally, seagrass meadow sediment accretion rates range from 0.13 to 9.9 mm yr−1 (Potouroglou et al. 2017), so 10 cm may represent accumulation over the past decade or longer. We collected 3 replicate 2.5-cm diameter cores (60-cc syringe) to 10 cm depth for sediment organic carbon analysis and an additional core for grain size analysis at each site, taking care to core between seagrass shoots to avoid collecting live seagrass tissue. The number of replicate cores for some sites was fewer than 3 because we discarded any cores that were < 10 cm depth. Sediment cores were sectioned into 3 intervals for processing: 0–3 cm, 3–6 cm, and 6–10 cm. We also collected 10 replicate 15-cm diameter cores to 20 cm depth to measure contemporary seagrass biomass at each site. Each replicate core was collected ~ 1 m apart for both sediment and biomass samples. All samples were kept on ice until transported to the laboratory where they were frozen (− 40 °C) until processed.

Sediment organic carbon content

We measured sediment sample bulk wet weight and dry weight after drying for 1–2 weeks at 60 °C to calculate bulk density (mass of dried soil divided by wet sample volume). Dried sediment samples were then homogenized using a mortar and pestle and plant tissue (both live and detrital) and shell fragments > 1 mm were removed from dried sample to isolate sediment organic matter from belowground biomass because the focus of this study was on the sediment component of carbon storage (Greiner et al. 2013; Oreska et al. 2017). While we mostly avoided sampling live plant tissue in sediment cores, pieces of live shoots were occasionally ensnared within the cores, especially at sites with more dense seagrass cover. Detrital seagrass tissue has potential to contribute to long-term carbon storage but was removed because retaining larger pieces of detrital tissue in the relatively small (2.5 cm-diameter) cores would have led to inflated percent organic carbon values. Removing shells and plant tissue from the dried bulk samples likely resulted in the underestimation of total organic carbon, however, the amount of dried material removed from each sample was minimal, typically < 0.1 g, and an order of magnitude smaller than the mean dry sample weight of ~ 8.5 g.

We indirectly measured percent total organic carbon by subtracting percent total inorganic carbon from percent total carbon. Total carbon was measured using a Carlo Erba NA1500 CNS Elemental Analyzer and total inorganic carbon was measured by coulometric titration using an acidification preparation device (AutoMate FX, Inc.) coupled to a CM5017 CO2 Coulometer (UIC, Inc.). We estimated carbon densities by multiplying percent total organic carbon by core section length (cm) and bulk density (g cm−3) and converted values to mg organic carbon (OC) cm−3. We calculated total core carbon as the summed carbon densities across all depth intervals as described in Fourqurean et al. (2014) for each replicate core and calculated the mean value for each site to generate mean surface carbon densities.

Grain size

We measured sediment grain size for each core depth interval for each site. Prior to analysis, samples were dried for 1–2 weeks at 60 °C and passed through a 2-mm sieve. We analyzed particle size distribution using a Beckman Coulter LS-13,320 multi-wave particle size analyzer equipped with PIDS technology at the University of Florida Environmental Pedology and Land Use Laboratory. To assess site-level variability in grain size, we analyzed a duplicate core for ~ 10% of samples and found the coefficient of variation for average particle diameter was < 15% across all duplicate cores. We calculated the percent of fine-grained sediment, hereafter referred to as percent fines, for the entire length of the cores (0–10 cm) as the depth-averaged sum of percent clay (< 2 μm) and percent silt (2–50 μm) grain size classes.

Historical seagrass data

We analyzed annual seagrass monitoring data (Florida Department of Environmental Protection, 2021) to calculate historical predictor variables for each site. We used all available seagrass cover data that preceded our sampling event in 2020 (13 years total within 2006–2019). Seagrass cover was determined using the Braun-Blanquet method, a visual cover estimation technique. The Florida DEP estimated cover using a slightly modified Braun-Blanquet scale where the lowest Braun-Blanquet score is 1 and spans values of 0.1–5% cover (Table S1 in Online Resource 1). The modified scale may be less precise in estimating small amounts of seagrass cover but still allowed us to compare differences in seagrass cover between sites and within sites over time.

Braun-Blanquet density values are ordinal, so we converted values to percent cover using a log-linear transformation (van der Maarel 2007) for ease of interpretation, though both transformed and untransformed Braun-Blanquet data perform well as proxies for continuous data in statistical analyses (Furman et al. 2018). For each site, we calculated the following historical seagrass cover variables from 4 replicate quadrats surveyed at each site annually across the monitoring period: mean total cover (includes all species), variability in cover (coefficient of variation), mean species richness, and mean relative dominance by species (species-specific cover divided by mean total cover).

Contemporary seagrass data

We thawed and rinsed seagrass samples and removed epiphytes by gently scraping leaves with a razor blade. After drying for ~ 1 week at 60 °C, we measured dry weight of biomass separated by species and by aboveground vs. belowground tissue. For each site, we calculated the following contemporary seagrass biomass variables from the 10 replicate biomass cores: mean total, mean above- and belowground biomass, and mean relative dominance by species (species-specific biomass divided by total biomass). Contemporary species richness was calculated as the number of unique species identified across all 10 replicate biomass cores.

Environmental variables

To account for environmental heterogeneity across our study area, we calculated a relative exposure index as described in Fonseca and Bell (1998) and measured elevation (as a proxy for depth) and distance to navigation channels (as a proxy for physical disturbance related to boating activity) for each site. The relative exposure index uses wind speed, wind direction, and fetch to estimate relative exposure related to wave energy. Wind data from 2018 to 2020 were obtained from the NOAA National Data Buoy Center (Station CDRF1; https://www.ndbc.noaa.gov/). Fetch was estimated using the R package ‘waver’ (Marchand and Gill 2018) by calculating site distance to shore from eight bearings based on the Florida Shoreline (1 to 40,000 Scale) dataset (available at https://geodata.myfwc.com/).

We used modelled elevation as a proxy for water depth to account for differences in light that may influence seagrass meadow distribution (Duarte 1991) and carbon cycling (Serrano et al. 2014). We measured elevation using the NOAA National Centers for Environmental Information Continuously Updated Digital Elevation Model (CUDEM)–1/9 Arc-Second Resolution Bathymetric-Topographic Tiles (available at: https://coast.noaa.gov). We measured site distance to navigational channels using the Coastal Maintained Channels in U.S. Waters dataset (Office of Coast Survey 2015). Portions of the navigation channels in Cedar Key were excluded from this dataset, so we georeferenced a more comprehensive map of channels (https://beta.marinenavigation.noaa.gov/gateway/) to the study area using the World Imagery Basemap in ArcGIS 10.8, and then digitized the channels to create a vector file with a more accurate representation of Cedar Key navigation channels.

Statistical analyses

All statistical analyses were performed in R version 4.2.0 (R Core Team 2022). We first tested for spatial autocorrelation of carbon densities to determine whether spatial weights would be necessary for further statistical analysis and found no evidence of significant spatial autocorrelation based on Moran’s I statistic (Moran’s I = − 0.04, p-value = 0.42) implemented using the ‘ape’ package (Paradis and Schliep 2019). We examined bivariate relationships between each predictor variable and carbon densities using simple linear regression with the ‘lm’ function. We checked assumptions of linearity, normality, and homoscedasticity of residuals by visually examining fitted values vs. residuals, scale-location plots, and performing Shapiro-Wilk tests using the ‘shapiro.test’ function, respectively. Species richness was log-transformed to meet model assumptions. We adjusted the threshold for statistical significance for bivariate analyses to increase the power to detect relationships with small effect sizes (r2 < 0.25) and a relatively small sample size (n = 25): we considered evidence for relationships to be strong, moderate, and weak at α < 0.01, 0.01 ≥ α < 0.05, and 0.05 ≥ α < 0.10, respectively (Mudge et al. 2012; Muff et al. 2022). While linear regression may oversimply the likely nonlinear and multivariate nature of factors driving carbon storage, the bivariate analyses were useful in assessing general patterns (i.e., direction and magnitude of effects) between potential driver variables and carbon densities.

We evaluated hypothesized causal relationships driving carbon storage with exploratory path analysis using the ‘lavaan’ package (Rosseel 2012). We selected a subset of the many potential predictor variables of carbon densities we measured–at least one variable from each of the variable types (i.e., sediment, environmental, historical seagrass, and contemporary seagrass variables)–to create a path model. First, we checked for collinearity between potential predictor variables within each variable type based on Pearson correlation coefficients calculated using the ‘Hmisc’ package (Harrell 2022). Within each variable type, there were significant correlations (α < 0.05) among most of the potential predictor variables (Table S2 in Online Resource 1). We therefore selected noncorrelated variable(s) from each type: for example, we selected contemporary total biomass as the contemporary seagrass variable, because total biomass was significantly correlated with contemporary above- and belowground biomass, species richness, and T. testudinum relative dominance (Table S2 in Online Resource 1).

The path model consisted of two exogenous variables (elevation and distance to channel), two historical endogenous variables (variability in cover and historical species richness), and three contemporary endogenous variables (total biomass, sediment percent fines, and carbon densities,) to test our hypothesis that environmental drivers would have direct effects and indirect effects mediated by historical drivers on contemporary features that influence carbon storage (Fig. 2). Data were centered and scaled using the ‘scale’ function to facilitate analysis across variables with different measurement units and scales. We iteratively tested variations of the model with additional direct links between predictor variables and carbon densities. We assessed model fit by examining the goodness-of-fit χ2 test statistic (χ2p > 0.05), comparative fit index (CFI ≥ 0.95), and standardized root mean squared residual (SRMR ≤ 0.10).

Fig. 2
figure 2

Causal diagram illustrating hypothesized relationships between historical and contemporary predictors of surface carbon densities. Variables names were abbreviated as follows: Distance to channel = distance to navigation channel; Historical variability = historical variability in cover (CV); Historical diversity = historical species richness. Solid lines indicate causal relationships and dashed lines indicate correlations that are either positive (black) or negative (red)

Full size image

Figures were generated in R version 4.2.0 (R Core Team 2022) using the following R packages: ‘DiagrammeR’ (Iannone 2022), ‘DiagrammeRsvg’ (Iannone 2016), ‘ggplot’ (Wickham 2016), ‘ggplotify’ (Yu 2021), ‘ggmap’ (Kahle and Wickham 2013), ‘ggsn’ (Santos Baquero 2019), ‘ggpubr’ (Kassambara 2020), ‘RColorBrewer’ (Neuwirth 2022), ‘rsvg’ (Ooms 2022), ‘semPlot’ (Epskamp 2022), and ‘tiff’ (Urbanek and Johnson 2022).

Results

Sediment organic carbon content and grain size

Sediment characteristics varied across sites (Table S3 in Online Resource 1). Dry bulk density ranged from 1.01 to 1.94 g cm−3, with a mean of 1.48 g cm−3, 95% CI [1.38, 1.59]. Total carbon spanned an order of magnitude from 0.10 to 2.37 total carbon by percent weight (mean = 0.92, 95% CI [0.64, 1.19]). Sediments from all sites contained low amounts of inorganic carbon, averaging 0.18 total inorganic carbon by percent weight, 95% CI [0.11, 0.24]. Therefore, organic carbon comprised the bulk fraction of the sediment total carbon at most sites. Carbon densities incorporated total organic carbon and dry bulk density measurements and varied considerably from 1.23 to 21.04 mg organic carbon (OC) cm−3, with a mean of 9.85 mg OC cm−3, 95% CI [7.53, 12.17], though high values for carbon densities were not spatially clustered (Fig. 1). Mean particle size and percent fines also varied across sites, with mean values of 315.53 μm, 95% CI [288.09, 342.97], and 4.25% volume, 95% CI [3.11, 5.39], respectively.

Historical and contemporary seagrass metrics

Both historical and contemporary seagrass metrics varied across sites (Table S3 in Online Resource 1). Historical seagrass mean total cover (log-linear transformed) and variability in cover (coefficient of variation) ranged from 22.94 to 154.82% and from 0.13 to 2.1, respectively. While up to five seagrass species occur in Cedar Key, historical species richness averaged over the entire monitoring period was low (mean = 1.4, 95% CI [1.26, 1.54]). Historically, T. testudinum had the highest relative dominance (mean = 0.52, 95% CI [0.38, 0.66]), followed by H. wrightii (mean = 0.34, 95% CI [0.18, 0.51]) and S. filiforme (mean = 0.15, 95% CI [0.06, 0.24]. Two sites were completely bare during the August 2020 sample collection, though they were historically vegetated. Contemporary seagrass biomass ranged from 0 (bare) to 752.55 g m−2. While the maximum contemporary species richness was 4, richness was lower at most sites (mean = 2.00, 95% CI [1.57, 2.43]). Contemporary relative dominance values were similar to patterns in the historical dataset where T. testudinum had the highest relative dominance (mean = 0.59, 95% CI [0.42, 0.77]), followed by H. wrightii (mean = 0.37, 95% CI [0.17, 0.56]) and S. filiforme (mean = 0.11, 95% CI [0.03, 0.19]). H. engelmannii was rarely observed in both the historical and contemporary datasets. R. maritima cover was recorded for a single site in 2006 and was not present in any of the contemporary biomass cores; historical and contemporary relative dominance of R. maritima were therefore excluded from analyses.

Environmental variables

Environmental characteristics varied across sites even though all sites spanned a relatively small area (Table S3 in Online Resource 1). The relative exposure index (unitless) ranged from 4.1 × 103 to 2.6 × 104, averaging 1.8 × 104, 95% CI [1.5 × 104, 2.1 × 104]. Elevation extracted from the digital elevation model ranged from − 2.95 to − 0.09 m, averaging − 1.02 m, 95% CI [1.40, − 0.64]. The distance to the nearest navigational channel ranged from 0.19 to 3.55 km, averaging 1.35 km, 95% CI [0.99, 1.71].

Bivariate relationships

We found statistically significant relationships (p-value < 0.10) relationships between sediment and seagrass variables and carbon densities but not between environmental variables (i.e., relative exposure, elevation, distance to channel) and carbon densities (Table 1). Percent fines was strongly related to carbon densities (Fig. 3). Both historical mean total cover and species richness were positively related to carbon densities while H. wrightii relative dominance and variability in cover were negatively related to carbon densities (Fig. 4). Contemporary H. wrightii relative dominance was also negatively related to carbon densities while contemporary T. testudinum relative dominance was positively related to carbon densities (Fig. 5). Contemporary S. filiforme relative dominance and species richness were both weakly positively related to carbon densities (Fig. 5).

Table 1 Comparison of effect sizes in simple linear regression models with surface carbon densities as the response variable. Relationships were considered statistically significant at α < 0.10 (bold text) 
Full size table
Fig. 3
figure 3

Scatterplot showing the relationship between sediment percent fines and surface carbon densities

Full size image
Fig. 4
figure 4

Scatterplots depicting a subset of relationships (p-value < 0.10) between surface carbon densities and historical predictor variables: a H. wrightii relative dominance, b cover, c variability in cover, and d species richness

Full size image
Fig. 5
figure 5

Scatterplots depicting a subset of relationships (p-value < 0.10) between surface carbon densities and contemporary predictor variables: a H. wrightii relative dominance, b S. filiforme relative dominance, c T. testudinum relative dominance and d species richness

Full size image

Path analysis

We found greatest support in terms of model fit for the hypothesized path model with an additional direct link from historical variability to carbon densities (Table 2; Fig. 6). Contrary to the hypothesized model (Fig. 2), the final, best-fitting model contained positive path coefficients between elevation and percent fines and distance from channel and historical species richness, and negative path coefficients between elevation and historical species richness and historical species richness and contemporary total biomass. Therefore, we found evidence that environmental drivers indirectly influence carbon densities via percent fines and historical species richness in unexpected ways and that historical variability influences carbon densities directly, or at least via mechanisms other than contemporary total biomass.

Table 2 Direct and total effects of each predictor variable on surface carbon densities from the path analysis
Full size table
Fig. 6
figure 6

Diagram of the best fitting path analysis model. Variables names were abbreviated as follows: Distance to channel = distance to navigation channel; Historical variability = historical variability in cover (CV); Historical richness = historical species richness. Solid lines indicate causal relationships and dashed lines indicate correlations that are either positive (black) or negative (red). Nodes are labeled with coefficients of determination (R2), edges are labeled with path coefficients, and edge thickness is scaled to reflect the relative strength of path coefficients

Full size image

Discussion

Seagrass meadows in Cedar Key, Florida generally stored less carbon compared to previous regional and global estimates, but surface carbon storage to 10 cm depth was spatially variable, ranging from 1.23 to 21.04 mg OC cm− 3 within the 25 km2 sampling area. Identifying drivers of variability observed at small spatial scales improves our ability to generate reliable global estimates of carbon storage and enhances our understanding of the role of seagrass meadows in the global carbon cycle. We found that historical measures of meadow cover were generally better predictors of carbon densities than contemporary measures of seagrass biomass, highlighting the importance of site history to ecosystem function. Species identity was also a significant driver of carbon storage: meadows dominated by the long-lived, large-bodied species, T. testudinum, stored the greatest amounts of carbon. Both historical and contemporary species diversity were weakly related to carbon densities and environmental variables had only indirect effects on storage. Our findings suggest meadow stability and late-successional species assemblages, and to a lesser degree species diversity, enhance seagrass sediment carbon storage at the meadow scale.

Comparison to global and regional carbon stocks

Carbon storage reported for different depth intervals is not directly comparable, so it is not surprising that mean carbon densities measured for 0–10 cm in our study were generally lower than densities reported for deeper cores. Mean percent total organic carbon (0.73%) was less than the global median of 1.8% in the top meter of sediments (Fourqurean et al. 2012) and mean carbon densities, 9.85 Mg organic carbon (OC) ha− 1 in the top 10 cm, were lower than the estimated global mean of 15.4 Mg OC ha− 1 in the top 20 cm (Kennedy et al. 2022). Densities were also lower than values reported for 0–25 cm in temperate Zostera marina meadows (27.2 Mg OC ha− 1; Röhr et al. 2018) and for 0-100 cm in tropical seagrass meadows (129.0 Mg OC ha-1; Alongi 2014). Regionally, densities were lower than average values reported for 0–20 cm from subtropical and tropical meadows from the Gulf of Mexico (25.7 Mg OC ha− 1; Thorhaug et al. 2017), but the maximum density from our study (21.04 Mg OC ha− 1) approached the maximum value reported for 0–20 cm from nearby sites up to 80 km south of Cedar Key (~ 20 mg OC cm− 3 from Barry et al. 2018). If scaled to 1-m stocks, average carbon densities from Cedar Key (9.85 Mg OC ha− 1 for 0–10 cm) would approach or exceed the reported global and regional averages. However, organic carbon content tends to decrease with depth (Fourqurean et al. 2012) so additional measurements would be needed to compare densities across varying depth intervals.

Maximum carbon densities from Cedar Key were greater than densities estimated from deeper cores in nearby sites influenced by a different river system further south along peninsular Florida’s west coast (Barry et al. 2018). Here, we suspect hydrologic controls, mainly river discharge, explain differences in densities. Seagrass meadows in this study and in Barry et al. (2018) are in open coastal systems, where local variation in carbon storage may be related to allochthonous carbon and sediment delivery rates (Ewers Lewis et al. 2018). The Suwanee River, with a mean annual discharge of 298.5 m3 s− 1 (Franklin et al. 1995), likely delivers more carbon and sediment to Cedar Key than do smaller rivers flowing into coastal meadows to the south; for example, the Waccasassa River has a mean annual discharge of only 8.0 m3 s− 1 (Florida Department of Environmental Protection 2005). The Suwanee and Waccasassa are partially and primarily spring-fed, respectively, so both deliver relatively small amounts of sediment to the coast compared to other river systems in the Gulf of Mexico, which limits the carbon storage potential in receiving meadows (Thorhaug et al. 2017). Watershed characteristics and delivery rates of fine-grained sediments may therefore represent broad-scale controls on carbon storage in this region.

As expected, we found a strong relationship between percent fines and carbon densities (Figs. 3 and 6) in congruence with past studies (Dahl et al. 2016; Howard et al. 2020; Krause et al. 2022). Serrano et al. (2016a) found that tight correlations between fine-grained sediment and sediment organic carbon in seagrass meadows typically indicate that seagrass-derived organic carbon contributes minimally to the sediment organic carbon pool. We therefore suspect that allochthonous organic matter, rather than seagrass-derived organic matter, comprises the bulk of carbon densities in Cedar Key, though organic carbon source analysis would be needed to confirm this hypothesis. Mean percent fines in our study (4.25%) was much lower than the global mean percent mud (18.93%) for 0–10 cm in seagrass sediments (Piñeiro-Juncal et al. 2022), which suggests a lower supply of fine-grained sediments and/or higher site erosional energy may explain the low mean carbon densities in Cedar Key relative to global averages. However, percent mud encompasses a slightly larger particle size class (< 63 μm) than does percent fines (< 50 μm), and particle size distributions assessed using different methods (e.g., sieving vs. laser diffraction) may not be directly comparable (Eshel et al. 2004). Percent fines is considered a driver of carbon storage because the greater particle surface area provides more adsorption sites for small organic particles, and the resulting mineral-associated organic carbon is largely protected from degradation (Mayer 1994; Miyajima et al. 2017). Percent fines may also be a covariate of carbon storage because fine-grained sediment accumulates in depositional environments where organic matter also accumulates (Mazarrasa et al. 2017; Novak et al. 2020). Whether percent fines is a causative or correlative agent of carbon densities, it is nonetheless a strong predictor of seagrass carbon storage.

Meadow persistence, composition, and diversity as drivers of carbon storage

While mean carbon densities were generally low in Cedar Key, there was significant variation across sites, and consistent with our hypotheses, historical seagrass features directly influenced carbon densities. Sites with persistent, high cover over time stored more carbon, whereas sites with higher contemporary biomass did not necessarily store more carbon (Table 1; Figs. 4 and 5). The path model emphasized the importance of historical variables in driving carbon storage, which aside from percent fines, had the greatest total effects on carbon densities (Table 2). Historical variability in cover had a direct negative effect on carbon densities beyond the indirect path via contemporary biomass (Fig. 6), suggesting that meadow history rather than standing biomass is the better predictor of carbon densities. In other Florida estuaries, Howard et al. (2020) similarly found a positive correlation between historical (two-year averaged) seagrass cover and carbon storage in Florida Bay, and Lebrasse et al. (2022) found meadow extent in St. Joseph’s Bay had not changed over a 30-year period and inferred that belowground carbon stores were also likely to be stable. Other studies have linked contemporary seagrass cover or biomass with carbon storage (Lavery et al. 2013; Belshe et al. 2018; Alemu et al. 2022), but historical drivers of seagrass biomass, such as long-term nutrient history, may underpin observed patterns (Armitage and Fourqurean 2016). Seagrass cover can change over intra- to interannual time scales (O’Brien et al. 2018), but sediment carbon accumulation occurs over longer time scales (Mateo et al. 1997; Orem et al. 1999; Mckee et al. 2007), so metrics that account for annual fluctuation in cover may better explain variation in carbon densities. It is unclear how long of a historical record is needed to evaluate meadow stability as a driver of sediment carbon densities and likely varies across environmental and ecological disturbance gradients (Snelgrove et al. 2014). In Cedar key, 14-year meadow cover history predicted surface (0–10 cm) carbon densities well, where meadows with stable seagrass cover confer greater carbon storage benefits than do meadows with more variable cover; however, it’s not clear whether this time scale is appropriate to predict carbon densities deeper than 10 cm.

We also detected species composition effects on carbon densities, presumably because the spatial extent of our study (25 km2) did not span pronounced geomorphic gradients that often overwhelm community drivers of carbon storage (Belshe et al. 2018). Sites with historically and contemporaneously high amounts of H. wrightii stored less carbon and sites with contemporaneously high amounts of T. testudinum–and S. filiforme to a lesser degree–stored more carbon (Figs. 4 and 5). T. testudinum, a large-bodied, slow-growing persistent species, may more strongly attenuate wave or current velocity and better stabilize sediments (Fonseca and Fisher 1986), thereby enhancing organic carbon deposition and storage in sediments. H. wrightii, a small-bodied, fast-growing pioneering species, may better withstand harsher hydrodynamic conditions in settings where fine sediments and organic carbon are unlikely to accumulate (Fonseca and Bell 1998).

H. wrightii and T. testudinum relative dominance may also be indicators of site disturbance history and ecological succession: contemporary relative dominance by T. testudinum was negatively correlated with variability in cover and both historical and contemporary relative dominance by H. wrightii were positively correlated with variability in cover (Table S2 in Online Resource 1). In this way, T. testudinum may indicate past meadow stability, whereas H. wrightii may indicate past meadow instability. Alemu et al. (2022) similarly found higher carbon stocks in Singaporean seagrass communities containing persistent seagrasses (Cymodocea species), but other studies from the Indo-Pacific have found higher carbon stocks associated with pioneering species in the genus Halophila (Lavery et al. 2013; Gillis et al. 2017). Our findings agreed with the global pattern of large-bodied, long-lived seagrass species driving meadow carbon storage (Kennedy et al. 2022; Piñeiro-Juncal et al. 2022) but species identity effects on carbon storage must be studied at fine (< 1 m) spatial scales to uncover the mechanisms generating this pattern.

The relationship between species richness and carbon storage was less clear: historical and contemporary species richness were both only weakly directly related to carbon densities (Table 1; Figs. 4 and 5) but historical species richness had a strong total effect on carbon storage in the path model, mainly via a positive effect on percent fines (Table 2; Fig. 6). Given the moderately strong relationships between species’ relative dominance and densities, we suspect the positive relationships between species richness and densities could be due to selection effects (Loreau and Hector 2001), where more diverse assemblages are more likely to contain the persistent species, T. testudinum. Competition between species may also drive the effects of species richness on surface carbon storage. For example, Krause et al. (2023) observed increased leaf lengths in mixed-species meadows in the Florida Keys where S. filiforme and T. testudinum competed for light; increased canopy height resulting from interspecific competition may enhance particle trapping and promote carbon storage. These findings suggest species identity as well as richness may be important to the ecosystem function of carbon storage at the site-level.

In our study, meadows contained a maximum of four species, so we might expect to find more pronounced species richness-carbon storage relationships in the comparatively species-rich meadows of the Indo-Pacific, where mixtures of species with different rooting depths better stabilize sediments than monocultures (Duarte et al. 1997). We also note that species richness is a simplified descriptor of species diversity that does not always reflect functional trait diversity (Schoolmaster et al. 2020), the purported mechanism driving biodiversity-ecosystem function relationships; therefore, fine-scale flume experiments or in situ hydrological measurements are needed to provide a mechanistic understanding of how seagrass species identity and trait diversity affect surface carbon storage via particle trapping and sediment stabilization.

Environmental drivers of carbon storage at the meadow-scale

We identified several site-level biological and sediment predictors of carbon storage, yet our linear regression results suggest environmental variables were less important drivers at the meadow scale (Table 1). We expected environmental variables, especially those related to hydrology (i.e., relative exposure, elevation) and anthropogenic disturbance (i.e., distance to navigation channels), to drive sediment deposition and carbon storage because water column temperature, salinity and nutrient levels are similar across nearshore meadows in this region (Sullivan et al. 2021). While we did not assess bulk sediment, sediment porewater, or seagrass tissue nutrient concentrations, we suspect seagrasses were nutrient-limited in Cedar Key, as is typical for tropical seagrasses in relatively unimpacted regions (Duarte and Chiscano 1999). There are also other potential hydrologic drivers of storage that we did not directly capture in this study, such as current speed (Fonseca and Bell 1998) and position within the meadow (i.e., edge vs. interior; Oreska et al. 2017).

We anticipated sites with greater exposure and higher elevation would experience greater wave energy and thus be less favorable to sediment deposition and sites near channels would have lower carbon densities due to sediment resuspension caused by boat traffic. However, the spatial scale relevant to boating impacts is difficult to identify. If boaters reliably utilize channels, we expect seagrass beds to experience highly localized impacts near channels, so we might not observe a trend between distance to channel and carbon densities for sites situated beyond 10–100 m from navigational channels. We did not find the hypothesized bivariate relationships between environmental variables and carbon densities at the meadow-scale, though the path analysis provided evidence that elevation and distance to navigation channels indirectly influenced carbon densities.

The path model highlights the indirect role of the environment on carbon storage, though the effects of elevation and distance to channels are difficult to interpret (Fig. 6). We expected higher-elevation sites closer to channels to host coarser sediments and more historically diverse seagrass assemblages with less consistent cover in response to mild disturbance. According to the intermediate disturbance hypothesis, moderate levels of environmental stress should favor greater diversity (Menge and Sutherland 1987). Our path model only partially supported this hypothesis, where distance to channel negatively affected variability, but elevation negatively influenced species richness and distance to the channel positively affected species richness. It is plausible that species richness decreases with elevation and increases with distance from channels because a particular species is better adapted to shallower, more disturbed environments, i.e., H. wrightii. Alternatively, elevation and distance to channels may not accurately represent wave energy and disturbance, respectively, because the elevation range was relatively small (− 0.1 to − 3.0 m) and boating activity commonly occurs outside of the navigation channels in Cedar Key. The path model was useful in exploring hypothesized relationships between several drivers that interact to affect carbon densities over different time scales but can and should be revised as our understanding of meadow-scale drivers of carbon storage improves.

The environmental variables we measured were likely too coarse to capture small-scale differences in hydrodynamic intensity mediated by seagrass cover, which closely correlates with sediment and carbon accretion rates (Zhu et al. 2022). Given the strong relationship between percent fines and carbon densities, meadow-facilitated organic matter trapping and burial may drive carbon storage in Cedar Key, where meadows have a limited sediment supply and likely accumulate allochthonous organic carbon at relatively slow rates. In other systems, seagrass productivity and species-specific traits, such as lignin content, may be more relevant drivers of carbon storage. Collectively, our results indicate that (unmeasured) site-level differences in hydrodynamic intensity, perhaps caused by or correlated to seagrass cover and composition, can drive variability in fine-grained sediment deposition and carbon storage at the meadow-scale. Therefore, we expect scaling up measurements from a single location over even a relatively small meadow area (~ 25 km2) to estimate carbon stocks may be inaccurate based on the high variability in carbon storage we and others found at this spatial scale (Oreska et al. 2017).

Conclusions

Species diversity-carbon storage relationships are likely scale and context dependent so that the strength and direction of relationships vary across spatial and temporal scales (Thompson et al. 2018; Qiu and Cardinale 2020). In terrestrial systems, plant species diversity often directly correlates with soil carbon storage, and this relationship is mainly mediated through on-site primary production and litter deposition (Chen et al. 2018; Zhou et al. 2019; Xu et al. 2020). In contrast, seagrass sediment carbon accrues over long times scales and storage is limited by the ability of a meadow to trap and retain fine-grained sediments and allochthonous organic matter. In seagrasses, biodiversity-ecosystem function theory may operate over longer time scales where historical species diversity promotes carbon storage directly through canopy particle trapping efficiency and indirectly via succession. Successional changes in seagrass species assemblage enable persistent species that may most efficiently capture and store carbon, such as T. testudinum, to establish in previously disturbed areas (Williams 1990). Our findings suggest seagrass species diversity may benefit carbon storage over smaller spatial and longer temporal scales. Meadows dominated by persistent species may contain the greatest carbon stores, but pioneering species can facilitate the establishment of persistent species. Preserving seagrass diversity ensures long-term processes such as succession sustain meadow development and ultimately promote sediment carbon storage.