An Integration of Numerical Modeling and Paleoenvironmental Analysis Reveals the Effects of Embankment Construction on Long‐Term Salt Marsh Accretion

There are still numerous uncertainties over the influence of anthropogenic interventions on salt marsh dynamics. This study uses the Ribble Estuary as a test case and an integrated approach of numerical modeling and paleoenvironmental analysis to investigate the contribution of embankment construction to long‐term marsh accretion. Accretion rates derived using optically stimulated luminescence dating (OSL) were combined with a multi‐proxy paleoenvironmental investigation on sediment cores extracted from the salt marsh, the mobile seafloor of the central Irish Sea and the river catchment area. These analyses provided a first evolutionary perspective on the Ribble Estuary preceding any management interventions. The paleoenvironmental analyses were then compared to simulations conducted using the hydrodynamic model Delft3D to investigate the effects of embankment construction on estuarine hydrodynamics and morphodynamics of the salt marsh over the period constrained by the OSL. The numerical simulations showed that embankments were responsible for an overall intensification of the ebb currents in the system which promoted sediment export. The paleoenvironmental analyses showed that the marsh has been accreting at a rate of 4.61 to 0.86 cm yr−1 over the last ca. 190 years and that the high sedimentation rate was caused by a naturally high rate of sediment supply. The model‐data integration showed that the effects of the embankment construction on sediment transport did not compromise the long‐term resilience of the salt marsh because of the high rates of sediment supply and the river dredging which enhanced the flood dominance of the tide near the tidal flat.

