Spatial patterns of organic matter content in the surface soil of the salt marshes of the Venice Lagoon (Italy)

Salt marshes are crucial eco-geomorphic features of tidal environments as they provide important ecological functions and deliver a wide range of ecosystem services. Being controlled by the interplay between hydrodynamics, geomorphology, and vegetation, the contribution of both organic matter (OM) and inorganic sediments drives salt marsh vertical accretion. This allows marshes to keep pace with relative sea level rise and likewise capture and store carbon, making them valuable allies in climate mitigation strategies. Thus, soil organic matter (SOM), i.e. the organic component of the soil, plays a key role within salt marsh environments, directly contributing to soil formation and supporting carbon storage. This study aims at inspecting spatial patterns of OM in surface salt marsh soils (top 20 cm), providing further insights into the physical and biological factors driving OM dynamics that affect salt marsh survival and carbon sink potential. Our results reveal two scales of variations in SOM content in marsh environments. At the marsh scale, OM variability is influenced by the interplay between surface elevation and changes in sediment supply linked with the distance from the marsh edge. At the system scale, OM content distribution is dominated by the gradient generated by marine and fluvial influence. The observed variations in SOM are explained by the combination of inorganic and organic input, preservation conditions, and sediment grain size. Our results highlight the importance of marshes as carbon sink environments, further emphasising that environmental conditions within a tidal system may generate strongly variable and site-specific carbon accumulation patterns, enhancing blue carbon assessment complexity.

1 Introduction

Salt marshes may be defined as halophytic grasslands or dwarf brushwoods encroaching sediment platforms bordering saline tidal water bodies (Beeftink, 1977; Silvestri and Marani, 2004). Salt marshes lie at the upper margin of the intertidal environment, and they are regularly flooded by the tides and occasionally by storm surges (Allen, 2000; Friedrichs and Perry, 2001; D'Alpaos and D'Alpaos, 2021; D'Alpaos et al., 2007; D'Alpaos, 2011). They occur worldwide, particularly in middle and high latitudes on low-energy coasts, from microtidal to macrotidal regimes (Allen and Pye, 1992). Coastal wetlands, situated at the interface between terrestrial and marine ecosystems, often lie in some of the world's most densely populated areas and deliver a wide range of ecosystem services that contribute to human wellbeing, as they are simultaneously some of the most economically relevant and vulnerable ecosystems on Earth (D'Alpaos and D'Alpaos, 2021; Costanza et al., 1997; de Groot et al., 2012). Salt marshes provide a critical habitat for numerous species, protect coastal areas by buffering floods and wave action, improve water quality, support commercial fisheries, provide recreational opportunities, and contribute to climate regulation by acting as carbon (C) sinks (e.g. Barbier et al., 2011; Boesch and Turner, 1984; Costanza et al., 1997; Lefeuvre et al., 2003; Temmerman et al., 2013; Fairchild et al., 2021; Möller et al., 2014).

