Towards unraveling total subsidence of a mega-delta – the potential of new PS InSAR data for the Mekong delta

Total subsidence in deltaic areas is the cumulative effect of a range of driving mechanisms, both natural and anthropogenic. The populous and low-lying Vietnamese Mekong delta is facing accelerating subsidence rates and effective mitigation strategies are urgently needed to save-guard the future sustainability of the delta. This paper gathers results from existing measurements and estimates of subsidence in the Mekong delta and presents new, delta-wide datasets of PSI observations of vertical velocity from 2014–2019. We describe the practical application of this new data in ongoing projects in Vietnam and outline a planned approach to determine depth-dependent subsidence rates, using this new dataset in combination with field surveys and physics-based numerical models, to advance towards improved quantitation of the contributions of individual subsidence mechanisms.

1 Introduction

1.1 Causes of land subsidence in deltas

Land subsidence in deltas can be caused by a range of drivers and processes. The actual subsidence that occurs and contributed to surface elevation change over time, is the sum of all processes. Land subsidence in deltas occurs because of natural processes, like the compaction of soft, unconsolidated sediments over time and with increased overburden load (e.g. new sediments or seasonal flooding), oxidation of organics, but also as a result of tectonic movements. Human activities can enhance the natural subsidence or trigger new processes, for example, by extracting fluids (decreasing pore-pressure), lowering of the surface water table or by adding weight to the surface by buildings or infrastructure. All these mechanisms act on different temporal and spatial scales, creating a complex subsidence response of an entire delta. To develop detailed adaptation and mitigation strategies for subsiding deltas subsidence measurements need to be untangled into different processes and drivers to gain understanding on where, how much and why individual subsidence drivers cause subsidence locally.

The populous Vietnamese Mekong delta, like many other deltas worldwide, is facing increased rates of subsidence during the past decades. Dedicated research focused on quantifying rates, drivers and processes have greatly increased insights and the availability of subsidence data in recent years. The Mekong delta is experiencing relatively high rates of natural compaction of its young, Holocene sediments (Zoccarato et al., 2018). On top of this, human activities and economic development during the past decades have led to intensification of agricultural practices and urbanization of the delta. Intensification of land use change from a natural to more human dominated land-use correlated to increasing rates of observed land subsidence (Minderhoud et al., 2018) as anthropogenic drivers, enhancing and triggering subsidence processes increased. One of the large drivers of the subsidence rates presently experienced, was the steady increase in groundwater extraction since the 90's from its thick, locally over 500 m, multi-aquifer system (Minderhoud et al., 2017). As new elevation data revealed that the deltas elevation relative to local sea level is even lower than previously assumed (on average the delta plain is elevated only ∼80 cm above local sea level) (Minderhoud et al., 2019), the accelerating subsidence rates strongly increase the vulnerability to sea-level rise-induced flooding, salinization (Eslami et al., 2019), coastal erosion and, ultimately, threaten the livelihood of 18 million delta inhabitants with permanent inundation (Minderhoud et al., 2020).

Table 1 Overview of the existing subsidence observations in the Mekong delta.

Quantifications of land subsidence are scarce in the Mekong delta and generally do not cover long measurement periods. In recent years, several direct and indirect measurements and derived quantifications on surface level lowering, the result of land subsidence, became available both from remote sensing sources (sensors on satellites) and (in)direct ground-based measurements (Table 1).

In this paper, we present new Persistent Scatterer Interferometry (PSI) data for the Mekong delta created by Gisat, through an EU's Copernicus Emergency Management Service – Risk & Recovery Mapping activation by the German International Cooperation Agency (GIZ) and Geological Survey (BGR). We explain the methodology, present the resulting InSAR-based subsidence estimates and discuss an example of practical application of the data in DEM projections as done in ongoing development projects. Finally we discuss the proposed approach to move forward towards interpreting the new PSI estimates to arrive at depth-dependent subsidence rates, which is the next step towards unraveling specific subsidence processes and mechanisms.

2 New subsidence estimates for the Mekong delta based on PS InSAR

Recently, new datasets of PSI-based subsidence estimates for the Mekong delta became available. This InSAR-based subsidence monitoring was carried out under the framework of Copernicus Emergency Management Service (© 2018, 2019 European Union), EMSN057: “Ground subsidence Analyses, Mekong Delta, Vietnam”, and EMSN062: “Assessing changes in ground subsidence rates, Mekong Delta, Vietnam”. EMSN057 focused on mapping areas of three cities in the Delta in high detail, i.e. density of detected persistent scatterer (PS) points: Ca Mau, Long Xuyen and Rach Gia (approximately 20×20 km large rectangles centered in the cities), while EMSN062 aimed at mapping of the whole Delta (more than 40 000 km²) at lower detail.

