Establishing sustainable sediment budgets is critical for climate-resilient mega-deltas

Many of the world’s major river deltas face a sustainability crisis, as they come under threat of increases in salinity and the extent of tidal zones forced by combinations of sea-level rise, changes in river discharge and channel geometry. The relative contribution of these factors to future increases in tidal extent remains unconstrained, with most prior work emphasising the role of climate-driven sea-level rise. Here we use new field data from the Mekong delta to measure variations of river discharge and changes of channel geometry, and project them into the future. We combine these with projections of future sea-level rise into a 2D hydrodynamic numerical model and quantify the influence of the different driving factors on future tidal extension into the delta. We show that within the next two decades, tidal extension into the Mekong delta will increase by up to 56 km due to channel deepening (92%), dominantly driven by anthropogenic sediment starvation. Furthermore, even under strong mitigation scenarios, sediment starvation still drives a long-term commitment to future tidal extension. Specifically, by 2098 eustatically rising sea-levels are predicted to contribute only modestly to the projected extension. These findings demonstrate the urgent need for policy makers to adopt evidence-based measures to reverse negative sediment budgets that drive tidal extension into sediment starved deltas.


Introduction
The world's deltas occupy less than 1% of the global land surface, but 4.5% of the global population lives on them (Edmonds et al 2020). The ecosystem services that deltas support provide a range of socio-economic functions (Szabo et al 2016) that underpin global food security (Smajgl et al 2015). However, many large deltas are being drowned (Syvitski et al 2009) due to rapid relative sea-level rise, which in 46 of the world's largest deltas averages rates of more than 6 mm yr −1 (Tessler et al 2018), presenting a major threat to these systems (Anthony et al 2015). Rapidly rising relative sea-levels are driven by a range of compound pressures including eustatic sea-level rise, accelerating subsidence and sediment starvation (Ericson et al 2006, Syvitski 2008, Vörösmarty et al 2009. Some studies estimate that by 2100 greenhouse gas-driven eustatic sea-level rise alone will increase the areas at risk of flooding in the world's deltas by more than 50% (Syvitski et al 2009), exacerbating saline intrusion into their sensitive ecosystems (Smajgl et al 2015, Ensign and Noe 2018). Meanwhile, delta subsidence, driven by compaction of sedimentary strata (Törnqvist et al 2008) is also being accelerated due to increased rates of ground water extraction (Zhang et al 2008b, Erban et al 2014 in many locations. Furthermore, many large deltas are also sediment starved (Dunn et al 2019) resulting from declines in the supply of fluvial sediment caused by upstream dam construction (Dai andLiu 2013, Kondolf et al 2014), rapid and accelerating sand mining (Bravard et al 2014, Hackney et al 2020, and the construction of flood embankments (Vörösmarty et al 2009, Auerbach et al 2015, Chapman et al 2017, Tamura et al 2020. Such sediment starvation potentially compromises the ability of deltas to offset rising sea levels through sediment deposition (Syvitski et al 2009, Giosan et al 2014. The combined influence of these multiple pressures of relative sea-level rise and changes in water and sediment discharges, driven by global warming (Lan et al 2016) and catchment engineering (Vörösmarty and Sahagian 2000), directly impacts the sensitive balance between fluvial and tidal forces that governs tidal extent (TE) in large deltas. TE is defined here as the region of the delta where flow discharge is dominated by tidal forcing. This is important because TE regulates saline intrusion (Eslami et al 2019), water and sediment routing (Nienhuis et al 2018), and channel stability (Zhang et al 2017). Disentangling the relative impacts of the distinct pressures influencing TE, and especially their future evolution, is pertinent to generating evidencebased policy that could underpin adaptation and mitigation strategies for exposed deltas.
We focus on the Mekong delta because it is classified by the Intergovernmental Panel on Climate Change (IPCC) as one of the three most vulnerable deltas globally (Parry et al 2007) and because it is representative of other large deltas not only in terms of the multiple drivers of change it is experiencing, but also in its response to these pressures. For example, similar to many other deltas (Dunn et al 2019), the supply of fluvial sediment to the Mekong has experienced a major decline, from 160 Mt yr −1 in the 1980's (Milliman and Meade 1983) to 87 Mt yr −1 in 2005 (Darby et al 2016). These declining rates of fluvial sediment supply have been attributed to the effects of shifting tropical cyclone tracks (Darby et al 2016) as well as extensive river impoundment and hydropower development across the Mekong's catchment (Kondolf et al 2014). Furthermore, many sections of the river and its delta have been subject to intensive sand mining activities, with previous studies documenting annual removal rates in excess of 50 Mt yr −1 (Bravard et al 2013) and a recent study suggesting that unregulated sand mining is also occurring (Jordan et al 2019). This produces a significant net negative sediment budget when compared to sand influx rates of 6.2 ± 2.0 Mt yr −1 at the delta apex (Hackney et al 2020 (Brunier et al 2014), which in the recent past has led to an intensification of saltwater intrusion (Eslami et al 2019), and reversing progradation of the channel mouths (Tamura et al 2020). Synchronous to this anthropogenically-driven channel bed level lowering, the delta is also experiencing: (a) a eustatic rising of sea-levels, forced by global warming, which for the Mekong delta region is projected to be up to 1 m by 2100 (RCP 8.5) (relative to 1998) (see supplementary figure 8 (available online at stacks.iop.org/ERL/16/064089/mmedia)), and (b) an anthropogenically accelerated sediment compaction which is projected to give cumulative subsidence of up to 1.4 m by 2100 (relative to 1990) (Minderhoud et al 2020). Here we use new field data and aligned numerical hydrodynamic modelling to develop the first analysis of the combined influence of future (to 2098) eustatic sea-level rise, river discharge and delta channel changes, the latter partitioned into sediment loss and land subsidence, on tidal extension in the Mekong delta.

Data and methods
We quantify historical bed-level lowering of the principal distributary channels of the Vietnamese Mekong Delta (VMD) by comparing successive bathymetric surveys for the years 1998, 2008 and 2018. We then propagate the identified trends of bed-level changes 20 years forward by assuming a 'business as usual' evolution of pressures and generate potential future alternate delta analogues. We develop a twodimensional hydrodynamic numerical model (based on the varying delta topographies) and undertake a series of numerical experiments (n = 105) to explore the response of delta hydraulics under a range of compound pressures, accounting for co-variations in: (a) eustatic sea-level rise; (b) river discharge, and; (c) channel bed lowering. We analyse annual time series of predicted hourly water discharges at 100 locations across the delta channel network to quantify tidal extension for each one of the scenarios investigated.

Historical bathymetric surveys
Historical bathymetric surveys (Brunier et al 2014) were obtained from the Mekong River Commission (www.mrcmekong.org). These surveys consist of point data for the years 1998 and 2008 and cover approximately 51% (432 km 2 ) of the delta's principal channels. They include depth measurements with a depth uncertainty of ±0.2 m for every 10 m which generates a mean vertical error of 0.19 and 0.21 m for the 1998 and 2008 datasets, respectively. They are referenced to the World Geodetic System 1984 (WGS84) and the Ha Tien 1960 mean sea level (MSL) datum. We adopted the Universal Transverse Mercator (UTM) coordinate system and projected the geographic coordinates into WGS84 UTM 48 N. We verified historic MRC data and procedures used for vertical referencing, and then converted the elevation values from Ha Tien 1960 to the Hon Dau 1992 MSL by subtracting the 0.167 m offset from each reported elevation. Hon Dau is the vertical datum currently used in Vietnam.

Contemporary bathymetric survey
Our contemporary bathymetry was obtained from a delta-wide survey conducted in May 2018. We used two vessels, equipped with identical single beam echo-sounding (SBES) systems consisting of a Garmin GT20-TM sonar transducer, fully submerged (∼0.5 m) below the water surface, linked to a global positioning system (GPS). Depending on channel depth, we altered the sonar pulse frequency between 800 and 455 kHz with beam angles of 1.6 • and 2.5 • , respectively. Data were recorded at a frequency of 1 Hz; vessel speed was kept below 20 km h −1 to minimise air entrainment and data noise. Garmin does not specify a measurement error for their system, hence we adopt the same vertical uncertainty as for the historical data, which generates a mean vertical error of 0.25 m for the 2018 dataset. We henceforth propagate the mean errors of the historical and 2018 surveys into our volumetric differencing and subsequent channel bed level lowering rates. A limited part of the principal delta channels (145 km 2 out of 854 km 2 ) had already been thoroughly surveyed in 2014 (supplementary figure 1). We surveyed these areas following a straight course along the channels to identify general trends of change within the 4 year period (see later). We then concentrated our efforts to extensively survey the remaining areas (709 km 2 out of 854 km 2 ), that had not been surveyed since 2008. To optimise channel coverage, our survey tracks followed oblique cross-sections crossing from bank to bank, thereby forming open triangles with their apexes spaced approximately two channel widths apart.
Acoustic recordings were converted into projected (WGS84 UTM 48N) coordinates and measurements of total water depth. The latter were transformed to elevations using hourly water level recordings from the Vietnamese hydrological agency's network of water level monitoring stations (supplementary figure 1), referenced to the Hon Dau MSL. We compared the bed elevations measured in 2014 with the corresponding ones measured in 2018 and found that within the 4 year period the landward section of the principal Mekong had been lowered by 2.9 m on average, whereas a limited section of the Bassac had been raised by 0.6 m on average. We adjusted the 2014 data accordingly and used them to augment our 2018 survey data.

Principal channel bed surfaces
We interpolated channel bed surfaces from the survey data and compared historical and contemporary channel morphology. To account for the anisotropy of channel bed elevations, due to flow directionality ( Surface differencing of contemporary and historical maps (see section 3 and supplementary figure 2) highlighted spatial trends for channel bedlevel lowering, with change being most prominent in the landward sections, whereas in areas very close to the shore almost no change was observed ( figure 1(a)). These patterns were summarised by extracting linear trends of laterally averaged channel change for the principal Mekong (section indicated as MM' in supplementary figure 2(c)) and Bassac (section indicated as BB' in supplementary figure 2(c)) (supplementary figure 3) channels. Trends of channel bed-level lowering were interpolated spatially into raster maps representing the historical trend of bed-level changes across the principal channel network (supplementary figures 4(a) and (b)).

Hydrodynamic model
A 2D hydrodynamic model (supplementary figure 5) was developed using the Danish Hydraulics Institute's MIKE21 FM modelling suite, using a combination of regular and flexible mesh elements (n = 244 000). Simulations were undertaken using a high-performance computing system, with a full year simulation requiring approximately 1.5 d to complete. Our modelling simulations included a hydrodynamic component, with channel bathymetry remaining constant during scenario runs but varied between scenarios as described in section 2.5.
The domain topography (

Scenarios
We designed a series of 105 simulation scenarios to explore the response of delta hydraulics under a range of compound pressures.
For eustatic sea-level rise (a), we used seven scenarios, based on projections from the Warming Acidification and Sea level Projector (WASP) (Goodwin et al 2017, Oppenheimer et al 2019 for the grid point closest to the Mekong delta, which are broadly consistent with the IPCC Fifth Assessment Report (AR5) projected range (Church et al 2013). We chose to assess 0.06, 0.1, 0.2 and 0.25 m of sea level rise (relative to 1998) by the years 2008, 2018, 2031 and 2038, respectively. These levels of sea-level rise correspond approximately to average projections across all RCPs for this period (supplementary figure 8). We also explored the effect of a 0.5 m sea-level rise which is predicted by all RCPs (but for different years) and that of a 1 m SLR, that reflects the upper end of the RCP 8.5 by 2100. Within the hydrodynamic model, we vertically offset the sea boundary (which is based on the 1998 tidal time series) by the specific mean sea-level rise scenarios.
The effects of river discharge variability (b) were assessed using three distinct river flow discharge data sets. Specifically, we used hourly observations of river discharge at Tan Chau and Chau Doc that covered the entire years 2013, 2010 and 2001 (supplementary figure 8). We selected these years because their mean annual discharge values are close to the median, lower and upper bounds, respectively, of the range of discharges observed during 2001-2018. Finally, to explore the effects of channel bed-level lowering (c) we repeated all simulations with five different bathymetries corresponding to historical, contemporary and three future projected conditions. The three future scenarios were developed to represent the potential evolution of channel geometry under different mitigation approaches. The first mitigation scenario (scenario M1) is a hypothetical scenario in which effective regulation is implemented immediately to maintain the delta channel geometry at its contemporary (i.e. 2018) state and as such matches the contemporary channel bed levels ( figure 1(a)). Scenario M2 corresponds to a hypothetical trajectory where mitigation is deferred and the delta channels continue to lower their beds (at their present rates) until 2028, when further bed-level lowering is halted. We produced M2 (

Anthropogenic channel bed-level lowering
Comparison of the interpolated channel bed surfaces for 1998, 2008 and 2018 reveals a substantial, accelerating, bed-level lowering of the delta channels. Specifically, our data show a mean bed-level lowering of 2.5 m (σ = 3.9 m) since 1998, with change being most prominent in the landward sections, whereas in areas very close to the shore almost no change was observed ( figure 1(a)). This is likely due to an increasing sediment influx from the coast (Ha et al 2018) and the increased sand mining activity that is particularly concentrated in the landward reaches of the delta. Comparisons at decadal intervals reveal that the bed levels of the common areas (51% of the principal delta channels) were lowered by 0.9 m (σ = 3.2 m) during 1998-2008, and by 1.6 m (σ = 2.9 m) during 2008-2018 (supplementary figure 2). Minderhoud et al (2020) have quantified a mean cumulative subsidence for the Mekong delta plain of 0.2 m for the period 1990-2018. Therefore, at least 2.3 out of the 2.5 m (92%) of bed-level lowering that we quantify here cannot be explained by subsidence and must be as a result of sediment loss. Accounting for this modest subsidence generates volumetric sediment loss estimates of approximately 341 × 10 6 ± 97 × 10 6 m 3 and 584 × 10 6 ± 112 × 10 6 m 3 which, assuming a bulk density of 1600 kg m −3 , corresponds to 546 ± 155 Mt and 934 ± 195 Mt for the 1998-2008 and 2008-2018 periods, respectively. Extrapolating the sediment losses found on the common areas across the entire principal channel network, suggests that these channels lost sediment at average rates of 108 ± 33 Mt yr −1 during 1998-2008 and 184 ± 39 Mt yr −1 during 2008-2018 ( figure 1(b)). Combining the observed trajectory of bed-level lowering with the 2018 survey, taking into consideration the spatial variability of quantified changes (see section 2 and supplementary figures 3 and 4), we project the future morphology of the delta channels (figure 1(c)) by assuming a 'business as usual' evolution of compound pressures, giving an estimated bed-level lowering of 5.3 (3.9) m by 2038 (2028), relative to the 1998 analogue.
The average GHG-driven sea-level rise is projected as 0.25 m in 2038, relative to 1998 regional sea levels (figure 1(d) blue line), across all RCP scenarios we consider here (supplementary figure 8). In contrast, the linear projected channel geometry for 2038 is 5.3 m lower than in 1998 ( figure 1(d), red line) suggesting that the channel beds in 2038 will be 5.55 m lower with respect to sea-level than they were in 1998 ( figure 1(d), black line). The cumulative subsidence for 2038 is projected to be 0.4 m by 2038 (Minderhoud et al 2020). Therefore, the vast majority of bed-level lowering in our scenarios for the recent past (2.3 of 2.5 m between 1998 and 2018) and future (2.6 of 2.8 m between 2018 and 2038) is attributable to anthropogenically-driven sediment loss ( figure 1(d), orange line) under the combined influence of reduced sediment supply from the catchment and excessive river-bed mining.

Anthropogenically-driven tidal extension
The hydrodynamic model was used to explore 105 scenarios (90 are shown in figure 2), spanning a wide range of compound pressures. Simulated water discharge values at 1 h intervals at 100 locations across the delta (supplementary figure 9) show strong seasonal variability, driven by the monsoon-governed variations in river discharge, and hourly fluctuations governed by the tidal forcing. We adopted the convention that, under complete tidal dominance the flow would vary equally between positive and negative (inversed flow) discharge values (50% Q > 0), resulting in a net zero discharge over a full year, whereas a 75% Q > 0 represents the zone where fluvial and tidal forces equalise. We calculate the ratio of positive and negative discharge at each one of the 100 study locations and used this ratio to map the transition between fluvial and tidal dominance across the delta. We quantify the TE for each scenario by measuring the mean distance between the location of equal tidal and fluvial dominance (75% Q > 0) and the channel mouths (n = 7). Henceforth we report values of TE corresponding to the median (2013)  The simulated average TE into the delta in 1998 was 41 km. Our simulations show that if the channel geometry had been maintained at 1998 conditions (HA1), the TE would have increased by only 1.3 km landward (TE = 42.3 km) during 1998-2018 and an additional 9.3 km (TE = 51.6 km) during 2018-2098, driven solely by the eustatically rising sealevels (figure 2, grey panels) predicted by the RCPs. However, for the simulations in which the bed-level lowering during 1998-2018 is included, the landward increase of TE was extended by an additional 13.2 km, so that its predicted location by 2018 was 54.2 km from the channel mouths. Should immediate (2018) action be taken to halt channel deepening (scenario M1), the TE is predicted to increase a further 7.1 km landward (TE = 61.3 km) by 2098 due to eustatic sea-level rise (figure 2, green panels). In the M1 scenario, the TE (TE = 61.3 km) predicted in 2098 reaches almost 10 km further landward than the HA1based equivalent (TE = 51.6 km). However, if channel degradation continues for another decade (scenario M2), the predicted TE in 2098 is 89.9 km, 35.7 km further inland than in 2018 (figure 2, orange panels). Moreover, if channel degradation is not halted until 2038 (scenario M3) the TE is predicted to more than These trajectories, showing the potential future tidal extension within the delta under different drivers of change, are summarised in figure 3. The VMD is subjected to two groups of pressures of different magnitude. GHG-driven sea-level rise forces a slow, progressive, tidal extension of ∼0.1 km yr −1 (figure 3, grey dotted line). However, channel bed-level lowering, dominated by sediment starvation during the last 20 years, forces a much more rapid (0.5 km yr −1 ) tidal extension landward (figure 3, black line). If the process of bed-level lowering continues unmitigated, the TE in 2098 is projected to be double that in 2018. Our simulations also show abrupt accelerations of the landward increase of the TE, interpreted here as being related to delta topology. Since the delta is increasingly bifurcated downstream (Edmonds and Slingerland 2008), and because channel widths increase seaward where tides become more dominant (Nienhuis et al 2018), both the number of distributary channels and the total channel cross-sectional area decrease landward (supplementary figure 10). Thus, as the TE stretches further landward, fewer and narrower channels exist to attenuate tidal energy. This mechanism creates a positive feedback which expands the TE further landward. Considering our historical trajectory (figure 3, black line) in relation to the delta topology (supplementary figure 10) suggests that channel narrowing further exacerbated the tidal extension between 2008 (TE = 44.5 km) and 2018 (TE = 54.2 km). Our future trajectories will also similarly accelerate should the TE become greater than 70 km, which falls within a transition zone (60-90 km) of acute reduction of total cross-sectional area and number of distributary channels (supplementary figure 10). This key distance of 70 km is also close to the major cities of Can Tho (c. 75 km) and Vinh Long (c. 80 km) and could be crossed by the tidal incursion in the near (before 2038) or medium-term (after 2064) future under scenarios M3 and M2, respectively. If immediate action to halt channel deepening is taken (scenario M1) the TE will remain relatively constrained (TE < 70 km), for the majority of seasonal fresh-water flux conditions, until the end of the century.

Discussion and conclusions
Should channel bed lowering continue unmitigated, the TE in the Mekong delta will rapidly increase by ∼56 km in the next two decades. This dramatic landward expansion of the TE risks causing severe socioeconomic impacts, via exacerbation of saline intrusion (Smajgl et al 2015, Eslami et al 2019). Moreover, tidal extension will also affect delta channel stability and flood risk because of changes in bifurcation function (Zhang et al 2017) and water and sediment routing (Nienhuis et al 2018). It is imperative, therefore, to identify correctly the dominant drivers of tidal extension in order to develop viable mitigation and adaptation strategies. Importantly, to date the dominant narrative has focused on the role of climate change in driving sea-level rise (Smajgl et al 2015, Ensign and Noe 2018) and land subsidence (Erban et al 2014, Minderhoud et al 2020 in driving tidal extension. Our results highlight that, although such factors do indeed contribute to tidal extension in the long term (i.e. by 2098), channel bed level lowering linked to cumulative anthropogenic sediment starvation, stemming from sand mining and damming, has a rapid and substantial influence that dominates all the other drivers of tidal extension combined. Critically, our findings demonstrate that present-day levels of channel bed level lowering stimulate such a powerful long-term commitment to tidal extension that, even if it is halted now, channel bed level lowering will continue to remain the net dominant driver of tidal extension until at least the end of the century. This means that attempts to mitigate against tidal extension must focus on reversing the negative sediment budgets that are driving channel bed level lowering. Our findings likely have widespread implications, as the Mekong is not the only large delta system that is sediment starved and which is highly vulnerable to tidal extension. Fluvial sediment reduction has been reported for many large rivers (Dunn et  . It is clear that attempts to maintain climate resilient deltas globally will fail unless they include efforts to establish, sustain and enhance positive sediment budgets into their fragile systems.

Data availability statement
The data that support the findings of this study are openly available at the following URL: https:// hydra.hull.ac.uk/resources/hull:17951. grant number: NE/P014704/1 and supported by the Viet Nam Ministry of Science and Technology, Project Code KC08.12/16-20 'Study on adverse impact of the river bed degradation of the lower Mekong River and proposed mitigation solution' . D R P additionally acknowledges funding received from the European Research Council under the European Union's Horizon 2020 research and innovation program (Grant Agreement 725955). G V thanks V Dung, M V Bao, H P Pham, T M Nhat, S Teasdale, J Johnson, and C Hackney for their help during the channel survey and C Collins for High Performance Computing support. The authors declare no conflicts of interest.