2 of 24 become colonized by vegetation (Reed, 1990). The resilience of salt marshes to sea-level rise is thought to be dependent upon sediment availability and rate of sediment supply to the estuarine accommodation space (FitzGerald et al., 2008;Ganju et al., 2017;Kirwan et al., 2010Kirwan et al., , 2016 and the effects of anthropogenic interventions for shoreline management on the hydrodynamic and morphodynamic feedbacks of the system (Lee et al., 2017;Li et al., 2020;Palmer et al., 2019;Pontee, 2014).
Salt marshes are generally resilient to sea-level rise when sediment supply and organogenic production are high enough to allow marsh accretion (Kirwan et al., 2010) or when salt marshes can migrate inland (Kirwan et al., 2016). Increasing suspended sediment concentrations promotes the long-term resilience of salt marshes by enhancing both vertical accretion (Temmermann et al., 2004) and lateral expansion (Willemsen et al., 2021). Willemsen et al. (2021) observed that, during mild weather, salt marshes that would experience retreat due to exposure to higher wave energy (between 0.10 and 0.15 m), are able to switch from a retreating to an expansional behaviour with an increase in sediment supply. However, the accumulation of such sediments on salt marsh surfaces depends on the age of the marsh, with young low marsh surfaces trapping sediments at a high rate until they reach an equilibrium level with mean high water, whereas high old marsh surfaces accrete at a slower rate comparable to the rate of sea-level rise (Temmermann et al., 2004). The minimum amount of suspended sediment concentrations to allow marsh survival under current rate of sea-level rise has been estimated to be ∼20 mg/L (Kirwan et al., 2010). Increasing rates of sea-level rise can compromise the stability of estuaries and salt marshes by causing marsh drowning and lateral retreat, as new accommodation space is created and, as a result, the amount of sediment inputs required to preserve marsh stability increases (FitzGerald et al., 2008;Ganju et al., 2017;Kirwan et al., 2010Kirwan et al., , 2016. Coastal management can significantly influence the accretion and expansion of salt marshes. Several studies around the world found that the presence of embankments alters the way sea-level rise affects the tidal signal by changing the resonance properties of the coastline (e.g., Carless et al., 2016;Holleman & Stacey, 2014;Pelling et al., 2013). These changes in the tidal signal affect sediment transport dynamics in estuaries and, consequently, sediment supply to coastal marshes (Guo et al., 2016;Palazzoli et al., 2020). It has also been observed that the construction of embankments can channelise sediments towards the estuarine accommodation space, inducing siltation of intertidal areas (Lee et al., 2017;Pontee, 2014).
Salt marshes in the Ribble Estuary, North-West England, have an extremely important environmental and economic value, as they are one of Europe's largest salt marsh systems (∼2,302 ha) and have been restored through managed realignment to provide flood protection (Tovey et al., 2009). Previous investigations undertaken in the estuary using bathymetric charts suggested that a very rapid accretion characterised the system in the 19th century, triggered by embankment construction started in 1810 AD, whereas a more natural accretion characterised the estuary in the 20th century following a stabilization of the system ( Van der Wal et al., 2002). Van der Wal et al. (2002) suggested that the construction of embankments caused an infilling of the estuary and accelerated the natural rate of accretion, concluding that anthropogenic interventions have outstripped those of changes in natural forcing on the morphological development in the area. There is no evidence of the direct effect of embankments on the hydrodynamics and morphodynamics of the system and there is no bathymetric data available for the period that precedes 1847, therefore no shift from a slower natural to a faster unusual accretion has been detected in the system through field data in the available literature ( Van der Wal et al., 2002). The sediment supply to the Ribble Estuary is thought to be predominantly marine (Lyons, 1997). The central Irish Sea is characterised by waves of mobile sediments that move up to 70 m per year (Van Landeghem et al., 2012). This high amount of mobile sediment has been linked to the legacy of sediments discharged through the melting and retreat of the Irish Sea Glacier (Scourse et al., 2021;Van Landeghem & Chiverrell, 2020) and Wright et al. (1971) suggested that these sediments are the major source of sediment supply to the Ribble estuary. Therefore, considering the great availability of sediment that can potentially reach the estuary, the infilling of the estuary and rapid accretion of the system might have alternatively been natural and been caused by a surplus of sediments available that infilled the estuary until the system reached equilibrium.
Despite of the numerous insightful studies on the subject, there are still many uncertainties over the effects of embankment construction on long-term salt marsh accretion with respect to natural sediment availability. This study uses the Ribble Estuary as a case study and an integrated approach of paleoenvironmental reconstruction and hydrodynamic modelling to investigate the effects that embankment construction undergoing in the estuary since 1810 AD had on the long-term accretion of Hesketh Out Marsh with respect to the high rates of natural 3 of 24 sediment supply. It specifically assesses how the embankments have affected the long-term sediment delivery to and accumulation on the salt marsh platform and whether they caused any changes in the long-term marsh accretion rate since the beginning of the formation of the tidal flat. The Ribble Estuary has been chosen as a test case because of the establishment of extensive tidal flats and marsh areas since the 19th century, high quantity of natural sediment availability and existing embankments. Hesketh Out Marsh has been chosen as a sampling site due to the environmental and economic importance of this site, as it has been managed realigned, after being reclaimed in 1980 AD, to provide more resilience to the coast against flooding (Tovey et al., 2009), making it essential to understand the factors that can influence the resilience of this site in the long term.
Sediment cores from Hesketh Out Marsh, located in the centre of the estuary, were recovered and analysed spanning the evolution from sand-flat to salt marsh, with time scales and accretion rates secured by optically stimulated luminescence dating (OSL). Sediment geochemistry (X-Ray fluorescence (XRF)) was analysed for the marsh cores, for samples collected from the Ribble catchment and for cores extracted from the adjacent seafloor of the Irish Sea to test a hypothesis that the sand-flat/marsh materials are predominantly Irish Sea in provenance. Changes in sediment accretion rates alongside down-core changes in particle size distribution, geochemistry and organic content were used to detect changes or shifts in sedimentation processes and depositional environment. A record of estuary management was then used to determine whether any change in accretion and sedimentation processes might be chronologically linked to the beginning of embankment construction. The hydrodynamic model Delft3D was employed to simulate sediment accumulation within the tidal flat-salt marsh system in scenarios that include (a) the existing embankments and (b) zero anthropogenic disturbance, to investigate the hydrodynamic and morphodynamic feedbacks that might cause changes in sedimentation when embankments are present, and further understand whether the presence of embankments might be responsible for a significant change in sedimentation. The evolutionary profile of the marsh was then compared to the rate of accretion, revealing whether any change in accretion rate can be connected to the natural evolution and stabilisation of the system and related changes in the mode of sediment supply.

Study Site
The Ribble Estuary is located on the Lancashire coast of North-West England (Figure 1a). The estuary is funnelshaped, tidally dominated and macro-tidal (Wakefield et al., 2011). It experiences moderate wave energy owing to the waves generated in the Irish Sea basin (Pye & Neal, 1994). Lyons (1997) suggested that sediments are introduced into the estuary mostly by tidal pumping, especially during high storm surges. It is thought that the formation of the extensive intertidal sand-silt flats and salt marshes lying along the riverbanks (Figure 1b) resulted from the combination of infilling of sandy sediments from the bed of the Irish Sea and deposition of the silt and clay-sized sediments coming from the River Ribble, with the sea influencing the sediment supply to the estuary in a much more significant quantity than the river ( Van der Wal et al., 2002). The accretion of the intertidal flat and salt marsh platform might have been further aided by the moderate wave climate characterizing the estuary insufficient to cause significant lateral erosion (Van der Wal et al., 2002). Glacial Isostatic Adjustment models show a long-term decrease in relative sea-level in North-West England for the past 2000 years (Shennan et al., 2018;Tooley, 1974), but historical tidal gauge records show that this was followed by an average increase in mean sea-level of 2 mm yr −1 since the beginning of the measurements in 1858 AD (PSMSL, 2019;Woodworth et al., 1999). Plater et al. (1993) also suggested that a secular rise in sea-level characterized the area, increasing the accommodation space. In 1810 AD, land reclamation of the estuary begun, and embankments were built along the northern bank of the estuary first, from Lytham to Preston, and along the southern bank later, north of Southport (Parsons et al., 2013d;Van der Wal et al., 2002). Between 1847 AD and 1910 AD the channel began to be dredged as the inner estuary was experiencing an infilling (Parsons et al., 2013d;Van der Wal et al., 2002). Hesketh Marsh and Hesketh New Marsh were reclaimed respectively in 1859 AD and 1883 AD for agricultural purposes (Figures 1c-1f), while Hesketh Out Marsh was only reclaimed in 1980 AD (Figure 1g), but a managed realignment project that restored Hesketh Out Marsh was conducted between 2007 AD and 2017 AD (Figure 1h; Tovey et al., 2009). 10.1029/2021JF006524 4 of 24

Paleoenvironmental Reconstruction
To perform the analysis for paleoenvironmental reconstruction, Hesketh Out Marsh was cored using a percussion corer with a gauge 1 m long and 7.5 cm in diameter. Two cores 3.20 m long (HOM19-1 and HOM19-2) were extracted from the middle marsh ( Figure 2h). The middle marsh was chosen as a sampling location because it is characterized by an accumulation rate slow enough to provide a record which more likely extends before embankment construction and is less likely characterized by sediment loss due to wave erosion compared to records that can be recovered by the lower marsh area (Berry & Plater, 1998). Furthermore, it was not enclosed by embankments before 1980 AD (Tovey et al., 2009), therefore it can provide a record of accretion extending more toward present day compared to the higher marsh (Hesketh New Marsh), which was enclosed in 1883 AD. HOM19-1 was cored in the sunlight and used for geochemical, particle size distribution and near-infra red analysis. HOM19-2 was cored using black sleeves to protect the sediments from sunlight exposure and was used for optically stimulated luminescence analysis. Stratigraphic analysis was performed on both cores to ascertain that they contain similar stratigraphic patterns and that the OSL chronology developed for the sediment profile can provide temporal context to the measurements of changing sediment composition (Figure 2i).

Sediment Provenance
Geochemistry is widely used as a proxy to reconstruct the provenance of marshland sediments (i.e., marine or riverine). For instance, it was employed by Plater et al. (2000) in the attempt to reconstruct sediment provenance in the Tees estuary, North-West England and by Hazermoshar et al. (2016) to characterize surficial sediments in Eynak Marsh, North of Iran. To investigate the provenance of the marsh sediments, the geochemistry of the salt marsh sediment was compared to that obtained for samples in the likely sediment source areas in the Ribble catchment and the eastern Irish Sea (Liverpool Bay -Sefton Coast). The HOM19-1 core was subsampled into a total of 320 samples at 10 mm continuous intervals. To characterize the river sediments, a set of 15 samples from five different locations (three per location) were collected from the surface sediments of the Ribble catchment; the locations included the lower (Figure 2g (Figure 2f) tributaries. Marine sediments are transported from the Irish basin into the estuary from south-west and, to a minor extent, from the southern sand dunes through erosion (Parsons et al., 2013d). Characterizing this onshore delivery of sediments to the estuary, samples were taken from coastal aeolian dunes near Southport (n = 3) and the upper 1 m of seafloor sediment ∼ 20 km offshore in 28-30 m water depths (Figures 2b and 2c). The seafloor samples came from 6 m (max length) vibrocores collected at 2 km spacing across the Irish Sea by a Fugro Geotechnical team during survey cruises for installing the Western High

of 24
Voltage -Direct Link marine cable (see Pearce et al., 2011). The top 1 m was targeted and subsampled at 0.1 m intervals progressing down three vibrocores (n = 30), with these depths from the more mobile surface sediments of the eastern Irish Sea.
Geochemical data were measured by XRF using a XEPOS 3 Energy-dispersive XRF for all the samples from the Hesketh Out Marsh core -divided into upper core (0-100 cm) and basal core (101-320 cm) -and probable source area materials from the Irish Sea, Southport sand dunes and Ribble catchment. All samples were lightly hand ground, pressed and then measured under a He atmosphere under combined Pd and Co excitation radiation and using a high resolution, low spectral interference silicon drift detector. The XEPOS 3 undergoes daily standardization procedure and has accuracies verified routinely using 18 certified reference materials (Boyle et al., 2015). Light elements were corrected for organic content, using loss-on-ignition (LOI) values measured by heating the samples at 105°C overnight to evaporate all moisture content and then igniting them in a furnace at 450°C for 4.5 hr to combust all organic matter (Boyle, 2000). Principal component analysis was performed using PAST3 (Hammer, 2019) to explore the geochemical compositions of the samples and to assess any association between the salt marsh materials and the offshore and fluvially sourced sediments. The parameters selected for this PCA were: Group 1 formed by silicates indicators (Si (mg/g)); Group 2 formed by salt water (Na (mg/g), S (mg/g)) and shell content indicators (Ca (mg/g), Sr (μg/g)); Group 3 formed by organic content indicators (LOI (%)), silt/ clay mineral indicators (K (mg/g), Al (mg/g), Fe (mg/g), Ti (μg/g), Rb (μg/g), Zr (μg/g), Mn (μg/g), Nb (μg/g)), pollutants (Mg (mg/g), P (mg/g), As (μg/g), Pb (μg/g), Zn (μg/g), Ni (μg/g), Cu (μg/g), V (μg/g), Cr (μg/g), Ga (μg/g), Ge (μg/g), Br (μg/g), Ba (μg/g), I (μg/g)) and the rare elements (Y (μg/g), La (μg/g), Ce (μg/g)). The division in groups was based on the association between parameters performed by the PCA (Figure 4) as well as the geochemical properties of the parameters (Boyle, 2000;Plater et al., 2000).

OSL Analysis of Accretion Rates
To date the Hesketh Out Marsh sequence and to quantify sediment accretion rates over time, OSL analysis was performed on coarse quartz grains (90-150) from HOM19-2. Sampling of the core was designed to characterize all visible changes in core stratigraphy; a total of 12 samples were extracted from the core (HOM19-2-1 to HOM19-2-12, see Figure 5 for down-core sampling profile). Quartz grains were used for the analysis because their signal resets rapidly in response to sunlight (Godfrey-Smith et al., 1988), thus they are likely well zeroed in young sediments (Ballarini et al., 2003). The dominant bedrock geology of the area, Triassic desert mud-, silt-and sand-stones, would also likely contribute to provide bright quartz grains with good luminescence signal (Fitzsimmons, 2011).
To isolate the coarse quartz grains for equivalent dose measurement, all samples were treated with a 10% (v/v) dilution of 37% HCl to dissolve carbonates and with a 10% (v/v) dilution of H 2 O 2 to remove organic content. Dry sieving was used to extract grains 90-150 μm in dimeter for samples HOM19-2-2 and HOM19-2-4 and 90-125 μm for the rest of the samples. Density separation using sodium polytungstate was performed to isolate the fraction of grains between 2.62 and 2.70 g cm −3 and remove feldspar and heavy metals grains. The quartz grains were then etched using 40% hydrofluoric acid to remove the outer portion of the grains contaminated by alpha particles and to remove any remaining feldspar grains. Another wash with a 10% dilution of HCl was then performed to remove soluble fluorides produced in the previous step.
Grains were mounted as 2 mm-diameter multiple-grain aliquots on to 9.8 mm-diameter aluminum discs and used for analysis. Forty-eight discs were analyzed for each sample except for HOM19-2-4, which only had enough quartz for 14 discs. No coarse quartz was recovered for sample HOM19-2-1. The luminescence measurements were made using an automated Risø TL/OSL DA-15 reader equipped with a 90 Sr/ 90 Y source (Bøtter-Jensen et al., 2003), the quartz grains were stimulated with a light source consisting of blue (470 nm) light emitting diodes and the luminescence signals were detected using a 5 mm thick U-340 filter and convex quartz lens beneath a photo-multiplier tube. The signal was recorded at 125°C for a total of 40 s, where the OSL signal was summed for the first 0.6-1.5 s and the background was calculated from 33.6 to 38.5 s of the decay curve. A 1.5% instrument measurement error was incorporated into the curve fitting. Individual aliquots were accepted or rejected from D e calculation according to their performance during the following screening criteria (considering the associated uncertainties): (a) the response to the test dose was greater than three sigma above the background; (b) the test dose uncertainty was <20%; (c) the recycling ratios and (d) OSL-IR depletion ratios were within the range 0.8-1.2, and (e) recuperation was <5% of the response from the largest regenerative dose (7-8 Gy).

10.1029/2021JF006524 7 of 24
A pre-heat plateau test of 21 aliquots from sample HOM19-6 was used to determine the preheat temperature used for OSL analysis. The results show that both the D e and recuperation values started increasing from 220°C; thus a pre-heat temperature of 200°C for 10 s and a cut-heat of 160°C for 0 s were used for the single aliquot regenerative dose (SAR) protocol (Murray & Wintle, 2000). A dose-recovery test was performed on 24 aliquots from sample HOM19-2-6 using a 4.15 Gy dose suggested that the SAR protocol was appropriate (ratio of 0.97 ± 0.06, and overdispersion value of 5 ± 0%). The D e distributions ( Figure S1 in Supporting Infomation S1) show asymmetries for all samples. Thus, the minimum age model (MAM; Galbraith & Laslett, 1993) was used for age calculation.
Dosimetry samples were prepared using 20-80 g of bulk sediment. The samples were dried at room temperature, weighted, ground to homogenize the material and stored for four weeks. A high-resolution gamma spectrometer was used to measure the radionuclide concentrations (K, U and Th) for 72 hr. The final water content was estimated calculating the weighted average of the mean water content during the period in which samples were transitioning from surface to depth and the mean water content during the period in which the samples were at depth. Both contents were estimated by calculating the weighted average of the water content during the period in which the samples were saturated (below mean water level) and the water content during the period in which the samples were at field water content (above mean water level). Following the model provided by Roberts and Plater (2005), for the period in which samples were transitioning from surface to depth, it was estimated that they were saturated 75% of the time, and for the period in which they were at depth, it was estimated that they were saturated 90% of the time. The saturated water content was estimated as 45%, while the field water content was measured by weighting the samples before and after drying them at room temperature. Cosmic dose was calculated using Prescott and Hutton (1994) equation that accounts for altitude, latitude, longitude and depth below ground level of each sample. Th, U and K concentrations were converted into beta and gamma dose-rates using Guerin et al. (2011) conversion factors and attenuated for grain size and etching (Guerin et al., 2012). The outer coating of the quartz grains was removed by hydrofluoric acid etching and so alpha dose-rates were negligible. The sum of the three beta dose rates and the three gamma dose rates was then calculated and corrected for water content and combined with the cosmic dose to determine the total dose-rate (Table S1 in Supporting Information S1). Dose-rates were calculated using the Dose Rate Age Calculator (Durcan et al., 2011). Th:U was calculated to check whether it was in equilibrium (3 < Th:U < 4) or in disequilibrium due to the presence of saltwater (Aitken, 1985). The age of each sample was calculated by dividing the D e by the dose rate and presented with total uncertainties (±1 sigma); the present age used is 2019 (Table S2 in Supporting Information S1).
All OSL ages were then modeled using the Bayesian package "OxCal" -version 4.4, which constructs age-depth sequences by integrating the IntCal20 calibration curve with statistical models (Bronk Ramsey, 1995, 2009Reimer et al., 2020). The OSL ages were entered in the model as a P-sequence, in which sediment deposition is modeled as a Poisson (random) process, and where a parameter (k) determines the extent to which sedimentation rates are allowed to vary. A uniformly distributed prior was used for k such that k 0 = 1, and log 10 (k∕k 0 ) ∼U (−2, 2); this allows k to vary between 0.01 and 100. The final probability distributions with ranges incorporating 68.2%, 95.4% and 99.7% of the total area of the distributions with the highest probability density were modeled using a Markov Chain Monte Carlo sampling routine. Outliers were handled as a student-T distribution using a general outlier model. A boundary was inserted at 120 cm, where the core exhibits a neat stratigraphic change and stops being sand dominated ( Figure 2i). The final accretion profile uses the probability weighted mean age and 95.4% uncertainties from the Bayesian age-depth model.

Particle Size Distribution Analysis
Particle size distributions (PSDs) can be used to infer whether sediment deposition is driven by accretion or infilling. According to the particle size distribution model developed by Rahman and Plater (2014), fine-skewed to near-symmetrical distributions characterized by well-sorted sand-sized sediment, typical of traction load delivered by the faster tidal flow velocities (i.e., during flood phase), are attributable to infilling; near-symmetrical distributions characterized by fine to very fine, poorly sorted silts, typical of the suspension load that settles during the turn of the tide (i.e., during ebb phase), is attributable to gradual accretion. PSDs were measured for HOM19-1 at 10 mm intervals using a Coulter LS 13 320 Single-Wavelength Laser Diffraction Particle Size Analyzer that determines the dimensions of individual particles 0.375-2000 μm. Subsamples were digested in 6% concentrated H 2 O 2 (50 ml of H 2 O 2 per 10 ml of sample) to remove any organic component, with the samples then dispersed and sonicated in Na 6 O 18 P 6 and analyzed under sonicating measurement conditions. The resulting PSDs are the average of three repeats after elimination of outliers. The Coulter LS 13 320 undergoes regular calibration checks using samples with known size distributions. End-member modeling analysis (EMMA) was conducted using the EMMAgeo R package to statistically derive the most recurrent modes that explain the variance in the PSDs down-core profile (Dietze et al., 2012). End-member modeling analysis has been successfully applied to sediment deposits in salt marsh settings by Clarke et al. (2014). Clarke et al. (2014) identified six main modal peaks occurring in consistent particle size ranges that account for more than 98.5% of the explained cumulative variance in PSDs, which agree with the six fixed particle size windows identified by Clarke et al. (2013) through the estimate of the particle size ranges deposited through settling, saltation and traction processes. Hence, the maximum number of end-members that the model could detect in the down-core PSD profile was set as seven. Particle size categories with several zero values were combined to enable end-member analysis of the entire data set. A robustness test was performed to check on each end-member ( Figure S2a in Supporting Infomation S1). The model showed that 86% of the variance in the PSDs present in the HOM19-1 core could be explained by four PSD end-members ( Figure S2b in Supporting Infomation S1). The loading of the end-members explaining most of the variance was plotted to observe what depositional processes (i.e., infilling or accretion) could be correlated with the PSDs in the core. The down-core profile of the end-members scores were then plotted to observe changes in distribution along the core and attribute whether possible changes in accretion over time detected by the OSL in HOM19-2 were attributable to infilling or gradual accretion.

Particle Size Distributions
To reveal whether any change in accretion rate can be connected to changes in the mode of sediment supply connected with the natural evolution and stabilization of the system, Plater et al. (2007) particle size interpretation was used to reconstruct the transition from sandflat to mudflat to salt marsh. Salt marsh environments are characterized by very slow currents and sediments are fine, poorly sorted, near-symmetrical distributed and platykurtic. Mudflat sediments are deposited with low energy conditions and characterized by fine, very poorly sorted, positively skewed to symmetrical distributed, meso-to platykurtic sediments. Sand flats sediments are deposited by high energy currents and are instead characterized by coarse, well-sorted to moderately sorted, positively skewed to symmetrical distributed, meso-to leptokurtic sediments. GRADISTAT 8.0v (Blott & Pye, 2001) was used for the statistical analysis of the PSDs measured for each sample of HOM19-1 to calculate D50, D90, standard deviation, skewness and kurtosis (geometric results).

X-Ray Fluorescence Analysis
To provide a more robust paleoenvironmental reconstruction, PSD statistical analysis was integrated with geochemical proxies. The geochemical composition of sediments is thought to be directly correlated with changes in stratigraphy and has been used in various case studies for the investigation of environmental transitions in marshland settings. Plater et al. (2007) used changes in the geochemical composition of the sediments as a proxy for stratigraphical changes to reconstruct the evolution of Romney Marsh/Dungeness Foreland depositional complex. Kolditz et al. (2012) implemented the same approach to study the environmental transitions in a marshland setting during sea-level rise in North-West Germany, and Hazermoshar et al. (2016) used it in a multi-proxy study that aimed to characterize surficial sediments in Eynak Marsh in North of Iran. The X-Ray Fluorescence (XRF) measurements for all HOM19-1 samples were used to produce the profiles of the following geochemical proxies: Si/K (silicates normalized toward mineral matter as sand indicator), Ca/K (calcium normalized toward mineral matter as shell content indicator), Fe/Mn (redox indicator), Rb (silt/clay mineral indicator); LOI was added to the list of geochemical proxies as organic content indicator (Plater et al., 2007).

Near-Infra Red Spectra Analysis
The XRF analysis was supplemented by near-infrared (NIR) spectra derived sediment composition. NIR spectra were measured by diffuse reflectance for all HOM19-1 samples using an integrating sphere on a Bruker MPA Fourier-Transform NIRS, with NIR spectra based on combining 24 scans collected at 8 cm −1 intervals across the range 3,595-12,500 cm −1 (Russell et al., 2019). NIR spectra were converted to first derivatives, calculated using a centrally weighted Savitzky-Golay smoothing algorithm minimizing spectral noise (Russell et al., 2019). The NIR spectra were interrogated using multiple regression of known end-member spectra onto the spectra derived for the unknown sediment samples (Russell et al., 2019). NIR spectra for the unknown samples are a conservative mixture of the spectra of component materials, and by using appropriate end-member spectra the multiple regression process can unmix the proportions of the components (Russell et al., 2019). The multiple regression analysis was undertaken in R using the (LM) function in R (R Core Team, 2013), using mineral matter (powdered quartz SiO 2 ), CaCO 3 shell material and organic matter (International Humic Substances Society (IHSS) humic/ fulvic acid standards); the output concentrations for these components are reported as weight percent ( Figure S3 in Supporting Infomation S1).

Cluster Analysis
A classical hierarchical clustering routine was run on all the components of the PSD, XRF and NIRS analysis using PAST3 (Hammer, 2019) to produce a dendrogram showing how data points can be clustered and statistically identify any change in the depth-profile. The unweighted pair-group average algorithm, which joins clusters based on the average distance between all members in the groups, with a stratigraphic constraint, was applied with a Euclidean similarity index.

Numerical Modeling
To investigate the accretion and evolution of the marsh in scenarios with no anthropogenic disturbance and with the addition of embankments, the FLOW module of the numerical finite-difference model Delft3D was used to simulate the hydrodynamics and sediment transport of the system. The model calculates non-steady flow and transport phenomena using Navier-Stokes and transport equations (DELFT Hydraulics, 2014). The model set-up for this study accounts for suspended-load, evaluated through the advection-diffusion equation, of multiple cohesive sediment fractions; erosion and deposition for cohesive sediments are computed using the Partheniades-Krone formulations (Partheniades, 1965). It also accounts for vertical diffusion of sediments due to turbulent mixing and sediment settling due to gravity.
The original design of the domain and the set-up and calibration of the model for this study were performed by Li et al. (2018Li et al. ( , 2019. The domain consists of a grid of 344 × 80 cells in the east to west and north to south direction (Figure 3a), and three equally spaced vertical layers. The cell size is approximately 20 × 20 m within the river, 70 × 70 m across the intertidal area (including the marsh platform) and 300 × 300 m in the outer estuary. Simulations with and without embankments presence were conducted using bathymetries from 1847 AD and 2008 AD  (Table 1; Van der Wal et al., 2002;Parsons et al., 2013d). The 1847 AD bathymetry is the oldest with data availability within the timeframe of the dates produced by the OSL analysis (Table S2 in Supporting Infomation S1) and it offers the opportunity to simulate the presence of an accommodation space before the formation of the intertidal flat and before the beginning of the channel dredging ( Figure 3b). The 2008 AD bathymetry was used to simulate an established intertidal flat and salt marsh platform and completed channel dredging (Figure 3c). Between 1847 AD and 1910 AD the channel began to be dredged as the inner estuary was experiencing an infilling, with dredge spoils disposed off the mouth of the estuary until 1,905 and offshore after then, and the depth of the channel increased approximately by 1 m when it was at its deepest (Parsons et al., 2013d;Van der Wal et al., 2002). It has been observed that channel dredging can modify tidal asymmetry (Van Maren et al., 2015;Zhu et al., 2015). To understand how channel dredging interfered with the effects of embankment construction on the sediment transport, a simulation was run using the 1847 AD bathymetry with a channel depth increased by 1 m (Figure 3d); this allowed to investigate the sole effect of channel dredging without considering the morphological evolution of the intertidal flat and the construction of further embankments that occurred during the period of the historical dredging. The bathymetric data for 1847 AD was obtained by digitizing the maps provided by Van der Wal et al. (2002) with reference Above Ordnance Datum Newlyn (AODN). Since the data provided by Van der Wal et al. (2002) covered a restricted area of the domain (53.59°N to 53.89°N, −3.13°E to −2.83°E), to provide a domain large enough to stabilize the model, this was integrated with present data obtained from the combination of two datasets: the bathymetry data for the open sea collected from EDINA DIGIMAP and the LiDAR data for the coastal regions downloaded from the Environment Agency's LiDAR data archive. The elevation of the inland terrain model and very offshore areas is, on average, of the same order of magnitude than what presented in Van der Wal et al. (2002). The choice of integrating the digitised bathymetry with present bathymetric data was based on limitations in terms of data availability for that period and the assumption that changes in the most offshore bathymetry and inland terrain elevation (where the marsh is already well formed) would affect the analyses and calculations performed on the salt marsh platform and within the intertidal area to a minor extent compared to the long-term changes of the river, inner and outer estuary bathymetry. The bathymetry influences the propagation of the tidal currents (e.g., Van Maren et al., 2015;Li et al., 2018;Palmer et al., 2019), therefore this choice could cause limitations such as an overall overestimation or underestimation of the net sediment transport within the estuary and of the final sediment budget of the salt marsh. The bathymetric data for 2008 AD was obtained directly by combining the bathymetry data for the open sea collected from EDINA DIGIMAP and the LiDAR data for the coastal regions downloaded from the Environment Agency's LiDAR data archive. Since the two datasets have different vertical reference levels, Low Astronomical Tide and AODN respectively, they were adjusted and referred to Mean Sea Level (MSL) following the spatially varying Vertical Offshore Reference Frame corrections provided by the UK Hydrographic Office prior to combining. The model is constrained within two open boundaries, one 20 km offshore and one across the River Ribble. The model was successfully calibrated by Li et al. (2018Li et al. ( , 2019 using the 2008 AD bathymetry and 2008 AD time-series for both river discharge and offshore water level. Data for the offshore boundary was provided by the Extended Area Continental Shelf Model fine grid (CS3X), which has approximately 12 km resolution, covers an area from 53.55°N to 53.92°N and from −3.33°E to −2.75°E. Data for the river boundary was collected from the National River Flow Archive and consists of a time series of daily averaged river discharge values; a constant discharge of 44 m 3 /s was estimated using the mean discharge for the simulated period. Each simulation was run for one year. To simulate past ocean boundary conditions, 10 tidal harmonics (M2, S2, N2, K2, O1, P1, Q1, K1, M4 and S1) were derived for the Ribble Estuary using the global tidal model GOT-e 4.10c (Ray, 1999;Zaron & Elipot, 2021). To ascertain that the tidal harmonics would allow for the modeling of realistic hydrodynamic for the simulated years, the water level at ocean boundary for the 2008 AD scenario forced with tidal harmonics was compared to the water level for the 2008 AD scenario forced with water level time-series ( Figure S4 in Supporting Information S1). Since the river contribution to sediment delivery in the intertidal area is minimal (Figure 4; Lyons, 1997;Van der Wal et al., 2002) and no time-series for river discharge are available for 1847 AD, the time-series for river discharge from 2008 AD were used to force the river boundary of the 1847 AD scenarios as well. The Ribble Estuary experiences moderate wave energy (Pye & Neal, 1994). However, for simplicity, the effects of wind waves have been neglected in this study to investigate exclusively the effects Figure 5. Age-depth profile of optically stimulated luminescence dating (OSL) ages modeled using Bayesian analysis and accretion rates calculated from the modeled ages. The change in color represents the observed change in core stratigraphy. The OSL ages were entered in the model as a P-sequence using a k parameter which varies between 0.01 and 100. The squared gaps, dark gray and light gray distributions incorporate respectively 68.2%, 95.4% and 99.7% of the total area of the distributions with the highest probability density. Outliers were handled as a student-T distribution using a general outlier model. The final accretion profile uses the probability weighted mean age and 95.4% uncertainties from the Bayesian age-depth model.

Table 1
List of Simulated Scenarios Indicating Name and Description of the Simulations of embankments on the sediment transport. A uniform and not erodible initial bed composition was applied, to avoid feedbacks between hydrodynamics and the initial bathymetry. This is because the focus of the study was the investigation of the redistribution of input sediments and deposition within the estuarine accommodation space with and without embankments, rather than the modeling of the morphological evolution of the estuary. The initial suspended sediment input and the sediment input at river and ocean boundaries were simulated using cohesive erodible sediments. The suspended sediment input at the river boundary was calculated as an average for the year 2008 in proportion to the river discharge, using the Wright-Parker formulation (without stratification) for the suspended load, which estimated a mean value of 0.002 Kg/m 3 . The input at the ocean boundary for the suspended load and the initial suspended sediment concentration in the domain were given using the value measured by Silva et al. (2016) in 2011 (0.001 Kg/m 3 ). To simulate the presence of embankments, thin dams were imposed, and the position of the dams was chosen based on the position of the embankments at the time of each bathymetric survey (Table 1; Van der Wal et al., 2002;Parsons et al., 2013d).
To quantify the contribution of embankment construction on the accretion of the tidal flat-salt marsh, the sediment budget of Hesketh Out Marsh platform was calculated at the end of the simulated year for each scenario by multiplying the cumulative deposition of each cell for the cell area and summing them (Ganju et al., 2015). Tidal asymmetry has been recognised as one of the main factors controlling the net import/export of sediments and the large-scale morphological evolution of estuaries Guo et al., 2016;Palazzoli et al., 2020) and there is evidence that the presence of embankments alters the way sea-level rise affects the resonance properties of the coastline, modifying the tidal signal and changing the sediment supply to coastal marshes (e.g., Carless et al., 2016;Holleman & Stacey, 2014;Pelling et al., 2013). To investigate the effect of embankment construction on the tidal signal, Fourier analysis was conducted for each scenario for an entire monthly tidal cycle. The MATLAB package T-TIDE was used to conduct the analysis (Pawlowicz et al., 2002). The distortion and asymmetry of the tidal signal were analysed using the principal lunar semi-diurnal constituent M2 and the largest shallow water constituent M4; the distortion was calculated using the ratio A 4−2 = AM 4 /AM 2 where A is the amplitude of the tidal height and the asymmetry was calculated using Δθ = 2θM 2 −θM 4 where θ is the phase of the tidal height (Blanton et al., 2002;Friedrichs & Aubrey, 1988). When Δθ is between 0° and 180° the flood phase dominates, whereas when it is between 180° and 360° the ebb phase dominates; the magnitude of A 4−2 is representative of the significance of the dominance (Friedrichs & Aubrey, 1988). The harmonic analysis was not performed on dry cells and on salt marsh or tidal flat cells intermittently covered by water. Figure 4 shows the first two components of the PCA performed on the XRF measurements of the sediments taken from Hesketh Out Marsh, the Irish Sea floor, Southport sand dunes and the Ribble catchment. Component 1 and component 2 summarize 72.5% of the variance in the data, respectively 63.8% and 8.7%. Figure 3a shows how the parameters are used by each component to separate the samples. Component 1 separates samples rich in salt water (Na and S), shells (Ca and Sr) and silicates (Si) from those rich in silt/clay mineral indicators (e.g., Zr, Ti, Rb, K, Al) and pollutants (e.g., As, Pb, Zn). Component 2 separates samples rich in silicates (Si) with and without shelly materials (Ca and Sr). The samples from Hesketh Out Marsh basal core overlap with the samples from the Irish Sea cores and one sample from the Southport sand dunes. Hesketh Out Marsh upper core is characterized by samples that in part cluster with the marine and basal core samples or form a second cluster by themselves. The samples from the Ribble catchment form a third cluster with two of the samples from the Southport sand dunes. This implies that there is a strong correlation between the geochemical composition of the Irish Sea floor sediments and the Hesketh Out Marsh basal core sediments, while most of the upper core does not show this association; conversely, the sediments collected from the Ribble catchment do not correlate with any of the other two clusters. This suggests that the first 200 cm of sediment deposited on the tidal flat-salt marsh platform was predominantly supplied by the Irish Sea floor, with this input from the seafloor becoming less significant in the sediment deposition captured within the top 100 cm of the core; the river instead provides minimal sediment input to the intertidal flat system. Figure 5 shows the age-depth distribution of the OSL ages modeled using Bayesian analysis and accretion rates calculated from the modeled ages of Hesketh Out Marsh. The ages between 300 and 125 cm overlap within uncertainties between approximately 1734-1806 AD and 1768-1849 AD, with a modeled accretion rate of 4.61 cm yr −1 . Between 100 and 0 cm, the ages range between 1768-1849 AD and 1926-1986 AD and the modeled accretion rate is 0.83 cm yr −1 .

Marsh Evolution and Sediment Provenance
The end-member analysis showed that 86% of the variance in the PSDs present in the HOM19-1 core could be explained by four PSD end-members (EM1, EM2, EM3 and EM4). These end members represent the recurring modes detected in the PSDs. Figure 6a shows that EM1 has a mixed near-symmetrical distribution dominantly clay to silt with a secondary very fine to fine sand mode, EM2 has a fine-skewed to near-symmetrical distribution dominantly very fine to fine sand, EM3 has a bimodal distribution with near-symmetrical clay to silt mode and near-symmetrical to fine skewed very fine to fine sand mode, and EM4 has a fine-skewed to near-symmetrical distribution of fine to coarse sand. Figure 6b shows that up to 300 cm EM4 and EM2 contribute equally to the PSDs of the core. Between 300 cm and to 150 cm the distribution becomes dominated by EM2 and punctuated by layers of EM4; EM1 and EM3 only give a negligible contribution in this section. Between 150 and 100 cm EM1 and EM2 contribute equally to the distribution and dominate this section of the core, however a minor contribution is also given by EM3. From 100 cm up to the top of the core EM1 is the dominant end-member, contributing to 70% of the distribution, while EM3 contributed for about 30%. This suggests that up to 300 cm the particle size distributions were indicative of a high energy environment dominated by infilling of coarser sand through load traction. Between 300 and 150 cm the environment became less energetic and dominated by infilling of finer sand occasionally layered by infilling of coarser sand and rarely punctuated by settling of silt-sized material. Between 150 and 100 cm the distributions reveal a progressive transition from an environment dominated by infilling through load traction into an environment dominated by accretion through settling. From 100 cm up to the top of the core the environment becomes dominated by settling of silt-sized particles layered with occasional infilling of fine sand.
This suggests that the rapid accretion that occurred between 1734-1806 AD and 1768-1849 AD was caused by an infilling of well-sorted, sand-sized sediment into the estuary accommodation. The slower accretion that characterizes the period comprised between 1768-1849AD and 1926-1986 AD, on the other hand, was triggered by the deposition of poorly sorted silts. Figure 7 shows that the core is dominated by coarse to very coarse sand up to 150 cm. From 150 to 100 cm fine sand becomes dominant. At 100 cm, grains become finer, ranging between fine sand and coarse silt. The grains are mostly well-sorted to moderately sorted up to 100 cm, while the top 100 cm of the core is mainly dominated by poorly to very poorly sorted grains. The grain size distribution is mostly positively skewed up to 150 cm, where near symmetrical distribution starts becoming more common up to 100 cm; from 100 cm to the top of Figure 6. End-members loading (a) and down-core profile of end-members scores (b) for HOM19-1. A fine-skewed to near-symmetrical distribution is typical of traction load delivered by the faster tidal flow velocities (i.e., during flood phase); a near-symmetrical distribution is typical of the suspension load that settles during the turn of the tide (i.e., during ebb phase).
the core the distribution becomes mostly near symmetrical. The core is mostly leptokurtic up to 150 cm. From 150 to 100 cm, it becomes mostly mesokurtic. From 100 cm to the top of the core it becomes mostly platykurtic. NIRS analysis suggests that mineral content constitutes more than 90% of the material collected from Hesketh Out Marsh, while shell and organic content together only constitute less than 0.2%. SiO 2 increases consistently up to 50 cm, where it reaches almost 100% of the total composition, and then decreases consistently up to 15 cm, where it drops sharply. CaCO 3 behaves in a specular way; it decreases up to 100 cm, where the quantity becomes negligible such that the model stops detecting it, except for a few short sections; it then increases again up to the top of the core; however, the latter increase is likely caused by the biological production of calcium carbonate by plants rather than by an increase in shell content. Humic/Fulvic acid is constant up to 50 cm, where it then starts increasing up to the top of the core. Rb stays constant up to 120 cm, where it starts decreasing up to the top of the core. Ca/K and Si/Al are constant up to 100 cm, where they start decreasing slightly up to 75 cm, and then remain constant up to the top of the core; Ca/K shows a slight increase toward the top which again could be due to biogenic production. Fe/Mn stays constant up to 140 cm, where it then starts decreasing up to 120 cm, after which it remains constant up to the top of the core. LOI is constant up to 100 cm, where it starts to increase up to the top of the core. The clustering analysis statistically shows that the core is divided into two main sections, one extending from the bottom of the core up to 120 cm, dominated by sand-sized, well-sorted grains, shelly material, and reduced conditions, and one extending from 120 cm to the top of the core, dominated by silt-sized, poorly sorted grains, and high mineral and organic content.
Overall, Figure 7 suggests that a sandflat environment dominates the core up to 150 cm, while a transition from sandflat to mudflat occurs between 150 and 70 cm. At 70 cm a mudflat environment becomes dominant up to 15 cm, where the transition to a fully vegetated marsh occurs. The most statistically significant change in the down-core composition occurs at c. a. 120 cm, where the sandflat is replaced by a mudflat. Figure 8 summarizes the results derived from the sediment budget calculations performed for the numerical simulations. Figure 8a shows that the sediment budget of Hesketh Out Marsh for present day decreases compared to the sediment budget for the 1847 AD bathymetry; this is true for both scenarios (i.e., with and without embankments). It also shows that, for all simulated bathymetries, the sediment budget decreases when embankments are added. Figure 8b shows that dredging of the channel causes an increase of 0.5% in the sediment budget when there are no embankments present and an increase of 1.25% when embankments are added.  Figure 9 shows the tidal analysis performed for each simulated bathymetry. The Ribble Estuary resulted to be overall flood dominated in all simulated years and for both scenarios (i.e., with and without embankments; Figures 9a and 9b). In a scenario with no embankments, the flood dominance increases for the present bathymetry with respect to the 1847 AD bathymetry (Figure 9a). When embankments are added, despite the estuary remaining flood dominated, there is an overall intensification of the ebb currents, particularly near the tidal flat ( Figure 9c). For the 2008 AD bathymetry, the intensification of the ebb currents is significant also in the outer estuary ( Figure 9c). This explains why the sediment budget decreases when embankments are added, since an intensification of ebb currents causes an increase in the sediment export out of the system. Figures 9d and 9e shows that the increase in ebb dominance follows an increase in the phase of the M 2 constituent and a decrease in the phase of the M 4 constituent. With channel dredging, the flood dominance near the tidal flat is enhanced in both scenarios (i.e., with and without embankments; Figures 9a and 9b). This explains the overall increase in the sediment budget caused by the channel dredging compared to the non-dredged bathymetry, as the enhanced flood dominance near the tidal flat promotes sediment import. When embankments are added, the intensification of the ebb currents near the tidal flat is lower compared to the non-dredged bathymetry (Figure 9c), which explains why the increase of sediment budget for the dredging scenarios compared to the no-dredging scenarios is greater when embankments are present.

Modeled Scenarios
Overall Figure 8 suggests that embankments promote sediment export and cause a decrease in the sediment budget of the salt marsh.

Sediment Provenance
The principal component analysis ( Figure 4) showed a correlation between the geochemistry of Hesketh Out Marsh sediments and the geochemistry of the Irish Sea floor sediments, indicating that the provenance of the marsh sediments is predominantly marine. Figures 6 and 7 show the causes that determined the strong correlation between the Irish Sea cores and the lower part of Hesket Out Marsh core, while the upper 0-1.25 m of salt marsh sediments do not show an association with the seabed sediments. The lower 1.25-3 m of sandflat deposits of the core are characterized by sediments entrained at the seafloor by flood currents and deposited directly in the estuarine accommodation space before a tidal flat was formed, while the upper part is characterized by fine suspended sediments deposited by ebb currents through tidal creeks when a tidal flat was already forming (Rahman & Plater, 2014). Our results agree with the interpretation of Lyons (1997) that sediments are delivered into the estuary mainly by tidal pumping. They further support previous hypotheses that abundant legacy of sediment left during the last glacial cycle (Scourse et al., 2021) and regularly remobilized through the central Irish Sea, often as sediment waves (Van Landeghem et al., 2012), has been the main source of sediment supply to the estuary (Wright et al., 1971).   Figure S5 in Supporting Infomation S1 for tidal distortion (A 4−2 ), M 2 amplitude (AM 2 ) and M 4 amplitude (AM 4 ) for the same scenarios.

Changes in Marsh Accretion
While some studies (e.g., Davis et al., 2010;Madsen et al., 2005;Madsen et al., 2007bMadsen et al., , 2007a have found well-bleached quartz in tidal flats and salt marsh environments, Madsen et al. (2009) found partially bleached quartz in three hurricane deposits made of coarse sandy particles rapidly deposited on Little Sippewissett Marsh. High levels of turbidity characterizing the water table during high energy events, which are a major source of sediment supply to the Ribble Estuary (Lyons, 1997), might have contributed to partial bleaching by attenuating the blue end of the daylight spectrum (Berger, 1990;Plater & Poolton, 1992). Moreover, the tidal inundation of Hesketh Out Marsh is semidiurnal, therefore sediment deposition can occur either at daytime or night-time (UKHO, 2019). This means that some grains experience shorter exposure to light if not none, further contributing to partial bleaching (Mauz et al., 2010). However, since the sediments analyzed are very young -between 67 ± 9 and 264 ± 26 years old -they might have also suffered incomplete resetting of prior stored charge which can lead to significant offsets in the equivalent dose (Murray & Olley, 2002). The asymmetry in the D e distribution could alternatively be caused by microdosimetry. It was indeed found that, in sandy environments, K-feldspar can affect quartz dose rate distribution and, as a consequence, the D e distribution skewness and dispersion of dose rate distributions increase as the number of potassium-rich feldspar grains is decreased relative to the number of quartz grains or if the quartz grain size increases, as this causes an increase in distance between the quartz grains and the source of potassium (Guerin et al., 2015). Nonetheless, regardless of the cause of the asymmetrical D e distribution, accurate OSL ages can be derived using small aliquots and the MAM (Galbraith & Laslett, 1993).
The Bayesian analysis of the OSL ages shows that between 1734-1806 AD and 1768-1849 AD there was a very rapid accumulation of sediments (ca. 4.61 cm yr −1 ), followed by a slower accretion (ca. 0.83 cm yr −1 ) that persisted until 1926-1986 AD ( Figure 5). The overlap of the OSL ages between 300 and 125 cm was likely to be due the resolution of the technique that is not able to differentiate individual ages in such fast deposition. However, the timeframe of the environmental reconstruction of the historical maps confirms the accuracy of the age intervals at a decadal scale. Indeed, according to the reconstruction of marsh evolution (Figure 7), the first period of fast accretion determined the formation of the sand flat substrate until 1768-1849 AD, when the rate of accretion decreased, and the sandflat progressed into a mudflat. Accordingly, the historical reconstruction shows that in 1840 AD Hesketh Out Marsh area was occupied by a sandflat (Figure 1b), while between 1890 AD and 1950 AD a mudflat was the dominant environment (Figures 1c-1e). The new OSL ages here extend this understanding by providing the first record of sediment dating and the longest record of accretion rates for the Ribble Estuary (Mamas et al., 1995) that then allow us to compare with hydrodynamical modeling to understanding the processes involved in salt marsh accretion. The numerical simulations show that the sediment budget of Hesketh Out Marsh decreases toward present day ( Figure 8a); this is true for both scenarios (i.e., with and without embankments). There is an agreement between the trends observed for the sediment budget of the 1847 AD and 2008 AD scenarios and the trends detected in the long-term accretion profile provided by the luminescence ages. The results of the two techniques are therefore comparable.
Historical records of the estuary show that in 1810 AD embankment construction started along the northern bank of the Ribble Estuary. Between 1810 AD and 1850 AD the whole northern bank, from Lytham to Preston, and part of the southern bank, north of Southport, was embanked, with the exception of the coastline surrounding Hesketh Out Marsh area (Van der Wal et al., 2002). The period where the majority of embankment construction occurred (1810-1850 AD) falls within the age intervals given by the OSL ages for the period of the fast sand flat accretion (from 1734-1806 AD to 1768-1849 AD). However, the OSL age intervals suggest that there is a possibility that the rapid sedimentation might have started before the beginning of the 19th century. Moreover, no shift was detected from a slower to a faster accretion rate in the lower part of the core, therefore there is no direct evidence that the embankment construction might have triggered the infilling that caused the rapid sedimentation.

Effects of Embankments on Marsh Sediment Budget
The calculation of the sediment budget for all the simulated bathymetries revealed that the final sediment budget for Hesketh Out Marsh is greater when no embankments are present (Figure 8a). The quantity of sediment entering and exiting an estuary during each tidal cycle is influenced by tidal asymmetry Guo et al., 2016;Palazzoli et al., 2020). A flood dominance (flood phase shorter and more intense than ebb phase) triggers a net landward sediment import as flood velocities can be sufficiently high to resuspend sediment, while an ebb dominance is responsible for a net export (Lanzoni & Seminara, 2002;Van Dongeren & De Vriend, 1994).
Tidal analysis for the Ribble Estuary (Figures 9a-9c) showed that, despite the estuary remaining flood dominated, the addition of embankments along the coastline affected the tidal propagation within the system by causing an overall intensification of ebb currents and weakening of flood currents. The Irish Sea, where the Ribble Estuary is located, is influenced by a phenomenon of resonance (Pickering et al., 2012). In San Francisco Bay (Holleman & Stacey, 2014), Yellow River Delta (Pelling et al., 2013) and the Patagonian Shelf (Carless et al., 2016), it has been observed that the introduction of embankments can affect the resonance properties of a basin, causing changes in the amplitude and phase of the M 2 and M 4 constituents, which result in changes in tidal distortion and asymmetry.
In the Ribble Estuary, this increase in ebb dominance when embankments are added, follows an increase in the phase of the M 2 constituent and a decrease in the phase of the M 4 constituent (Figures 9d and 9e).
When the channel is dredged, the areas of the channel adjacent to the tidal flat become flood dominated (Figures 9a and 9b). It has been observed that an enhancement of flood dominance following channel dredging is caused by a tidal amplification which increases estuarine circulation (Van Maren et al., 2015;Zhu et al., 2015).
Other estuaries, such as the Ems Estuary, located between the Netherlands and Germany, have showed similar behavior (Van Maren et al., 2015). In the Ribble Estuary, this resulted in a lower intensification of ebb currents when embankments are added (Figure 9c), hence attenuating the effects of the embankments on the sediment budget ( Figure 8b).
Overall, our numerical modeling shows that embankments enhanced sediment export out of the system (Figures 8  and 9). Embankments begun to be built in 1810 AD; however, the oldest bathymetry available is from 1847 AD (Parsons et al., 2013d;Van der Wal et al., 2002). Therefore, we cannot directly discuss changes in sediment transport caused by the beginning of the embankment construction. Luminescence derived accretion rates ( Figure 5), in line with historical maps (Figure 1) and later bathymetric studies (Van der Wal et al., 2002), show a continuous infilling of the estuary and progressive formation of the tidal flat over the period studied. Accordingly, tidal analysis shows that in 1847 AD the estuary was overall flood dominated despite the presence of the embankments (Figure 9b). For the 2008 AD bathymetry, with no embankments present, the flood dominance of the Ribble Estuary is greater than it is for the 1847 AD (Figure 9a). Both in the Ribble Estuary (Li et al., 2018;Pannozzo et al., 2021a) and in other sites (e.g., Donatelli et al., 2018;Palmer et al., 2019), it has been observed that an increase in accommodation space promotes ebb dominance; the strengthening of flood dominance for the 2008 AD no-embankment scenario could be therefore caused by more extensive intertidal flats present in the 2008 AD bathymetry, responsible for a decrease in accommodation space. A similar behavior has been observed in the Tamar Estuary (Australia), where infilling of the accommodation space causes an increase in flood dominance (Palmer et al., 2019). However, when embankments are added, the 2008 AD scenario is characterized by a weakening of the flood currents compared to the 1847 AD scenario which extend also to the outer estuary ( Figure 9b). This indicates that further land reclamation and the addition of more embankments reduced the overall flood dominance of the estuary toward present day, which explains the sharper drop in sediment budget when embankments are added to the 2008 AD bathymetry compared to when they are added to the 1847 AD bathymetry. It further suggests that in 1810 AD the estuary might have been characterized by greater accommodation space and might already been flood dominated, but the beginning of embankment construction might have already started reducing the natural strengthening of flood dominance driven by the infilling of the accommodation space.

Natural Accretion of the Marsh
The OSL age intervals suggested that there is a possibility that the rapid sedimentation might have started before the beginning of the 19th century, preceding the start of the embankment construction ( Figure 5). Moreover, no shift was detected from a slower to a faster accretion rate in the lower part of the core, showing no direct evidence that the embankment construction might have triggered the rapid sedimentation ( Figure 5). The sediment budget and tidal analysis suggest that in the Ribble Estuary, embankments were responsible for an overall intensification of the ebb currents, promoting sediment export and causing a decrease in the sediment budget of the tidal flatsalt marsh platform (Figures 8 and 9). Therefore, embankment construction was unlikely to be responsible for an increase in sedimentation. This suggests that the infilling of the estuary responsible for the high sedimentation rate must have had a different cause.
The particle size distribution analysis ( Figure 6) suggests that the rapid accretion, detected by the Bayesian analysis of the OSL ages and occurred between 1734-1806 AD and 1768-1849 AD, was mainly caused by infilling of the estuary with sand-sized sediments during flood tide (Rahman & Plater, 2014) and the reconstruction of Hesketh Out Marsh evolution ( Figure 6) suggests that in this period a sandflat was the dominant environment (Plater et al., 2007). In the first period (bottom 300 cm of the core) the Hesketh Out Marsh area was a highly energetic environment that favored deposition of coarser sand into the estuary; this was then followed by a transition (between 300 and 150 cm of the core) in which the environment became less energetic and dominated by infilling of finer sand ( Figure 6). Glacial Isostatic Adjustment models show a long-term decrease in relative sea-level in North-West England for the past 2000 years (Shennan et al., 2018;Tooley, 1974). Historical tidal gauge records showed that this was followed by an average increase in mean sea-level of 2 mm yr −1 since the beginning of the measurements in 1858 (PSMSL, 2019;Woodworth et al., 1999), and Plater et al. (1993) suggested that a secular rise in sea-level characterized the area, increasing the accommodation space. The secular rise in sea-level of the past centuries and consequent increase in accommodation space, might have led to the formation of the estuary, favoring an infilling of the area with marine sediments. The central Irish Sea is characterized by waves of mobile sediments that move up to 70 m per year (Van Landeghem et al., 2012). This high amount of mobile sediment is a legacy of materials left by the former Irish Sea Glacier that is mobilized on the seafloor in the form of sediment waves and sand sheets (Van Landeghem et al., 2012). Wright et al. (1971) suggested that this sediment was the main source of sediment supply to the estuary, and this is confirmed by the analysis of sediment provenance (Figure 4) that shows a strong correlation between the floor sediment of the Irish Sea and the sediment of the initial Hesketh Out Marsh sandflat platform. The progressive increase of the accommodation space and the overall flood dominance of the estuary (Figure 9a) might have facilitated the deposition of sea-bed sediments transported from the Irish Sea floor into the estuary during flood tide (Plater et al., 1993), therefore explaining such a rapid deposition of sediment into the area. As the platform started increasing in elevation, the deposition of coarse sand was replaced by the deposition of finer sand, which was responsible for the formation of the sandflat (Reed, 1990). The slower accretion that characterized the period comprised between 1768-1849 AD and 1926-1986 AD, on the other hand, was caused by the deposition of poorly sorted silts ( Figure 6, Rahman & Plater, 2014). The reconstruction of Hesketh Out Marsh evolution (Figure 7) suggests that the accretion rate decreases with the progression from a sandflat to a mudflat (Plater et al., 2007) and the particle size distribution analysis ( Figure 6) suggests that between 150 and 100 cm a transition occurred from an environment dominated by infilling through load traction into an environment dominated by accretion through settling, while in the top 100 cm of the core the environment becomes dominated by settling of silt-sized particles layered with occasional infilling of fine sand (Rahman & Plater, 2014). This mode of deposition is typical of mudflat-salt marsh systems. When mudflat platforms reach an elevation which is above sea-level, tidal creeks form and become the main mean through which marine sediments can reach the marsh platform and deposit; they transport silt-sized sediments to the marsh areas of higher elevation during flood tide and deposit them during ebb tide (Reed et al., 1999). The sandy layers that punctuate the upper 100 cm of the core are instead likely to be caused by infilling of fine sand through load traction during occasional flooding of the marsh platform owed to spring tide high waters or storm surges (Pannozzo et al., 2021a(Pannozzo et al., , 2021b. Accordingly, the numerical simulations show that the sediment budget of Hesketh Out Marsh decreases toward present day both for embankments and no-embankments scenarios ( Figure 8a). The 1847 AD bathymetry is characterized by a greater accommodation space preceding the formation of the intertidal flat, whereas the 2008 AD bathymetry was characterized by an established intertidal flat and salt marsh platform. Therefore, it is likely that the increase in elevation of the sandflat platform, following the infilling of the estuary accommodation space, and the conversion into a mudflat-salt marsh platform, was responsible for the decrease in accretion rate detected by OSL around 1768-1849 AD. Hence, the rapid accretion of the Ribble Estuary tidal flat-salt marsh system seems to be a natural evolution of the system that adapted to the high availability of sediments.

Long-Term Marsh Resilience to Embankment Construction
The numerical simulations showed that the estuary has been flood dominated through the investigated period, but embankments were responsible for an overall intensification of the ebb currents in the system which promoted sediment export. The paleoenvironmental analysis showed that Hesketh Out Marsh has been accreting at a rate of 4.61 to 0.86 cm yr −1 over the last ca. 190 years and that the estuary has experienced continuous infilling over the period studied caused by a natural high rate of sediment supply. Furthermore, the modeling showed that the channel dredging reduced the intensification of the ebb currents caused by the construction of the embankments. The combination of these analyses suggests that the effects of the embankment construction on sediment transport did not compromise the long-term resilience of the salt marsh because of the high rates of sediment supply and undergoing river dredging that enhanced flood dominance near the tidal flat.
Part of Hesketh Out Marsh was enclosed by embankments in 1980 AD and managed realigned between 2007 AD and 2017 AD (Tovey et al., 2009). Managed realigned sites have been found to be resilient to rates of sea-level rise below 5 mm yr −1 (Orr et al., 2003). The resilience of managed realigned salt marshes has been found to be highly dependent upon the breach design and size of the site (Kiesel et al., 2020(Kiesel et al., , 2022 as well as sediment supply, size of sediment supplied to the marsh and tidal range (Xu et al., 2022). It was also observed that the vegetation cover at different stages of marsh maturity and related ability to trap sediments in restored salt marshes are comparable with those of forming marshes (Xu et al., 2022). However, for Hesketh Out Marsh, the effects on the reclamation and subsequent managed realignment on the marsh accretion cannot be discussed in relation to the accretion rates provided by the luminescence analysis, due to the low resolution of the technique that was only able to provide with ages ranging between 1734-1806 AD and 1926-1986 AD with a minimum error of ±30 years.
Other studies have shown that embankment construction can trigger higher sedimentation rates in certain areas of an estuary. This was observed, for instance, by Pontee (2014) in the Humber Estuary and Lee et al. (2017) in Chesapeake Bay and Delaware Bay, who concluded that, by channelizing sediment, embankments can favor the delivery of sediment to a specific area of an estuary accommodation space, inducing siltation of intertidal areas and aiding vertical accretion of tidal flat-salt marsh systems. This suggests that the effects of embankments on sediment deposition can change from case to case depending on the hydrodynamic and morphodynamic feedbacks characteristic of the site and that results for single basin studies cannot be generalized. Along the North-West coast of England, other estuaries of great economic importance and at high flood risk, experience tidal dynamics and rate of sediment supply similar to the Ribble Estuary (Parsons et al., 2013c). Liverpool Bay estuaries -that is, Dee Estuary and Mersey Estuary -and Morecambe Bay estuaries are all tidally dominated, macrotidal, flood dominated and receive sediments from the central Irish Sea that accumulate in intertidal areas forming extensive tidal flats and salt marsh platforms (Parsons et al., 2013a(Parsons et al., , 2013b(Parsons et al., , 2013c. They were also all affected by land reclamation and embankment construction over the past c. a. 200 years (Parsons et al., 2013a(Parsons et al., , 2013b(Parsons et al., , 2013c. Hence, the findings of this study could be used to inform about the effects that embankment construction might have on the sediment transport dynamics of these estuaries in relation to sediment availability.

Conclusions
This study used Hesketh Out Marsh, in the Ribble Estuary, and an integrated approach of paleoenvironmental analysis and numerical modeling to investigate the contribution of embankment construction to long-term salt marsh accretion. OSL analysis showed that Hesketh Out Marsh has been accreting at a rate of 4.61 to 0.86 cm yr −1 , and that this rapid accretion responsible for the infilling of the estuary might have started before the beginning of embankment construction in 1810 AD, which is in contrast with previous hypotheses. Moreover, no shift from a slower to a more rapid accretion was noticeable in the trend that could link the beginning of the rapid accretion to the beginning of the embankment construction. Furthermore, the numerical simulations showed that in the Ribble Estuary the addition of embankments is responsible for an overall intensification of the ebb currents, which enhances sediment export. This confirms that embankments cannot be considered responsible for an increase in sedimentation. XRF analysis showed that sediment supply to the Ribble Estuary is predominantly from the Irish Sea and, in conjunction with PSD and NIRS analyses, it showed that the changes in accretion rate detected by OSL were connected to the natural evolution and stabilization of the system, supporting the hypothesis that the high sediment supply to the Ribble Estuary is likely to be connected to the high availability of mobile sediment in the central Irish Sea that was discharged by the Irish Sea Glacier during the last deglaciation. The rapid sedimentation is, therefore, more likely to be attributed to the natural sediment supply to the estuary. The following decrease in sedimentation is instead attributable to the natural evolution of the mudflat-salt marsh system that reached equilibrium. Despite the numerical modeling showing that embankments promote sediment export, the long-term accretion of the salt marsh was not compromised by the effects of the embankment construction. The numerical modeling showed that the estuary has been flood dominated throughout the studied period and that the infilling of intertidal flats, combined with undergoing river dredging, further enhanced this flood dominance. It is likely that a combination of overall flood dominance, natural high rates of sediment supply, infilling of the accommodation space with this sediment and undergoing river dredging, contributed to the formation and accretion of the salt marsh despite the negative impact of the embankments on the marsh sediment budget.