2.1 Methodology

Land subsidence was mapped using time series of Sentinel-1 SLC images from descending track 18. Time series covered period from October 2014 to January 2019 or September 2018 in case of EMSN057 respectively. In total more than 180 images were acquired, stacked and analyzed by PS InSAR (Persistent Scatterer Interferometry) technique. Images were first co-registered to optimal master image using TanDEM-X DEM and precise orbital information to achieve sub-pixel accuracy. In the next step atmospheric delay was estimated as the low-pass component of the phase residuals from displacement estimations on set of points with the lowest amplitude dispersion. In case of EMSN057, every city was processed separately and reference points were selected individually for each processed area. Therefore, their stability and comparability of absolute figures of estimated subsidence rates between areas and with EMSN062 could not be guaranteed. In EMSN062, the reference point was selected in the hilly area at the border with Cambodia, whose rock outcrops are supposed to be stable in contrast to the most of the Delta. In this case, each data swath was originally processed separately but due to deteriorating PS's quality at the swath edges each swath was divided into two sub-areas. Therefore, six sub-areas were processed independently and results were merged at the end of processing assuring each sub-area was referenced to the common stable reference point. Due to size of the Delta, density of resulting PS points could not be as high as for 3 cities and their number was reduced by application of thresholding considering amplitude stability and coherence. Standard PS InSAR (Ferretti et al., 2001) algorithm implemented in the SARProZ© was applied on each point in order to estimate the displacement rate (average annual velocity), its standard deviation, displacement time series and other attributes. The rate was estimated using linear displacement model assuming that characters of displacement trends are mostly linear in time and that their magnitude is significantly smaller than 2.8 cm per acquisition frequency, i.e. half of SAR sensor's wavelength per 12 d until September 2016 and per 6 d afterwards. Displacements were measured in satellite's line of sight (LOS), which ranged from 30–45^∘ from the vertical direction. Subsidence rate was recalculated from LOS to vertical plane direction under assumption that most of PS points in the Delta should be subject to subsidence while not possessing significant displacement component in horizontal direction. Each point was accompanied by attributes associated to quality of detected variables and showing the trend of detected ground displacement in time. As there had been many construction activities in the Delta during 2014–2019 period, each point was supplied with additional information about the start and the end date of the reflector estimated using evolution of detected radar intensity signal.

2.2 PSI-estimates of subsidence rates

The detected movement of the persistent scatterers was used to estimate vertical velocity of the objects (Fig. 1). Generally, there is higher PS point density in urbanized area, while considerably lower density in rural or natural landscapes. Initial general findings are in line with patterns and outcomes identified by previous studies (Erban et al., 2014; Minderhoud et al., 2018) though quantitative assessment still needs to be conducted: subsidence hot spots are located mostly in built-up areas with average annual subsidence rates ranging from 20–50 mm yr⁻¹ and agriculture areas are being affected to lesser extent with rates amounting to 0–20 mm yr⁻¹.

Figure 1 PSI-estimated average subsidence velocity in the Mekong Delta between 2014 and 2019. Based on Copernicus Emergency Management Service product (© 2019 European Union), [EMSN062] Mekong Delta, Vietnam: Average Subsidence Velocity map.

The results of the two datasets each serve their own purpose in application. The EMSN062 results facilitate comparative assessment of subsidence patterns at regional level with high reliability. They reveal areas more prone to subsidence and may be assessed in conjunction with built-up, land use or groundwater extraction information for the whole Delta. However, point density is insufficient for analysis at local level, e.g. for individual city. As demonstrated for 3 cities this gap is bridged by results from localized and detailed analysis from EMSN057. As shown in Fig. 2, some small and many medium-size buildings were covered by one more multiple detected PS points. High variability of subsidence rates was detected at very local (building-to-building) level in all analyzed areas. Such micro-level variability are likely caused by subsoil compaction due to different load contributions of infrastructure and buildings depending on their size, built-up material, foundations and also the object's age.

Figure 2 Estimated average subsidence velocity rate of movement in the area of Long Xuyen hospital. Subsidence rates are assigned to individual buildings. Outcomes of PS InSAR analysis were confirmed by field validation. Based on products of Copernicus Emergency Management Service (© 2018 European Union), [EMSN057] Long Xuyen, Vietnam: Average Subsidence Velocity map. © Google Maps.