However, salt marsh ecosystems are exceptionally vulnerable to the effects of climate changes and increasing human pressure such as accelerating relative sea level rise, declining sediment supply, and coastal land reclamation (D'Alpaos et al., 2011; Gedan et al., 2009; Kirwan et al., 2010; Marani et al., 2007; Morris et al., 2002; Mudd et al., 2004; Mudd, 2011; Temmerman et al., 2005; Silvestri et al., 2018; Breda et al., 2022), having lost between 25 % and 50 % of their global historical coverage (Mcowen et al., 2017).

Salt marshes are intrinsically dynamic environments, and their evolution importantly depends on complex feedbacks between physical and biological components. The interplay between hydrodynamics, geomorphology, and vegetation controls marsh vertical accretion, which has allowed marshes to keep pace with relative sea level rise over thousands of years (D'Alpaos et al., 2007; Marani et al., 2007; Fagherazzi et al., 2012; Mudd et al., 2009, 2004; Allen, 2000; Brückner et al., 2019). Vertical accretion in tidal marshes is driven by the contribution of both organic matter (OM) and inorganic sediments (Mudd et al., 2009; Fagherazzi et al., 2012; Nyman et al., 2006; Neubauer, 2008). Tides, wind waves, and storms are the main mechanisms for sediment transport, deposition, or mobilisation on the marsh platform, shaping marsh morphology, which in turn influences water depth and flow velocity. Halophytic plants, spatially organised in characteristic patches related to soil elevation and hydroperiod (Bertness and Ellison, 1987; Pennings and Callaway, 1992; Silvestri et al., 2005; Silvestri and Marani, 2004), contribute considerably to marsh accretion, influencing the water flow, enhancing sediment settling, damping sediment resuspension, and producing biomass, which results in OM accumulation (D'Alpaos et al., 2007; Marani et al., 2006; Mudd et al., 2009). In some cases, OM contribution to marsh surface vertical accretion can be much greater than that of sediment deposition (Ewers Lewis et al., 2020; Morris et al., 2016; Allen et al., 2021). The organic material that helps to build marsh elevation is likely a combination of in situ production of belowground root tissue inserted into the sediments (Craft et al., 1993; Day et al., 1999) and autochthonous or allochthonous organic materials that are deposited over the surface (Nyman et al., 2006; Mudd et al., 2009; Ewers Lewis et al., 2019; Mueller et al., 2019). Surface litter produced during the annual cycles of plant growth and decay settles on the ground and is trapped within the inorganic sediments deposited. A variable proportion of the salt marsh soil organic matter (SOM) (Drexler et al., 2020; Middelburg et al., 1997; Mueller et al., 2019) has an allochthonous source, deriving from suspended particulate organic matter, and is often adsorbed on mineral matter and on estuarine and marine phytoplankton, microphytobenthos, and non-local macrophyte litter carried to the marsh surface by waves and tides. Furthermore, tidal flooding inhibits microbial aerobic activity and slows down decomposition, fostering C accumulation in marsh soils (Keuskamp et al., 2013; Mueller et al., 2018; Kirwan et al., 2014; Morris et al., 2016). As a result of these dynamics, the C captured through plant photosynthesis is buried and preserved as soil organic carbon (SOC) and may be locked away from the atmosphere over centennial to millennial timescales (Perillo et al., 2009; Duarte et al., 2005). This process allows salt marsh environments to act as carbon sinks, serving as natural or artificial reservoirs that accumulate and store carbon-containing compounds, thereby helping to offset the effects of greenhouse gas emissions on the Earth's climate (Watson et al., 2000). The C sink function of vegetated coastal ecosystems, including salt marshes, mangrove forests, and seagrass meadows, has been increasingly recognised in recent years, and the term “blue carbon” was coined to indicate the C sequestered in these ecosystems, with a potential role in climate change mitigation (Chmura et al., 2003; Duarte et al., 2005; McLeod et al., 2011; Macreadie et al., 2019; Nellemann et al., 2009).

Recent models of intertidal system dynamics predict marsh evolution incorporating both physical and ecological processes (e.g. D'Alpaos et al., 2007; Kirwan and Murray, 2007; Marani et al., 2007, 2013; Morris et al., 2002; Mudd et al., 2004, 2009). However, to refine current representations of OM accumulation, a deeper understanding is needed of crucial biogeomorphological factors driving OM dynamics in tidal environments. Drivers of SOM variability in salt marshes are innately difficult to identify, as SOM content is the result of many complex processes acting simultaneously. The study of the interplay between these processes is complicated by the fact that they occur over a wide range of spatial scales (Marani et al., 2006; Ewers Lewis et al., 2020).

While, at the global scale, climatic parameters appear to be important drivers for C sequestration and storage (Chmura et al., 2003), at the system scale, the position within the gradient generated by marine and fluvial influence was observed to be a key predictor of organic content, with the environmental conditions resulting from the interplay between physical and biological factors (Kelleway et al., 2016; Van de Broek et al., 2016). Geomorphological factors and vegetation characteristics appear to be important drivers of soil organic content variability when considering finer spatial scales (e.g. Chen et al., 2016; Leonard et al., 2002; Roner et al., 2016). For instance, elevation and hydroperiod importantly affect vegetation characteristics, organic matter supply, and sediment deposition (Chmura et al., 2003) other than microbial community and organic matter preservation conditions (Kirwan et al., 2013; Marani et al., 2006; Mudd et al., 2009; Yousefi Lalimi et al., 2018). The latter are also influenced by sediment type and grain size, which are additionally associated with different concentrations of suspended particulate organic matter (Kelleway et al., 2016; Van de Broek et al., 2016; Saintilan et al., 2013). Therefore, sediment grain size was pointed out as a key predictor of organic content by Ford et al. (2019) and Kelleway et al. (2016), whereas Van de Broek et al. (2016) reported variations in suspended sediment concentration and stable particulate OC as the main drivers of observed differences in organic storage. The vegetation type was also observed to significantly affect soil organic content at different scales (Ewers Lewis et al., 2020; Saintilan et al., 2013). Above- and belowground production greatly varies across and within plant species, both in terms of quantity and quality, the latter determining different degrees of decay resistance (Scarton et al., 2002; Stagg et al., 2018). In addition, Van de Broek et al. (2018) show that SOC from allochthonous sources (e.g. particulate organic matter, phytoplankton, microphytobenthos, and non-local macrophyte litter) may be a main component of SOC stock preserved in marsh sediments, and this finding could significantly impact organic carbon sequestration assessments (Van de Broek et al., 2018).

Here we aim to investigate spatial patterns of OM in salt marsh soils, providing further insights into the physical and biological factors driving OM dynamics, which affect salt marsh survival and C sink potential. Toward these goals, we analysed soil organic content and SOC stock for the surface soil layer (0–20 cm) in 10 salt marshes of the Venice Lagoon from 60 sediment cores together with different variables including soil and morphological and vegetation characteristics. The choice to analyse surface marsh soil for assessing spatial patterns of soil organic content is driven by the need to capture the layer most directly influenced by current environmental variables. By contributing to a better understanding of OM distribution and its controls in tidal marshes, our study has the potential to improve spatially explicit biomorphological models of salt marsh evolution, allowing the evaluation of marsh resilience and associated C storage as an ecosystem service. Improved representation of SOM dynamics in biomorphodynamic models will provide critical information for management and conservation strategies, favouring the identification of optimal locations and actions to direct protection and restoration efforts.

Figure 1 Location of the study areas within the Venice Lagoon (Italy). Insets show the position of the 30 m long transects (orange lines) within each study area located perpendicularly to the main marsh edge: MI, CA, PA, SA, SF, SE, VB, CO, FO, and CV (Map data: © 2022 Google-Landsat/Copernicus, Maxar Technologies).

2 Material and methods

2.1 Study site

The Venice Lagoon, situated in the northwestern Adriatic Sea (Italy), is the largest brackish tidal waterbody within the Mediterranean Sea. It is a shallow, microtidal, back-barrier system, characterised by a semidiurnal tidal regime with an average range of 1.0 m and a mean water depth over tidal flats of about 1.5 m, and is connected to the Adriatic through three inlets known as Lido, Malamocco, and Chioggia (Fig. 1).

The Venice Lagoon is also highly influenced by meteorological forcings, with the southeasterly Sirocco wind responsible for severe storm surges flooding the city of Venice, northeasterly Bora wind responsible for large wind waves (Tommasini et al., 2019). The present-day morphology of the Venice Lagoon is strongly affected by the intense human pressure occurring over the last centuries, going from the diversion of the main rivers, which are now debouching out of the lagoon, to the more recent construction of jetties at the inlets and the excavation of large navigable channels, to the rash groundwater exploitation in the 20th century. The above interventions have seriously modified sediment supply and dynamics, local hydrodynamics and morphology, and rates of relative sea level rise, favouring sediment starvation and the deepening of the lagoon (Zecchin et al., 2009; Brambati et al., 2003). The recent activation of storm surge barriers, designed to prevent flooding of the city of Venice, was proved to further alter the lagoon hydrodynamics (Mel et al., 2021) contributing to the deepening of the tidal flats and reducing salt marsh sedimentation (Tognin et al., 2021, 2022). Salt marsh areas of the Venice Lagoon have dramatically shrunk in the last centuries, with a decrease in extension from about 180 km² in 1811 to about 43 km² in 2002 (Carniello et al., 2009; D'Alpaos, 2010; Tommasini et al., 2019; Finotello et al., 2023). In addition, historical river diversions have significantly reduced freshwater inputs into the lagoon, thereby impacting water salinity and vegetation characteristics. Freshwater inputs currently flow into the lagoon through 12 main tributaries distributed along the landward boundary of the lagoon, with a mean annual contribution of about 35 m³ s⁻¹ and a peak discharge of 344 m³ s⁻¹ (Zuliani et al., 2005). Spatial and temporal variability in salinity in the Venice Lagoon is additionally influenced by groundwater inputs (Gieskes et al., 2013). Estimates of the volume of underground freshwaters entering the lagoon floor vary widely, from 15 % of total freshwater flow to more than 100 % (Zirino et al., 2014). Salinity levels in the Venice Lagoon today vary from approximately 20 PSU at the northeastern mainland edge to about 34–35 PSU at the three inlets (Zirino et al., 2014).

The study sites are located in 10 salt marshes of the Venice Lagoon, at variable distances from the inlets (Fig. 1). In the northern lagoon, the Sant'Erasmo (SE), San Felice (SF), and Saline (SA) salt marshes border large tidal channels departing from the Lido inlet. These marshes host exclusively halophytic vegetation, mainly constituted by Salicornia veneta, Limonium narbonense, and Sarcocornia fruticosa, associated with Spartina maritima, Juncus maritimus, Puccinellia palustris, Inula crithmoides, Suaeda maritima, and Arthrocnemum macrostachyum (Silvestri, 2000; Silvestri et al., 2005; Yang et al., 2020). In contrast, the Pagliaga (PA) and Campalto (CA) sites are located at the lagoon–mainland boundary in the northern lagoon, with PA being close to the estuary of the spring water Dese River which debouches into the lagoon, carrying an average 65 m³ s⁻¹ of freshwater discharge and negligible sediment supply. Freshwater inputs at the PA site maintain a near-freshwater plant community dominated by Phragmites australis, together with the halophytic species Juncus maritimus, Salicornia veneta, and Halimione portulacoides (Silvestri, 2000). The CA marsh, whose main edge faces a shallow tidal flat exposed to Sirocco wind, hosts halophytic species dominated by Limonium narbonense, associated with Sarcocornia fruticosa, Spartina maritima, Salicornia veneta, and the scarce Suaeda maritima, Triglochin maritima, and Juncus gerardii. In the southern part of the Venice Lagoon, the Mira (MI) marsh is close to the landward boundary, and it is characterised by the presence of Phragmites australis on the marsh edge, together with the halophytic species Juncus maritimus, Puccinellia palustris, Limonium narbonense, and Sarcocornia fruticosa. Canale Virgilio (CV) and Fossei (FO) are located within the marsh belt in front of the Malamocco and Chioggia inlets and host halophytic vegetation consisting mainly of Limonium narbonense, Sarcocornia fruticosa, and Puccinellia palustris, together with Juncus maritimus, Spartina maritima, Salicornia veneta, and Suaeda maritima. The Conche (CO) salt marsh fringes the mainland and faces the wide subtidal flat that occupies the central southern Venice Lagoon, being exposed to Bora wind. CO hosts halophytic species dominated by Sarcocornia fruticosa, Suaeda maritima, Inula crithmoides, and Halimione portulacoides. Finally, the Valle di Brenta (VB) salt marsh, at the southern end of the lagoon, is separated from the mainland by a channel and is characterised by the dominance of Juncus maritimus, Sarcocornia fruticosa, Limonium narbonense, and Puccinellia palustris, with the presence of Triglochin maritima. The marsh areas under study were already mapped by Sebastiano Alberti in 1611, and, considering an overall accretion rate of about 2.0–3.0 mm yr⁻¹, the study deposits were probably accumulated over the past century (Tommasini et al., 2019).

2.2 Experimental design and sample analysis

Six sediment cores were collected along transects from the marsh edge to the inner part of each of the 10 salt marshes for a total of 60 cores. The selected sites are distributed across the Venice Lagoon in order to represent the different environmental conditions typical of the system (Fig. 1). Transects typically start on the marsh edge facing a channel and therefore show a tide-dominated sedimentation. In contrast, the CA and CO transects face tidal flats and are exposed to energetic wind waves. At each site, cores were collected at an increasing distance from the marsh edge (0, 2.5, 5, 10, 20, and 30 m) so as to represent morphological and vegetation zonation. The position of each core was surveyed using a GPS (Leica GS16) in RTK mode, and vegetation characteristics were recorded within a 1 m × 1 m quadrat by the Braun-Blanquet method (Braun-Blanquet, 1964) and registering species cover percentages (Table S3 in the Supplement). The Shannon diversity index was used to measure the diversity of species in each 1 m × 1 m quadrat.

Soil samples were taken every 5 cm up to the depth of 20 cm (0, 5, 10, 15, and 20 cm) from each core, collected through a Eijkelkamp auger, 3 cm diameter corer (Howard et al., 2014), and subsamples were oven-dried at 60 °C for 48 h or to constant weight and prepared for different analyses, including OM and C content, soil density, and grain size distribution. The percent of OM of each sample was determined through a loss-on-ignition (LOI) procedure by combusting ∼ 2 g of ground sediment in a muffle furnace at 375 °C for 16 h of ground samples (Ball, 1964; Frangipane et al., 2009; Roner et al., 2016). The difference in weight between pre- and post-treatment provided the OM content as a weight percentage. The percent of organic C was calculated from LOI using the following local conversion equation built on a subsample analysed with an elemental analyser (Puppin, 2023; Puppin et al., 2023b):

$$\begin{array}{}\text{(1)}& \mathrm{OC}=(\mathrm{0.37}\pm \mathrm{0.01})\mathrm{LOI}+(\mathrm{0.0024}\pm \mathrm{0.0004}){\mathrm{LOI}}^{\mathrm{2}}.\end{array}$$

Sediment dry bulk density (DBD) was calculated from the water content as the difference in weight between wet and dry samples and between organic and inorganic fractions, according to Kolker et al. (2009), by assuming a water density of 1.02 and a mineral and organic sediment density of ρ_i = 2.6 g m⁻³ and ρ_o = 1.2 g m⁻³, respectively, knowing the inorganic sediment and OM percentage from LOI analysis (Kolker et al., 2009). The product of the percent of OC and DBD determined soil carbon density (SCD). Soil carbon stocks in the top 20 cm of the soil were calculated based on Howard et al. (2014), summing the C content in each soil interval along the core in order to better represent OM content variability within the depth interval. Inorganic particle size distribution analyses were performed, after the removal of the OM through a treatment with 35 % hydrogen peroxide (H₂O₂) for 36 h, using laser granulometry (Mastersizer 2000 version 5.40, Malvern Instruments).

As to the analysis of OM content, we considered sediments up to the depth of 20 cm because we wanted to investigate the soil layer directly impacted by current vegetation. Indeed, according to Trumbore (2009), in most vegetated ecosystems the majority of underground plant biomass and microbial activity exist within the top 20 cm of soils. Conversely, in order to analyse sediment characteristics with reference to current hydrodynamics and morphological dynamics, we limited the reference depth for grain size analysis to 5 cm, considering the accreting behaviour of salt marshes.

Values of water salinity for each study area were calculated as the mean of the last 5 years of quarterly measurements from the CTD multiparameter probes of the Veneto Region Environmental Protection Agency (ARPAV, 2024) (Table S5).

2.3 Statistical analysis

Statistical analyses were performed using MATLAB R2021a. Analyses included the use of the Kendall test, the Kruskal–Wallis H test, and post hoc median comparisons. The Kendall test is a non-parametric, rank-based method, aimed to determine whether there exists a monotonic relationship between two variables of interest. The value of the τ coefficient ranges from 1 to −1, indicating a positive or a negative association, respectively. The Kruskal–Wallis test is a non-parametric version of the classical one-way ANOVA using ranks of data to compute the chi-square statistic and compare the medians of the groups of data to determine if the samples come from the same population. When the Kruskal–Wallis test showed a significant difference between groups, a multiple comparison test was used to determine which pairs of means were significantly different. For the statistical tests, we included all available data in the analysis without attempting to identify and eliminate possible outliers.

Figure 2 Distribution of surface sediment variables analysed and the surface elevation profile along the transects in the northern lagoon: organic matter content (LOI %; mean value and standard deviation in top 20 cm), dry bulk density (g cm⁻³; mean value and standard deviation in top 20 cm), vegetation cover (%), grain size distribution (D₅₀ µm and sand–silt–clay percentage; mean value and standard deviation in top 5 cm), and surface elevation ($\mathrm{m}\mathrm{a}.\mathrm{m}.\mathrm{s}.\mathrm{l}.$). The dashed line represents the linear regression between the percent of organic matter values.

Figure 3 Distribution of surface sediment variables analysed and the surface elevation profile along the transects in the southern lagoon: organic matter content (LOI %; mean value and standard deviation in top 20 cm), dry bulk density (g cm⁻³; mean value and standard deviation in top 20 cm), vegetation cover (%), grain size distribution (D₅₀ µm and sand–silt–clay percentage; mean value and standard deviation in top 5 cm), and surface elevation ($\mathrm{m}\mathrm{a}.\mathrm{m}.\mathrm{s}.\mathrm{l}.$). The dashed line represents the linear regression between the percent of organic matter values.

3 Results

3.1 Distribution of surface sediment variables along the transects

We examined the distribution of the surface sediment OM content, DBD, grain size distribution, and the vegetation cover along the surface elevation profile of study transects in the northern (Fig. 2) and southern (Fig. 3) lagoon.

The topographic GPS field surveys (Figs. 2 and 3d, g, j, m, and p) show that marsh surface elevations along the transects range between about 0.20 and 0.60 $\mathrm{m}\mathrm{a}.\mathrm{m}.\mathrm{s}.\mathrm{l}.$ CO is the highest site, with an average elevation along the transect of 0.55 $\mathrm{m}\mathrm{a}.\mathrm{m}.\mathrm{s}.\mathrm{l}.$, while FO is the lowest site, with an average elevation of 0.24 $\mathrm{m}\mathrm{a}.\mathrm{m}.\mathrm{s}.\mathrm{l}.$ Surface elevations are generally higher near the marsh edge and adjacent to the channel levees, and they decrease toward the inner marsh. Some sites show a pronounced levee morphology (i.e. SF, SA, and FO), while at other sites the latter is less distinct due to irregular profiles (i.e. MI, CV, and VB) or flattened profiles (i.e. PA and SE). The CA and CO transects, which face a tidal flat and are exposed to wave energy, are characterised by higher surface elevations toward the inner marsh, with a convex shape profile (Tognin, 2022).

The vegetation cover is generally high (never below 50 %), with no detectable trends along the transects.

OM content in surface soils (mean value in top 20 cm) ranges between 4 % and 30 % (green lines in Figs. 2 and 3b, e, h, k, and n) and shows important differences both within and among study areas. CA, FO, PA, VB, CV, and MI show generally higher organic content, with average OM content values along the transects ranging between 15 % and 21 %, than CO, SF, SE, and SA, with average OM content values along the transects ranging between 5 % and 7 %. Trends in OM content along the transects suggest an increase toward the inner marsh in most of the study sites (e.g. CA, PA, SF, SE, CV, CO, and VB). The DBD in surface soil (grey lines in Figs. 2 and. 3b, e, h, k, and n) shows a distribution mirrored to that of OM content, with values ranging between 0.28 and 1.15 g cm⁻³, displaying a more or less marked decrease toward the inner marsh in most of the study sites. The highest values of DBD are found at the CO, SE, SF, and SA sites, with average values of 1.02, 0.81, 0.94, and 0.80 g cm⁻³, respectively.

The grain size of the surface inorganic sediment fraction (Figs. 2 and 3c, f, i, l, and o) ranges from clay to very fine sand, with a D₅₀ (mean value in top 5 cm) ranging between 7 and 28 µm, with a mean value of 14 µm. Finer median grain size values are found at the PA site (average D₅₀ value of 5 µm), which is characterised by higher percentages of clay, whereas coarser median grain size values are found at the SF and CO sites (average D₅₀ value of 19 and 20 µm, respectively), with the latter being characterised by some shell layers on the marsh surface. As to the distribution of the median grain size, D₅₀, along the transects, we observe coarser grains along the marsh edges (i.e. along the tidal channel levees) at SF, SA, and VB. CA and CO median grain sizes show a smooth trend with coarser particles in the central part of the transects. The SE transect shows slightly coarsening sediments toward the inner marsh, bordered by a tidal flat.

Figure 4 Effects of the distance from the marsh edge on OM, DBD, and grain size: organic matter content (LOI %; mean value in top 20 cm) (a), dry bulk density (g cm⁻³; mean value in the top 20 cm) (b), and sediment median grain size distribution (D₅₀ µm; mean value in the top 5 cm) (c). In the box plots the central mark indicates the median, and the bottom and top edges of the box indicate the 25th and 75th percentiles, respectively. The whiskers extending to the most extreme data points are not regarded as outliers. Swarm plots show single values (circles).

Figure 5 Vegetation distribution and influence on SOM. Species cumulative percent cover at each station along the transect (a). Organic matter content in marsh surface soils (LOI %; mean value in top 20 cm) as a function of dominant species (c). Organic matter content in marsh surface soils (LOI %; mean value in top 20 cm) (d) and OC stock in the top 20 cm (e) as a function of the Shannon diversity index H in the 1 m × 1 m plot around core sites. Limonium narbonense roots exposed by marsh edge erosion at the SE site (b) and Phragmites australis roots and stems within PA site soil (f). Panel (c) includes mean values (horizontal-line markers) and standard deviations. Swarm plots show single values (closed grey circles). Panels (d) and (e) include the median and standard deviation of binned values. Open circles represent single values, and bold lines represent their linear regressions (OM = 9.91 +H4.43, R² = 0.0679, p value = 0.0444; C stock (20 cm) = 0.39 +H0.08, R² = 0.0982, p value = 0.0147). INU: Inula crithmoides. AST: Aster tripolium. LIM: Limonium narbonense. SAL: Salicornia veneta. SAR: Sarcocornia fruticosa. SPA: Spartina spp. HAL: Halimione portulacoides. SUA: Suaeda maritima. PUC: Puccinellia palustris. JUN: Juncus maritimus. PHR: Phragmites australis. TRI: Triglochin maritima.

3.2 Dependence of sediment variables on the distance from the marsh edge

Considering the data from all our study areas, OM content in surface soils was significantly positively correlated with the distance from the marsh edge (Kendall's tau test; τ = 0.1551, p value = 0.0002) (Fig. 4a), whereas the DBD was significantly negatively correlated with the distance from the marsh edge (Kendall's tau test; τ = −0.2302, p value < 0.0001) (Fig. 4b). A significant negative correlation was also found between the median sediment grain size (D₅₀) on the marsh surface and the distance from the marsh edge (Kendall's tau test; τ = −0.0924, p value = 0.0269) (Fig. 4c). The distance from the marsh edge was negatively correlated to surface elevation in our study case (Kendall's tau test; τ = −0.1204, p value = 0.0040); thus, as expected, DBD and D₅₀ were positively correlated to marsh surface elevations (τ = 0.1387 and 0.1503, p value = 0.0003 and 0.0001, respectively). Nevertheless, no significant relationship was found between the OM content and surface elevation.

Figure 6 Surface sediment variables at different study sites: organic content (LOI %; mean value in top 20 cm) and water salinity at each study area (a), dry bulk density (g cm⁻³; mean value in top 20 cm) (b), OC stock in the top 20 cm (g cm⁻²) and SOC density (g cm⁻³; mean value in top 20 cm) (c), grain size distribution (D₅₀ µm; mean value in top 5 cm) (d). In the box plots the central mark indicates the median, and the bottom and top edges of the box indicate the 25th and 75th percentiles, respectively. The whiskers extending to the most extreme data points are not regarded as outliers. Swarm plots show single values (circles). Different letters above the box plots indicate significant differences on the basis of the Kruskal–Wallis test and post hoc multiple comparison test (when two groups do not show a significant difference, they have at least one letter in common).

3.3 Vegetation distribution and influence on SOM

The vegetation surveys revealed that Sarcocornia fruticosa, Limonium narbonense, Puccinellia palustris, and Juncus maritimus are among the most abundant species in terms of percent cover along the analysed transects (Fig. 5a). Sarcocornia fruticosa exhibits an increasing trend in cover with distance from the marsh edge, peaking at approximately 5 m from the edge before declining. Conversely, the cover of Limonium narbonense and Puccinellia palustris increases towards the inner marsh. The cover of Juncus maritimus varies along the transect. At the lagoon scale, Phragmites australis, typically found in brackish water environments, was observed only at two sites, PA and MI (Fig. S1 in the Supplement). Additionally, Puccinellia palustris was more abundant at sites located in the southern part of the lagoon (Fig. S1). We found no relationship between OM content and vegetation cover. When considering the dominant species, OM content of marsh surface soil displayed a wide variability (Fig. 5c). Higher mean organic content was observed in the presence of Limonium narbonense (Fig. 5b) and Phragmites australis (Fig. 5f) as dominant species. Although highly variable, organic content showed a significant positive correlation with vegetation species diversity (Kendall's tau test for Shannon diversity index H and OC stock in top 20 cm: τ = 0.1931, p value = 0.0300. Mean SCD in top 20 cm: τ = 0.1818, p value = 0.0412) (Fig. 5d and e).

3.4 Variability in sediment characteristics among study sites

As to the variability in marsh soil characteristics in different areas, the Kruskal–Wallis analysis indicated that there was a significant difference between OM content at the 10 sites (p value = < 0.0001), and the post hoc test confirmed that CO, SA, SE, and SF have significantly lower mean OM content in surface soil (Fig. 6a). Study sites were characterised by different mean water salinity ranging between 24.3 and 32.4 PSU, which revealed a significant negative relationship with OM content (Kendall's tau test; $\mathit{\tau }=-\mathrm{0.3524}$, p value = 0.0002). CO, SA, SE, and SF also show significantly higher values of DBD (Kruskal–Wallis; p value = < 0.0001) (Fig. 6b).

Mean SOC density in the top 20 cm (product of the percent of OC, which is function of the OM and DBD determined) ranged between 0.014 and 0.026 g cm⁻³ with lower values at SA, SF, and SE (Fig. 6b). Carbon stock in the top 20 cm ranged between 0.21 and 0.66 g cm⁻², with lower values at SA, SF, and SE (Fig. 6c). CO, which showed a mean organic matter percentage similar to SA, SF, and SE, displays, in contrast, higher values of SCD and C stock compared to SA, SF, and SE, due to particularly high DBD values at this site.

Median sediment grain size (D₅₀) on the marsh surface was not equal at different study sites (Kruskal–Wallis; p value = < 0.0001), although it did not show clear distribution patterns (Fig. 6d).

4 Discussion

The sources of OM content in salt marsh soils are influenced both by local and non-local processes. Firstly, OM is the result of in situ production of belowground (root, rhizome, and tuber tissue) (Craft et al., 1993; Rybczyk and Cahoon, 2002) and aboveground biomass, thus directly depending on the local primary production. OM content is also affected by the accumulation of organic material produced at other sites (Nyman et al., 2006; Mudd et al., 2009; Ewer Lewis, 2019; Mueller et al., 2019), which is transported eventually deposited on the marsh surface by hydrodynamic processes (i.e. tides and waves) acting at larger spatial scales . Both autochthonous and allochthonous organic materials, once part of marsh soil, are also affected by decomposition resulting from local topographic, sedimentological, and environmental conditions (Chen et al., 2016). Due to the intricate interaction of these local and non-local dynamics, two main spatial scales in OM variations can be identified: the marsh scale (metres to tens of metres) and the system scale (ranging from kilometres, encompassing the entire lagoon or estuary). In the following discussion, we will first examine evidence of OM variations in our results at the marsh scale by considering trends along the transects. Subsequently, we will shift our focus to the system scale by comparing results among study sites.

At the marsh scale, our results show that OM content in surface soils displays a significant trend with the distance from the marsh edge, with OM content generally increasing towards the inner marsh (Fig. 4a), consistently with previous findings (e.g. Chen et al., 2016; Leonard et al., 2002; Roner et al., 2016). This overall trend in OM with distance from the edge is influenced by various processes acting at the marsh scale, among which the interaction between sediment delivery and local topography plays a pre-eminent role. Suspended material is primarily delivered onto the marsh platform through inundation by overbank flow along tidal channels (Bayliss-Smith et al., 1979) and by apical flow at creek heads (Torres and Styles, 2007). As soon as the flow reaches the vegetated marsh platform, current velocities and turbulent energy rapidly decrease (D'Alpaos et al., 2007; Mudd et al., 2010), thus promoting the deposition of more abundant sediments during the flooding initial phase in close proximity of the marsh edge (Christiansen et al., 2000; Roner et al., 2016). As a result of the larger deposition close to the edge, inner marshes generally present slightly lower elevation than that of the marsh margin (Figs. 2 and 3). Lower elevations promote the persistence of an anaerobic environment, which slows down OM decomposition by reducing microbial respiration (Halupa and Howes, 1995; Kirwan et al., 2013; Puppin et al., 2023a; Roner et al., 2016).

Progressive energy dissipation over the vegetated marsh platform also promotes selective material settling. Coarser, denser inorganic sediment is mainly deposited near the marsh edge, while the inner marsh receives a higher proportion of finer sediments, as supported by the observed grain size distribution along the transect (Fig. 4c), and less dense organic material (Leonard et al., 2002; Miller et al., 2022). In addition to the already larger proportion of organic material settling in the inner marsh, the supply of finer inorganic sediment (Fig. 4c) may also promote conditions favourable for OM preservation. This can occur due to the reduced oxygen exchange resulting from the lower porosity and drainage capacity of finer sediments, as well as their greater potential for protecting C from decay through organic mineral interactions and the formation of micro- or macro-aggregates (Kelleway et al., 2016).

Considering all available data along transects allows us to capture global trends and average out variabilities related to local conditions. However, analysing site-specific trends can offer an interesting perspective on drivers locally affecting OM dynamics. Overall, 7 out of 10 analysed transects show an increasing trend in OM with the distance from the marsh edge (CA, PA, SF, SE, CV, CO, and VB). At the SA, MI, and FO study sites the trend in OM content deviates from the average behaviour, and this can be attributed to the effects of local variability. At the SA site, observed OM content is very low (Fig. 3n), and its distribution pattern may be masked by the intrinsic variability in the measurements. At the MI site, Phragmites australis grows on the marsh edge, providing a contribution of OM that outcompetes that provided by the halophytic vegetation of the inner marsh (Fig. S1). At the FO site, we observed abundant beach-cast seagrass wracks on the marsh edge, which are likely transported from the extensive seagrass meadows located on the tidal flats adjacent to the FO area (Fig. 1f). These wracks may serve as an additional source of OM, locally influencing the OM trend with the distance from the marsh margin.

Conversely, in some transects where the organic matter trend aligns with the average increasing trend with distance from the marsh margin, we observe particular cases in the behaviour of grain size and topographic variables, which occasionally deviate from the average trend. For example, at the inner end of the SE transect, we can observe an unexpected slight increase both in sediment grain size (Fig. 2l) and in topographic elevation (Fig. 2m). This may be related to the presence of a tidal flat at the inner border of the SE marsh (Fig. 1j), which can represent an additional source of sediment supply. Two notable exceptions in terms of marsh topography are also represented by the CA and CO marshes, which face a tidal flat (Fig. 1b and h) and are exposed to energetic wind waves. As a consequence of the influence of wind waves, the elevation profiles exhibit a subtly convex shape, with the disappearance of the raised margin and a slight inward shift (i.e. between 5 and 10 m from the edge) of the maximum elevation along with the locations of higher median grain size values. However, these variations do not appear to affect the overall increasing trend of OM content towards the inner marsh at the SE, CA, and CO sites.

We also observe variations and patterns of OM at the system scale. The position within the gradient generated by marine and fluvial influence was previously observed to be a key predictor of OM content (e.g. Van de Broek et al., 2016; Ewers Lewis et al., 2020; Kelleway et al., 2016; Macreadie et al., 2017). At the lagoon scale, our results show lower mean OM content in surface soil in areas which are directly affected by marine influence, being closer to the inlets or along the main channels branching from them (CO, SE, SF, and SA) (Figs. 1 and 6a). Conversely, higher OM contents are observed at sites closer to the mainland, e.g. CA, MI, and PA. In agreement with this result, a significant negative relationship is observed between OM content and water salinity. However, it is unlikely that salinity as such directly controls soil organic content, as previous observations suggest an inverse relationship between soil salinity and decomposition (Hemminga et al., 1991; Wang et al., 2019). The effect of salinity on decomposition may be overcome by the co-occurring effects of other factors acting at different positions within the lagoon, such as vegetation characteristics, hydrodynamic conditions, sediment supply, and freshwater inputs. OM increase at less-saline sites is likely minimally related to the supply of already-stabilised organic suspended particles from terrestrial sources as suggested for other study areas (e.g. Van de Broek et al., 2016, 2018; Gorham et al., 2021; Omengo et al., 2016) because, in the Venice Lagoon, after historical river diversions, the fluvial supply of organic and inorganic material dramatically decreased. However, residual freshwater inputs, especially in terms of groundwater, can still locally reduce salinity levels and consequently affect vegetation characteristics, with usually increasing macrophyte biomass at lower salinity values (Hansen et al., 2017; Van de Broek et al., 2016).

A relationship between OM and grain size can also be observed at the lagoon scale. Higher values of median grain size are observed at sites closer to the inlet (i.e. SF and SE) or adjacent to first-order channels connected to them (SA) (Fig. 6d), whereas higher fractions of fine sediments are observed at the lagoon–mainland boundary near the Dese River mouth (PA) (Fig. 2). This pattern is consistent with the general grain size gradient observed in the Venice Lagoon, reflecting the typical pattern of decreasing hydrodynamic energy conditions from the inlets to the landward shore (Zonta et al., 2018). Considering the landward decreasing grain size gradient observed within the lagoon, enhanced C preservation capacity of fine sediments (Kelleway et al., 2016) may have a role in the observed organic content pattern. Furthermore, we may hypothesise that at sites where inorganic sediment inputs are greater, where hydrodynamic energy is higher, the soil organic fraction is proportionally lower.

Another potentially important factor controlling OM content is vegetation type, as different species exhibit varying biomass production and decomposition resistance (e.g. Van de Broek et al., 2016; Ewers Lewis et al., 2020; Ford et al., 2019; Saintilan et al., 2013; Yuan et al., 2020).

While our data do not allow for a full statistical analysis of the relationship between OM content and vegetation type, we can derive some interesting qualitative observations. We find no discernible trend in vegetation cover along the transects nor a significant relationship between OM content and vegetation cover.

Considering the relationship between SOM and dominant species, we observe that higher OM percentages in surface soil are not necessarily associated with dominant species having greater aboveground biomass (as reported in the literature; see Table S4). We speculate that the lack of a relationship between aboveground biomass and OM may be due to the continuous transport and mixing of locally produced litter by marsh flooding, weakening the effect of local aboveground biomass production. Moreover, previous studies indicate a major impact from belowground biomass, which, inserted into the sediments, directly contributes to OM content (Craft et al., 1993; Rybczyk and Cahoon, 2002). The highest mean organic contents are observed in the presence of Limonium narbonense and Phragmites australis as dominant species (Fig. 5c). Phragmites australis, characterised by high aboveground biomass and even higher belowground biomass, forms a dense and deep network of leathery stems, roots, and rhizomes  (Fig. 5f) (Moore et al., 2012; Scarton et al., 2002). Limonium narbonense, despite low aboveground biomass (Table S4), produces massive woody roots (Fig. 5b) (Brooks et al., 2021). The belowground biomass of these plants may importantly contribute to SOM content quantitatively and qualitatively, as belowground litter decomposition was observed to decline with increasing lignin content (Stagg et al., 2018; Puppin et al., 2023a). Interestingly, these species are more abundant in locations associated with higher values of SOM, such as the inner marsh for Limonium narbonense and brackish areas for Phragmites australis (Figs. 5a and S1). However, it was not feasible to directly measure or estimate the relationship between aboveground and belowground biomass and their effect on SOM content in this study.

In addition, soil organic content shows a significant positive correlation with vegetation species diversity (Fig. 5d and e) in agreement with Ford et al. (2019). Ford et al. (2016) found that plant species richness was one of the most important explanatory variables of root biomass, and Xu et al. (2020) suggested that species richness may increase biomass productivity due to multiple mechanisms including competition reduction, niche complementarity, selection effects, and biotic and abiotic facilitation.

Overall, the interplay of these dynamics at both marsh and system scales also impacts sediment bulk density. We observe higher soil densities along marsh edges and at sites such as CO, SE, SF, and SA, where organic content is generally lower and coarser sediments predominate (Figs. 4b and 6b). This observation aligns with previous research indicating that soil density is influenced by both organic matter content and grain size. Indeed, several studies have reported a significant and negative correlation between SOM and bulk density (e.g. Holmquist et al., 2018; Morris et al., 2016), while sand content has been shown to correlate positively with soil bulk density (Tanveera et al., 2016).

As a consequence of the variability in soil organic content and density both at the marsh and at the system scale, SOC density and C stock show significant variations across and within different salt marshes, enhancing the complexity of blue C assessment. Based on our estimates of mean SOC density in the top 20 cm, and considering an expected accretion rate of about 0.3 cm yr⁻¹ for salt marshes in equilibrium with relative sea level rise (RSLR) (Day et al., 1998), the average C accumulation rate for the salt marshes in the Venice Lagoon is estimated to be approximately equal to 69 $\mathrm{g}\mathrm{C}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$. Our results are consistent with the mean C accumulation rate in the Australian tidal marshes of 54.52 $\mathrm{g}\mathrm{C}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$ calculated by Macreadie et al. (2017) from 323 soil cores to the depth of 30 cm all around Australia and using a mean accretion rate of 0.21 cm yr⁻¹. However, we observed that, under the same accretion rate, estimated SOC from our study may result in C accumulation rates varying up to 50 % from one place to another because of the two-scale variability highlighted by our measurements. This underscores the need for careful consideration of local variability when assessing blue carbon sequestration and storage potential in wetland environments.

5 Conclusions

Overall, our data reveal two scales of variations in sedimentary OM content in salt marsh soils. One occurs at the marsh scale and is influenced by the interplay between surface elevation and sediment supply dynamics along the distance from the marsh edge. The other occurs at the system scale, dominated by the gradient generated by marine and fluvial influence.

The observations on morphological, sedimentological, and vegetational features emphasise the dynamic feedbacks between hydrological dynamics, sediment supply, surface elevation, and vegetation characteristics and allowed us to relate them to SOM spatial patterns.

Our results show that OM content generally follows an inverse relationship with the distance from the marsh edge and, broadly, is greatest in salt marshes subjected to freshwater inputs and characterised by the presence of fine-grained sediments. Although the interplay between physical and biological factors makes it difficult to unravel the direct cause–effect relationship, variations in allochthonous inputs, sediment grain size, and/or preservation condition may explain the observed variation in SOM. Our data on local vegetation do not allow us to detect plant community characteristics controlling overall variations in SOM. Hence, we highlight the need for further analyses of quantitative and qualitative characteristics of different tidal marsh vegetation types. However, species diversity was observed to be an important factor influencing SOM content.

Our data on the spatial distribution of SOM provide additional insights on organic matter dynamics in salt marshes and inform model representations of SOM accumulation, with the possibility to help improve the ability for biogeomorphological models to describe marsh response to the effects of climate change and anthropogenic perturbation. Moreover, our analyses emphasise the potential complications introduced by local variability in assessing the blue carbon sequestration and storage potential in wetland environments. This underscores the importance of carefully considering both the local variability and the factors influencing it. Finally, our findings may inform conservation strategies and restoration interventions, providing information on conditions promoting OM storage and preservation, such as the maintenance or recovery of freshwater inputs, the supply of finer sediments, and the enhancement of vegetation diversity.