2.3 Vertical velocity of individual buildings or infrastructure

The high spatial resolution provided by the Sentinel-1 derived InSAR data allows the identification of vertical motion rates of objects on the ground. Horizontal inaccuracies in the dataset make it difficult to attribute specific PSI points to individual houses in dense urban settlements. It is however possible to correlate specific PSI point to larger buildings or free standing objects, like high-voltage power pylons in a rice paddy. Data from the EMSN057 dataset is especially suitable to identify individual objects because of the higher PSI density and provides, for the first time in the Mekong delta, to compare velocity difference between different, adjacent objects which provides relevant information on depth-dependent subsidence (Fig. 3). The relation between objects acting as PSI points and actual subsidence of the surrounding delta surface does not correlate one-on-one as objects may sink faster or slower, depending on structural weight and foundation depth.

Figure 3 Differential subsidence of individual objects and their relation between foundation depth and weight. Combining this relation with the PSI-estimated subsidence rates can be used to investigate and quantify depth-dependent subsidence rates.

Figure 4 PSI-estimated average subsidence velocity of An Giang University buildings in Long Xuyen. Subsidence rates are assigned to individual buildings. The founded buildings (green outline) subside at a speed of ∼0.5 mm yr⁻¹ while the surrounding area subsides at a rate of ∼20 mm yr⁻¹, resulting in about 150–200 mm in ten years since construction (inset photo). Outcomes of PS InSAR analysis were confirmed by field validation. Based on products of Copernicus Emergency Management Service (© 2018 European Union), [EMSN057] Long Xuyen, Vietnam: Average Subsidence Velocity map.

Preliminary, qualitative in-situ field observations show that roads, parking spaces and similar infrastructure have the same subsidence rates as nearby ground. Small buildings with few floors show mostly similar rates of sinking as roads while most bigger, taller buildings tend to have significantly lower subsidence rates, likely related to a deeper foundation depth. The University of An Giang in the city of Long Xuyen, constructed ten years ago, is an example of a building complex with significantly lower subsidence rates compared to nearby roads and open spaces. The university buildings move downward with ∼0.5 mm yr⁻¹ while the roads subside with ∼20 mm yr⁻¹ (Fig. 4).

2.4 Practical application of new PSI data to create elevation projections for urban flood analysis

The PSI data is used within the framework of the ongoing GIZ project: “Mekong Urban Flood Resilience and Drainage Programme” to create projections of elevation in three cities of the Mekong delta (Rach Gia, Long Xuyen, Ca Mau) to investigate the effect of subsidence on future urban flooding and drainage. The PSI-based subsidence estimates were projected, assuming a gradual reduction in rates of 1 % yr⁻¹ (assuming settlement decreases over time with increasing compaction) and subtracting from the Digital Elevation Model (DEM) of the cities. Although interpolation of PSI points estimate leads to uncertainties, as individual points do not necessarily represent subsidence rates of its surrounding area, this approach does provide a reasonably well first estimate of future elevation.

Figure 5 Elevation map of Long Xuyen city in An Giang province in 2018, based on leveling survey in 2010, corrected to 2018 using additional survey data and new PSI estimates. And elevation projected for 2050, following substraction of interpolated PSI-estimated subsidence rates. The projected DEM is used to evaluate future urban flood resilience and drainage capacity of the city.

GIZ used projected DEMs for the layout of the drainage systems for Ca Mau, Long Xuyen and Rach Gia with the aim of facilitating smooth operations at least until the middle of the century. For the city of Ca Mau it appears unlikely that purely gravity driven drainage systems will be viable. Within the next two decades relative sea-level rise will necessitate the protection of the city with a ring of dykes and the use of pumps.

3 Approach towards unravelling the contribution of individual subsidence drivers in the Mekong delta

The new PSI-based estimates of subsidence of individual objects in the Mekong delta enable a detailed and, literately, in-depth analysis aimed to quantify depth-dependent subsidence rates, which is the next step towards quantitative unraveling measured subsidence rate for different driving mechanisms. Planned activities consist of measuring height differences in-situ during field (validation of PSI data) and correlating the results with the building's age, foundation depth, estimated load, surrounding land use, (ground) water dynamics and geological setting. Through this depth-dependant subsidence rates can be correlated to (hydro)geological and lithological record and potentially upscaled to the entire delta. These results will also provide valuable new insights to improve numerical models and therefore our predictive capabilities of process-based subsidence.